斷層錯動數(shù)值模擬中的幾個關(guān)鍵影響因素

張炳焜

(四川省建筑科學(xué)研究院，四川成都610081)

針對斷層錯動類問題的有限元數(shù)值模擬，有幾個關(guān)鍵問題值得注意，全面的考慮這些問題無疑能使計算結(jié)果更為準(zhǔn)確。本項(xiàng)研究分別考慮了以下幾個影響因素：準(zhǔn)靜態(tài)法、網(wǎng)格尺寸、應(yīng)變軟化和材料阻尼。為了能夠說明這四個影響因素的重要性，本文建立簡單均一上覆砂土層二維有限元模型，深度H=20 m，寬度L=80 m，土層內(nèi)摩擦角φ=35°，剪脹角為ψ=15°，彈性模量E=60 MPa，泊松比ν=0.3，不考慮內(nèi)聚力的影響。以各影響因素為變量展開計算，模擬逆斷層60°，垂直位錯量d=1 m時上覆土層應(yīng)變與位移變形情況。

1 準(zhǔn)靜態(tài)分析

準(zhǔn)靜態(tài)法是一種用靜力學(xué)方法近似解決動力學(xué)問題的簡易方法，該方法能在有限程度上反映荷載的動力特性，它物理概念清晰，與全面考慮結(jié)構(gòu)物動力相互作用的分析方法相比，計算方法較為簡單，計算工作量很小、參數(shù)易于確定，并積累了豐富的使用經(jīng)驗(yàn)，易于設(shè)計工程師所接受。在進(jìn)行數(shù)值模擬時，有兩個方面值得注意，一是在對地震導(dǎo)致的地面震動問題進(jìn)行分析時，動力分析法是必不可少的。而另一方面，對于斷層錯動問題，若是綜合考慮斷層錯動產(chǎn)生的波動和位移這兩種因素的影響，則模擬結(jié)果會缺乏針對性。為了能夠有針對性的說明斷層錯動時上覆土層的響應(yīng)情況，選擇準(zhǔn)靜態(tài)法進(jìn)行模擬。ABAQUS/EXPLICIT模塊常常在解決準(zhǔn)靜態(tài)類型問題是非常有效的。

本文對斷層垂直錯動速度分別為v=1 m/s，v=0.5 m/s，v=0.1 m/s和v=0.05 m/s的情況進(jìn)行了計算和對比，從上覆土層塑性應(yīng)變云圖(圖1)可以發(fā)現(xiàn)，當(dāng)v=1 m/s時，塑性應(yīng)變帶未曾貫通上覆土層，而是較為集中的分布在了上覆土下部，由于高速的下部斷層錯動，導(dǎo)致上覆土層下部，即靠近斷層處的土體受到高速的撞擊而產(chǎn)生局部的應(yīng)變，該區(qū)域的應(yīng)變集中，消耗了本應(yīng)傳遞至上覆土上部的應(yīng)變能量。隨著錯動速度的降低，當(dāng)v=0.1 m/s時，塑性應(yīng)變帶較為完整的貫通了上覆土層。根據(jù)以上分析得出結(jié)論：相對高速的斷層錯動速度會導(dǎo)致上覆土下部應(yīng)變的集中，并且地面傾斜值低于相對低速的情況。而當(dāng)斷層錯動速度減小到一定程度后，若是繼續(xù)降低斷層錯動速度，覆土的塑性應(yīng)變和地面位移將不再發(fā)生明顯的變化，反而增加了計算所需的時間。因此，本文推薦采用的斷層錯動速度為0.1 m/s。

(a)v=1 m/s

(b)v=0.5 m/s

(c)v=0.1 m/s

(d)v=0.05 m/s圖1 塑性應(yīng)變云圖

2 網(wǎng)格尺寸

