基于周邊收斂的初始地應力場反分析

傅鶴林，李亮，鄭剛強，張加兵

(中南大學 土木工程學院，湖南 長沙 410075)

?

基于周邊收斂的初始地應力場反分析

傅鶴林，李亮，鄭剛強，張加兵

(中南大學 土木工程學院，湖南 長沙 410075)

摘要:基于隧道開挖導致的周邊收斂量測值，結合FLAC3D正分析得出的一些位移基礎量，通過分解、模擬自重應力場及邊界構造應力場，建立周邊收斂位移的線性方程組，并采用最小二乘法矛盾方程，然后反推出初始地應力場。將反演分析得出的初始地應力場的結果代入FLAC3D進行數值模擬計算，正算該典型斷面在開挖完成以后，下臺階收斂測線的最終凈空收斂值，并將其與現場監控量測得到的實測數據進行對比分析。分析表明兩者基本吻合，驗證了反演分析的合理性。

關鍵詞：反分析；地應力場；數值模擬；位移反分析法

初始地應力場是影響巖土工程穩定性的重要因素之一，也是地震預測和巖體工程設計、施工的重要依據之一。采用的初始地應力場是否可靠，將直接影響到工程設計與施工的可靠性與安全性[1]。初始地應力場主要由自重應力場和構造應力場2個部分組成，自重應力場由地心引力(自重)引起，而構造應力場由地質構造運動引起。目前，原始地應力場一般采用黎曼法、水壓破裂法來測定。但由于工程地質條件復雜，上述測試方法不僅工作量、技術難度大，費用昂貴；而且測試結果的離散度也大，測點過少也難對應力場做出全面地判斷，難以推廣應用[2]。反分析是指利用反映系統力學行為的某些物理量推算該系統的巖體初始參數或初始地應力場。在巖土工程領域內，則被稱為反分析法。其中反映系統力學行為的某些現場實測物理量，被稱為反分析法的基礎信息。根據反分析基礎時所利用的基礎信息不同，反分析法可分為應力反分析法、位移反分析法和混合反分析法。位移反分析法是通過測量開挖硐室周邊或圍巖內的相對位移和絕對位移來反推巖性參數及初始地應力場的。本文以彈塑性修正后的現場隧道周邊收斂量測值為基礎，結合FLAC3D正分析出的一些位移基礎量，依據彈性力學疊加原理，建立了周邊收斂位移的線性方程組，求解初始構造應力的大小，從而反分析出初始應力場。

1基本原理與基本假設

1.1基本原題

地下工程的開挖導致圍巖產生應變，發生應變處必然伴隨著位移。一般來說，地下工程開挖產生的絕對位移如拱頂沉降，可用水準儀、全站儀等易測得；產生的相對位移可用收斂計、鉆孔伸長儀等易測得。根據測點某方向的位移等于各種荷載引起的同點位移之和建立方程組，并采用最小二乘法求解矛盾方程，反推出初始地應力場。

假設在線彈性受力狀態下，某點在自重應力和構造應力作用下引起量測方向的位移分別為dσ和dt，測點在量測方向的總位移為dc，則有：

dc=dσ+dt

(1)

(2)

(3)

設相互獨立的實測測點位移值共有n個。將每個測點的實測位移值記為DK，則根據最小二乘法原理，可將n個測點的總誤差寫為

(4)

(i,j=1, 2;i≤j)

(5)

求出σij后，利用相應的地應力計算公式可求出初始地應力分量大小。

1.2絕對位移與相對位移轉換

隧道周邊收斂量測是一種沿測線方向的相對位移的量測，而利用數值分析軟件正算得到的是監測點的絕對位移。因此，有必要對2種位移進行相互轉換。

設兩相鄰點沿量測方向的位移為Um，位移量測基線與x軸的夾角θ，則量測方向相對位移與絕對位移的關系，由圖1可知：

(6)

圖1　兩相鄰點間的位移分量關系圖Fig.1 Displacement components relationship between two adjacent points

1.3初始地應力場的彈塑性校正

