李勇輝，孫 勇，朱伯文，唐江濤

(1.貴州大學(xué) 國土資源部喀斯特環(huán)境與地質(zhì)災(zāi)害重點(diǎn)實(shí)驗(yàn)室，貴陽 5500033； 2.成都理工大學(xué) 地質(zhì)災(zāi)害與防治與地質(zhì)環(huán)境保護(hù)國家重點(diǎn)實(shí)驗(yàn)室，成都 610059)

西南地區(qū)主要是山區(qū)，在修建公路鐵路的過程中會(huì)遭遇大量的邊坡。而邊坡的穩(wěn)定性分析是重中之重，首先要分析該邊坡的穩(wěn)定與否，如若需要進(jìn)行開挖，則需要對(duì)該開挖方案進(jìn)行穩(wěn)定性分析，從而判斷該方案的可靠度與可行性。近年來，數(shù)值模擬在邊坡穩(wěn)定性分析中得到大量的運(yùn)用，其中FLAC3D軟件因其模擬結(jié)果與實(shí)際情況更加相符，從而得到更廣泛的應(yīng)用。

1 FLAC3D原理

FLAC3D是由美國itasca公司開發(fā)的一種模擬仿真計(jì)算軟件，是對(duì)二維有限差分程序FLAC2D的擴(kuò)展，能夠?qū)ν馏w、巖石和其它材料進(jìn)行三維結(jié)構(gòu)的受力特性模擬和塑性流動(dòng)分析。通過對(duì)三維網(wǎng)格中的多面體單元進(jìn)行不斷調(diào)整來擬合實(shí)際的結(jié)構(gòu)單元。單元材料主要有線性和非線性兩種本構(gòu)模型，在外力作用下，當(dāng)材料達(dá)到屈服流動(dòng)后，網(wǎng)格可以產(chǎn)生相應(yīng)的變形和移動(dòng)(大變形模式)。FLAC3D主要采用顯式拉格朗日算法和混合-離散分區(qū)技術(shù)，能夠十分準(zhǔn)確模擬材料的塑性流動(dòng)和破壞。因?yàn)椴挥眯纬蓜偠染仃?，所以該軟件用較小的內(nèi)存空間就可以求解相當(dāng)大范圍的三維問題。FLAC3D主要有三大類模型：開挖模型、彈性模型和塑性模型。FLAC3D網(wǎng)格中的每個(gè)區(qū)域都可以賦予不同的材料模型，并且可以指定材料參數(shù)的統(tǒng)計(jì)分布以及變化梯度。而且模型中還包含有節(jié)理單元，也稱之為界面單元，可以模擬兩種以及多種材料界面，不同材料性質(zhì)的間斷特性。節(jié)理允許發(fā)生滑動(dòng)或者分離，因此能夠用于模擬巖體中的斷層、節(jié)理和摩擦邊界。FLAC3D中的網(wǎng)格生成器gen，可以通過匹配、連接在網(wǎng)格生成器中生成局部網(wǎng)格，從而方便地生成模擬所需要的三維結(jié)構(gòu)網(wǎng)格。還可以自動(dòng)產(chǎn)生交叉結(jié)構(gòu)網(wǎng)格(比如說相交的巷道)，三維網(wǎng)格由整體坐標(biāo)系X、Y、Z系統(tǒng)所確定，不同于FLAC程序是由行列方式確定。這就提供了比較靈活的產(chǎn)生和定義三維空間參數(shù)。FLAC3D計(jì)算步驟為建立網(wǎng)格，劃分網(wǎng)格，選擇本構(gòu)模型，賦參數(shù)，選擇初始邊界條件達(dá)到原始平衡，模擬開挖監(jiān)測坡內(nèi)響應(yīng)[1-2]。

2 研究區(qū)三維模型建立

2.1 模型的建立

