網格精度與對稱性對深部巖體分區破裂化的影響

摘 要：【目的】深部巖體工程開挖卸荷后，圍巖普遍存在分區破裂化現象，容易導致隧洞大變形、大范圍塌方、涌水等工程事故，利用數值模擬軟件可以較好地重現該現象，但數值模擬網格質量會對模擬結果產生較大影響，為保證模擬結果的相對可靠，有必要對其進行研究分析。【方法】基于深部巖體應變軟化模型，建立深埋多隧洞數值模型，研究網格對稱性與網格精度對圍巖分區破裂形態、破裂帶數目、破裂范圍的影響。【結果】結果表明，網格的對稱性和精度會對圍巖分區破裂化模擬結果產生較大影響，非對稱網格及粗糙網格會使分區破裂數值模擬結果產生較大誤差。【結論】網格劃分過程中要注意網格的對稱性以及單元尺寸的選取，本研究對類似研究具有一定的參考價值。

關鍵詞：深部巖體；應變軟化；圍巖分區破裂化；網格對稱性；網格精度

中圖分類號：TU45" " "文獻標志碼：A" " 文章編號：1003-5168（2024）20-0061-05

DOI：10.19968/j.cnki.hnkj.1003-5168.2024.20.012

Effects of Grid Accuracy and Symmetry on Zonal Fracture of Deep Rock Mass

Abstract： [Purposes] After excavation and unloading of deep rock mass engineering， there is a common phenomenon of zonal fracture in the surrounding rock， which can easily lead to engineering accidents such as large deformation of the tunnel， large-scale collapse， and water inflow. The use of numerical simulation software can better reproduce this phenomenon， but numerical simulation mesh quality will have a greater impact on simulation results. In order to ensure the relative reliability of the simulation results， it is necessary to study and analyze the numerical simulation mesh quality. [Methods] A numerical model of deeply buried multi-tunnels was established based on the strain softening model to study the effects of grid symmetry and grid accuracy on the zonal fracture morphology， number and range of surrounding rock fracture zones. [Findings] The results show that the symmetry and accuracy of the grid will have a great impact on the simulation results of zonal fracture of the surrounding rock. Asymmetric grids and rough grids will cause large errors in the numerical simulation results of zonal fracture. [Conclusions] Therefore， attention should be paid to the symmetry of the mesh and the selection of the unit size during the meshing process. This study has certain reference value for subsequent scholars to study this phenomenon through numerical simulation.

Keywords： deep rock mass; strain softening; zonal fracture of surrounding rock; grid symmetry; grid accuracy

0 引言

近年來，地表層資源逐漸枯竭，因此深部地下工程的規劃和開發受到學界與相關領域的廣泛關注。但由于深部巖體工程獨特的工程地質環境致使圍巖的變形和破壞呈現明顯的非線性特征，分區破裂化現象［1-2］就是其中典型變形特征之一。該現象是指在深部巖體開挖硐室或隧道時，其兩側和工作面前的圍巖中產生膨脹帶和壓縮帶，破裂區和未破裂區交替出現的結構性破壞現象，容易引起圍巖大變形和失穩，造成人員生命財產的重大損失。

數值模擬方法在模擬復雜地質環境下的深部巖體圍巖分區破裂化現象具有良好的適用性。唐春安等［3］利用RFPA對巷道圍巖中出現的環狀間隔破裂現象進行模擬分析，結果表明深埋隧洞圍巖在高應力狀態下發生分區破裂化的概率更大。李樹忱等［4］結合數值模擬軟件研究分區破裂化的衍生過程，并計算巷道圍巖發生破裂的狀況，結果顯示破裂區的數量及范圍信息與現場勘探情況具有高度一致性。王紅英等［5］采用FLAC模擬軟件深入探討分區破裂化現象，得出時間和空間是影響分區破裂化的重要因素。陳旭光等［6］采取輪廓法研究深部巖石裂紋的發展軌跡，將其結果導入有限元軟件中進行模擬計算，能夠較好地重現圍巖分區破裂化現象。

在模擬深部巖體開挖圍巖分區破裂化的過程中，計算模型的網格精度與對稱性對模擬結果有較大影響，本文結合深部巖體峰后應變軟化模型，借助FLAC3D有限差分軟件，以錦屏二級水電站引水隧洞為背景建立數值計算模型，研究數值模擬中網格的精度與對稱性對圍巖分區破裂化的影響，為同類研究提供了參考。