根據上述方法獲得的地應力能夠滿足平衡條件，但由于采用彈性假設，彈塑性條件下，并不能保證反分析所得到的初始地應力場內各點的應力狀態都能滿足塑性屈服條件，直接應用這樣的地應力場進行分析會對結構施加多余的作用和荷載[3]。因此，有必要對反分析得到的初始地應力場進行彈塑性校正。目前，國內外對基于實測地應力線性回歸的初始地應力場進行彈塑性校正研究較多，但對基于相對位移的初始地應力場進行彈塑性校正研究較少。本文為研究的方便，將初始地應力場的彈塑性校正，主要通過對現場監測數據的彈塑性校正來體現，即將彈塑性校正后的監測數據，作為反分析的基本量。根據文[4]和文[5]指出，綜合考慮巖土體變形的“時間效應”和“空間效應”，以及塑性變形等因素，認這實測位移值僅是真正位移釋放值的60% ～ 80%。

1.4基本假設

本文的研究將基于以下幾點假設：

1)地層視為各向同性的線彈性模型；

2)區域內原始應力場不變；

3)研究問題為平面應變問題；

4)初始地應力僅由自重應力和構造應力組成，不考慮地下水的作用；

5)研究的是小范圍區域地應力場，假定地區邊界構造應力為均布分布，自重應力場為沿豎向線性分布；

6)隧道開挖前后，測線之間的夾角不變。

2計算模型及方案

2.1計算斷面概況

壁板坡隧道是滬昆客專最長隧道，全長為14 788 m。隧道進口位于貴州省盤縣紅果鎮上紙廠村，出口位于云南省富源縣后所鎮三丘田村。隧道最大埋深約為735 m，中心里程為DK985+061，起訖里程DK977+667～DK992+455。本文計算斷面靠近斷層破碎帶，圍巖級別為V級，里程為DK981+900，斷面走向為北10°～30°西，與附近斷層走向基本正交。巖性主要為厚層狀凝灰巖與玄武巖互層。附近斷層為一逆斷層，其走向為北10°～40°東。設計支護方案中采用復合式襯砌，初期支護由系統錨桿、單層鋼筋網、噴射混凝土組成，二次襯砌為模筑混凝土，其中初襯采用厚度為24 cm的C30噴射混凝土，二次襯砌采用C35混凝土，厚度為55 cm。

2.2模型的建立及計算參數

本文針對厚層狀凝灰巖與玄武巖互層體，玄武巖與凝灰巖厚度比為2∶1，玄武巖在上，凝灰巖在下。采用FLAC3D數值分析軟件建立地應力反分析平面應變模型，用于初始地應力場。根據計算需要，經過一定程度簡化后建立的計算模型如圖2。圍巖實體單元和襯砌結構單元本構模型均采用理想彈性模型。

圖2　計算模型網格劃分Fig.2 Calculation model

根據室內試驗結果并參考壁板坡隧道地質資料，確定各巖層和支護結構物理力學參數如表1。

表1巖層和支護結構物理力學參數

Table 1 Mechanical parameters of rock mass and supporting structure

巖層密度/(kg·m-3)變形模量/GPa泊松比凝灰巖23002.8790.25玄武巖25807.5120.20初期支護2400300.30

2.3模型邊界條件及實現

由于假設初始地應力僅由自重應力和構造應力組成，故可將計算域內的初始地應力場視為自重應力場和邊界構造應力場的線性疊加。通過分解、模擬自重應力場及邊界構造應力場，最后得出初始地應力場。鑒于本文將現場實測的隧道收斂視為一種平面內的變形體現，因此，可認為圍巖變形僅由重力和豎直平面內的構造應力引起。通常，自重應力場采用巖體實測密度和重力加速度模擬；而豎直平面內的構造應力場由以下3個子構造因素共同作用形成：1)X方向水平擠壓構造運動引起的應力場；2)Z方向豎直擠壓構造運動引起的應力場；3)豎直平面內剪切變形構造運動引起的應力場。

對于構造運動場的模擬比較困難，這主要是因為其邊界條件的復雜性。模型的不同邊界條件處理對反演結果具有一定的影響。構造運動場的邊界條件分為應力邊界條件和位移邊界條件。根據實際情況需要和模擬的方便，本文采用應力邊界條件模擬構造應力場的邊界條件，構造應力為均勻分布。

綜上所述，本模型的二維初始地應力場邊界模式可表示為[6]：

1)自重應力場的模擬，如圖3(a)；2)構造應力場通過在邊界上施加的單位應力p來體現，如圖3(b)，圖3(c)和圖3(d)，但反映構造運動作用力的最終值決定于p與相應的線性解系數或回歸系數的乘積。