對有限元模擬計算而言，除了材料參數(shù)的選取以及模型的建立外，模型的網(wǎng)格劃分也是非常關(guān)鍵的，網(wǎng)格劃分的好壞直接影響計算的精度。斷層錯動情況的網(wǎng)格尺寸是個關(guān)鍵因素，與其他模擬不同，斷層錯動伴隨著較大的位移和應(yīng)變的發(fā)生，若是模型網(wǎng)格尺寸較大，則可能導(dǎo)致塑性應(yīng)變區(qū)的寬度發(fā)生變化，即破裂帶可能變寬，而影響了計算精度，對后續(xù)的研究造成誤導(dǎo)。本文有限元模型中部斷層區(qū)域是研究的重點(diǎn)，在該區(qū)域內(nèi)的網(wǎng)格尺寸更應(yīng)該進(jìn)行嚴(yán)格的控制，力求達(dá)到均勻和精細(xì)。

本文對模型中部網(wǎng)格尺寸分別為L=10 %H，L=7.5 %H，L=5 %H和L=2.5 %H的情況進(jìn)行了計算和對比，H為覆土深度。從上覆土層塑性應(yīng)變云圖(圖2)可以發(fā)現(xiàn)，隨著網(wǎng)格尺寸的減小，塑性應(yīng)變帶的寬度逐漸減小，當(dāng)L=2.5 %H時，塑性應(yīng)變帶寬度達(dá)到合理水平。因此，本文推薦網(wǎng)格尺寸L=2.5%H。

(a)L=10 %H

(b)L=7.5 %H

(c)L=5 %H

(d)L=2.5 %H圖2 塑性應(yīng)變云圖

3 應(yīng)變軟化

應(yīng)變軟化是指材料試件經(jīng)一次或多次加載和卸載后，進(jìn)一步變形所需的應(yīng)力比原來的要小，即應(yīng)變后材料變軟的現(xiàn)象。在斷層發(fā)生錯動后，涉及到較大的上覆土體變形和應(yīng)變，土體強(qiáng)度的降低是應(yīng)該考慮的，若是不考慮此因素的影響而進(jìn)行數(shù)值模擬，得出的模擬結(jié)果與實(shí)際可能會產(chǎn)生一定的差別，之前學(xué)者們對此得出了較好的對比結(jié)果[1-3](Roth et al.1981,1982;Loukidis 1999)，針對路堤邊坡穩(wěn)定性而言，亦有類似結(jié)論[4-5](Potts et al.1990,1997)。針對摩爾庫倫模型，內(nèi)摩擦角、剪脹角和內(nèi)聚力是該模型最主要的強(qiáng)度參數(shù)，將上覆土層的塑性剪應(yīng)變e作為控制參數(shù)，考慮內(nèi)摩擦角、剪脹角和內(nèi)聚力隨塑性剪應(yīng)變的增加而數(shù)值逐漸降低，當(dāng)塑性剪應(yīng)變達(dá)到某一個臨界值后，內(nèi)摩擦角、剪脹角和內(nèi)聚力將不再發(fā)生變化，保持其殘余值。此處，內(nèi)摩擦角和剪脹角的峰值分別記為φmax和ψmax，殘余值分別記為φres和ψres，塑性剪應(yīng)變臨界值記為e0。則內(nèi)摩擦角和剪脹角與塑性剪應(yīng)變e之間的關(guān)系為：

(1)

(2)

本文分別對考慮應(yīng)變軟化和未考慮應(yīng)變軟化的情況進(jìn)行數(shù)值模擬，當(dāng)考慮應(yīng)變軟化時，內(nèi)摩擦角和剪脹角的殘余值分別為φres=25°，ψres=0°，塑性剪應(yīng)變臨界值e0=0.1。

從以上兩種情況的塑性應(yīng)變云圖可以看出(圖3)，二者的塑性應(yīng)變帶均貫通了上覆土層，并且分布較為均勻，但是，當(dāng)未考慮應(yīng)變軟化時，破裂帶和水平面的夾角明顯較小。從某種意義上降低了斷層的傾角。為了能夠充分的估計斷層錯動對上覆土層的影響，考慮應(yīng)變軟化模型是合理的，有必要的。

(a)考慮應(yīng)變軟化

(b)未考慮應(yīng)變軟化圖3 塑性應(yīng)變云圖

4 材料阻尼