正確的地質(zhì)模型是計(jì)算的基礎(chǔ)，在計(jì)算中不可能也無必要將所有細(xì)節(jié)考慮進(jìn)去，因此對(duì)工程地質(zhì)條件的深入認(rèn)識(shí)與抽象是建立合理地質(zhì)概化模型的前提和關(guān)鍵。研究區(qū)邊坡順擬建公路走向約121 m，垂直公路走向約151 m，最大垂直開挖深度約28 m。邊坡上堆積的殘坡積物較薄，基本為一巖質(zhì)邊坡，巖性主要為白云巖、泥質(zhì)白云巖、泥巖等三大類，根據(jù)現(xiàn)場調(diào)查，研究區(qū)邊坡并無大型斷層，結(jié)構(gòu)面面等，因此不考慮斷層及結(jié)構(gòu)面等的影響，為模型長度取121 m，寬度取151 m，底部高程為800 m，參照壩址區(qū)工程地質(zhì)平面圖、剖面圖以及邊坡巖級(jí)分帶圖等，建立研究區(qū)數(shù)值模型。FLAC3D數(shù)值模型建立方法大體上有兩種，一種是通過前處理較為方便的其他軟件建立數(shù)值模型并劃分網(wǎng)格，將網(wǎng)格節(jié)點(diǎn)和單元信息導(dǎo)出，經(jīng)處理后再導(dǎo)入FLAC3D中[3]；另一種是直接在FLAC3D程序中采用堆砌法生成模型,本文采用第一種方法，利用cad-surfer-ansys-flac聯(lián)合建模方法，建立其數(shù)值模型。見圖1。

圖1 研究區(qū)三維模型Fig.1 Three dimensional model of research area

2.2 模型本構(gòu)關(guān)系的選取

數(shù)值計(jì)算中，將巖土體看作理想彈塑性材料，屈服準(zhǔn)則采用內(nèi)置的Mohr-Coulomb準(zhǔn)則，即：

(1)

式中：NΦ=(1+sinφ)(1-sinφ)；σ1為最大主應(yīng)力；σ3為最小主應(yīng)力；φ為內(nèi)摩擦角；C為內(nèi)聚力。

若fs<0，表示巖體發(fā)生剪性屈服；fs≥0，表示巖體未發(fā)生剪性屈服。

當(dāng)法向應(yīng)力為張應(yīng)力時(shí)，則超出了莫爾-庫侖準(zhǔn)則的力學(xué)有效性范圍，此時(shí)要求最小主應(yīng)力不得超過巖體的抗拉強(qiáng)度，否則巖體將出現(xiàn)張性屈服，其判別公式為[4]：

ft=σ3-σt

(2)

若ft≥0，表示巖體張性屈服；若ft<0，表示巖體未發(fā)生張性屈服。

對(duì)于FLAC-3D數(shù)值計(jì)算所要求的參數(shù)體積模量(K)和剪切模量(G)可以通過下式進(jìn)行換算[5]：

(3)

(4)

式中：E為彈性模量；γ為泊松比。

實(shí)際上，在FLAC3D中使用M-C本構(gòu)模型時(shí)也可直接輸入彈模和泊松比，程序可實(shí)現(xiàn)自動(dòng)換算。

2.3 模型邊界的選取

FLAC3D提供兩種應(yīng)力邊界條件，一種是位移(速度)邊界條件，一種是應(yīng)力邊界條件。應(yīng)力邊界一般用于需考慮原巖應(yīng)力的情況下。對(duì)模型底部采用完全約束，約束其3個(gè)方向的自由度；而對(duì)模型的前后左右均只作單向約束，即僅約束水平方向速度而不約束其豎向速度。除此之外，由于上部房屋的影響，因此需要考慮施加均布載荷，就需要使用應(yīng)力邊界荷載條件。這樣，巖土體在生成自重應(yīng)力階段可發(fā)生自由沉降，符合實(shí)際情況。

3 邊坡穩(wěn)定性分析

3.1 自然斜坡穩(wěn)定性分析

為了解自然斜坡在天然條件下的穩(wěn)定性情況，考慮上部房屋荷載情況下斜坡在自身重力作用下的穩(wěn)定性，分別進(jìn)行斜坡應(yīng)力分析、斜坡變形分析以及斜坡穩(wěn)定性分析

3.1.1 斜坡應(yīng)力分析

斜坡應(yīng)力分析主要從斜坡的重力場以及最大最小主應(yīng)力場3個(gè)方面進(jìn)行，斜坡重力場模擬圖見圖2-圖7。

由圖2-圖3可知，坡體應(yīng)力量值范圍在-2.04～0.13 MPa，坡體應(yīng)力分布均勻，自上而下逐步增加，應(yīng)力最大值分布在斜坡地面，坡面由于施加有房屋荷載，所以在一定范圍內(nèi)，坡表的應(yīng)力值較大，極少數(shù)地方出現(xiàn)拉應(yīng)力。

圖2 坡體自重應(yīng)力云圖Fig.2 The slope self-weight stress nephogram

圖3 坡表自重應(yīng)力云圖Fig.3 The slope surface self-weight stress nephogram