(a)自重；(b)X方向均勻擠壓構造；(c)Z方向均勻擠壓構造；(d)XZ平面均勻剪切構造圖3　二維初始地應力場邊界模式Fig.3 Two-dimensional initial stress field boundary mode

計算中具體的實現過程為：

1)自重應力場：模型左右兩側邊界約束X方向位移，模型前后邊界約束Y方向位移，底部約束3個方向位移，頂面為自由面。模型頂部根據埋深施加相應的面力代替頂部的自重應力，本文中頂部邊界施加豎向均布荷載P=-2.74MPa。

2)構造應力場：①X方向單位構造應力作用時，模型前后邊界約束Y方向位移，底部約束三個方向位移，頂部面為自由面；②Z方向單位構造應力作用時，模型左右兩側邊界分別約束X方向位移，模型前后邊界約束Y方向位移，底部和頂面為自由面；③XZ平面內剪切應力的作用時，模型前后邊界約束Y方向位移，其余面為自由面。

2.4測線布置及現場監控量測結果

根據現場監控量測的實際情況，在DK981+900斷面上臺階開挖完成之后，進行監測點的埋設，布置了3條周邊收斂測線(見圖4)。監測測線布置完成后，及時利用數顯收斂儀進行第1次的量測。經過一段時間的持續監控量測，凈空收斂數據最后趨于穩定，3條測線的量測相對位移值見表2。根據上文1.3的內容指出，需對現場實測相對位移進行修正，本文采用70%的修正系數修正實測位移值，得修正后的相對位移值。

圖4　位移收斂量測布置圖Fig.4 Displacement convergence measurement plan

測線編號距離/m角度/(°)實測位移值修正后位移值AB77-6410.826844-4.01-5.71BC77-1107.49970-8.12-11.59CA110-647.4997136-4.56-6.52

備注：伸長為+，縮短為-。

3地應力求解

3.1各子應力場的模擬結果

利用FLAC3D分別模擬了自重、單位水平構造力應力、單位豎直構造應力和單位豎直平面內剪切應力作用下上臺階開挖的情況，圍巖變形結果如下。

3.1.1重力作用下

在重力作用下模型中各點X方向和Z方向位移見圖5。

(a)X向位移分布云圖;(b)Z向位移分布云圖單位：m圖5　自重作用下圍巖位移分布Fig.5 Displacement distribution of surrounding rock under gravity action

利用fish語句提取對應測線測點的位移見表3。

表3　測點位移統計表

3.1.2單位水平應力(1 MPa)作用下

在單位水平應力作用下模型中各點X方向和Z方向位移見圖6。

利用fish語句提取對應測線測點的位移見表4。

(a)X向位移分布云圖;(b)Z向位移分布云圖單位：m圖6　單位水平應力作用下圍巖位移分布Fig.6 Displacement distribution of surrounding rock under horizontal stress

ID號位移/mmX方向Z方向640.0010.881771.2920.159110-1.2900.158

3.1.3單位豎直應力(1 MPa)作用下

在單位豎直應力作用下模型中各點X方向和Z方向位移見圖7。

(a)X向位移分布云圖;(b)Z向位移分布云圖單位：m圖7　單位豎直應力作用下圍巖位移分布Fig.7 Displacement distribution of surrounding rock under vertical stress

利用fish語句提取對應測線測點的位移見表5。

表5　測點位移統計表

3.1.4單位剪應力(1 MPa)作用下

在單位剪應力作用下模型中各點X方向和Z方向位移見圖8。

(a)X向位移分布云圖;(b)Z向位移分布云圖單位：m圖8　單位剪應力作用下圍巖位移分布Fig.8 Displacement distribution of surrounding rock under shear stress

利用fish語句提取對應測線測點的位移見表6。

表6　測點位移統計表

根據絕對位移與相對位移之間的轉換，可得出各應力邊界條件下沿測線方向的相對位移如表7。

表7　測線方向相對位移

3.2構造地應力值的確定

根據1.1節可知，構造地應力值可通過建立如下方程組進行求解：

(7)

求解上述方程組(7)可得：

σh=5.78MPa；σv=2.25MPa；τ=0.15MPa