對于準(zhǔn)靜態(tài)問題而言，雖然該方法大大的降低了動力波動對模型的影響，但是材料的阻尼仍然是不能忽視的問題，本文考慮Rayleigh阻尼的影響。

在有限元分析中，當(dāng)采用振型分解法時，即把多自由度結(jié)構(gòu)的復(fù)雜振動分解為按各個振型的獨(dú)立振動的疊加，則常假定阻尼矩陣[C]滿足振型的正交條件，即:

(3)

瑞利(Rayleigh)曾提出阻尼矩陣采用質(zhì)量矩陣與剛度矩陣的線性組合，必滿足正交條件，即：

[C]=α[M]+β[K]

(4)

阻尼矩陣用轉(zhuǎn)置振型矩陣及振型矩陣前后相乘后，將得: 如令α=0，即C=β[K]，此時β=2ξn/ωn

ξn=(α/ωn+βωn)/2

(5)

如令β=0，即C=α[M]，此時α=2ξnωn

可見，當(dāng)阻尼矩陣僅與質(zhì)量矩陣有關(guān)時，阻尼比與振動頻率成反比，高階振型的阻尼就非常??；同樣，當(dāng)阻尼矩陣僅與剛度矩陣有關(guān)時，阻尼比正比于頻率，高階振型的阻尼非常大，這就使結(jié)構(gòu)的高階振型對結(jié)構(gòu)動力反應(yīng)的影響減弱甚至消除了。

系數(shù)α、β可由給定的第i、第j階振型的阻尼比ξi、ξj反算得到：

(6)

通常結(jié)構(gòu)振動分析中廣泛采用的是簡單的Rayleigh阻尼模型。

在工程中常采用簡化求法，即：

α=ξ1·ω1β=ξ1/ω1

(7)

本文分別對考慮材料阻尼和未考慮材料阻尼的情況進(jìn)行數(shù)值模擬，從塑性應(yīng)變云圖可以看出(圖4)，當(dāng)考慮材料阻尼時，塑性應(yīng)變帶貫通了上覆土層，并且分布較為均勻，而未考慮材料阻尼時，破裂帶除了貫通上覆土層外，在破裂帶左側(cè)地面出現(xiàn)了一個局部的塑性應(yīng)變區(qū)域，上覆土層破裂的不規(guī)律。因此，雖然采用準(zhǔn)靜態(tài)法進(jìn)行數(shù)值模擬，但是材料的阻尼仍是不容忽視的關(guān)鍵因素之一。

(a)考慮材料阻尼

(b)未考慮材料阻尼圖4 塑性應(yīng)變云圖

5 結(jié)論

通過以上計算分析，得出以下結(jié)論：對于地震引發(fā)的地面震動問題，動力分析方法是切實(shí)有效的。不過對于下覆斷層錯動類問題，為了能夠較有針對性的對問題進(jìn)行分析，建議采用準(zhǔn)靜態(tài)法進(jìn)行研究，并綜合考慮有限元網(wǎng)格尺寸、土體應(yīng)變軟化特性以及材料阻尼等影響因素。避免對上覆土層及結(jié)構(gòu)的影響估計不足，使分析結(jié)果與實(shí)際貼切，為日后的工程設(shè)計和研究提供幫助。

[1] Roth WH, Scott RF, Austin I (1981). Centrifugemodelling of fault propagation through alluvial soils.Geophys Res Lett 8(6):561-564

[2] Roth WH, Sweet J, Goodman RE (1982) Numerical and physical modelling of flexural Slip phenomena and potential for fault movement. Rock Mech 12:27-46 (Suppl.)

[3] Loukidis D (1999). Active fault propagation through soil. Diploma Thesis, School of Civil Engineering, National Technical University, Athens

[4] Potts DM, Dounias GT, Vaughan PR (1990). Finite element analysis of progressive failure of Carsington Embankment. Géotechnique 40(1):79-101

[5] Potts DM, Kovacevic N, Vaughan PR (1997). Delayed collapse of cut slopes in stiff clay. Géotechnique 47(5):953-982