Hardin-Drnevich本構(gòu)模型在ABAQUS中的實(shí)現(xiàn)

歐陽(yáng)君 林 飛

（湖南省水利水電勘測(cè)設(shè)計(jì)研究總院 長(zhǎng)沙市 410007）

引 言

ABAQUS軟件具有強(qiáng)大的計(jì)算功能和廣泛模擬性能，特別是對(duì)強(qiáng)非線性和不連續(xù)非線性問(wèn)題的求解具有明顯的優(yōu)勢(shì)[1]。利用ABAQUS軟件進(jìn)行土石壩地震反應(yīng)分析時(shí)，要涉及到土體的動(dòng)力本構(gòu)模型的選擇問(wèn)題。ABAQUS軟件中材料的本構(gòu)模型都是建立在傳統(tǒng)的塑性力學(xué)基礎(chǔ)上的，然而傳統(tǒng)的彈塑性理論不能充分反映土料的動(dòng)力變形機(jī)制[2]，同時(shí)，軟件中對(duì)材料的阻尼的定義采用Rayleigh法和直接模態(tài)阻尼法，按上述方法確定的材料阻尼隨著材料變形的變大而減小，而實(shí)際土體材料的阻尼比隨著其變形的變大而增大，且變化幅度約為0.2～30[3]。

本文針對(duì)ABAQUS軟件在土石壩地震動(dòng)力反應(yīng)分析中的不足，利用其提供的二次開(kāi)發(fā)用戶(hù)子程序接口，完成了反映土體動(dòng)力變形特性的Hardin-Drnevich[4]本構(gòu)模型的開(kāi)發(fā)工作。選擇了新疆某土石壩作為研究對(duì)象，用鄧肯-張E-B模型完成壩體在穩(wěn)定滲流期的應(yīng)力與變形計(jì)算[5]，再輸入1985年烏恰7.4地震中喀什水電站6.8級(jí)余震的實(shí)測(cè)地震加速度時(shí)程曲線[6]進(jìn)行計(jì)算。計(jì)算結(jié)果驗(yàn)證了本文模型的準(zhǔn)確性和實(shí)用性，為土石壩動(dòng)力反應(yīng)分析材料模型的選取提供有益的參考。

1 Hardin-Drnevich本構(gòu)模型

目前用于土石壩地震反應(yīng)分析的本構(gòu)模型大致可分為兩大類(lèi)：一類(lèi)是等效線性粘彈性模型，另一類(lèi)是真非線性粘彈塑性模型。真非線性粘彈塑性模型根據(jù)加載歷史、屈服條件、硬化規(guī)律、流動(dòng)法則等建立，在理論上較為合理。由于真非線性粘彈塑性模型比較復(fù)雜，雖然相比等效粘彈性模型有一些優(yōu)點(diǎn)[7]，但因其計(jì)算參數(shù)較難確定，方程建立及求解工作量大、所需時(shí)間長(zhǎng)，受應(yīng)力路徑的影響較大，且目前缺乏合理的計(jì)算模型，故工程上的實(shí)際應(yīng)用較少。等效線性粘彈性模型盡管存在一些缺陷，但概念明確、參數(shù)易確定、應(yīng)用方便、且計(jì)算結(jié)果有一定的精度，補(bǔ)充一些相關(guān)的計(jì)算模式，包括殘余變形模式，就能夠全面分析地震反應(yīng)，而且在參數(shù)的確定和應(yīng)用方面積累了較豐富的試驗(yàn)資料和工程經(jīng)驗(yàn)，故被工程界所接受，實(shí)用性強(qiáng)。

