地下工程的巖石力學數(shù)值模擬綜述自20世紀70年代以來，數(shù)值模擬方法作為一種科學有效而又快速簡便的分析方法，被逐步引人到交通、水工、采礦、建筑等地下工程的力學分析、穩(wěn)定性評價、方案比較等工作中。目前，這種方法已作為一種常用的手段在實際工程中大量采用。1數(shù)值模擬方法選擇現(xiàn)今可用于對巖體工程結(jié)構(gòu)進行力學分析的數(shù)值方法多種多樣，每一種方法有其針對性和特點，對一個具體的問題用數(shù)值模擬方法進行分析時，應選擇一種最適合該問題的方法進行研究。數(shù)值模擬方法的選擇，取決于研究對象即巖體工程結(jié)構(gòu)的巖石力學性質(zhì)和數(shù)值模擬的目的。嚴格地講，巖石除具有彈性性質(zhì)外，還具有塑性性質(zhì)和粘性性質(zhì)，只不過在特定情況下，某種性質(zhì)占主導地位而已。在巖體工程實踐中，硬巖及應力水平不甚高的中硬巖，其力學性質(zhì)主要呈現(xiàn)為彈性或彈塑性；高應力環(huán)境下的軟巖，其力學性質(zhì)主要呈現(xiàn)為塑性或粘塑性；對于服務(wù)時間較長的地下工程，巖石極軟或軟且應力水平很高，則在計算分析中不能忽視巖石的流變性質(zhì)。為了達到了解整個巖體工程系統(tǒng)的應力及變形規(guī)律的目的，各種數(shù)值方法均可采用，但以彈、塑性有限單元法或拉格朗日元法最為適宜。這兩種方法的單元劃分靈活，計算所需參數(shù)較少且易獲得，軟件也易于得到，成本較低。局部工程結(jié)構(gòu)的應力及變形分析，若巖石中硬以上，則各種方法均可采用；若巖石軟弱，則宜采用能進行大變形分析的拉格朗日元法；若巖體可能發(fā)生非連續(xù)破壞，則宜采用離散單元法、非連續(xù)變形法或流形元法或大變形分析的拉格朗日元法；若巖體可能發(fā)生非連續(xù)破壞，則宜采用離散單元法、非連續(xù)變形法或流形元法。要實現(xiàn)了解巖體結(jié)構(gòu)的破壞特征及動態(tài)破壞過程的目的，則只能采用離散單元法、非連續(xù)變形法或流形元法，因為這些方法就是針對巖體介質(zhì)的非連續(xù)性而提出的。此外，對于一個具體問題，是進行平面分析或是進行三維分析，也需做出恰當選擇。嚴格地講，所有的問題都是三維的，但如果采用平面分析既能達到目的，計算結(jié)果誤差也不大，為了降低費用和快速方便起見，則以采用平面分析為宜。反之，則應采用三維程序代碼進行計算分析。表1給出了幾種數(shù)值模擬方法所依據(jù)的基本原理、求解方式、離散化方法及其適用條件，可供選擇模擬方法時參考。表1 幾種數(shù)值模擬方法的對比數(shù)值模擬方法基本原理求解方法離散方式適用條件有限單元法最小勢能原理解方程組全區(qū)域劃分單元巖石中硬以上，小變形，巖體不會發(fā)生非連續(xù)性破壞如滑動、轉(zhuǎn)動、分離等邊界單元法Betti互等定理解方程組邊界上劃分單元同上離散單元法牛頓運動定律顯式差分按結(jié)構(gòu)弱面分布特征劃分單元巖石中硬以上，低應力水平，大變形，巖體沿弱面發(fā)生非連續(xù)性破壞非連續(xù)變形法最小勢能原理解方程組按主要結(jié)構(gòu)弱面實際情況劃分單元大變形，巖體發(fā)生非連續(xù)性破壞數(shù)值流形法最小勢能原理解方程組全區(qū)域劃分單元中硬以上巖體的連續(xù)或非連續(xù)變形拉各朗日元法牛頓運動定律顯式差分全區(qū)域劃分單元巖石軟弱，大變形，巖體的破壞以變形為主2 計算模型的建立數(shù)值模擬方法確定以后，建立一個能充分符合原型幾何特征和物理性質(zhì)的計算模型至關(guān)重要。建立的模型與實際原型越接近，則計算結(jié)果越能代表原型的實際情況。但是，實際工程的幾何形狀及巖體的物理力學性質(zhì)是極其復雜的，要使建立的計算模型與工程原型在幾何和物理性質(zhì)上完全相同幾乎是不可能的，因此，應根據(jù)工程的實際情況和研究目的對原型進行適當?shù)墓こ烫幚恚员氵M行計算分析。21 幾何形狀若需對較大范圍內(nèi)的地下工程進行力學分析，可按照工程施工設(shè)計建立數(shù)值模擬分析所需的幾何模型。如要對具體部位巖體的應力及變形進行分析，則應按照開挖形成的實際幾何形狀或統(tǒng)計平均情況建立模型。22 計算范圍計算范圍原則上按照開挖(采)活動所引起的圍巖應力重新分布的顯著影響范圍確定。對于隧道型或洞室型開挖工程，理論分析表明，開挖活動的顯著影響范圍為最小開挖尺寸的35倍。對較大范圍的地下工程進行力學分析時，其顯著影響范圍取決于工程的賦存深度、開挖規(guī)模與賦存深度之比等。一般情況下，計算范圍向上取至地表，向下取開挖空間較小尺寸的O510倍即可，水平方向的影響范圍自開挖邊界向外推最大開挖深度的距離。23 邊界條件在地下工程巖石力學數(shù)值模擬計算分析中，計算范圍通常取至開挖影響范圍之外，除地表為應力邊界(一般為零應力邊界)外，下邊界可為全約束或鉛垂方向約束邊界。前后左右邊界視情況不同而進行相應處理：若僅考慮巖體的自重作用，則前后左右邊界可為垂直于邊界平面方向上約束或全約束；若有水平構(gòu)造應力存在，且計算分析所用代碼不能對地應力初始化，則前后左右邊界應為應力邊界。24 初始狀態(tài)初始狀態(tài)包括原始地應力狀態(tài)、大型地質(zhì)構(gòu)造面(巖層層面、巖性界面、斷層破碎帶、軟弱夾層、裂隙等)、地下水賦存及其徑流狀態(tài)等。對一個具體的地下工程而言，其所處區(qū)域的原始地應力可以認為其大小和方向是恒定不變的，三個主應力的性質(zhì)為壓應力。大量的實測結(jié)果表明，鉛垂應力分量為原始地應力場的主應力之一，且其量值取決于巖體自重，即。原始地應力的水平主應力方向可以通過分析該地區(qū)長期以來的地質(zhì)構(gòu)造運動過程、地質(zhì)構(gòu)造特征以及該地區(qū)地形地貌特征來確定。原始地應力的大小可用實測方法獲得，如鉆孔應變法、孔徑變形法和水壓致裂法等，也可用采用反分析方法。前者應用較為成熟，但測量過程復雜，成本較高，不宜大量采用。后者方法簡單，成本低廉，便于大范圍長期進行。值得指出的是，測量得到的原始地應力主方向，尚需與地質(zhì)分析方法得到的主應力方向進行比較和綜合分析，才能得出正確結(jié)論。大型地質(zhì)構(gòu)造面中的巖層層面和巖性界面如兩側(cè)巖石結(jié)合良好，即巖層層面和巖性界面具有與兩側(cè)巖石相同數(shù)量級的C(粘聚力)、(內(nèi)摩擦角)值，則在計算分析中可不考慮這些構(gòu)造面的影響，但對于影響巖體應力分布和巖體穩(wěn)定性的斷層破碎帶、軟弱夾層和大型裂隙則必須予以考慮。斷層破碎帶和厚度較大的軟弱夾層可按具有相應物理力學性質(zhì)的巖層處理，厚度較小的軟弱夾層和大型裂隙用節(jié)理單元進行處理。在水體下(如江、河、湖等)的地下工程，若水對采礦有顯著影響，則需在探明巖體滲透性的條件下，用偶合力學分析和滲流分析的程序代碼進行計算分析，以便考慮地下水的靜、動壓對巖體應力及位移的影響。25 模型離散化模型離散化主要有單元劃分多少和如何劃分單元兩個問題。對于將巖體看作連續(xù)介質(zhì)的計算分析方法如有限單元法和拉格朗日元法，劃分多少個單元，取決于計算范圍大小、計算精度要求、計算時間(與計算成本成正比)。計算范圍越大，計算精度要求越高，則需劃分的單元也應越多，反之亦反。試驗表明，計算范圍一定時，單元越多，計算精度越高，但單元數(shù)達到某一數(shù)量后，隨著單元數(shù)的增加，計算精度的提高越來越小，即對特定的問題，存在一個適宜的劃分單元數(shù)范圍。對于將巖體看作非連續(xù)介質(zhì)的計算分析方法，如離散單元法，劃分多少單元取決于計算范圍大小，計算精度要求、巖體結(jié)構(gòu)特征、工程結(jié)構(gòu)單元尺寸和計算時間。與有限單元法不同的是：靠近開挖空間的巖體劃分的單元大小，與工程結(jié)構(gòu)尺寸相比，應不會從根本上改變實際巖體的工程特征，并盡量與實際巖體的巖塊大小一致。劃分單元的疏密變化的原則是：應力分布變化劇烈的部位，單元劃分應相對密一些；應力分布均勻的部位，單元劃分相對疏一些。26 支護結(jié)構(gòu)的作用這里所謂的支護特指錨桿、錨索和噴射混凝土支護。因目前所用的大多數(shù)軟件只能解彈性或彈塑性問題，即認為巖體為彈性材料或彈塑性材料，每一步開挖或開采后，巖體的變形及破壞是瞬時完成的。所以，若按照施工順序施加支護結(jié)構(gòu)，因該步開挖引起的巖體變形已完成，施加的支護對本次開挖無效，這顯然低估了支護的作用；若在開挖或開采前施加支護，則對于開挖或開采后再施加支護的實際情況而言，無疑又夸大了支護的作用。實際情況是：由于巖體具有滯后和流變性質(zhì)，開挖引起的變形不是瞬時完成，而是持續(xù)相當一段時間，達到穩(wěn)定的時間取決于巖體性質(zhì)和載荷大小等。一般情況下，施加支護時，變形仍未完成，支護結(jié)構(gòu)對本次開挖或開采仍起一定作用，而作用大小取決于達到穩(wěn)定的時間和施加支護的時刻。施加支護時已釋放的載荷或變形，通常是根據(jù)經(jīng)驗取一適當?shù)南禂?shù)，因各人的認識和經(jīng)驗不同，這一系數(shù)的取值主觀性較大，不便于在實際工程中采用。為此，(1)通過觀測開挖工程巖體變形，得到巖體變形隨時間的變化規(guī)律；(2)按照支護的施工時間，由變形一時間曲線確定支護時巖體變形已完成的比例；(3)假設(shè)支護時開挖等效荷載釋放比例與巖體變形完成比例相等(因為巖體變形與載荷成正比)，計算支護時已釋放的載荷；(4)將本次開挖支護分為兩步施工(第二步計算時撤除第一步的殘余載荷)進行計算分析即可。27 巖體性質(zhì)及力學指標用有限單元法或拉格朗日元進行計算分析時，對于硬度很大、彈性性質(zhì)十分明顯的巖體，可以進行彈性計算分析。但在一般情況下，巖體既有彈性性質(zhì)，同時也有塑性性質(zhì)和粘性性質(zhì)，對于埋深不大的中硬以上巖體，將巖體視為彈塑性體進行計算分析已能滿足工程需要。巖體力學性質(zhì)參數(shù)指標的選取，關(guān)系到計算結(jié)果是否可信，必須慎重對待。實驗室試驗獲得的巖石力學性質(zhì)指標一般不能直接用于計算分析，必須進行工程處理。巖體的強度參數(shù)可用格吉方法和費申科方法等進行弱化，變形參數(shù)(彈性模量)可用霍克布朗方法進行處理，反分析方法獲得的巖體力學參數(shù)一般具有較高的可信度。3 計算結(jié)果分析用數(shù)值模擬方法對地下工程進行巖石力學分析的目的主要有：(1)研究開挖過程中巖體應力及變形變化規(guī)律，開挖對其他重要工程、地表重要構(gòu)筑物的影響及影響程度，開挖引起地表沉降的規(guī)律，巖體加固的作用和作用機理等；(2)根據(jù)計算獲得的巖體應力大小及分布狀態(tài)，判斷圍巖及重要工程構(gòu)元的穩(wěn)定狀態(tài)或破壞情況；(3)依據(jù)巖體應力分布狀態(tài)和破壞情況，進行開挖方案比較，優(yōu)化開挖工藝和工程結(jié)構(gòu)參數(shù)及支護參數(shù)等。支護對巖體的加固作用和作用機理可以通過分析加固前后圍巖的應力狀態(tài)變化、關(guān)鍵部位應力及變形大小的變化來表示。例如用錨桿加固后，原來巖體中的拉應力區(qū)或塑性破壞區(qū)減小或完全消失，或者拉應力轉(zhuǎn)變?yōu)閴簯Γ敯宄两碉@著減小，則說明錨桿加固發(fā)揮了應有作用。如將巖體看作彈性材料進行分析，則有限元或拉格朗日元法計算結(jié)果并不能直接給出巖體的穩(wěn)定狀態(tài)和破壞情況，在此情況下，一方面要看圍巖體中是否存在拉應力，拉應力大小是否達到或超過巖體的抗拉強度，如是則表明巖體發(fā)生破壞；另一方面，對關(guān)鍵部位或工程結(jié)構(gòu)構(gòu)元的壓應力區(qū)用恰當?shù)钠茐臏蕜t(如HoekBrown破壞準則、Griffith破壞準則等)判斷其是否發(fā)生破壞。作彈塑性分析時，有的應用軟件還能給出巖體的塑性破壞區(qū)域，但值得注意的是，巖體發(fā)生塑性破壞，并不一定意味著工程結(jié)構(gòu)構(gòu)元完全失去效用。用離散單元法進行分析時，除作上述分析外，還可依據(jù)巖體沿節(jié)理裂隙滑動、開裂和冒落情況進行分析。用數(shù)值模擬的計算結(jié)果進行方案比較，是數(shù)值模擬方法最擅長的方面。通常是比較方案之間的若干參考量來選擇最優(yōu)或較優(yōu)方案，這些參考量包括：巖體中重要部位的拉應力區(qū)及拉應力大??；巖體破壞區(qū)大??；關(guān)鍵部位的變形量；關(guān)鍵部位的應力分布狀態(tài)。在工程量和施工成本基本相同的條件下，巖體中拉應力區(qū)及拉應力越小，巖體破壞區(qū)越小，巖體變形量越小，應力分布狀態(tài)越有利于巖體穩(wěn)定，則方案越優(yōu)。方案之間參考量矛盾時，需進行綜合分析確定最優(yōu)方案。