說明：上述結果正號表示構造應力方向與程序加載的單位應力方向相同，力的性質也相同。由于程序中加載的單位正應力為壓應力，剪應力對微元體任意一點取矩為逆時針。根據材料力學中對于應力的正負號規定為：正應力受拉為正，受壓為負；剪應力以對微元體內任意一點取矩為順時針者為正[7]。故模型邊界的實際構造應力為：

sh=-σh=-5.78MPa;

sv=-σv=-2.25MPa;τ=-0.15MPa

3.3初始地應力分量的確定

3.3.1等效側壓力系數的確定

對于各向同性巖體，僅自重作用下，巖體側壓力系數公式為：

(8)

式中，μ為巖體泊松比。

由于本文研究的是凝灰巖與玄武巖互層巖體，其泊松比并不能單一的采用凝灰巖或玄武巖的泊松比，而應采用等效泊松比代替。根據楊春和等[8]研究，平面應變問題下，僅自重作用時，巖體等效泊松比計算公式為：

μ12=αAμA+αBμB

(9)

式中：μ12表示巖體等效泊松比；αA和αB表示A和B 2種巖石的含量百分比；μA和μB表示A和B 2種巖石的泊松比；EA和EB表示2種巖石的壓縮彈性模量。

結合室內試驗結果，經計算得等效泊松比為μ12=0.22，故僅自重作用時，巖體等效側壓力系數為：λ=μ/(1-μ)。

3.3.2初始地應力分量值的確定

對于平面問題，假設Z方向為重力作用方向，X方向為水平方向。則任一點的應力狀態可表示為σ=(σx,σz,τxz)，其中σx，σz和τxz為初始地應力的各個應力分量。根據地下工程中地應力場構成的表達式和邊界應力的插值，可將監測斷面的初始地應力分量近似表示為：

(10)

3.4反演結果驗證

將上述反演出來的結果代入FLAC3D進行數值模擬計算，正算該典型斷面在開挖完成以后，下臺階收斂測線(如圖4中DE位置處)的最終凈空收斂值，并將其與現場監控量測得到的實測數據(圖中的實測數據已進行過修正)進行對比分析。分析表明兩者基本吻合(如表8)，但仍存在一定誤差。

分析其原因可能主要有以下2個方面：1)地應力反演時，構造應力的邊界條件采用應力邊界條件，且為均布荷載，這與實際初始構造地應力有一定出入。2)對實測數據的修正時，修正系數依據的是文獻[4]和文獻[5]給出的建議值，但對于不同的地下工程，由于其所處地質條件和圍巖情況的不同，對相應測量數據的修正應該視具體情況而定。

表8　下臺階測線凈空收斂位移對比

此外，根據反演出來的初始地應力分量值及材料力學的最大主應力求解公式，可知最大主應力的作用方向與水平X軸夾角約為6°，最大主應力方向位于計算斷面內，方向為北60°～80°西。由于計算斷面位于一逆斷層附近，根據1.4節第(2)條假設認為計算斷面的構造應力場與斷層區域的構造應力相同。該逆斷層的走向為北10°～40°東，根據赤平極射投影[10]和逆斷層判別主應力[11]的方法，分析出該區域構造應力場方向為東南-西北擠壓，垂直斷層走向的水平應力應是比較大的，即北50°～80°東(或南50°～80°西)方向為最大主應力方向。這與初始地應力場反演的結果基本相同。

4結論

1)以彈塑性修正后的現場周邊收斂量測數據為基礎，結合彈性力學疊加原理，將初始地應力場假設為3個均布子構造應力場和沿豎向線性分布的自重應力場的迭加場，利用FLAC3D有限差分軟件進行分析計算，得出了計算斷面位置處的構造應力分量的值；

2)將計算得出的構造應力分量的值代入到FLAC3D，正算該典型斷面在開挖完成以后，下臺階收斂測線的最終凈空收斂值，并將其與現場監控量測得到的實測數據對比，驗證了初始地應力場反演的結果的合理性和準確性，為其他類似反演分析提供一定參考。

參考文獻：

[1] 汪波, 何川, 吳德興,等. 蒼嶺特長公路隧道地應力場反演分析[J]. 巖土力學, 2012, 33(2): 628-634.

WANG Bo, HE Chuan, WU Dexing, et al. Inverse analysis of in-situ stress field of Gangling super-long highway tunnel[J]. Soil Mechanice, 2012, 33(2): 628-634.