等效線性粘彈性模型把土視為粘彈性體，采用等效剪切模量Geq和等效阻尼比λeq這兩個(gè)參數(shù)來(lái)反映土的動(dòng)應(yīng)力動(dòng)應(yīng)變關(guān)系的兩個(gè)基本特征：非線性與滯后性，并且將剪切模量與阻尼比均表示為動(dòng)應(yīng)變幅的函數(shù)，即 Geq=G（γ）和 λeq=λ（γ）。 這種模型的關(guān)鍵是要確定最大動(dòng)剪切模量Gmax與平均有效主應(yīng)力σ0′的關(guān)系以及剪切模量和阻尼比與動(dòng)剪應(yīng)變幅之間的關(guān)系式，國(guó)內(nèi)外許多學(xué)者對(duì)此進(jìn)行了研究，其中以Hardin-Drnevich模型[4]最有代表性、應(yīng)用最廣。

Hardin-Drnevich模型假定應(yīng)力-應(yīng)變關(guān)系的骨架曲線為一條雙曲線，當(dāng)動(dòng)剪應(yīng)變?chǔ)?∞時(shí)，曲線以靜力極限剪應(yīng)力τmax為漸近線。當(dāng)γ=0時(shí)，曲線的切線斜率為最大剪切模量Gmax，見(jiàn)圖1。將此曲線的縱坐標(biāo)轉(zhuǎn)變?yōu)?γ/τ，繪制成 γ/τ～γ 的關(guān)系線，見(jiàn)圖 2。 此線為一條直線，在縱軸的截距為 1/Gmax，斜率為 1/τmax，有

等效線性剪切模量為：

式中γ為動(dòng)剪應(yīng)變幅；γr為參考剪應(yīng)變，定義γr=τmax/Gmax；Gmax為最大動(dòng)剪切模量；τmax為最大剪應(yīng)變。根據(jù)土石料動(dòng)力試驗(yàn)資料，最大剪切模量Gmax可以表示為：

式中 σ0′為平均有效主應(yīng)力，可由 σ′0=[σ′1+（1+K0）σ′3]/3估算。K和n為實(shí)驗(yàn)參數(shù)。τmax由下式可以求得：

式中 σ′1、σ′3分別是最大和最小有效主應(yīng)力，c′、φ′為有效內(nèi)聚力和內(nèi)摩擦角。

圖1 Hardin-Drnevich模型的主干線

圖2 Hardin-Drnevich模型中的γ/τ～γ關(guān)系線

圖3為Hardin-Drnevich模型的滯回環(huán)，陰影部分的面積為滯回曲線面積的一半,主干線在原點(diǎn)的切線斜率為Gmax，開(kāi)始卸載時(shí)曲線的切線坡度亦為Gmax。假定滯回環(huán)的面積AL為三角形面積Aabc的百分?jǐn)?shù)，有：

能量損耗系數(shù)的定義為 η=△w/（2πw），△w 為應(yīng)變一周內(nèi)所損耗的能量，即滯回環(huán)的面積，w為應(yīng)力應(yīng)變一周內(nèi)物體內(nèi)部積累的最大彈性變形能，即三角形adc的面積。根據(jù)結(jié)構(gòu)動(dòng)力學(xué)知識(shí)η=2λ，加上圖示的幾何關(guān)系和上述假定，得到等效阻尼比為：

當(dāng)γ接近無(wú)窮大時(shí)，Geq=0，定義λmax=2K1/π為最大阻尼比，同時(shí)根據(jù)式（2）和式（6）可以得出等效阻尼比對(duì)動(dòng)剪應(yīng)變幅的關(guān)系式：

圖3 Hardin-Drnevich滯回環(huán)模型

由于各種土的性質(zhì)差異很大，試驗(yàn)表明各類(lèi)土的應(yīng)力應(yīng)變關(guān)系并不完全符合雙曲線的形狀，故需要加以適當(dāng)?shù)男拚?，所以Hardin[8]（1972年）建議將雙曲線模型修正為如下形式：

其中

式中γh為剪應(yīng)變，a、b為試驗(yàn)參數(shù)，與土的振動(dòng)次數(shù)、振動(dòng)頻率、土單元的平均有效主應(yīng)力有關(guān)。

