Bouc-Wen土體動力模型的閾值應(yīng)變研究

魯麗雪，遲世春

（大連理工大學(xué) 海岸與近海工程國家重點實驗室，遼寧 大連 116024）

1 引 言

土體的動力模型即動荷載下土體的應(yīng)力-應(yīng)變關(guān)系是表征土體動力學(xué)特性的基本關(guān)系，是土動力學(xué)研究的中心問題之一，也是分析土動力反應(yīng)等一系列工程問題的重要基礎(chǔ)。土體的動力特性主要表現(xiàn)出非線性、滯后性、變形累積以及強度與剛度退化特性。到目前為止，所提出土體動力模型大多不能令人滿意。Bouc-Wen土體動力微分模型由Bouc[1]1967 年提出，Wen[2]和Baber等學(xué)者[3－4]對其進行了改進，使之能夠模擬土體的強度和剛度退化特性。Bouc-Wen模型既包含非線性阻尼，又包含非線性剛度，可通過合理選擇參數(shù)得到不同形狀的滯回圈，用以描述土體在動力荷載作用下非線性、滯后性、變形積累以及強度與剛度退化等特性。Bouc-Wen土體動力模型也能很好地反映剛度和阻尼隨剪應(yīng)變的變化，適于描述多種類型土的循環(huán)動力特性。該模型還無需處理加載拐點，方便編程，具有良好的通用性。

圖1為GDS動三軸試驗得到的摻礫土心墻料滯回圈與Bouc-Wen模型滯回曲線的比較，可以看出，試驗點與理論模型吻合良好。

國外學(xué)者已就Bouc-Wen模型的熱力學(xué)可容許性開展大量研究工作，指出盡管Bouc-Wen模型獨立于內(nèi)時理論發(fā)展起來，但其仍屬于基于內(nèi)時概念的內(nèi)部變量理論類模型[5－6]。由于某些應(yīng)力路徑引起負的能量耗散而可能違反熱力學(xué)第二定律[7]。不過，單個滯回圈自動閉合而不發(fā)生應(yīng)力或應(yīng)變滑移，因此服從Drucker（或 Ilyushin）公設(shè)[8－9]。

Collins等[10]基于熱力學(xué)基本原理建立了土體本構(gòu)關(guān)系模型的方法，采用量綱分析，給出耗散函數(shù)的可能形式，對受力變形過程中的能量耗散機制等問題也進行了初步探討。郭曉霞[11]、遲世春[12]等從熱力學(xué)基本定律出發(fā)，對土體Hardin-Drnevich模型及 Ramberg-Osgood模型的增量耗散函數(shù)進行了構(gòu)造，但由于采用了Masing準則假設(shè)，并未考慮實際阻尼比隨剪應(yīng)變的變化等。因此，有必要結(jié)合壩料的實際阻尼研究土石料動力變形機制。

本文以Bouc-Wen模型為研究對象，探討其適用條件以及參數(shù)影響，從熱力學(xué)基本定律出發(fā)，構(gòu)造壩料的自由能函數(shù)和耗散函數(shù)，研究其耗散特征及動力變形機制，繪制其耗散應(yīng)力空間和真實應(yīng)力空間屈服曲線，探討其屈服面的變化規(guī)律。

圖1 理論滯回曲線與試驗結(jié)果比較Fig.1 Comparison between theoretically predicted hysteretic loops and test results

2 Bouc-Wen模型耗散函數(shù)的構(gòu)造

2.1 Bouc-Wen土體動力模型

Bouc-Wen土的動力模型的應(yīng)力-應(yīng)變關(guān)系為

式中：τ、γ分別為切應(yīng)力和切應(yīng)變；Gmax為最大剪切模量；α為屈服比，控制屈服后剪切剛度，取值范圍(0,1)，主要影響滯回圈的面積。α越大，滯回圈越瘦，α對曲線的影響見圖2。

圖2 參數(shù)α 對滯回曲線的影響Fig.2 Influence of parameter α on hysteretic loop

參數(shù)ζ=ζ(t)是遲滯變量，用來控制土的非線性，由下列微分方程確定：

式中：A、g、b、n均為模型形狀參數(shù)。

由式（1）、（2）可知，Bouc-Wen模型中加、卸載拐點由遲滯變量ζ和應(yīng)變速率d/dtγ控制，無需特殊處理，這樣可大大簡化程序的編寫。

