巖土工程數(shù)值計(jì)算方法

1、巖土工程數(shù)值計(jì)算方法報(bào)告學(xué)院：土木與環(huán)境工程學(xué)院姓名：xxxxxx學(xué)號(hào)：xxxxxxxx三維有限差分穩(wěn)定性分析一、 FLAC3D 基本原理FLAC(Fast Lagrangian Analysis of Continua)是由美國(guó)Itasca 咨詢公司研究開發(fā)的顯式有限差分程序，可用于工程力學(xué)計(jì)算，模擬巖石、土等材料的力學(xué)行為。由于其采用了顯式拉格朗日算法及混合離散劃分單元技術(shù)，使得該程序能較好地模擬地質(zhì)材料在達(dá)到強(qiáng)度極限或屈服極限時(shí)發(fā)生的破壞和塑性流動(dòng)，分析漸進(jìn)破壞和失穩(wěn)，特別適用于模擬大變形。材料通過(guò)單元和區(qū)域表示，根據(jù)計(jì)算對(duì)象的形狀構(gòu)成相應(yīng)的網(wǎng)格。每個(gè)單元在外載和邊界約束條件下，按照給出

2、的本構(gòu)關(guān)系產(chǎn)生力學(xué)響應(yīng)。FLAC 軟件主要是為巖土工程穩(wěn)定性分析開發(fā)的巖石力學(xué)計(jì)算程序，它包括了反映地質(zhì)材料力學(xué)效應(yīng)的特殊計(jì)算功能，能夠計(jì)算地質(zhì)類材料的高度非線性(包括應(yīng)變硬化/軟化)、不可逆剪切破壞和壓密、粘彈(蠕變)、空隙介質(zhì)的應(yīng)力滲流耦合及動(dòng)力學(xué)行為等。FLAC 提供了多種材料本構(gòu)模型：各向同性彈性模型、橫觀各向同性彈性模型、摩爾-庫(kù)侖塑性模型、應(yīng)變硬化/軟化塑性模型、德魯克-普拉格塑性模型、遍布節(jié)理模型、雙屈服塑性模型、霍克-布朗模型、空單元模型等。另外，程序設(shè)有界面單元，可以模擬斷層、節(jié)理和摩擦邊界的滑動(dòng)、張開和閉和行為。支護(hù)結(jié)構(gòu)，如砌襯、錨桿、支架等與圍巖的相互作用也可以在FLAC

3、 中進(jìn)行模擬。同時(shí)，用戶可根據(jù)自己的需要在FLAC 中創(chuàng)建自己的本構(gòu)模型，進(jìn)行各種特殊修正和補(bǔ)充。FLAC 采用顯式算法來(lái)獲得模型全部運(yùn)動(dòng)方程的時(shí)間步長(zhǎng)解，從而可以追蹤材料的漸進(jìn)破壞和跨落，這對(duì)研究開采的時(shí)間效應(yīng)和空間效應(yīng)是非常重要的。此外，程序允許輸入多種材料類型，亦可在計(jì)算過(guò)程中改變某個(gè)局部的材料參數(shù)，增強(qiáng)了程序使用的靈活性，用來(lái)提供采動(dòng)區(qū)域的跨落過(guò)程和開采中的充填過(guò)程。FLAC 具有強(qiáng)大的后處理功能，用戶可以直接在屏幕上繪制圖形，或以文件形式創(chuàng)建和輸出打印多種形式的圖形。使用者還可以根據(jù)需要，將若干個(gè)變量合并在同一幅圖形中進(jìn)行研究分析。FLAC3D 的計(jì)算過(guò)程如圖1.1 所示，在每個(gè)時(shí)步

4、內(nèi)，首先，根據(jù)高斯定律，由節(jié)點(diǎn)速度求出新的應(yīng)變速度；根據(jù)本構(gòu)方程計(jì)算出各時(shí)步內(nèi)單元新的應(yīng)力；然后，由運(yùn)動(dòng)方程計(jì)算新的結(jié)點(diǎn)速度和位移。在各個(gè)時(shí)步內(nèi)進(jìn)行循環(huán)計(jì)算，當(dāng)計(jì)算達(dá)到平衡狀態(tài)時(shí)，不平衡力趨近于0，如果有塑性流動(dòng)產(chǎn)生，則不平衡力趨近某一數(shù)值。圖1.1 FLAC3D循環(huán)求解過(guò)程二、節(jié)理化模型節(jié)理化模型是各向異性塑性模型，它包括了包含在摩爾庫(kù)侖體內(nèi)特殊方向上的弱面。根據(jù)應(yīng)力狀態(tài)、弱面走向以及模型體和弱面的材料特性的不同，屈服可能發(fā)生在模型體內(nèi)，或者發(fā)生在弱面上，或者在兩個(gè)部位同時(shí)發(fā)生。這種模型在FLAC中的實(shí)現(xiàn)方法是首先判別總體破壞，同時(shí)應(yīng)用到和FLAC中摩爾庫(kù)侖模型中相同的相關(guān)塑性修正，然后對(duì)更