由圖4-圖5可知，邊坡最大主應(yīng)力的量值范圍為-0.85～0.58 MPa,坡體的最大主應(yīng)力位于底面，與重力方向的應(yīng)力相似，坡體中大主應(yīng)力分布均勻，由上往下逐漸增大，坡體中無應(yīng)力集中現(xiàn)象。在坡表最大位移云圖中可以知道，位于右側(cè)的局部地區(qū)出現(xiàn)了拉應(yīng)力。主要是地形起伏大以及層面間的接觸面處。

圖4坡體最大主應(yīng)力云圖Fig.4 The slope max principal stress nephogram

圖5 坡表最大主應(yīng)力云圖Fig.5 The slope surface max principal stress nephogram

由圖6-圖7可知，斜坡的最小主應(yīng)力的量值范圍為-2.04～0.09 MPa,坡體中最小主應(yīng)力分布與重力分布類似，由上到下逐漸增大，均勻分布。坡表的最小主應(yīng)力云圖可知，只有作用均布荷載的地方應(yīng)力出現(xiàn)一點(diǎn)異常，其他都正常。

圖6 坡體最小主應(yīng)力云圖Fig.6 The slope mini principal stress nephogram

圖7 坡表最小主應(yīng)力云圖Fig.7 The slope surface mini principal stress nephogram

3.1.2 斜坡變形分析

坡體變形一定程度上可以反映邊坡的穩(wěn)定狀態(tài)，斜坡的變形分析主要從X方向以及Y方向的位移入手，模擬結(jié)果圖見圖8-圖9。

圖8 斜坡X方向位移云圖Fig.8 X direction displacement nephogram of slope

圖9 斜坡Y方向位移云圖Fig.9 Y direction displacement nephogram of slope

圖8-圖9為數(shù)值模擬結(jié)果圖，從作圖可以知道，坡表X方向的最大位移為2.2 mm，主要集中不同巖性相交位置處，坡表Y方向最大位移為0.6 mm，主要集中在施加均布荷載的前段。自然斜坡的變形量非常小

3.1.3 斜坡穩(wěn)定性分析

塑性區(qū)以及剪應(yīng)變的增量可以反映邊坡的穩(wěn)定狀態(tài)，一般來說，邊坡塑性區(qū)越大，深度越深，則對(duì)邊坡的不離程度也相即增大，對(duì)模擬結(jié)果進(jìn)行切片，分別沿斜坡橫向以及縱向切片，得到橫向以及縱向的最大剪應(yīng)變速率。見圖10、圖11。

圖10 斜坡橫向剪應(yīng)變速率云圖Fig.10 Transversal shear strain rate nephogram of slope

圖11 斜坡豎向剪應(yīng)變速率云圖Fig.11 Longitudinal shear strain rate nephogram of slope

由圖10、圖11可知，橫向以及豎向的剪應(yīng)變速率都非常小，幾乎可以忽略不計(jì)，說明斜坡處于一個(gè)穩(wěn)定的狀態(tài)。

3.2 開挖條件下邊坡穩(wěn)定性研究

開挖方案為沿著開挖線垂直開挖至指定高程843.0 m處，最大的開挖深度局部可以達(dá)到28.0 m，開挖將會(huì)形成陡立邊坡，模擬開挖之后邊坡的穩(wěn)定性以及變形。

首先在自然斜坡達(dá)到上述自重平衡的條件下，對(duì)原本存在的位移進(jìn)行清零處理，以便可以清楚掌握在開挖情況下邊坡產(chǎn)生的位移情況。

3.2.1 邊坡應(yīng)力分析

開挖邊坡主要從最大主應(yīng)力、最小主應(yīng)力進(jìn)行初步分析，分析是否存在應(yīng)力集中部位，模擬結(jié)果見圖12-圖16。

由圖12-圖13可以看出，開挖情況下，坡體的最大主應(yīng)力總體分布較均勻，自上而下逐步增加，最大主應(yīng)力的量值范圍為0.055～-0.85 MPa，但在開挖的坡腳出現(xiàn)了應(yīng)力集中的情況，應(yīng)力的大小約為-0.2 MPa。此外，在坡表作用于均布荷載的地方，應(yīng)力也出現(xiàn)“異常”，其他地方都應(yīng)力分布形式都為正常情況。

圖12 坡體最大應(yīng)力分布圖Fig.12 Max stress distribution in slope

圖13 坡表最大應(yīng)力分布圖Fig.13 Max stress distribution in slope surface

圖14 坡體最小應(yīng)力分布圖Fig.14 Mini stress distribution in slope

圖15 坡表最小應(yīng)力分布圖Fig.15 Mini stress distribution in slope surface