A控制0ζ=時滯回圈的斜率，影響滯回圈的高度和面積。描述土的滯回特性時，A取值(0,1]，其單獨變化時對滯回圈的影響見圖3。

圖3 參數(shù)A對滯回曲線的影響Fig.3 Influence of parameter A on hysteretic loop

g控制dγ/dt變號時曲線的剛度。為確保能量耗散值為正，應(yīng)滿足 g＞0，在此基礎(chǔ)上具有物理意義的5種滯回圈分別是：①g+b＞0和g?b＞0（見圖4）；②g+b＞0和g?b＜0（見圖5）；③ g+b＞0和g?b=0（見圖6）；④g+b = 0 和g?b＞0（見圖7）；⑤g+b＜0和g?b＞0（見圖8、9）。圖2～9中 n = 1，Gmax=600000 MPa。

由圖可知，當 g+b＜0 時滯回圈內(nèi)凹，g+b＞0時滯回圈外凸，g+b = 0 時滯回圈呈擬線性。

圖4 參數(shù)g對滯回曲線的影響Fig.4 Influence of parameter g on hysteretic loop

圖5 參數(shù)b對滯回曲線的影響Fig.5 Influence of parameter b on hysteretic loop

圖6 參數(shù)g+b對滯回曲線的影響(g-b=0,A=0.40)Fig.6 Influence of parameter g+b on hysteretic loop(g-b=0,A=0.40)

圖7 參數(shù)g-b對滯回曲線的影響(g+b=0,A=0.50)Fig.7 Influence of parameter g-b on hysteretic loop(g+b=0,A=0.50)

圖8 參數(shù)g+b對滯回曲線的影響(g-b=40000,A=0.15)Fig.8 Influence of parameter g+b on hysteretic loop(g-b=40000,A=0.15)

圖9 參數(shù)g-b對滯回曲線的影響(g+b=20000,A=0.10)Fig.9 Influence of parameter g-b on hysteretic loop(g+b=20000,A=0.10)

n為屈服指數(shù)，控制加載拐點初始加載階段線性至非線性轉(zhuǎn)換時的屈服前后滯回圈光滑程度，即屈服尖銳程度系數(shù)，n越小，過渡越平滑，其取值范圍為0至無窮，n＞10時，近似理想彈塑性；n趨向無窮時，應(yīng)力-應(yīng)變曲線趨向雙線性。n對單調(diào)加載曲線的影響見圖10。

圖10 參數(shù)n對滯回曲線的影響Fig.10 Influence of parameter n on hysteretic loop

2.2 模型耗散函數(shù)的構(gòu)造

為確保模型不違背熱力學(xué)定律，本文基于熱力學(xué)基本定律構(gòu)造了Bouc-Wen土體動力模型的增量耗散函數(shù)表達式。

假定塑性中心移動沿骨架曲線建立耗散函數(shù)，采用Ziegler正交條件，得到耗散應(yīng)力空間中屈服函數(shù)的表達式。引入耗散應(yīng)力及真實應(yīng)力間的差別項，即遷移應(yīng)力，得到真實空間中的屈服函數(shù)。由于Bouc-Wen模型的表達式是微分形式，為簡化表達式，取n=1構(gòu)造耗散函數(shù)。

滯回曲線統(tǒng)一表達為

骨架曲線為

式（3）與第1象限中骨架曲線走向相同的再加荷曲線g、γ1取負值，與第3象限中骨架曲線走向相同的卸荷曲線b取負值，與第3象限中骨架曲線走向相同的再加荷曲線g、γ1、b取負值；式（4）對于第3象限中骨架曲線g、b取負值，γ0為剪應(yīng)變幅值，γ1為滯回曲線與橫軸的交點，可由反對稱性得出，下同。

綜合公式（3）、（4）可得到塑性中心移動為即時骨架曲線時的耗散函數(shù)：

真實應(yīng)力表達式為

式中：xbl、ρbl分別為耗散應(yīng)力和遷移應(yīng)力，遷移應(yīng)力，其公式為

3 屈服曲線的繪制

結(jié)合某心墻壩工程選擇粗堆石料、細堆石料及反濾料Ⅰ的動力試驗成果，使用上述方法探討其屈服面的變化規(guī)律。粗堆石料為角礫巖與花崗巖的混合物，細堆石料為弱風化以下花崗巖開挖料加工而成，二者采用大三軸動力試驗機進行動力試驗；反濾料I則為弱風化以下花崗巖加工而成，采用共振柱試驗機及中型三軸儀進行動力試驗。

3.1 最大動剪切模量

