條形荷載下梯度非均勻非飽和土的動(dòng)力響應(yīng)分析1)

施力維 馬 強(qiáng) ,?, 舒進(jìn)輝

* (青海大學(xué)土木工程學(xué)院,西寧 810016)

? (青海省建筑節(jié)能材料與工程安全重點(diǎn)實(shí)驗(yàn)室,西寧 810016)

引言

動(dòng)荷載作用下地基動(dòng)力響應(yīng)問(wèn)題的研究一直是土木工程界一個(gè)重要的研究課題.諸多學(xué)者對(duì)動(dòng)荷載作用下均質(zhì)單相彈性地基[1-2]、兩相飽和土地基[3-4]和三相非飽和土地基[5-10]的動(dòng)力響應(yīng)問(wèn)題進(jìn)行了一系列研究,有關(guān)動(dòng)荷載作用下均質(zhì)土地基動(dòng)力響應(yīng)的研究已有很多重要的成果.

然而需要注意的是,上述對(duì)土體動(dòng)力響應(yīng)問(wèn)題的研究絕大部分工作是基于土體具有統(tǒng)計(jì)各向同性且均勻等假設(shè)為前提條件的.而在自然環(huán)境下的土體絕大多數(shù)都為固-液-氣三相組成的非飽和多孔介質(zhì),且由于沉積年代的不同表現(xiàn)出不同的物理力學(xué)性質(zhì),土體沿深度方向具有非均勻性.因此采用層狀非飽和土土地基模型更加符合實(shí)際工程情況.目前關(guān)于單相層狀非均勻介質(zhì)的計(jì)算模型較多,主要有薄層單元法[11-12]、傳遞矩陣法[13-15]、剛度矩陣法[16-18]、傳遞和反射矩陣方法[19-21]等.除了上述方法,由于回傳射線矩陣法(reverberation-ray matrix method,RRMM)推導(dǎo)的矩陣方程中不包含正指數(shù)項(xiàng),因此較好地防止了在高速、高頻動(dòng)載荷作用下因土層厚度或參數(shù)相差較大而出現(xiàn)的數(shù)值計(jì)算問(wèn)題,被廣泛運(yùn)用于層狀模型計(jì)算中[22].文獻(xiàn)[23]擴(kuò)展和改進(jìn)了回傳射線矩陣的方法,用于分析各向異性層狀彈性介質(zhì)中自由波的傳播,證明了回傳射線矩陣法的高精度.柳偉等[24]將回傳射線矩陣法推廣至樁土系統(tǒng)的振動(dòng)分析中,分析了外露長(zhǎng)度、埋置深度、樁端約束情況對(duì)埋置結(jié)構(gòu)自振特性的影響.文獻(xiàn)[25]將回傳射線矩陣法推廣到地基梁自振特性的研究中,對(duì)不同邊界條件下的自振頻率、衰減系數(shù)及模態(tài)函數(shù)進(jìn)行了分析.隨后文獻(xiàn)[26]建立了在動(dòng)載荷下梯度非均勻飽和土體的回傳射線矩陣法計(jì)算列式,提出回傳射線矩陣法對(duì)非均勻兩相飽和土動(dòng)力響應(yīng)分析也具有很好的適用性.Ma 等[27]基于回傳射線矩陣法研究了飽和多孔功能梯度材料對(duì)飽和土地基中移動(dòng)荷載引起的振動(dòng)的隔振效果.

綜上可知目前關(guān)于層狀彈性地基和飽和土地基波動(dòng)響應(yīng)問(wèn)題的研究較多,有關(guān)層狀非飽和土地基動(dòng)力響應(yīng)問(wèn)題的研究還鮮有報(bào)道.因此,鑒于非飽和土是土體在自然界更為普遍的一種存在狀態(tài),考慮土體的非均勻性,本文基于多孔介質(zhì)混合物理論,采用回傳射線矩陣法研究了條形荷載作用下梯度非均勻非飽和土地基的動(dòng)力響應(yīng)問(wèn)題.通過(guò)傅里葉積分變換和Helmholtz 矢量分解原理,結(jié)合邊界條件,建立了在條形荷載作用下非均勻非飽和土體的計(jì)算列式.分析了非飽和土物理力學(xué)性質(zhì)在深度方向呈梯度變化時(shí)應(yīng)力、位移以及孔隙壓力等物理量的變化規(guī)律.

1 非飽和土控制方程及通解

1.1 非飽和土控制方程

基于非飽和土多孔介質(zhì)混合物理論,三相非飽和土體動(dòng)力控制方程如下[8,10]