2 Hardin-Drnevich模型的程序?qū)崿F(xiàn)

在進(jìn)行土石壩的動(dòng)力反應(yīng)時(shí)，土料的應(yīng)力-應(yīng)變關(guān)系通常用八面體上的應(yīng)力-應(yīng)變關(guān)系表示[9]，設(shè)八面體上剪應(yīng)力為τ，八面體上的剪應(yīng)變?yōu)棣茫瑒t剪切變形模量為：

八面體上的剪應(yīng)變?yōu)椋?/p>

式中J′為應(yīng)變偏量的不變量，用下式計(jì)算：

式中：εx，εy，εz，γxy，γyz和 γzx分別為各單元的正應(yīng)變和剪應(yīng)變。

根據(jù)彈性力學(xué)原理，應(yīng)力應(yīng)變關(guān)系為：

其中：

式中δij稱(chēng)為克羅納克爾，ν為土的泊松比。

將式（13）寫(xiě)成增量形式為：

采用隱式積分求解非線性動(dòng)力問(wèn)題時(shí)，雅克比（Jacobi）矩陣的定義將直接影響到積分格式的收斂速度，上述模型的雅克比矩陣在二維條件下定義如下：

在ABAQUS軟件中的隱式積分模塊Standard和顯式積分模塊Explicit中[10]，分別采用外掛子程序的辦法實(shí)現(xiàn)Hardin-Drnevich模型的二次開(kāi)發(fā)。

3 算例分析

3.1 計(jì)算模型

以新疆某位于基巖上的均質(zhì)壩為例進(jìn)行數(shù)值計(jì)算。壩高44.0 m，壩頂寬8.0 m，上游壩坡坡比為1∶2.5，下游壩坡坡比為1:2.75。為獲得動(dòng)力計(jì)算所需的初始應(yīng)力場(chǎng)，靜力計(jì)算采用Duncan-chang（鄧肯-張）E-B模型，模型參數(shù)如表1所示。本文模型進(jìn)行動(dòng)力時(shí)計(jì)算參數(shù)見(jiàn)表2?；鶐r面運(yùn)動(dòng)輸入1985年烏恰7.4地震中喀什水電站6.8級(jí)強(qiáng)余震的實(shí)測(cè)地震加速度時(shí)程曲線，歷時(shí)17.80 s，步長(zhǎng)0.01 s。地震動(dòng)水平向加速度幅值為2.43 m/s2，豎向加速度幅值為2.07 m/s2。水平和豎向地震波時(shí)程曲線見(jiàn)圖4。

表1 靜力模型計(jì)算參數(shù)

表2 動(dòng)力模型計(jì)算參數(shù)

圖4 喀什地震加速度-時(shí)程曲線

3.2 計(jì)算結(jié)果與分析

壩體絕對(duì)加速度響應(yīng)是大壩對(duì)地震動(dòng)力作用最直接的反應(yīng)，為了驗(yàn)證本文子程序的實(shí)用性和準(zhǔn)確性，選取壩頂節(jié)點(diǎn)，分析其峰值加速度和放大倍數(shù)，并繪制其絕對(duì)加速度響應(yīng)時(shí)程時(shí)程曲線見(jiàn)圖5。

圖5 壩頂節(jié)點(diǎn)絕對(duì)加速度響應(yīng)時(shí)程曲線

