韓偉歌 崔振東 涂新斌 金永軍 王彥兵

(①中國科學(xué)院地質(zhì)與地球物理研究所， 中國科學(xué)院頁巖氣與地質(zhì)工程重點(diǎn)實(shí)驗(yàn)室， 北京 100029， 中國)

(②石家莊鐵道大學(xué)， 省部共建交通工程結(jié)構(gòu)力學(xué)行為與系統(tǒng)安全國家重點(diǎn)實(shí)驗(yàn)室， 石家莊 050043， 中國)

(③中國科學(xué)院地球科學(xué)研究院， 北京 100029， 中國)

(④中國科學(xué)院大學(xué)， 地球與行星科學(xué)學(xué)院， 北京 100049， 中國)

(⑤國網(wǎng)經(jīng)濟(jì)技術(shù)研究院有限公司， 北京 102209， 中國)

0 引 言

隨著層狀巖體工程建設(shè)的不斷增多，由于軟弱層理面破裂而導(dǎo)致的工程施工事故也日益凸顯。層狀復(fù)合巖體由于其富含層理面極易發(fā)生彎曲、拉伸及剪切破壞，對(duì)隧道圍巖穩(wěn)定性有重大影響。因此，揭示層狀巖體隧道圍巖層理變形破裂特征對(duì)隧道安全施工及穩(wěn)定運(yùn)營具有重要意義。

為保證層狀巖體工程的安全施工運(yùn)營，國內(nèi)外學(xué)者采用眾多手段開展了大量研究。已有學(xué)者從理論上推導(dǎo)了層狀巖體本構(gòu)模型，為后續(xù)含層理物理模型試驗(yàn)及數(shù)值模擬研究奠定了基礎(chǔ)(劉卡丁等， 2002； 張玉軍等， 2002； 韓昌瑞等， 2008)。室內(nèi)隧道圍巖物理模型試驗(yàn)常采用相似材料進(jìn)行模擬，蘇士龍等(2020)采用該方法研究了層狀巖體巷道在不同支護(hù)條件下的變形破裂特征。而隧道圍巖在錨桿支護(hù)下的力學(xué)變化行為也被郝登云等(2017)在現(xiàn)場(chǎng)進(jìn)行了原位觀測(cè)研究。沙鵬等(2020)通過室內(nèi)層狀巖體試驗(yàn)對(duì)隧道層狀圍巖BQ分級(jí)方法進(jìn)行了改進(jìn)。室內(nèi)及現(xiàn)場(chǎng)模型試驗(yàn)工作量大且成本較高，因此數(shù)值模擬方法成為層狀巖體圍巖穩(wěn)定性研究的重要手段。通過數(shù)值模擬結(jié)果和力學(xué)分析相結(jié)合，可有效分析層狀巖體變形破裂特征，并可根據(jù)研究結(jié)果提供支護(hù)建議(魯海峰等， 2014； Tang et al.，2001)。王惠風(fēng)等(2020)利用離散元方法研究了不同層理厚度對(duì)巷道圍巖變形破壞的影響，得到層理厚度對(duì)圍巖破裂范圍有較大影響。此外，不同層理面強(qiáng)度對(duì)層狀巖體破裂特征的影響也被研究(韓偉歌等， 2019)。同時(shí)，不同層理面傾角也對(duì)層狀圍巖的變形有較大影響。柳厚祥等(2014)研究了大傾角層狀圍巖在不同支護(hù)條件下的穩(wěn)定性特征，而在不同層理角度下錨桿作用位置對(duì)圍巖約束支護(hù)效果也有較大差異(彭焱森等， 2013)。劉紅兵(2013)采用數(shù)值模擬手段研究了層理面角度產(chǎn)生的偏壓作用，并提出圍巖最大應(yīng)力梯度與層理面具有相同的方向。目前，層理面在層狀圍巖變形破裂過程中的角色定位還需進(jìn)一步細(xì)化分析，并且當(dāng)前對(duì)層狀巖體隧道的數(shù)值模擬研究主要采用有限元法和離散元方法，針對(duì)兩種方法各自的特點(diǎn)，常采用有限元方法進(jìn)行圍巖應(yīng)力分析，而離散元方法用于圍巖破裂研究。

因此，本文基于ABAQUS有限元軟件，通過對(duì)模型嵌入0厚度Cohesive單元方法建立天然層理面，從而在有限元中實(shí)現(xiàn)離散思想。通過該方法研究了層狀巖體隧道開挖后層理面在不同角度下的變形破裂特征，并結(jié)合Python和Matlab編程手段在有限元里實(shí)現(xiàn)了聲發(fā)射模擬，最終結(jié)合數(shù)值模擬結(jié)果及聲發(fā)射數(shù)據(jù)提出了采用三維地質(zhì)建模技術(shù)進(jìn)行層狀圍巖隧道監(jiān)測(cè)的建議方法。