式中,pw和pa分別為孔隙水壓力和孔隙氣壓力;ρ=(1-n)ρs+nSrρw+n(1-Sr)ρa(bǔ),其中 ρ 代表非飽和土介質(zhì)的總密度,ρm(m=s,w,a)分別代表固、液、氣三相的密度,Sr代表飽和度,n代表孔隙率;頂標(biāo)·表示對(duì)時(shí)間t求導(dǎo),wi和vi分別代表孔隙水和孔隙氣沿i方向相對(duì)于土骨架ui的位移.χ代表有效應(yīng)力系數(shù),受含水量等因素影響,在本文中根據(jù)文獻(xiàn)[28]將χ的值等同于Sr的值;λ和 μ代表Lame 常數(shù);a=1-Kb/Ks,Kb和Ks分別代表土骨架和土顆粒的壓縮模量,其中Kb=λ+2μ/3 ;Kw和Ka分別代表孔隙水和孔隙氣的體積壓縮模量;e=?·u代表土骨架的體積應(yīng)變,u代表土骨架位移張量,w和v代表孔隙水和孔隙氣的位移張量;系數(shù)A11～A24表達(dá)式詳見(jiàn)文獻(xiàn)[10].g代表重力加速度,kw和ka分別代表孔隙水和孔隙氣的滲透系數(shù).根據(jù)Fredlund[29]理論,水和空氣的滲透系數(shù)可表示為

式中,ηw,ηa分別代表孔隙水和孔隙氣的黏性系數(shù);krw,kra分別代表孔隙水和孔隙氣的相對(duì)滲透系數(shù);κ代表土的固有滲透率.相對(duì)滲透系數(shù)可由土-水特征曲線導(dǎo)出,本文采用V-G 模型[30]可得

式中,Se代表有效飽和度;α,m,d代表擬合參數(shù),m=1-1/d;Se=(Sr-Sw0)/(1-Sw0),Sw0代表水的束縛飽和度;pc=pa-pw為基質(zhì)吸力.

1.2 方程的求解

在簡(jiǎn)諧振動(dòng)中,所有變量可表示為f=的形式,其中表示變量f的幅值(簡(jiǎn)化起見(jiàn),下文公式推導(dǎo)略去～號(hào)),ω 為圓頻率.對(duì)空間變量x進(jìn)行傅里葉變換及其逆變換

式中,ξ 為x方向上的波數(shù);f是空間域變量,是相應(yīng)頻率域變量;i為虛數(shù),t代表時(shí)間.

式(1a)和(1b)可改寫(xiě)為

對(duì)式(5)兩邊取散度后代入式(1c)和(1d)中可得

式中b11=baρa(bǔ)A11,b12=baρa(bǔ)A12,b13=baρa(bǔ)A13+ω2ρa(bǔ),b21=bwρwA21,b22=bwρwA22,b23=bwρwA23+ω2ρw.

將式(5)代入到式(1e)中可得

式中,b1=-aχ-ω2/bw,b2=-a(1-χ)-ω2/ba,b3=ρω2+ω4ρw/bw+ω4ρa(bǔ)/ba.

同樣,對(duì)式(5)兩邊取散度后代入到式(7)并聯(lián)立式(6)可得

對(duì)方程式(6)和式(8)中的空間變量z進(jìn)行傅里葉變換

求解常微分方程組(9)可得通解如下

將式(7)進(jìn)行傅里葉變換可重寫(xiě)為

結(jié)合式(10)和式(11)可得兩個(gè)二階線性非齊次微分方程,求解可得其通解為

對(duì)土骨架體積應(yīng)變e=?·u進(jìn)行傅里葉變換并聯(lián)立式(10c)和式(12)可得

將式 εij=代入到非飽和土應(yīng)力-應(yīng)變關(guān)系式 σij=λeδij+2μεij-δijap中,式中P為等效孔隙流體壓力,表達(dá)為p=χpw+(1-χ)pa.通過(guò)對(duì)應(yīng)力-應(yīng)變關(guān)系式和式(5)進(jìn)行傅里葉變換后可得

結(jié)合方程式(10)和式(12),最后可得應(yīng)力、孔隙水和孔隙氣的相對(duì)位移在傅里葉變換域的通解為

方程式(10a)、式(10b)、式(12)和式(15)為非飽和土體動(dòng)力響應(yīng)問(wèn)題的應(yīng)力、位移和孔壓等物理量在頻率域的解答.

2 回傳射線矩陣法求解

本文考慮如圖1 所示的條形簡(jiǎn)諧荷載作用下梯度非均勻非飽和土地基的動(dòng)力響應(yīng)模型.非飽和土厚度為H,表面受到條形簡(jiǎn)諧荷載,其中荷載幅值為q0,分布長(zhǎng)度為 2l.

圖1 條形荷載作用下梯度非均勻非飽和土地基示意圖Fig.1 Physical model of graded non-homogeneous unsaturated soil under strip load