由圖4和圖5可見(jiàn)壩頂節(jié)點(diǎn)水平向和豎向加速度幅值與地震加速度幅值出現(xiàn)的時(shí)間相比有一些滯后。計(jì)算所得地震加速度反應(yīng)峰值和壩頂加速度放大倍數(shù)列于表3，由此可見(jiàn)，考慮豎向和水平向地震加速度時(shí)的水平向和豎向放大倍數(shù)分別為3.11和1.38，這和現(xiàn)行《水工建筑抗震設(shè)計(jì)規(guī)范》（SL 203-97）中的5.1.3條所推測(cè)設(shè)定值相吻合：設(shè)計(jì)烈度7度時(shí)壩頂加速度放大倍數(shù)為3.11與規(guī)范所規(guī)定值3.0相近。以上結(jié)論表明本文基于ABAQUS軟件的二次開(kāi)發(fā)平臺(tái)，利用其子程序UMAT開(kāi)發(fā)的Hardin-Drnevich動(dòng)力本構(gòu)模型在土石壩地震反應(yīng)分析中的研究是可行的。

表3 地震加速度反應(yīng)峰值及放大倍數(shù)

4 結(jié) 語(yǔ)

基于ABAQUS有限元軟件的二次開(kāi)發(fā)平臺(tái)，編制了Hardin-Drnevich模型的外接子程序，對(duì)該子程序的實(shí)用性和準(zhǔn)確性進(jìn)行了試驗(yàn)和理論驗(yàn)證，說(shuō)明編制的子程序能較好的模擬土體的動(dòng)力本構(gòu)特征；選取新疆某土石壩作為研究對(duì)象，基于ABAQUS軟件強(qiáng)大的非線性有限元分析環(huán)境，對(duì)其進(jìn)行地震反應(yīng)分析，給出了壩頂節(jié)點(diǎn)加速度時(shí)程曲線；計(jì)算結(jié)果表明：基于ABAQUS軟件二次開(kāi)發(fā)的Hardin-Drnevich動(dòng)力本構(gòu)模型是可行的，對(duì)土石壩地震反應(yīng)的計(jì)算結(jié)果也符合相應(yīng)規(guī)律；同時(shí)，利用ABAQUS軟件強(qiáng)大的非線性求解平臺(tái)，進(jìn)行土石壩的動(dòng)力反應(yīng)分析，并利用其用戶(hù)子程序UMAT，修改不同類(lèi)型的Hardin-Drnevich本構(gòu)模型也易于實(shí)現(xiàn)。

[1]徐遠(yuǎn)杰，王光琪.在ABAQUS中開(kāi)發(fā)實(shí)現(xiàn)Duncan-Chang本構(gòu)模型[J].巖土力學(xué)，2004，25（7）：1032-1036.

[2]鄭穎人，沈珠江，龔曉南.巖土塑性力學(xué)原理[M].北京:中國(guó)建筑工業(yè)出版社，2002.

[3]陳國(guó)興，謝君斐，張克緒.土的動(dòng)模量和阻尼比的經(jīng)驗(yàn)估計(jì)[J].地震工程與工程振動(dòng)，1995，15（1）：73-84.

[4]顧淦臣.土石壩地震工程 [M].南京:河海大學(xué)出版社，1989.

[5]歐陽(yáng)君,徐千軍,嚴(yán)新軍,等.基于ABAQUS軟件的土石壩應(yīng)力應(yīng)變分析[J].水資源與水工程學(xué)報(bào)，2009，20（6）：112-115.

[6]新疆維吾爾地震局.1985年新疆烏恰7.4級(jí)地震強(qiáng)余震觀測(cè)報(bào)告[M].北京:地震出版社，1987.

[7]趙劍明,汪聞韶,常亞屏,陳寧.高面板壩三維真非線性地震反應(yīng)分析方法及模型試驗(yàn)驗(yàn)證[J].水利學(xué)報(bào)，2003，（9）:12-18.

[8]Hardin B O,Drnevich V P.Shear modulus and damping in soils design equations and curves [J].Journal of Soil Mechanics and Foundation,ASCE,1972,98（7）:603-642.

[9]朱伯芳.有限元法原理與應(yīng)用[M].北京:中國(guó)水利水電出版社，1998.

[10]莊茁,張帆,岑松，等.ABAQUS非線性有限元分析與實(shí)例[M].北京:科學(xué)出版社，2005.