1 數(shù)值計(jì)算原理

本文采用黏聚力模型(Cohesive zone model)構(gòu)建層理面，該模型可解決裂尖奇異性問題，被Dugdale(1960)提出后，在復(fù)合材料領(lǐng)域得到廣泛發(fā)展(沈珉等， 2018)。Strom et al. (2014)采用該模型研究了層理面存在時(shí)的裂紋擴(kuò)展行為，證實(shí)了模型的準(zhǔn)確性。

1.1 黏性單元本構(gòu)模型

Cohesive單元采用基于牽引分離規(guī)律的線彈性本構(gòu)模型，損傷開始前應(yīng)力-應(yīng)變滿足線彈性關(guān)系，即：

(1)

式中：t、ε分別表示應(yīng)力、應(yīng)變； 下標(biāo)n、s、t分別表示法向及兩個(gè)切向方向；K為單元?jiǎng)偠取?/p>

1.2 裂紋擴(kuò)展準(zhǔn)則

1.2.1 損傷起始準(zhǔn)則

本文采用最大主應(yīng)力準(zhǔn)則作為裂縫起裂準(zhǔn)則：

(2)

(3)

由于軟件中負(fù)應(yīng)力表示壓縮，所以單元在純壓縮狀態(tài)下不會(huì)產(chǎn)生起始損傷。

1.2.2 損傷演化準(zhǔn)則

本文以有效位移來表征損傷變量D：

(4)

1.3 聲發(fā)射模擬方法

結(jié)合Python語言和Matlab編程程序?qū)BAQUS模擬結(jié)果進(jìn)行處理，在有限元里實(shí)現(xiàn)了聲發(fā)射模擬技術(shù)。首先，提取所有損傷單元的節(jié)點(diǎn)坐標(biāo)及單元耗散能，從而確定聲發(fā)射事件定位坐標(biāo)點(diǎn)和聲發(fā)射能量； 其次，統(tǒng)計(jì)不同時(shí)間下的損傷單元數(shù)，以此作為聲發(fā)射事件數(shù)(Han et al.，2020)。最后，獲取損傷單元參數(shù)MMIXDME以此判斷裂紋破裂類型。MMIXDME參數(shù)定義如下：

(5)

式中：Gn、Gs和Gt分別表示3種破裂類型的斷裂能；MMIXDME為0時(shí)表示純拉破裂，為1時(shí)為純剪破裂，數(shù)值在兩者之間時(shí)則表示為拉剪混合破裂。

2 數(shù)值模型

以某公路隧道為研究對(duì)象，隧道埋深100m，寬度和高度分別為10.86m和7.19m，隧道斷面形狀為仰拱和三心圓(牛澤林等， 2020)。為準(zhǔn)確繪制隧道斷面，利用AutoCAD軟件建立隧道斷面圖，將DFX格式的斷面圖文件導(dǎo)入ABAQUS模擬軟件，生成隧道斷面模型。為減少邊界約束的影響，數(shù)值模型左右邊界據(jù)隧道中心線距離不小于3～5倍洞徑，隧道底部據(jù)模型底部邊界不小于3倍隧道高度。因此，本次數(shù)值模型寬度設(shè)置為60m，模型底部邊界離隧道25m，隧道頂部距地面30m，其余70m埋深通過模型頂面施加載荷來實(shí)現(xiàn)。為研究層理面在隧道開挖后的變形破裂特征，采用黏聚力方法在數(shù)值模型中嵌入不同角度的0厚度Cohesive單元作為層理面，層理面間距設(shè)置為2m，層理面角度分別設(shè)置為0°、15°、30°、45°、60°、75°和90°(圖1)。

圖1 數(shù)值模型

根據(jù)《鐵路隧道設(shè)計(jì)規(guī)范》確定圍巖彈性模量E為1.2GPa，泊松比μ為0.4，容重為22kN·m-3(朱正國等， 2019)。設(shè)置層理面抗拉強(qiáng)度為0.4MPa，抗剪強(qiáng)度為0.8MPa。采用三節(jié)點(diǎn)線性平面應(yīng)變?nèi)切螁卧獙?duì)模型進(jìn)行網(wǎng)格劃分，模型全局布種尺寸為0.7。設(shè)置邊界條件約束隧道模型兩側(cè)水平方向位移，模型底部約束水平和豎直兩個(gè)方向位移，對(duì)模型施加體力并在模型頂部施加上覆土層的重力載荷。定義Geostatic和Static， General兩個(gè)分析步，并設(shè)置初始地應(yīng)力場(chǎng)進(jìn)行地應(yīng)力平衡，其中側(cè)壓力系數(shù)設(shè)置為0.75(朱正國等， 2017)。在Static， General分析步時(shí)，采用生死單元功能，實(shí)現(xiàn)隧道的開挖模擬。