1 數值計算模型

1.1 應變軟化模型

深部地下工程中，在高地應力條件下，巖體在開挖卸荷過程中表現出顯著的軟巖大塑性變形和應變軟化特性［7-8］，因此實際巖體更符合應變軟化模型，所以擬采用應變軟化模型來描述其本構關系，如圖1所示。

應變軟化模型基于摩爾-庫侖本構模型提出，其特點是與拉力流動法則相關聯，而不與剪切流動法則關聯。該模型的屈服函數、勢函數、塑性流動法則以及應力修正與摩爾-庫侖模型完全一致，主要的區別在于應變軟化模型中，黏聚力、內摩擦角、剪脹角等在塑性屈服開始后會發生變化，而在摩爾-庫侖模型中，這些參數都保持為常數。

式中：[ep1]、[ep2]、[ep3]分別對應第一、第二、第三主應變上的塑性應變量。

基于摩爾-庫倫應變軟化模型，FLAC3D采用塑性參數eps計算塑性剪切應變，并且在軟化過程中控制強度參數的變化，其增量形式定義為式（2）。

1.2 計算模型與參數

1.2.1 計算模型與邊界條件。考慮到邊界效應，計算模型左右邊界及上下邊界與隧洞間距均大于5倍洞徑，隧洞編號從左到右依次為1#、2#、3#、4#，采用六面體單元劃分，共剖分為58 638個單元，117 186個節點，如圖2所示。近隧洞單元大小尺寸為d/96（d為隧洞直徑），隧洞周邊單元網格均對稱，隧洞斷面形式為馬蹄形，如圖3所示。本研究對模型四周位移施加法向約束，對其底部位移全約束；模型施加的應力邊界條件依據現場監測的地應力反演數據［9］，考慮自重應力場的作用，分別在[x]、[y]、[z]方向施加49.81 、51.68 、58.09 MPa的構造應力，模擬深部地下高應力環境。

1.2.2 計算參數。根據眾多學者在高應力下巖體應變軟化的試驗結果可知［10-13］，巖體應變軟化后的內摩擦角會衰減至其原有強度的60%～85%，黏聚力下降至其原有強度的0%～70%，因此本研究采用的巖體內摩擦角和剪脹角均降至70%，黏聚力降至40%，采用相對收斂標準進行計算，模型初始力學參數見表1。

1.3 計算步驟

本研究旨在分析數值模擬中網格精度與對稱性在高應力狀態下對圍巖分區破裂形態、破裂帶數目以及破裂范圍的影響，因此分以下兩種情況進行計算分析。

①情況一。建立非對稱多隧洞網格模型（見圖4）與對稱多隧洞網格模型，對計算結果進行對比分析。

②情況二。取近隧洞單元尺寸大小分別為d/64、d/96、d/128的3種不同多隧洞對稱網格模型進行計算對比分析。

2 模擬結果及分析

2.1 網格對稱性對圍巖分區破裂化的影響

對稱和非對稱網格下圍巖分區破裂形態如圖5和圖6所示。由圖5與圖6可知，當模型接近平衡狀態，應力完全釋放時，對稱網格與非對稱網格對計算結果有較大影響。從圍巖分區破裂形態來看，在對稱網格計算下的每個隧洞圍巖分區破裂化形態上下左右呈一定的鏡像關系，分區破裂化程度相差不大，而在非對稱網格計算下的1#與4#隧洞左右兩幫圍巖分區破裂化表現出明顯的非鏡像關系，2#與3#隧洞圍巖分區破裂化呈一定的鏡像關系，且1#與4#隧洞的圍巖分區破裂化程度明顯弱于2#與3#隧洞。從圍巖分區破裂帶數目來看，在對稱網格下，1#隧洞到4#隧洞的破裂帶數目相差不大，都在30條左右，而在非對稱網格下，1#與4#隧洞的分區破裂帶數目在20條左右，而2#與3#隧洞分區破裂帶數目在33條左右，相差較大。從圍巖分區破裂范圍來看，經測量，在對稱網格計算下的1#～4#隧洞破裂帶最大徑向距離（最外緣分區破裂帶到隧洞中心的距離）為15.4～15.6 m，而在非對稱網格計算下，1#與4#隧洞破裂帶最大徑向距離為14.6 m左右，2#與3#隧洞為15.6 m左右。

2.2 網格精度對圍巖分區破裂化的影響

