隧道力學(xué)特征及數(shù)值模擬方法

1、2隧道力學(xué)特征及數(shù)值模擬方法2.1隧道幵挖生成的圍巖二次應(yīng)力場特征巖體在開挖前處于初始應(yīng)力狀態(tài)，初始應(yīng)力主要是由于巖體的自重和地質(zhì)構(gòu)造所引起的。在巖體進行開挖后改變了巖體的初始應(yīng)力狀態(tài)，使巖體中的應(yīng)力狀態(tài)重新分布，引起巖體變形甚至破壞。在這個時間工程中，地層應(yīng)力是連續(xù)變化的，特別地，洞室開挖后在未加支護的情況下，地層應(yīng)力所達到的新的相對平衡稱為圍巖的二次應(yīng)力狀態(tài)。一般來說，二次應(yīng)力場是三維場。在隧道施工過程中，橫向的二次應(yīng)力作用使得洞周圍巖的應(yīng)力狀態(tài)和變形狀態(tài)發(fā)生了顯著的變化，可將洞周圍巖從周邊開始逐漸向深部分為4個區(qū)域：（1）松動區(qū)由于施工擾動（例如施工爆破），區(qū)內(nèi)巖體被裂隙切割，越靠近洞室

2、周圍越嚴重，其內(nèi)聚力趨近于零，內(nèi)摩擦角也有所降低，強度明顯削弱，基本無承載能力，在重力的作用下，產(chǎn)生作用在支護上的松動壓力。（2）塑性強化區(qū)這一區(qū)域是圍巖產(chǎn)生變形的根源。隧道開挖后破壞了地層的原狀力線，在洞體四周產(chǎn)生了很高的應(yīng)力集中，此時，該處只存在切向應(yīng)力和指向隧道中心的徑向不平衡力，切向應(yīng)力由承載拱承擔(dān)，而對于徑向應(yīng)力，毛洞是無法承擔(dān)的，所以要釋放（在有支護的情況下一部分被初期支護承擔(dān)）。這就造成了洞體開挖后四周的圍巖向隧道中心發(fā)生位移，周邊的徑向應(yīng)力逐漸趨向零，而切向應(yīng)力隨著徑向位移而增大。這一應(yīng)力狀態(tài)的變化導(dǎo)致巖體從初始的二軸3）（這里只考察平面應(yīng)力狀態(tài)）受壓狀態(tài)轉(zhuǎn)變?yōu)閱屋S受壓狀態(tài)，使

3、得這一區(qū)域圍巖處于非常不利的受力狀態(tài)，當(dāng)這一應(yīng)力狀態(tài)超過巖體的強度極限時，洞室周圍出現(xiàn)了塑性區(qū)域或者破壞區(qū)域，產(chǎn)生塑性變形。如果洞室周圍塑性區(qū)域擴展不大，隨著徑向位移的出現(xiàn)，地層塑性區(qū)域達到穩(wěn)定的平衡狀態(tài)，圍巖沒有達到承載能力的極限值；但是如果塑性區(qū)域繼續(xù)擴展，則必須采取支護措施約束地層運動，才能保持洞室圍巖處于穩(wěn)定狀態(tài)，這時為了阻止地層運動，就顯出塑性變形壓力。彈性各點變形區(qū)域這一區(qū)域內(nèi)巖體在二次應(yīng)力作用下仍處于彈性變形狀態(tài)，的應(yīng)力都超過原巖的應(yīng)力，應(yīng)力解除后能恢復(fù)到原巖應(yīng)力狀態(tài)。其次，開挖面前方地層對已開挖區(qū)域的在開挖過程中不同橫斷面上的二次應(yīng)力分布和變形是不同的，這種開挖面效應(yīng)使得在開挖