坡體最小主應(yīng)力圖的分布形式與坡體最大主應(yīng)力分布形式類似，應(yīng)力分布總體均勻，自上而下應(yīng)力逐漸增大，應(yīng)力分布的量值為0.008～-2.054 MPa，最小主應(yīng)力的最大值分布在坡地，最小值分布于坡表面，在開挖的坡腳段也出現(xiàn)了應(yīng)力集中的現(xiàn)象，作用于均布載荷的上部應(yīng)力“異?！保吰缕渌课欢继幱谡顩r。

從圖16可以看出，邊坡的最大剪應(yīng)力分布在開挖段坡腳部分，剪應(yīng)力的量值范圍為0.005～0.682 MPa，圖16紅色區(qū)域就是最大剪應(yīng)力分布區(qū)域，值的大小約為0.67 MPa。綜合以上分析可知，在開挖坡腳部分出現(xiàn)應(yīng)力集中的情況，可能對(duì)邊坡的穩(wěn)定性造成不小的危害，需要進(jìn)行進(jìn)一步的分析論證。

圖16 邊坡最大剪應(yīng)力分布圖Fig.16 Max shear stress distribution of slope

3.2.2 邊坡變形分析

對(duì)邊坡的變形分析主要從開挖之后坡體在X以及Y方向的位移及其位移值，布設(shè)的監(jiān)測點(diǎn)進(jìn)行分析。見圖18-圖19。

從邊坡X方向的位移云圖可知，X方向最大位移量值方位為-13.6～12.5 mm，最大位移主要是層面接觸方向以及均布荷載作用的邊緣，如圖左上角紅色區(qū)域?yàn)榫己奢d作用邊緣處。在Y方向的位移云圖中可知，位移的量值范圍為3.2～-93.2 mm，最大位移集在相對(duì)巖體力學(xué)性質(zhì)軟弱處，在邊坡的右側(cè)，由于其開挖量比較小，邊坡Y方向的位移值比較小，呈現(xiàn)mm級(jí)別的，但是在邊坡左側(cè)由于上部房屋的作用以及開挖高度大的影響，邊坡Y方向的位移呈cm級(jí)別，最大位移為9.3 cm，邊坡處于不穩(wěn)定狀態(tài)。

3.2.3 邊坡穩(wěn)定性分析

塑性區(qū)以及剪應(yīng)變的增量可以反映邊坡的穩(wěn)定狀態(tài)，一般來說，邊坡塑性區(qū)越大，深度越深，則對(duì)邊坡的不離程度也相繼增大，對(duì)模擬結(jié)果進(jìn)行切片，分別沿斜坡橫向以及縱向切片，得到橫向以及縱向的最大剪應(yīng)變速率。見圖19-圖22。

圖17 邊坡X方向位移云圖Fig.17 X direction displacement nephogram of slope

圖18 邊坡Y方向位移云圖Fig.18 Y direction displacement nephogram of slope

由圖19-圖22可以看到，X方向的最大剪切應(yīng)變處于開挖坡腳線附近，其量值范圍為1.00×10-6～8.48×10-6，Y方向的最大剪切應(yīng)變位于坡表，量值范圍為1.00×10-7～9.50×10-7。從Y方向塑性區(qū)的分布云圖可以知道，開挖線后約20 m范圍內(nèi)巖體處于受剪切的狀態(tài)，受剪區(qū)與開挖線的夾角大約為45°。在X方向塑性區(qū)分布云圖可知，在高程為850.0 m以上的巖體基本處于受剪切狀態(tài)，表面巖體受拉張作用。

圖19 X方向最大剪應(yīng)變云圖Fig.19 X direction max shear strain nephogram

圖20 Y方向最大剪應(yīng)變云圖Fig.20 Y direction max shear strain nephogram

圖21 X方向塑性區(qū)分布云圖Fig.21 X direction plastic zone distribution nephogram

圖22 Y方向塑性區(qū)分布云圖Fig.22 Y direction plastic zone distribution nephogram

4 結(jié) 論

由FLAC3D數(shù)值模擬軟件做出的以上兩種分別在自然工況和開挖工況的穩(wěn)定性分析，得到以下結(jié)論：

1) 自然條件下斜坡處于一個(gè)穩(wěn)定的狀態(tài)。

2) 在該方案開挖條件下，邊坡的位移量值非常大，邊坡開挖，房屋的穩(wěn)定性會(huì)受到極大的影響，容易造成房屋的開裂，建筑物發(fā)生拉裂縫，開挖線后20 m以及高程850.0 m以上的巖體都會(huì)因開挖受到非常大的影響，在此種方案下，邊坡處于不穩(wěn)定狀態(tài)，極易產(chǎn)生失穩(wěn)變形破壞。