對(duì)于非均勻非飽和土,由于土體的非均勻性導(dǎo)致土體動(dòng)力問(wèn)題的控制方程為變系數(shù)偏微分方程,一般情況下很難獲得解析解.由于回傳射線矩陣法(RRMM)對(duì)非均勻單相彈性介質(zhì)和非均勻飽和多孔介質(zhì)都具有很好的適用性,因此本文嘗試采用RRMM 法對(duì)非均勻非飽和多孔介質(zhì)動(dòng)力控制方程進(jìn)行求解.為此,通常將非均勻材料在材料物理力學(xué)性質(zhì)變化方向上被簡(jiǎn)化為多個(gè)均勻材料薄層,層數(shù)N的數(shù)量取決問(wèn)題所需精度.回傳射線矩陣法的關(guān)鍵在于建立一對(duì)對(duì)偶坐標(biāo)系,圖2 給出了第i層的對(duì)偶坐標(biāo)系,(xi(i+1),zi(i+1))和 (xi(i-1),zi(i-1))中上標(biāo)i表示局部坐標(biāo)系位于第i層,上標(biāo) (i+1)和 (i-1)分別表示指向第i+1和i-1 層,通過(guò)建立局部坐標(biāo)系,可以有效避免計(jì)算過(guò)程中的數(shù)值問(wèn)題.

圖2 局部坐標(biāo)示意圖Fig.2 Dual local coordinates at the interface

2.1 散射關(guān)系

在局部坐標(biāo)系 (xi(i-1),0)下,式(10a)、式(10b)、式(12)以及式(15)可改寫(xiě)為

式中X,Mn,Nn,Ln,Hn的表達(dá)式見(jiàn)附錄A.

節(jié)點(diǎn)i處的應(yīng)力、位移以及孔壓在對(duì)偶坐標(biāo)系下的連續(xù)性條件為

利用式(16)和式(17)整理后得

式中入射波和出射波的波幅向量分別為節(jié)點(diǎn)i處的局部散射矩陣為Si,散射矩陣中各元素詳見(jiàn)附錄B.

考慮非飽和土地基表面受條形均布荷載作用,底面固定的邊界條件為

在z=0 處

在z=H處

利用式(16)、式(19)和式(20),整理后可得

式中

矩陣S0和SN中各元素詳見(jiàn)附錄C.

合并方程式(18)、式(21)和式(22)可得

方程(23)中提供的 8N個(gè)方程不滿足向量A和B中 16N個(gè)未知量解的要求,因此須再尋找一組方程.

2.2 相位關(guān)系和回傳射線矩陣

將式(16d)代入式(24)中可得

引入新矢量

結(jié)合式(23)、式(25)和式(26)可得

式中P為整體相位矩陣,定義為

U為整體置換矩陣,定義為

將式(27)代入式(23),最終可得

式中R=SPU稱(chēng)為回傳射線矩陣,I為單位矩陣.

通過(guò)式(30)求出C和D后,代入式(16),就可以獲得非飽和土地基中任一點(diǎn)處的應(yīng)力、位移和孔壓等物理量,再通過(guò)傅里葉逆變換即可求得相應(yīng)的各物理量在空間域上的表示.

3 數(shù)值算例

3.1 驗(yàn)證

為了驗(yàn)證本文算法的正確性,將本文非均勻非飽和土地基退化為均勻非飽和土地基,并與Shi[10]求得的條形荷載下非飽和土地基的解答進(jìn)行比較.驗(yàn)證時(shí)取計(jì)算參數(shù)為[10]:H=20 m ,μs=1.91×107Pa,ρs=2700 kg/m3,Sr=0.8,υ=0.2,n=0.6,Ks=3.6×1010Pa,k=1×10-12m2,q0=1 kPa,ω=1 rad/s,γ=0 .圖3 給出了本文數(shù)值解與文獻(xiàn)解的計(jì)算結(jié)果,從圖中可以看出,本文計(jì)算結(jié)果與文獻(xiàn)結(jié)果[10]比較接近,從而驗(yàn)證了本文算法的有效性.

圖3 本文數(shù)值解與文獻(xiàn)解的比較Fig.3 Comparison between present work and analytical

3.2 數(shù)值結(jié)果和討論

考慮沿深度方向非均勻非飽和土層的物理力學(xué)參數(shù)按冪函數(shù)變化,則有

式中,G(x3)為在x3處的 λ,μ,n,ρ,k,υ,Sr等物理力學(xué)參數(shù);γ 為梯度因子,γ 值的不同代表非飽和土體的非均勻程度不同.可以看出,當(dāng) γ=0 時(shí),非飽和土物理力學(xué)參數(shù)為土層表面處的材料參數(shù);當(dāng)γ→∞時(shí),非飽和土物理力學(xué)參數(shù)趨向于土層底面處的材料參數(shù).