3 數(shù)值模擬結(jié)果

分別獲取隧道開挖后不同層理面角度工況下的數(shù)值模擬結(jié)果，分析層理面的變形破裂特征。

3.1 Mises應(yīng)力云圖

首先，提取隧道開挖后不同層理面角度下的Mises應(yīng)力分布云圖，并與不含層理面模型進(jìn)行對(duì)比分析，如圖2所示。

圖2 不同層理面角度下Mises應(yīng)力云圖

由應(yīng)力分布云圖可以發(fā)現(xiàn)，不含層理面時(shí)Mises應(yīng)力呈對(duì)稱光滑分布，應(yīng)力集中部位為隧道拱腳處。而層理面的存在直接決定了隧道開挖后的應(yīng)力重分布特征。層理面阻礙了重分布應(yīng)力的傳播，應(yīng)力集中受限于某些層理面內(nèi)，應(yīng)力重分布區(qū)域減小，尤其是層理面角度為0°時(shí)這一現(xiàn)象最為明顯。應(yīng)力在穿越層理面時(shí)呈現(xiàn)非連續(xù)分布，層理面效應(yīng)明顯，這也導(dǎo)致應(yīng)力集中區(qū)的轉(zhuǎn)移，應(yīng)力集中區(qū)受層理面的影響主要集中在離隧道邊界最近的層理面處。而層理面角度為0°和90°時(shí)，由于模型的對(duì)稱性，應(yīng)力云圖也呈左右對(duì)稱分布。其中層理面角度為75°和90°時(shí)，應(yīng)力分布特征與無層理面時(shí)最接近，但層理面的影響同樣不可忽略。

3.2 水平位移云圖及層理面變形特征

層理面對(duì)應(yīng)力分布的影響必然會(huì)帶來圍巖位移的變化，提取不同層理面角度下的水平位移分布云圖，同時(shí)獲取與之相對(duì)應(yīng)的層理面變形分布圖(圖3)，為了更直觀觀測(cè)層理面變形情況，本文將層理面變形放大100倍。

圖3 水平方向位移云圖及層理面變形特征

由圖3可以發(fā)現(xiàn)，層理面對(duì)水平位移云圖分布有較大影響，主要的位移變形均發(fā)生在層理面附近，位移云圖有較為明顯的順層分布特征，其將直接決定圍巖的穩(wěn)定性。盡管不同層理面角度下的模型主要受到豎向重力載荷的影響，但從層理面變形特征上可以發(fā)現(xiàn)，層理面大變形表現(xiàn)在垂直層理面角度方向(夏彬偉等， 2012)，且變形影響范圍較廣，越靠近隧道開挖面變形越明顯。層理面角度為0°時(shí)，可發(fā)現(xiàn)較為明顯的拱頂沉降和底部隆起現(xiàn)象，且沉降和隆起變形的影響范圍較大，隨著層理面角度的增大，這一影響范圍逐漸減小。

進(jìn)一步提取不同層理面角度下的最大二次分布應(yīng)力值及最大水平變形位移值，如圖4所示。

圖4 隧道圍巖最大位移及應(yīng)力

由圖4可以發(fā)現(xiàn)，當(dāng)無層理面存在時(shí)，最大Mises應(yīng)力僅為5.29×10-3MPa，而最大水平變形位移也僅為6.23×10-3mm，兩者與有層理時(shí)相差3個(gè)數(shù)量級(jí)。因此，層理面的存在對(duì)隧道圍巖的穩(wěn)定性有著不可忽視的影響。圍巖變形程度與層理面角度呈非線性關(guān)系，與其他層理面角度相比，層理面角度為75°和90°時(shí)，有較小的應(yīng)力和變形，結(jié)合上述應(yīng)力云圖分布可知這兩種角度下的應(yīng)力分布與無層理面時(shí)最接近，圍巖變形最不明顯。而層理面角度為60°時(shí)則有較大的變形，圍巖最不穩(wěn)定，其次為30°和0°。

3.3 層理面破裂形態(tài)及聲發(fā)射特征