4、面前方和后方一定范圍內(nèi)的二次應(yīng)力場呈三維分布狀態(tài)，這種空間效應(yīng)可以通過隧道開挖后的位移場反映出來（如圖2.1）omTi矗道開挖后的位移揚（4）原巖狀態(tài)區(qū)巖體未受影響仍處于原巖狀態(tài)。新奧法認為圍巖壓力類型可以相互轉(zhuǎn)化，即變形壓力可以轉(zhuǎn)化為松動壓力。支護結(jié)構(gòu)對圍巖起到約束作用，如果支護強度不足，圍巖位移過大進入松弛變形階段，則形成了松弛壓力。松弛壓力的實質(zhì)是部分圍巖脫離母體后作用于支護結(jié)構(gòu)上的重力。圍巖中常分布著大量的裂隙和結(jié)構(gòu)面，而這些薄弱面的剪切強度遠低于巖石塊體強度，而且不能承受拉力的作用，因此巖體力學(xué)性質(zhì)往往表現(xiàn)為各向異性，在特定方向上的強度遠低于巖石的塊體強度。當(dāng)結(jié)構(gòu)面和隧道臨空面形成不

5、利組失去支撐的圍巖塊體脫離母體（圖22），其重力作用于支護結(jié)構(gòu)上形成松弛壓力，這時巖壓力將不在取決于圍巖的特征曲線。圖2.2圍巖塊體塌落是松動壓力的主要來源2.2隧道施工過程中圍巖對支護作用力的探討隧道施工時，由于承載拱效應(yīng)，原始的地層應(yīng)力并不是全部轉(zhuǎn)化為作用在結(jié)構(gòu)上的荷載，即使在隧道建成若干年后，作用在隧道襯砌上的壓力仍然小于初始應(yīng)力。這已被大量的工程測試資料所驗證。對于這種現(xiàn)象，比較權(quán)威的解釋是1515王夢恕地下工程淺埋暗挖技術(shù)通論M安徽：安徽教育出版社，2005隧道開挖后洞室周圍地層應(yīng)力的釋放，隧道的拱形形狀以及地層內(nèi)部摩擦力等導(dǎo)致承載拱發(fā)揮作用，周圍地層應(yīng)力進行重分布產(chǎn)生兩種變化，即一

6、部分被釋放，另一部分想深部和其他方向轉(zhuǎn)移。當(dāng)施作襯砌支護后，地層應(yīng)力的釋放過程受到擬制，一部分釋放荷載作用在襯砌結(jié)構(gòu)上，這部分荷載的大小正是我們需要了解的作用于襯砌結(jié)構(gòu)上的壓力。圍巖的松動塌方與提供支護的時機有尖，如果支護愈早，提供的抗力就愈大（圖23）巖就能穩(wěn)定。反之力作用下會松動塌落，所以要維持圍巖穩(wěn)定，既要維持巖的極限平衡，還要維持松動區(qū)內(nèi)滑移體支護遲，提供的支護抗力愈小，不足以維持圍巖的穩(wěn)定，松動區(qū)中的巖體在重AdvancingtunnelRasticzone圖2.3與距掌子面位置相矢的支護壓力1616劉波，韓彥輝.FLAC原理、實例與應(yīng)用指南M北京：人民交通出版社，2005所需的支護

7、抗力小于維持圍巖極限平衡狀態(tài)所需的支護抗力，那么只需要松動區(qū)還保持在極限平衡狀態(tài)之中，松動區(qū)內(nèi)滑移體就不會松動塌落。反之，則會松動塌落。由此，我們可以把維持松動區(qū)內(nèi)滑移體平衡所需的抗力等于維持極限平衡狀態(tài)的抗力，作為圍巖出現(xiàn)松動塌落和確定Pimin條件。要確定最佳支護結(jié)構(gòu)或最佳支護時間，必須確定最小圍巖壓力Pimin。最小圍巖壓力和圍巖允許最大位移兩者是等價的（圖2.4）。目前，無論是確定最小圍巖壓力還是確定圍巖允許最大位移都沒有好的計算方法。對于一1的情況，可以用下面方法估算1414何滿潮，黃潤秋王金安等.工程地質(zhì)數(shù)值方法M北京：科學(xué)出版杜，2006:AneesMJdEddnsqulKxIc