[2] 莫海鴻. 某地下廠房初始地應力場的反演分析[J].華南理工大學學報(自然科學版), 1995, 23(1): 159-164.

MO Haihong. Inverse analysis of in-situ stress field of an underground powerhouse [J]. Journal of South China University of Technology(Social Science Edition), 1995, 23(1): 159-164.

[3] 楊強, 劉福深, 任繼承. 三維初始地應力場的多尺度彈塑性校正[J]. 水力發電學報, 2007, 26(6): 24-29.

YANG Qiang, LIU Fushen, REN Jichen. Multi-scale elasto-plastic calibration of the regression analysis of 3D initial crustal stress field[J]. Journal of Hydroelectric Engineering, 2007, 26(6): 24-29.

[4] 王建宇. 隧道開挖引起的圍巖變形過程及巖體變形模量的反分析[J]. 土木工程學報,1986,19(2):70-80.

WANG Jianyu. Deformation process in surrounding rock induced by tunnel excavation and back analysis of deformation parameters[J]. China Civil Engineering Journal, 1986, 19(2): 70-80.

[5] 劉一波, 屠厚澤, 陸文興. 現場量測位移反分析法及在地下工程中的應用[J]. 地質與勘探,1997, 33(3): 48-52.

LIU Yibo, TU Houze, LU Wenxing. Back analysis of measuring displacements on the spot and application to underground engineering[J]. Geology and Prospecting, 1997,33(3):48-52.

[6] 郭懷志, 馬啟超, 薛璽成, 等. 巖體初始應力場的分析方法[J]. 巖土工程學報, 1983, 5(3): 64-75.

GUO Huaizhi, MA Qichao, XUE Xicheng, et al. The analytical method of the initial stress field for rock masses[J]. Chinese Journal of Geotechnical Engineering, 1983, 5(3): 64-75.

[7] 劉慶潭. 材料力學[M]. 北京: 機械工業出版社, 2003.

LIU Qingtan. Mechanics of materials[M]. Beijing: China Machine Press, 2003.

[8] 楊春和, 李銀平. 互層鹽巖體的Cosserat介質擴展本構模型[J]. 巖石力學與工程學報,2005,24(23):4226-4232.

YANG Chunhe, LI Yinping. Expanded cosserat medium constitutive model for laminated salt rock [J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(23): 4226-4232.

[9] 彭文斌. FLAC3D實用教程[M]. 北京：機械工業出版社, 2007.

PENG Wenbing. FLAC3D practical tutorial [M]. Beijing: Mechanical Industry Press,2007.

[10] 鐘林. 赤平投影在工程地質勘探中的應用研究[D].長沙: 中南大學, 2009.

ZHONG Lin. Research on the stereographic projection in geology engineering drilling and its application[D]. Changsha: Central South University, 2009.

[11] 尤哲敏. 大斷裂區深埋隧道地應力特征及圍巖穩定性分析[D]. 武漢: 中國地質大學, 2013.

YOU Zhemin. Analysis of characteristics of in situ stress and surrounding rock stability in deeply buried tunnel with large faults[D].Wuhan: China University of Geosciences, 2013.

(編輯陽麗霞)

Inverse analysis on initial stress field based on peripheral displacement convergence

FU Helin, LI Liang, ZHENG Gangqiang, ZHANG Jiabing

(School of Civil Engineering, Central South University University, Changsha 410075, China)

Abstract:Based on peripheral convergence measurement value caused by tunnel excavation and basic displacement proposed by FLAC3D, system of linear equations of peripheral displacement convergence was established in this paper by disassembled and simulated self-weight stress field. By virtue of least square method, the inconsistent equation was solved, and the inital stress field was then datermined. Plugging the inverse analysis results into FLAC3D, the clearance convergence value in the bottom bench convergence survey line after the excavation of the typical sections is analyzed and the result was compared with the test data. The result show that the measured value basically agrees well with calculated value, which verifies the rationality of the inverse analysis.

Key words:inverse analysis; stress field; numerical simulation; displacement anti-analytic method

中圖分類號：TU 452

文獻標志碼：A

文章編號：1672-7029(2016)02-0301-08

通訊作者：傅鶴林(1965-)，男，江西高安人，教授，博士，從事隧道、巖土工程等研究；E-mail: 517336864@qq.com

基金項目：國家自然科學基金資助項目(51008308，50878213)

收稿日期：2015-06-16