采用上述聲發(fā)射模擬手段進(jìn)一步提取隧道開挖后的層理面破裂特征及聲發(fā)射定位圖，如圖5所示。其中聲發(fā)射定位圖中聲發(fā)射事件點(diǎn)大小表示聲發(fā)射能量的大小，顏色則表示聲發(fā)射破裂類型，不同層理面角度下層理面破裂均表現(xiàn)為剪切破裂類型，這是由于隧道開挖后引起的不協(xié)調(diào)變形導(dǎo)致層理面剪切滑移破裂。

圖5 層理面破裂特征及聲發(fā)射定位圖

由圖5a層理面破裂特征可以發(fā)現(xiàn)，越靠近隧道開挖面層理面的破裂范圍越廣，層理面從隧道開挖面起裂，逐漸順著層理面向兩側(cè)擴(kuò)展，層理面與隧道開挖面的距離越大則該裂紋擴(kuò)展距離越短，當(dāng)二次重分布應(yīng)力不足于誘使層理面破裂時(shí)，裂紋止裂，此時(shí)層理面也許還在隧道開挖影響范圍內(nèi)，但層理面強(qiáng)度已足以抵擋二次應(yīng)力的擾動(dòng)影響。此外，在遠(yuǎn)離隧道開挖面區(qū)域?qū)永砻嫫屏淹耆俺尸F(xiàn)斷續(xù)破裂特征，這一現(xiàn)象表明巖石破裂可呈現(xiàn)跳躍式斷裂，這與崔振東在微觀層面觀測(cè)到的實(shí)驗(yàn)現(xiàn)象一致(崔振東等， 2018)。除此之外，根據(jù)圖5中聲發(fā)射定位圖可以發(fā)現(xiàn)，層理面破裂聲發(fā)射能量從0°至90°大致呈“n型”分布，層理面角度為60°和45°時(shí)有較強(qiáng)聲發(fā)射能量，而層理面角度為90°時(shí)，聲發(fā)射能量較小，并且越遠(yuǎn)離隧道開挖面聲發(fā)射能量越小，聲發(fā)射能量的大小可間接表征層狀圍巖的變形破壞程度。因此，從聲發(fā)射定位圖中也可發(fā)現(xiàn)不同層理面角度下圍巖的主要破裂區(qū)域分布特征。

進(jìn)一步提取層理面在整個(gè)破裂過程中的累計(jì)聲發(fā)射事件數(shù)和累計(jì)聲發(fā)射能量，得到圖6。由圖6可知，層理面角度為30°和45°時(shí)，累計(jì)聲發(fā)射事件數(shù)最多，表明隧道圍巖破裂單元最多，破裂范圍最廣。而層理面角度為60°時(shí)，累計(jì)聲發(fā)射能量最大，此時(shí)圍巖變形破壞最明顯，但其影響范圍并不大。層理面角度為90°和75°時(shí)，累計(jì)聲發(fā)射事件數(shù)和累計(jì)聲發(fā)射能量均較小，此時(shí)的層理圍巖受隧道開挖的影響最小，隧道圍巖相對(duì)更加穩(wěn)定，其次為15°和0°。

圖6 聲發(fā)射事件數(shù)及能量

3.4 層理角度對(duì)隧洞不同部位位移影響規(guī)律

根據(jù)上述層理圍巖變形及破裂結(jié)果，可設(shè)置多個(gè)隧道斷面監(jiān)測(cè)點(diǎn)進(jìn)行位移變形監(jiān)測(cè)，監(jiān)測(cè)點(diǎn)分布如圖7所示。

圖7 監(jiān)測(cè)點(diǎn)分布

提取不同層理面角度下6個(gè)監(jiān)測(cè)點(diǎn)的豎向位移，如圖8所示。由于數(shù)值模型主要受到自上而下的重力載荷的影響，因此監(jiān)測(cè)點(diǎn)的豎向位移自上而下呈現(xiàn)越來越小的變形趨勢(shì)，而監(jiān)測(cè)點(diǎn)6處于隧道底部，受隧道底部隆起的影響有較大豎向變形。層理面角度為0°和90°時(shí)，模型水平對(duì)稱，因此，監(jiān)測(cè)點(diǎn)2和3與監(jiān)測(cè)點(diǎn)4和5分別對(duì)應(yīng)相同的豎向位移。而在其他層理面角度下由于層理面的存在導(dǎo)致的變形不協(xié)調(diào)，即使監(jiān)測(cè)點(diǎn)呈對(duì)稱分布，則兩者監(jiān)測(cè)到的變形量也不相同，層理面上部的變形量往往大于下部。隧道頂部的沉降量在層理面角度為60°時(shí)最大為0.0197m，其次分別為45°、30°、0°、15°、75°，角度90°時(shí)最小為0.0182m。隧道底部隆起變形也是60°最大為0.0154m，其次為0°、45°、30°、15°、75°，角度90°時(shí)最小為0.0147m。