8、imfrrwvrdradialdiBpfaaQorrwrrty*圖2.4圍巖特征曲線和支護特征曲線16佝劉波，韓彥輝.FLAC原理、實例與應(yīng)用指南M北京：人民交通出版社，2005由巖體力學(xué)可知，在的情況下，圍巖松動區(qū)內(nèi)的滑裂面為一對對數(shù)螺線。假設(shè)松動區(qū)內(nèi)強度已大大下降，可以認為滑移巖體己無自承作用，以致于松動區(qū)內(nèi)滑移體的全部重量都要由支護抗力Pimin(2.1)來承受由此有:Pimin=G考慮到實際情況，真正作用在支護結(jié)構(gòu)上的壓力應(yīng)當(dāng)是重力與變形壓力的疊加，則式(2應(yīng)該寫為：2G(2.2)式中：G為滑移體的重量；b為滑移體的地寬?；企w的重量可以近似取下式：(2.3)G二*b(Rm.x-f0)

9、2為巖體容重式中：R為與Ri相應(yīng)的允許最大的松動取半徑;RmaxRmin帶入式(2.2)得pimln(2.4)按塑性區(qū)半徑公式即可得：P+cctgsin半YlP2sinjRmax=ro；_L(Rmin+濁(2.5)J+sinJ計算R喻時，采用的C值應(yīng)該再降低Pimin的大小主要取決于松動區(qū)半徑Rmax，當(dāng)原巖應(yīng)力愈大，C,值愈低和C,值損失愈多時，則入縱和Pimin就愈大。此外，還與巖體構(gòu)造狀況，施工爆破情況、外界條件等有矢，這些都會影響圍巖C值的降低。由于隧道圍巖力學(xué)性質(zhì)的復(fù)雜性及不確定性，要純粹利用解析的方法準(zhǔn)確計算巖對支護的壓力是很困難的，但是可以通過數(shù)值模擬及實測數(shù)據(jù)預(yù)測圍巖對支護的作

10、用。來三維力學(xué)分析結(jié)果及國內(nèi)外大量的實測數(shù)據(jù)均表明，原始地層應(yīng)力的釋放率與地面沉降和拱頂下沉之比有很好的一致1414何滿潮，黃潤秋王金安等工程地質(zhì)數(shù)值方法M北京：科學(xué)出版杜，2006,即=u地面/u拱頂(2.6)而/“(2.7)式中匚v為實際作用在結(jié)構(gòu)上的壓力，匚。為原始地層應(yīng)力于是.二心：匚=(U地面/U拱頂)對于拱部：:V(u地詁u拱種(2.9)對于邊墻：jv(U地詁u拱皿)koh(2.10)對于底部：W地顯U拱頂)H(2.11)DL上式中，w為襯砌結(jié)構(gòu)的自重；H為覆蓋土厚度；D為隧道結(jié)構(gòu)跨度；i為隧道側(cè)向任意一點力地面的距離，L為閉合長度。對于作用在結(jié)構(gòu)上的圍巖壓力與覆跨比的矢系，由上面

11、各式可以看出，作用在結(jié)構(gòu)上的圍巖壓力與荷載釋放率6密切相矢，而芳即為地面沉降與拱頂下沉之比因此只要研究U與H/D的尖系，就可以弄清楚作用在結(jié)構(gòu)上的圍巖壓力與覆跨比的矢系了。日本的島田隆夫通過試驗得到以下公式:(2.12)式中，H/D為覆跨比，、:是與施工方法有矢的常數(shù)。Atkinson等人根據(jù)模型試驗得到下式；U地面/U拱頂=1.0-:H/D(2.13)將2.12式和2.13式代入2.6式和2.11式就可以算得隧道在開挖過程中的圍巖釋放率和圍巖作用在支護上的作用力。2.3數(shù)值模擬方法隨著計算機技術(shù)的迅速發(fā)展，借助數(shù)值模擬方法與計算機圖形、圖象技術(shù)、可視化技術(shù)相結(jié)合，對地下工程的開挖步驟、支護工

12、藝的工程性態(tài)(穩(wěn)定和變形)模擬和過程再現(xiàn)已成為現(xiàn)實，這種工程性態(tài)模擬和過程的圖形顯示技術(shù)被稱為計算機仿真技術(shù)。在數(shù)值模擬中，不僅僅可以模擬巖體中的斷層、節(jié)理、裂隙等地質(zhì)結(jié)構(gòu)面，而且還可以模擬分布開挖、支護等施工過程，揭示不同施工步驟、施工工藝以及不同支護條件下的應(yīng)力與位移。尤其采用不連續(xù)分析技術(shù)，模擬和顯示洞室開挖過程，洞頂及邊墻有些部位巖塊失穩(wěn)而下落或滑移，為支護設(shè)計提過可靠依據(jù)。這些都是解析方法難以實現(xiàn)的。本文利用巖土工程界最為常用的有限單元法和有限差分法，對廈門海底隧道不同施工段進行了施工力學(xué)行為模擬分析。2.3.1彈塑性有限元分析巖體彈塑性本構(gòu)矢系是巖石主要非線性問題之一。巖石的彈塑性

13、是指巖石材料的應(yīng)力應(yīng)變矢系在屈服之前呈線性矢系，當(dāng)應(yīng)力達到屈服應(yīng)力時，應(yīng)力應(yīng)變矢系就變?yōu)榉蔷€性。由于彈塑性模型中應(yīng)變不僅依賴于受載的應(yīng)力狀態(tài)，而且與加載路徑有矢，因此一般彈塑性本構(gòu)矢系不能用應(yīng)力應(yīng)變?nèi)渴赶禍?zhǔn)確描述，只能用能反映與加載路徑有矢的應(yīng)力應(yīng)變增量矢系描述。在巖石非線性本構(gòu)矢系有限元分析中，一般采用初應(yīng)力法和初應(yīng)變法求解非線性平衡方程組。初應(yīng)力法是將荷載以微小增量形式逐級加在模型上，每加一級荷載增量：dF./,就會產(chǎn)生相應(yīng)的位移增量dj?應(yīng)變增量ch和應(yīng)力增量|對于具有初應(yīng)力的彈塑性應(yīng)力應(yīng)變本構(gòu)矢系可以寫J。成：-bptd:2.14):d;。：-bpid(2.15)其中Dpi為塑T生矩