5、新的應(yīng)力在弱面上產(chǎn)生的破壞進(jìn)行分析，同時(shí)對(duì)這些應(yīng)力分別進(jìn)行進(jìn)一步的校正。弱面內(nèi)的破壞準(zhǔn)則存在于包含了拉應(yīng)力路徑的摩爾庫(kù)侖屈服條件的局部形式中，與局部剪切流動(dòng)法則不相關(guān)聯(lián)而與局部拉應(yīng)力流動(dòng)法則相關(guān)聯(lián)。圖1.2顯示了整體坐標(biāo)系（x，y）和局部坐標(biāo)系(x，y)下存在于摩爾庫(kù)侖體內(nèi)的軟弱面。圖1.2 整體坐標(biāo)系中沿角方向的軟弱面為了簡(jiǎn)化這部分的符號(hào)，定義ij，它對(duì)應(yīng)于由各階段總體破壞的塑性修正的應(yīng)用引起的應(yīng)力分量，這種總應(yīng)力分解成局部應(yīng)力后可表示為：式中：節(jié)理夾角（從x坐標(biāo)軸逆時(shí)針?lè)较蛩闫穑?。依照約定，由表示弱面上切向引力分量的大小，相應(yīng)的應(yīng)變變量為，可以看到：有了這兩個(gè)符號(hào)，彈性增量的法則的局部表達(dá)

6、式可以表示為：式中：1=K4G3 ，2=K2G3 ，上標(biāo)e代表“彈性部分”。軟弱面的破壞準(zhǔn)則可以在圖1.3所示的，平面內(nèi)表示出來(lái)。圖1.3 FLAC中的軟弱面破壞準(zhǔn)則將摩爾庫(kù)侖破壞準(zhǔn)則定義為0，從點(diǎn)A到點(diǎn)B的局部破壞包絡(luò)線定義為： 將拉應(yīng)力破壞準(zhǔn)則定義為0，從點(diǎn)B到點(diǎn)C的局部破壞包絡(luò)線定義為： 式中：，分別為弱面的摩擦角、粘聚力和抗拉強(qiáng)度。對(duì)于摩擦角不為零的弱面，抗拉強(qiáng)度的最大值定義如下：剪切和張拉勢(shì)函數(shù)和對(duì)應(yīng)于不相關(guān)聯(lián)的流動(dòng)法則，剪脹角對(duì)應(yīng)于相關(guān)聯(lián)的流動(dòng)法則，它們分別表示為：破壞準(zhǔn)則邊界附近的流動(dòng)法則用FLAC中摩爾庫(kù)侖模型中所述的方法定義。這里，由函數(shù)h(，)=0表示(，)平面中0和0曲線

7、的對(duì)角線（圖1.4），此函數(shù)形式為：圖1.4 節(jié)理化模型中用于定義弱面流動(dòng)法則的區(qū)域其中和是兩個(gè)常量，定義如下：彈性假設(shè)和破壞準(zhǔn)則不一致，分別由(，)平面中位于1區(qū)或2區(qū)（分別對(duì)應(yīng)于h=0區(qū)域內(nèi)-或+區(qū)域）。如果位于1區(qū)，說(shuō)明平面內(nèi)是剪切破壞，應(yīng)用由勢(shì)函數(shù)確定的流動(dòng)法則，應(yīng)力點(diǎn)回歸到0的曲線上。如果位于2區(qū)，說(shuō)明是局部拉應(yīng)力破壞，應(yīng)用由勢(shì)函數(shù)確定的流動(dòng)法則，應(yīng)力點(diǎn)回歸到0的曲線上。首先考慮平面內(nèi)的剪切破壞，流動(dòng)法則如下：其中上標(biāo)表示弱面上同破壞有關(guān)的塑性部分，是待定參數(shù)，應(yīng)用式，進(jìn)行偏微分法后，上式可寫為：最終，由于局部應(yīng)力修正而分解到總體坐標(biāo)軸下得到的弱面上剪切破壞的總應(yīng)力修正式為：這些修正

8、被添加到應(yīng)力分量上， 包括了總體破壞的應(yīng)力修正，即這些應(yīng)力修正會(huì)計(jì)算出該時(shí)步新的應(yīng)力狀態(tài)。現(xiàn)在考慮弱面上的拉應(yīng)力破壞，這種情況下，流動(dòng)法則的形式為：這里是待定的參數(shù)，應(yīng)用式，進(jìn)行偏微分后，上式可寫為：應(yīng)用如上的推理，可以得到：分解到整體坐標(biāo)軸下以后，應(yīng)力修正變成了： 在大應(yīng)變模式下，為考慮由于剛體轉(zhuǎn)動(dòng)和變形引起的轉(zhuǎn)動(dòng)，對(duì)軟弱面的方位角進(jìn)行調(diào)整，修正值由一個(gè)區(qū)域內(nèi)所有角度的平均值得出，其表達(dá)式如下：其中： 且以弧度形式表示。三、采場(chǎng)開采模擬分析鞍千礦業(yè)有限公司許東溝采場(chǎng)設(shè)計(jì)長(zhǎng)1740m、寬280490m，最高海拔高度為+242.3m。設(shè)計(jì)邊坡階段高度12m，并段高度24m，階段坡面角上盤（東幫）