圖8 不同監(jiān)測(cè)點(diǎn)豎向位移

進(jìn)一步提取不同層理面角度下不同監(jiān)測(cè)點(diǎn)的水平位移，如圖9所示。由圖9可以發(fā)現(xiàn)，不同層理面角度下隧道斷面不同位置處的水平位移有較大差異。層理面角度為0°時(shí)，最大水平位移發(fā)生在右側(cè)拱腳處； 層理面角度為15°、30°和45°時(shí)，最大水平位移發(fā)生在左拱肩處； 層理面角度為60°和75°時(shí)，最大水平位移發(fā)生在右拱肩處； 層理面角度為90°時(shí)，最大水平位移發(fā)生在兩側(cè)拱肩處。

圖9 不同監(jiān)測(cè)點(diǎn)水平位移

針對(duì)不同層理面角度下的隧道圍巖變形特征，可進(jìn)行針對(duì)性的圍巖變形監(jiān)測(cè)。光纖監(jiān)測(cè)技術(shù)已在圍巖變形監(jiān)測(cè)中得到了廣泛應(yīng)用(韓賀鳴等， 2019； 張平松等， 2019； 朱鴻鵠等， 2020)，而目前得益于三維地質(zhì)建模技術(shù)的快速發(fā)展，三維激光掃描點(diǎn)云數(shù)據(jù)建模技術(shù)由于其精度高、速度快、抗干擾性強(qiáng)等特點(diǎn)可為隧道變形監(jiān)測(cè)提供更多技術(shù)手段(馬自軍等， 2020)。點(diǎn)云掃描可快速、高效獲取隧道斷面數(shù)據(jù)，可大大降低建模成本與作業(yè)時(shí)間，將其作用于隧道圍巖監(jiān)測(cè)具有獨(dú)特優(yōu)勢(shì)。但想要獲取精確的斷面點(diǎn)云模型往往需要生成數(shù)以千萬計(jì)的點(diǎn)云數(shù)據(jù)，而如此龐大的數(shù)據(jù)量對(duì)軟硬件設(shè)備將有更高要求。為解決這一問題，可引入分級(jí)重采樣思想，結(jié)合本文不同層理面角度下的圍巖變形特征可知，不同層理角度工況下的圍巖變形區(qū)域有較大差異，若采用點(diǎn)云數(shù)據(jù)進(jìn)行多次全斷面的精細(xì)化建模無疑是增加了無效工作量，無變形區(qū)域圍巖只需粗略點(diǎn)云建模即可滿足工程需求。因此，可僅對(duì)圍巖易變形區(qū)進(jìn)行局部點(diǎn)云數(shù)據(jù)細(xì)化處理，通過實(shí)時(shí)對(duì)比易變形區(qū)域點(diǎn)云數(shù)據(jù)差異實(shí)現(xiàn)隧道圍巖的精細(xì)化變形監(jiān)測(cè)。

4 結(jié) 論

(1)層理面的存在對(duì)二次應(yīng)力重分布具有一定阻礙作用，使得應(yīng)力分布呈現(xiàn)非連續(xù)特征，應(yīng)力集中區(qū)域局限于某些層理面內(nèi)，層理面角度為0°時(shí)應(yīng)力分布區(qū)域較小，最大應(yīng)力集中在離隧道邊界最近的層理面處。

(2)層理面直接決定隧道圍巖穩(wěn)定性，圍巖最大位移變形均發(fā)生在層理面附近，層理面變形表現(xiàn)在垂直層理面角度方向，越靠近隧道開挖面變形越明顯。

(3)層理面角度為60°時(shí)圍巖層理變形最大，而層理角度為30°時(shí)隧道圍巖層理面破裂范圍最廣，且層理面破裂均為剪切破裂類型，同時(shí)越靠近隧道開挖面層理破裂聲發(fā)射能量越大。

(4)不同層理角度圍巖表現(xiàn)出不同的易變形區(qū)，針對(duì)具體隧道圍巖的變形特征可采用三維掃描點(diǎn)云數(shù)據(jù)進(jìn)行分級(jí)重采樣，實(shí)現(xiàn)區(qū)域化精細(xì)實(shí)時(shí)建模，從而實(shí)現(xiàn)高效、便捷的隧道圍巖變形監(jiān)測(cè)。