14、陣，它與加載前的應(yīng)力水平有矢，而與應(yīng)力增量無矢。初應(yīng)力法是通過對；0打的處理將應(yīng)力修正到正確的水平上，初應(yīng)力心亠？不僅與加載增量前應(yīng)力水平有矢，還與本級所加荷載增量引起的應(yīng)變增量;有矢。量形式平衡方程為KoId二加；:,dFi(2.16)式2.16中，K。為線彈性計算中總剛度矩陣；】為校正荷載項，由式2.17決定:dF1-7IB:bipId:dv(2.17)由于OF?隨位移變化而變化，所以計算時必須進行迭代求解。初應(yīng)力法求解按照以下步驟實現(xiàn)：(1) 把全部荷載劃分成若干個增量，在每一級增量段內(nèi)，按照增量彈塑性平衡方程進行求解(2) 計算各單元的應(yīng)力增量及當(dāng)前應(yīng)力仏j=bId3j丿叭上dj匕片二

15、匕L+的片(2.1下標(biāo)i表示第i級荷載增量；j表示第j次迭代。(3)根據(jù)巖石的屈服準(zhǔn)則，由各單元應(yīng)力判斷單元是否屈服，對于塑性單元，計算應(yīng)力修正項并修正應(yīng)力(2.19)(4)塑性單元通過修正項s計算等效節(jié)點力，所有塑性單元的等效結(jié)點力疊加構(gòu)成總的修正荷載矢量IdFi=XtdCTpj.dV(5)在修正荷載作用下進行下次迭代運算，此時基本方程為(2.21)Kd.j二dFi重復(fù)進行(Q(5)步計算直至所有的塑性單元達到收斂精度要求。增量計在進行下一步的荷載算重新施加下一級荷載增量fdFj，重復(fù)計算(1)(5)步，直至計算完畢。通過累加各級載荷作用的計算結(jié)果，求得總位移和總應(yīng)力匚門。一般初應(yīng)力法的收斂

16、速度比較緩慢，因此通常采用常剛度和變剛度法相結(jié)合的方法加速收斂。在ANSYS中巖土工程問題的分析一般使用Drucker-Prager屈服模型(ffl2.5)，其等效應(yīng)力表達式為：ce=3TIMIS(2.22)其中：二；X，z/3為平均應(yīng)力或靜水壓；S；為偏應(yīng)力；1為材料常數(shù)；M為mises屈服準(zhǔn)則中的M-a圖25德魯克一普拉格屈服準(zhǔn)則上面的屈服準(zhǔn)則是一種修正的Mises屈服準(zhǔn)則，他考慮了靜水應(yīng)力分量的影響，靜水應(yīng)力(側(cè)向壓力)越高，則屈服強度越大，材料常數(shù)1的表達式如下：p2sin$*3(3_sin)(223)屈服準(zhǔn)則的表達式如下：6Ccos(2.24)y3(3-sin)最后的屈服準(zhǔn)則表達式為

17、：(2.25)2-3-2有限差分法分析概述有限差分法是較早用于求解給定初值和(或)邊值微分方程組的數(shù)值方法。在有限差分法中，基本方程組和邊界條件(一般為微分方程)近似地改用差分方程(代數(shù)方程)來表示，及由空間離散點處的場變量(應(yīng)力、位移)的代數(shù)表達式代替。這些變量在單元內(nèi)是非確定的，從而把求解微分方程的問題轉(zhuǎn)化成求解代數(shù)方程的問題。而有限元是將連續(xù)求解域離散成一組有限個相互連接的單元體，利用每一個單元內(nèi)假設(shè)的近似函數(shù)來分片求解區(qū)域上待求的未知場函數(shù)。單元內(nèi)的近似函數(shù)通常由未知場函數(shù)及其導(dǎo)數(shù)在各節(jié)點的數(shù)值和插值函數(shù)來表示，這樣有限元分析中的未知函數(shù)及其導(dǎo)數(shù)在各個節(jié)點上的數(shù)值就成為新的未知量，一經(jīng)

18、求出這些未知量，就可以通過插值函數(shù)計算出各單元內(nèi)場函數(shù)的近似值，進而得到整個求解域上的近似解。有限差分法和有限元法都產(chǎn)生一組待解方程組。盡管這些方程是通過不同方式推導(dǎo)出來的，但是兩者產(chǎn)生的方程是一致的。另外，有限元程序通常要將單元矩陣組合成大型整體剛度矩陣，而有限差分則無需如此，因為它相對高效的在每個計算步重新生成有限差分方程。在有限元法中，常采用隱式、矩陣求解方法，而有限差分法則通常采用顯式、時間遞步法解代數(shù)方程。FLAC采用拉格朗日分析方法，由于它不需要形成整體剛度矩陣，大變形計算時在每個計算步都很容易修正坐標(biāo)。位移增量施加到坐標(biāo)上以致網(wǎng)格隨著材料發(fā)生移動和變形，這就是所謂的拉格朗日法；若

19、材料移動和變形是相對于固定的網(wǎng)格，就是與此對應(yīng)的歐拉法。(2)有限差分法基本方程在彈性體上用相隔等間距h并平行于坐標(biāo)軸的兩組平行線劃分成網(wǎng)格(圖2.6)設(shè)f=fX,y為彈性體內(nèi)某個連續(xù)函數(shù)它可能是某個應(yīng)力分量或位移分量、也可能就是應(yīng)力函數(shù)、溫度、滲流等。這些函數(shù)，在平行于某軸的一根格線上，坐標(biāo)的只隨變化而改變。在臨近節(jié)點0處，函數(shù)f可以展開為泰勒級數(shù):1-。匚在節(jié)點3和節(jié)點1，X分別等于Xo其代入式2.26,得f3二fo-hSilr2gfi=toh丿o21入_X3)+_Jo(X_xj+(2.26)h和X。等于0h即Xo-h分別-h和h,將(2.27)(2.28)假定h充分小，因而可以不計它的三

20、次幕及更高次幕的各項，則式可以簡化為:聯(lián)立求解(2.29)(2.30)，可得到差分公式:ef、(3刃2同樣，可以得到1f2-f4f2(2.34)(2.29)(2.30)(2.31)(2.32)(2.33)2hf2J-2bh2公式(2.31)(2.34)是基本差分公式，通過這些公式可以推導(dǎo)出其他的差分公式例如利用公式(2.31)(2.33),可以導(dǎo)出混合二階導(dǎo)數(shù)的差分公式(2.35)(cf；,cfcfIrr|r用同樣的方法，由公式(2.32)(2.33)可以導(dǎo)出四階導(dǎo)數(shù)的差分公式。(3)平面的有限差分拉格朗日法數(shù)值原理對于平面問題，將具體的計算對象劃分成四邊形網(wǎng)格域的有限差分網(wǎng)格，每個圖2.7有

21、限差分單元劃分示意圖三角形網(wǎng)格域的有限差分公式可以用高斯散度定理的廣義形式推導(dǎo)得出:(3.36)一工dAs“ifdSA其中：為繞閉合面積邊界積分，為對應(yīng)表面S的單位法向量，f為標(biāo)量、snjSf矢量或張量，人為位置矢量，ds為微量弧長，A表示對整個面積A積分在面積A上定義f的梯度平均值為：(3-37)汙AA丄dAA;X,將式(2.36)代入上式可得：對于一個三角形的網(wǎng)格域式(2.38)的有限差分形式為:/邑-一(3.39)丿一cx(/As其中，As為三角形網(wǎng)格域的邊長，求和是對該三角形的三個邊進行；邊的平.；v取均值。平面問題有限差分法是基于物體運動與平衡的基本規(guī)律。最簡單的例子是物體質(zhì)量(m)

22、、加速度(du/dt)與施加力F隨時間變化的矢系，牛頓定律描述的運動方程為：(2.40)dum=FdtF(t)圖2.8物體對隨時間變化作用力的響應(yīng)當(dāng)幾個力同時作用于該物體時，如果加速度趨于零，即aF=0(對所有作用力求和)。式(2.40)也表示該系統(tǒng)處于靜力平衡狀態(tài)。對于連續(xù)固體介質(zhì)，式(2.40河以寫成如下廣義形式:pdu_卩(2.41)dt:Xj晶為重力加速度(體式中，丁為物體的質(zhì)量密度，t為時間，薯為坐標(biāo)矢量分量，力)分量，5為應(yīng)力張量分量禾Ij用式(2.39)，將f替換成網(wǎng)格域每邊平均速度矢量，這樣，網(wǎng)格域的應(yīng)變速率ej可以用網(wǎng)格點速度的形式表述:e1邑;u，J2仏&6(243河以求岀

23、應(yīng)變張量的所有分量根據(jù)力學(xué)本構(gòu)定律，可以由應(yīng)變速率張量獲得新的應(yīng)力張量：G=M(G,eq,k)(2.44)式中，M.表示本構(gòu)定律的函數(shù)形式；k為歷史參數(shù)，取決于特殊本構(gòu)矢系；表示“由替換”。在FLAC程序中，單元應(yīng)變率是計算各主要參數(shù)的紐帶由于非線性應(yīng)力應(yīng)變本構(gòu)矢系不具有唯一性，一般用增量形式表示，當(dāng)已知單元舊的應(yīng)力張量和應(yīng)變速率(應(yīng)變增量)時，可以通過式(244)確定新的應(yīng)力張量。例如各向同性最簡單的彈性本構(gòu)矢系張量差分形式如下：2I一i+v勺一Gekk2Gej=t(2.45)GJ：CTK3式中At為時間步，G,K分別是剪切模量和體積模量。在一個時步內(nèi)，單元的有限轉(zhuǎn)動對單元應(yīng)力張量有一定的影

24、響。對于固定參照系，此轉(zhuǎn)動使應(yīng)力分量有如(2.46)其中，打=(2.47)在大變形計算中先通過式(2.46)進行應(yīng)力校正，然后利用式(2.45)或(2.44)計算等前時步的應(yīng)力。一旦計算出作用在網(wǎng)格域上的應(yīng)力后，就可以確定作用到每個網(wǎng)格點上的等價力。在每個三角形子網(wǎng)格域中的應(yīng)力如同在三角形邊上的作用力，每個作用力等價于作用在相應(yīng)邊端點上的兩個相等的力。每個角點受到兩個力的作用，分別來自各相鄰的邊。因此：斤二丄6(n(S(1)+nfS(2)由于每個四邊形網(wǎng)格域有兩組兩個三角形域，在每組中對每個角點處相遇的三角形節(jié)點力求和，然后將來在這兩組的力進行平均，得到作用在該四邊形網(wǎng)格點上的的力。在每個網(wǎng)格點處，對于所有圍繞該網(wǎng)格點四邊形的力求和汗，得到作用于該節(jié)點的純粹節(jié)點力矢量。該矢量包括所有施加的荷載作用以及重力引起的體FjAg.mg(2.49)其中mg是凝聚在網(wǎng)格點處的質(zhì)量，定義為連接該