作為本文解答的具體應(yīng)用,取土層厚度H=40 m,層數(shù)N=40 ,均布荷載幅值q0=1 kPa,頻率ω=1 rad/s,l=1 m .非飽和土地基的物理力學(xué)參數(shù)見(jiàn)表1所示[31].簡(jiǎn)化起見(jiàn),假設(shè)沿深度方向各物理力學(xué)參數(shù)按同樣的梯度因子變化.

表1 非均勻地基的物理力學(xué)參數(shù)Table 1 Physico-mechanical properties of non-homogeneous foundation

為了研究梯度因子對(duì)非均勻非飽和土地基動(dòng)力響應(yīng)的影響規(guī)律,圖4 給出了梯度因子 γ 對(duì)豎向位移沿深度z的變化曲線.從圖中可以看出,隨著梯度因子 γ 的增大,豎向位移幅值剛開(kāi)始隨之顯著減小,但隨著 γ 的繼續(xù)增大豎向位移降低幅度開(kāi)始變小.這是由于梯度因子 γ 的增大,使得抵抗變形能力強(qiáng)的基底材料增多,隨著梯度因子 γ 的繼續(xù)增大導(dǎo)致土的整體剛度較大,因此豎向位移幅值增量較小.

圖4 梯度因子 γ 對(duì)豎向位移沿深度z 變化的影響曲線Fig.4 Influence of the gradient factor γ on vertical displacement variations against with depth z

圖5 和圖6 給出了梯度因子 γ 對(duì)孔隙水壓和孔隙氣壓沿深度z的變化曲線.從圖5 中可以看出,隨著梯度因子 γ 的增大,孔隙水壓的幅值先增大后減小,并且土體非均勻程度越高,孔隙水壓的幅值越大,當(dāng)非飽和土體趨于均質(zhì)土?xí)r(γ=0,50),孔隙水壓幅值大小相接近.從圖6 中可以發(fā)現(xiàn),隨著梯度因子γ的增大,孔隙氣壓幅值不斷增大,并且其沿深度方向的振動(dòng)頻率增大,波峰值不斷靠近地表處附近.

圖5 梯度因子 γ 對(duì)孔隙水壓沿深度z 變化的影響曲線Fig.5 Influence of the gradient factor γ on pore water pressure variations against with depth z

圖6 梯度因子 γ 對(duì)孔隙氣壓沿深度z 變化的影響曲線Fig.6 Influence of the gradient factor γ on pore air pressure variations against with depth z

圖7 給出了梯度因子 γ 對(duì)正應(yīng)力沿深度z的變化曲線.圖中可以看出,正應(yīng)力幅值隨著梯度因子γ的增大先增大后減小,并且土體非均勻程度越高,正應(yīng)力的幅值越大,當(dāng)非飽和土體趨于均質(zhì)土?xí)r(γ=0,50),正應(yīng)力的幅值大小相接近,這與孔隙水壓幅值大小受土體非均勻程度的影響相同.綜合圖4～圖7還可以發(fā)現(xiàn)梯度因子 γ 顯著地改變了豎向位移、正應(yīng)力、孔隙水壓和孔隙氣壓沿深度z方向上的振動(dòng)模態(tài).

圖7 梯度因子 γ 對(duì)正應(yīng)力沿深度z 變化的影響曲線Fig.7 Influence of the gradient factor γ on normal stress variations against with horizontal depth z

4 結(jié)論

本文基于多孔介質(zhì)混合物理論,研究了條形荷載作用下梯度非均勻非飽和土地基的動(dòng)力響應(yīng)問(wèn)題,建立了非均勻非飽和土地基在條形荷載作用下回傳射線矩陣法的計(jì)算列式,分析了非飽和土體物理力學(xué)性質(zhì)沿深度方向按梯度變化時(shí)應(yīng)力、位移以及孔隙壓力等物理量的變化規(guī)律.結(jié)果表明:(1)土體非均勻性對(duì)非飽和土的動(dòng)力響應(yīng)具有顯著影響,梯度因子改變了豎向位移、正應(yīng)力和孔隙壓力在深度方向上的振動(dòng)模態(tài),其中孔隙氣壓在深度方向的振動(dòng)頻率隨著梯度因子的增加而不斷增大,峰值不斷靠近地表處附近；(2)豎向位移隨梯度因子的增大不斷減小,正應(yīng)力和孔隙壓力隨著梯度因子的增大先增大后減小,并且土體非均勻程度越高,其相應(yīng)幅值越大;而孔隙氣壓隨梯度因子的增大不斷增大.

附錄A

附錄B

附錄C

附錄D