最大動剪切模量[13]由下式確定：

式中：Gmax為最大動剪切模量，σm為土體所受的初始平均靜應(yīng)力；pa為大氣壓力；k（與Gmax同量綱）、n1（無量綱）為試驗常數(shù)，均由動力試驗確定，筑壩粗堆石料、細堆石料及反濾料Ⅰ的模量系數(shù)k分別為2455.6、1303.9、976.0，指數(shù) n1分別為0.6004、0.6109、0.4980。

壩料動剪模量比隨動剪應(yīng)變和阻尼比隨動剪應(yīng)變的變化見圖11～13。

圖12 細堆石料對應(yīng)的模量衰減和阻尼增長曲線Fig.12 Corresponding modulus decline and damping growth curves of fine rockfill materials

圖13 反濾料Ⅰ對應(yīng)的模量衰減和阻尼增長曲線Fig.13 Corresponding modulus decline and damping growth curves of the filter materialⅠ

3.2 Bouc-Wen模型參數(shù)辨識

根據(jù)動力試驗提供的動剪切模量 Gs、等效阻尼比ξeq與動剪應(yīng)變[14]的關(guān)系（見圖11～13），用遺傳算法對Bouc-Wen模型進行參數(shù)辨識。

遺傳算法中選擇算法采用比例選擇方法，交叉采用非均勻算數(shù)交叉的方法，變異采用均勻變異，種群數(shù)量為80，迭代次數(shù)為1000，雜交概率為0.8，變異概率為0.1。由于參數(shù)n僅影響遲滯曲線的光滑度，取 n = 1，因此需辨識參數(shù)有α、A、g、b。

從模型參數(shù)識別的結(jié)果來看，隨著動力水平的提高，阻尼增加，滯回圈變高、變胖，α、A、g?b、g+b 均逐漸降低。圖14給出半對數(shù)坐標下粗堆石料對應(yīng)的模型參數(shù)與動剪應(yīng)變的變化規(guī)律。

圖14 參數(shù)α、A、g-b、g+b隨剪應(yīng)變的變化規(guī)律Fig.14 Variation trend of parameters α，A，g-b，g+b with the change of shear strain

在半對數(shù)坐標下擬合α-γ、A-γ、gb?-γ及g b+-γ曲線，得

3.3 屈服曲線的繪制

以粗堆石料為例，結(jié)合模量衰減和阻尼比增長建立Bouc-Wen模型，討論p-q（單位為MPa）的應(yīng)力空間中塑性中心的移動為骨架曲線時的屈服曲線的變化規(guī)律。隨著γm不斷增大，繪制屈服曲線如圖15～17所示，其中（γm，τm）為滯回曲線的頂點。

圖15 γm= 5×10-6時的屈服曲線Fig.15 Yield curves when γm= 5×10-6

圖16 γm= 4×10-5時的屈服曲線Fig.16 Yield curves when γm= 4×10-5

圖17 γm= 8×10-4時的屈服曲線Fig.17 Yield curves when γm= 8×10-4

圖中分別繪制了γm=5× 10?6、4 × 10?5和8×10?4時試樣的屈服曲線。由圖可見，應(yīng)變水平為5× 10?6時，不同γ/γm的屈服函數(shù)幾乎為一條直線，滯回圈應(yīng)變達到 4 × 10?5時，屈服面為斜率不同的直線。不同γ/γm的屈服面在p-q應(yīng)力空間為同一直線說明土體屈服時土體顆粒之間的常摩擦系數(shù)起控制作用。因此，試樣滯回圈應(yīng)變在 5 × 10?6～4 × 10?5之間存在一個剪切屈服曲線的直線斜率，不再保持一致的界限應(yīng)變，稱之為第1閾值應(yīng)變[15]，記為γt1。隨著應(yīng)變的增大，且達到 8 × 10?4時，屈服面已明顯變彎。剪切屈服面彎曲代表控制堆石體試樣塑性變形的機制除摩擦外，還有結(jié)構(gòu)變化等機制。因此，試樣滯回圈應(yīng)變在 4 × 10?5～8 × 10?4之間存在剪切屈服曲線由直變彎的界限應(yīng)變，稱之為第2閾值應(yīng)變，記為γt2[15]。