數值模擬是把實際模型劃分為多個微小單元計算，當單元接近無限小時便接近實際模型，但同時單元越小，計算時間也會越長，因此需要選取合適的單元尺寸來模擬實際工程。從圖5與圖6可以看出，當進行全斷面爆破開挖時，四個硐室的破裂化形態沿中軸線呈鏡像關系，因此以下研究只截取左側兩個硐室的結果進行分析。選取近隧洞單元尺寸大小分別為d/64、d/96、d/128的對稱網格進行計算，結果如圖7所示。從圍巖分區破裂形態與破裂帶條數來看，隨著網格精度的增加，圍巖分區破裂形態表現出明顯的不同，破裂帶數目呈先增大后減小的趨勢，網格越粗糙，其破裂帶數目越少，圍巖分區破裂化現象越不明顯，誤差越大。從圍巖分區破裂范圍來看，d/64網格模型的破裂帶最大徑向距離均在14 m左右，而d/96與d/128網格模型的破裂帶最大徑向距離均在15.6 m左右，相差幅度不大，選用d/96的單元尺寸已滿足實際計算精度需求。

3 結論

本研究基于深部巖體應變軟化模型，以錦屏二級引水隧洞工程為研究背景，分析網格精度與對稱性對深埋多隧洞圍巖分區破裂化數值模擬的結果影響，主要得出以下結論。

①在非對稱網格處，其圍巖分區破裂形態也不具對稱性，與對稱網格區域內圍巖分區破裂形態呈鏡像對稱關系明顯不同。在非對稱網格處的分區破裂帶數目與破裂范圍對比對稱網格區域顯著減小。

②網格越粗糙，其圍巖破裂形態越不明顯，破裂帶數目與破裂范圍相對網格精度較高的數值計算模型相差較大。同時當網格精度達到一定程度時已滿足計算需求，過高的網格精度會浪費大量的計算時間。

因此在構建數值計算模型時，應重視網格的對稱性以及單元劃分時選擇合適的尺寸，以確保模擬結果能夠準確反映實際巖體破壞行為，避免模擬結果存在較大誤差以及耗費過多的計算時間。

參考文獻：

［1］王明洋，陳昊祥，李杰，等.深部巷道分區破裂化計算理論與實測對比研究［J］.巖石力學與工程學報， 2018， 37（10）： 2209-2218.

［2］王猛，牛譽賀，于永江，等.主應力演化影響下的深部巷道圍巖變形破壞特征試驗研究［J］.巖土工程學報， 2016， 38（2）： 237-244.

［3］唐春安，張永彬.巖體間隔破裂機制及演化規律初探［J］.巖石力學與工程學報， 2008， 27（7）： 1362-1369.

［4］李樹忱，錢七虎，張敦福，等.深埋隧道開挖過程動態及破裂形態分析［J］.巖石力學與工程學報，2009，28（10）：2104-2112.

［5］王紅英，張強，張玉軍，等.深部巷道圍巖分區破裂化數值模擬［J］.煤炭學報，2010，35（4）：535-540.

［6］陳旭光，張強勇，李術才，等.基于擴展有限元的深部巖體分區破裂化現象初步數值模擬［J］.巖土力學，2013， 34（11）： 3291-3298.

［7］王水林，鄭宏，劉泉聲，等. 應變軟化巖體分析原理及其應用［J］. 巖土力學，2014，35（3）：609-622，630.

［8］BROWN E T， BRAY J W， SANTARELLI F J. Influence of stress-dependent elastic moduli on stresses and strains around axisymmetric boreholes［J］.Rock Mechanics and Rock Engineering， 1989， 22： 189-203.

［9］孫飛躍，劉希亮，郭佳奇，等.巖爆預測評估方法的動力數值分析［J］.應用力學學報，2022，39（1）：26-34.

［10］高富強，康紅普，林健.深部巷道圍巖分區破裂化數值模擬［J］.煤炭學報，2010，35（1）：21-25.

［11］王學濱，潘一山，伍小林.不同強度巖石中開挖圓形巷道的局部化過程模擬［J］.防災減災工程學報，2010，30（2）：123-129.

［12］蘇永華，鄭璇.深部節理巖體分區破裂化機制數值研究［J］.公路交通科技，2013，30（1）：94-101.

［13］許磊，郭帥， DAVIDE E，等.深部不同斷面巷道分區破裂形態與圍巖結構控制［J］.巖土工程學報，2023，45（4）：720-729.