9、65o，下盤（西幫）55o，運(yùn)輸平臺(tái)、清掃平臺(tái)、安全平臺(tái)寬度分別取18m、8m、8m。最終邊坡角：上盤（東幫）45o55o、下盤（西幫）40o50o。礦層走向145°165°，傾角很陡，一般都超過(guò)80°。露天采場(chǎng)現(xiàn)狀最低開采標(biāo)高為+96m，+144m水平以上的平臺(tái)已經(jīng)靠幫到界，形成邊坡垂直高度為98m的高陡邊坡，開采現(xiàn)狀圖見(jiàn)圖1.5。局部邊坡+168m至+240m設(shè)計(jì)邊坡角為55°，由于開采初期沒(méi)有設(shè)計(jì)，實(shí)際邊坡境界+204m水平以上的部分已經(jīng)超挖，所以現(xiàn)邊坡角度為50°，小于設(shè)計(jì)的邊坡55°的角度。表1.1 FLAC數(shù)值模擬計(jì)算參數(shù)

10、表巖性名稱密度(g/cm3)粘聚力(kPa)內(nèi)摩擦角(°)體積模量(GPa)剪切模量(GPa)抗拉強(qiáng)度(kPa)綠泥石英片巖2.80400329.005.40200綠泥石英片巖片理20030遼河群千枚巖3.35450349.505.70500含鐵石英巖2.701000408.336.25500圖1.5 開采現(xiàn)狀圖圖1.6 開采最終設(shè)計(jì)圖采場(chǎng)邊坡分7步FLAC模擬開挖，模擬開采至72米、48米、 24米、 0米、-24米和-48米等四個(gè)水平的最大主應(yīng)力1等值線圖、最小主應(yīng)力3等值線圖、最大剪應(yīng)變max增量等值線圖、水平向位移ux等值線圖、垂直向位移uy等值線圖。 a開挖到72m時(shí)塑性區(qū)

11、 b開挖到48m時(shí)塑性區(qū) c開挖到24m時(shí)塑性區(qū) d開挖到0m時(shí)塑性區(qū) e開挖到-24m時(shí)塑性區(qū) f開挖到-48m時(shí)塑性區(qū)圖1.7 礦區(qū)各開挖步塑性區(qū)變化情況（1）礦山分6步模擬開采后，上盤邊坡192米96米水平間的巖體在臺(tái)階開挖的瞬間出現(xiàn)了巖體單元的塑性屈服和弱面破壞，開采擾動(dòng)塑性區(qū)寬度約30米，邊坡應(yīng)力場(chǎng)重新調(diào)整后，屈服單元重回彈性狀態(tài)；下盤邊坡巖體則僅在臺(tái)階表層出現(xiàn)塑性屈服狀態(tài)。 a b c d e f圖1.8 各開挖步主應(yīng)力變化情況 a b c d e f圖1.9 各開挖步最小主應(yīng)力變化圖（2）由圖1.8和圖1.9各開挖部最大和最小主應(yīng)力變化云圖可知，受礦山開采影響，邊坡巖體進(jìn)行應(yīng)力場(chǎng)

12、重新調(diào)整，坡面處應(yīng)力完全釋放，各剖面的最大主應(yīng)力1和最小主應(yīng)力3從坡面向坡體內(nèi)逐漸增加，并過(guò)渡到原巖應(yīng)力狀態(tài)。坡面附近最小主應(yīng)力3總體上為00.5MPa，個(gè)別平臺(tái)坡頂處出現(xiàn)了拉應(yīng)力，如168米水平、48米水平、24米水平和24米水平。 a b c d e f圖1.10 各開挖步水平向位移變化圖 a b c d e f圖1.11 各開挖步垂直位移變化情況（3）圖1.10和圖1.11為各開挖步水平和垂直位移云圖，分析水平向位移ux等值線圖和垂直向位移uy等值線圖表明，邊坡巖體質(zhì)點(diǎn)主要為彈性卸荷回彈變形，量級(jí)較小，位移較大區(qū)域位于邊坡中下部，邊坡中下部比較危險(xiǎn)需要重點(diǎn)考慮一些防護(hù)措施。 a b c d e f圖1.12 各開挖步剪應(yīng)變?cè)隽浚?）圖1.12為各開挖部剪應(yīng)變?cè)隽吭茍D，從最大剪應(yīng)力max（即主應(yīng)力差值13）等值線