求出 3種壩料的閾值應(yīng)變，并將其繪制于圖11～13，得到第1、第2壩料閾值應(yīng)變在Bouc-Wen模型對應(yīng)的模量衰減及阻尼衰減情況。從粗堆石料及細堆石料的第1閾值應(yīng)變所對應(yīng)的歸一化動剪切模量看，第1閾值應(yīng)變對應(yīng)的 G/Gmax在0.93～0.99之間，阻尼比為1%～4%。第2閾值應(yīng)變對應(yīng)的割線模量與最大動剪切模量之比范圍為 0.54～0.79，與傳統(tǒng)意義上以孔壓開始升高或體積開始變化為標準定義的門檻應(yīng)變相當。

4 結(jié) 語

本文探討了Bouc-Wen土體動力微分模型的適用條件以及參數(shù)影響。從熱力學(xué)基本定律出發(fā)，克服Hardin-Drnevich模型及Ramberg-Osgood模型未考慮實際阻尼比隨剪應(yīng)變變化的局限，結(jié)合模量衰減和阻尼增長給出基于Bouc-Wen微分模型的增量耗散函數(shù)，研究了筑壩土石料的動力耗散特征及動力變形機制，得到筑壩土石料的閾值應(yīng)變。其中第二閾值應(yīng)變與傳統(tǒng)試驗方法得到的門檻應(yīng)變相當。

[1]BOUC R.Forced vibration of mechanical systems with hysteresis[C]//Proceedings of the Conference on Nonlinear Oscillators.Prague: [s.n.],1967: 315.

[2]WEN Y K.Method for random vibration of hysteretic systems[J].Journal of Engineering Mechanics Division,ASCE,1976.102: 249－263.

[3]BABER T T,WEN Y K.Random vibration of hysteretic systems[J].Journal of Engineering Mechanics Division,ASCE,1981,107(6): 1069－1086.

[4]BABER T T,NOORI M N.Random vibration of degrading,pinching systems[J].Journal of Engineering Mechanics,ASCE,1985,111(8): 1011－1027.

[5]CASCIATI F.Stochastic dynamics of hysteretic media[J].Structural Safety,1989,6(2－4): 259－269.

[6]SIVASELVAN M V,REINHORN A M.Hysteretic models for deteriorating inelastic structures[J].Journal of Engineering Mechanics,ASCE,2000,126(6): 633－640.

[7]CHARALAMPAKIS A E,KOUMOUSIS V K.A Bouc-Wen model compatible with plasticity postulates[J].Journal of Sound and Vibration,2009,322(4－5): 954－968.

[8]CAPECCHI D,DE FELICE G.Hysteretic systems with internal variables[J].Journal of Engineering Mechanics,2001,127(9): 891－898.

[9]ERLICHER S,BURSI O S.Bouc-Wen-Type models with stiffness degradation: Thermodynamic analysis and applications[J].Journal of Engineering Mechanics,2008,134(10): 843－855.

[10]COLLINS I F,HOULSBY G T.Application of thermomechanical principles to the modeling of geomaterials[C]//Proceedings of the Royal Society of London,Series A: Mathematical,Physical and Engineering Sciences.London: [s.n.],1997: 1975－2001.

[11]遲世春,郭曉霞,楊峻,等.土的動力 Hardin-Drnevich模型小應(yīng)變特性及其閾值應(yīng)變研究[J].巖土工程學(xué)報,2008,30(2): 243－249.CHI Shi-chun,GUO Xiao-xia,YANG jun,et al.Small strain characteristics and threshold strain of dynamic Hardin-Drnevich model for soils[J].Rock and Soil Mechanics,2008,30(2): 243－249.

[12]郭曉霞,遲世春,林皋.基于熱力學(xué)定律的土體動力Hardin-Drnevich模型再認識[J].巖土力學(xué),2008,29(9):2335－2340.GUO Xiao-xia,CHI Shi-chun,LIN Gao.Recognition of dynamic Hardin-Draevich model for soils based on generalized thermodynamics[J].Rock and Soil Mechanics,2008,29(9): 2335－2340.

[13]謝定義.土動力學(xué)[M].西安: 西安交通大學(xué)出版社,1988.

[14]王志良,韓清宇.黏彈塑性土層地震反應(yīng)的波動分析法[J].地震工程與工程振動,1981,1(1): 117－136.WANG Zhi-liang,HAN Qing-yu.Sticky elastoplastic seismic response of soil layer wave analysis method[J].Earthquake Engineering and Engineering Vibration,1981,1(1): 117－136.

[15]郭曉霞.熱力學(xué)方法在土體本構(gòu)模型中的應(yīng)用研究[D].大連: 大連理工大學(xué),2009.