基于三維網格模型重構的巖溶地質構造穩定性數值分析

王啟明 車愛蘭 周晗旭

摘要 由于巖溶地質的離散性大，傳統建模方法并不能真實反映其空間分布狀態，數值分析結果的可靠性受到限制。為精確分析無錫市某超高層項目基礎落位區巖溶地質構造穩定性，本文基于三維網格模型重構技術，建立該區域巖溶地質構造的三維有限元分析模型，分別分析了溶洞處治前及注漿加固后的穩定性，結果表明：溶洞未處治時，溶洞周邊巖體應力較大且存在應力集中，不均勻沉降現象明顯。對溶洞進行注漿加固處治后，溶洞周邊巖體應力集中現象減弱，模型頂面最大沉降量下降了54.55%，不均勻沉降趨勢大大減小。通過上述分析，不僅對該超高層項目基礎的設計和溶洞的處治具有重要指導意義，也為類似工程的數值分析提供參考。

關 鍵 詞 巖溶地質;三維網格模型重構;注漿加固;數值分析

中圖分類號 TU47? ? ?文獻標志碼 A

Abstract Because of the large discreteness of unfavorable geology， the traditional modeling method cannot truly reflect its spatial distribution， the reliability of the numerical analysis is limited. In order to accurately analyze the stability of karst geological structure in the foundation area of a super-high-rise project in Wuxi， a 3D finite analysis model of karst geological structure is established based on 3D mesh model reconstruction technology， and the stability of karst cave is analyzed before and after grouting reinforcement. The results show that when the cave is not treated， the rock mass around the cave has large stress concentration and obvious non-uniform settlement. After the grouting reinforcement， the stress concentration of rock mass around karst cave is weakened， the maximum settlement on top of the model is reduced by 54.55%， and the non-uniform settlement trend is greatly reduced. The above analysis not only has an important guiding significance for the foundation design of the super-high-rise project and the treatment of karst caves， but also provides a reference for the numerical analysis of similar projects.

Key words karst caves; 3D mesh model reconstruction; grouting reinforcement; numerical analysis

0 引言

我國幅員遼闊、地質復雜，巖溶地區分布廣泛，大量在建大型工程位于巖溶發育地區，在巖溶地區進行基礎工程設計與施工時，下伏溶洞對基礎工程穩定性有至關重要的影響[1]。溶洞發育程度、溶洞跨度、頂板厚度、圍巖特性等與溶洞上方的地基承載力密切相關[2]。地基滑動和坍塌、承載力不足、不均勻沉降等是溶洞存在常常引起的地基變形問題，影響基礎的穩定性[3]。因此，施工前必須對巖溶地質的情況做好仔細的研究和認真的準備，充分考慮其對施工造成的各種不利影響以及應當采取的主要技術措施，以確保施工過程順利、安全。否則，將會埋下嚴重的安全隱患，極有可能造成工期延誤、重大經濟損失和人員傷亡[4-5]。

對于巖溶地質構造穩定性進行準確的分析及評價，為工程問題的及時預測和適時決策提供有力支持，已經成為巖土工程界熱點研究課題，大量學者從室內試驗、數值分析、理論研究等多方面對巖溶地質展開研究，取得了豐富的成果。張永杰等[6]以湖南某高速公路線巖溶地質發育段為基礎，分別從定量和定性評價指標對巖溶區公路路基穩定性進行二級模糊綜合評判，使公路路基穩定性評價更具合理性。張聰、陽軍生[7]等對新研制的環保可控膏狀注漿材料（TGM）進行室內試驗，測試了TGM材料的基本力學性能，并應用于工程實際，取得了良好的效果。近年來，數值分析方法在地下巖溶研究中得到了越來越廣泛的應用，潘東東等[8]通過建立巖溶隧道承壓隱伏溶洞突水模型進行數值分析，研究了地下溶洞突水對隧道開挖的影響，趙明階[9]針對層間空隙巖溶建立相應數值模擬模型，分析了圍巖在開挖過程中的應力、位移和涌水情況。上述數值分析方法根據單一巖溶形態建立數值模型，研究特定巖溶地質狀態的數學和力學特性，取得了良好的效果。而實際地下巖溶發育極不均勻且復雜多變，因此根據某種特定巖溶狀態建立模型進行數值分析，并不能如實反映地下巖溶地質的空間分布狀態，計算結果的可靠性受到限制，故難以在工程中應用。

因此，為精確分析無錫市某超高層項目基礎落位區的巖溶地質構造穩定性，對該區域進行彈性波CT探測，利用三維網格模型重構技術建立該區域巖溶地質構造三維計算模型，通過數值模擬分析巖溶地質構造的應力和沉降值，對該項目基礎的設計與溶洞處治有重要指導意義。

1 彈性波CT探測現場試驗

1.1 彈性波CT工作原理

地球物理勘探最先應用于油氣資源開發，隨著物理學、數學、電子技術、計算機技術的發展，在不良地質探測領域取得明顯的研究進展，比較成熟的運用于不良地質探測方面的方法有高密度電阻率法、探地雷達法、彈性波CT法、電磁波CT法等幾種，其中彈性波CT法作為一種領先的地球物理勘探技術，具有探測效率高、探測精度高、抗干擾能力強等優勢，在巖溶地質探測中取得良好的應用效果。

利用彈性波對被探測對象一定厚度的層面進行掃描，在某一鉆孔中的激發點激發時，分別在其他鉆孔不同接收點進行接收，掃描方式有水平同步掃描、斜同步掃描、定點發射扇形接收掃描。圖1為彈性波CT射線網格示意圖。

1.2 工程概況

無錫市擬建某超高層建筑，初步勘察資料顯示，該地區地層巖性中上部為第四系全新統、更新統沉積物，主要成分為粘性土及少量粉土，下部基巖為三疊系下青龍（T1x）灰巖，灰-肉紅色。現場鉆孔取芯，芯樣如圖2所示，在基巖面以下，標高-78.00～-88.00 m范圍內，發育有空溶洞和半填充溶洞。根據《巖土工程勘察規范》對該項目場地進行巖溶地質彈性波CT探測試驗。其中，在基礎落位區域共布置鉆孔25個（ZK1～ZK25），深度約90 m，鉆孔孔徑不小于70 mm，孔口至孔底均有PVC管以保證鉆孔暢通。鉆孔布置如圖3所示。

1.3 數據采集

于2018年1月25日至2月10日，對基礎落位區域進行彈性波CT數據采集。彈性波CT采集設備主要有數據記錄系統、井下水聲檢波器線陣、大功率電火花震源（如圖4），各部分工作原理如下。

數據記錄系統：由高精度寬頻帶數據采集模塊（Goed數字地震儀）和筆記本電腦構成，負責對換能器線陣接受到的信號進行模數轉換和數字記錄。

井下水聲檢波器線陣：由12組換能器（傳感器）按0.5 m間隔排列構成12道檢波器線狀陣列用于在井孔中接受彈性波。水聲檢波器由寬頻帶高靈敏度壓電陶瓷元器件和放大電路構成，把作用在器件上的壓力波轉換成電信號。

大功率電火花震源：由電極高壓放電形成電弧，電弧瞬間將電極周圍的水體加熱成蒸汽，蒸汽膨脹并爆炸在地層中激發出彈性波。

具體采集參數設定如下：接收道數為12道，道間距為0.5 m，激發點間距為0.5 m，每組激發點數為12個，完成一組激發后，12道接收道上移3 m。采樣間隔20.833 μs，持續時間0.15 s。在25個鉆孔之間，采集72組彈性波CT剖面數據。采集到的波形例圖5所示。

1.3 數據處理及結果

由于數據采集現場含有大量的施工噪音、雜波、多次波，而且地層介質對彈性波也有不同程度的吸收，所以需要對數據采取一定的處理，最大可能消除干擾波，使波速構造圖像能最有效地展現地下特征。具體處理工作如下。

1）數據預處理。現場記錄數據中并不包含檢波器和激發點（震源）的位置信息。這些信息一般記錄在操作記錄中，為了便于數據的一體化管理，需要把這些數據添加到數據文件的文件名中，使數據和采集信息合為一個整體。

2）波形處理。首先進行濾波、視窗和疊加處理抑制噪音，提高信噪比。除此之外，利用數字濾波技術再一次抑制噪音，提高數據質量，以至于首波更為明顯。

3）首波走時提取。運用自行編制軟件自動讀取每個波形的首波運行時間，并自動進行對應于位置信息的記錄。

4）波速計算與巖溶分布評價。獲取首波走時數據后，采用射線追蹤算法、最佳路徑算法等方法，通過計算機處理，建立井間介質的二維圖像，依據特定插值方法將一系列上述二維圖像轉化為三維空間分布圖像，將地下溶洞的分析結果繪制在三維地層空間內。

以剖面ZK6-ZK19為例，對現場測試結果進行地質解釋。測試范圍標高-77～-88 m，剖面寬度9 m。剖面內測點間平均縱波速度[Vpa=2 760]m/s，最小波速[Vmin=16 266]m/s，最大波速[Vmax=3 399]m/s。剖面巖層中存在低速區，巖溶較為發育，多為黏土或碎石充填。圖6為剖面CT測試分析成果圖。

2 巖溶地質三維網格模型重構

2.1 建立三維可視化模型

根據彈性波CT探測結果，對圖2虛線標識范圍內6個鉆孔間探測區域進行三維網格模型重構，首先按照厘定的勘測線依次對彈性波CT切片進行三維空間展布，如圖7所示。

在彈性波CT法探測中，CT切片相交將探測區域劃分成若干個五面體，若直接進行Kriging插值，每個五面體內離散數據點的插值對象不僅包括五面體的3個側面，還包括其他彈性波CT切片，模型的精確度大大降低。在插值過程中，以彈性波CT切片及其交線為邊界進行區塊化Kriging插值，獲得三維可視化模型。

彈性波CT切片在三維規范化空間中的區域范圍是8.00 m×8.00 m×11.00 m的六面體，根據目標區域大小及所要重構模型的網格單元規模，將三軸向上離散數據點的間隔值設置為0.25 m×0.25 m×0.25 m，則生成37×73×45的規則離散數據點陣。對離散數據點陣進行區塊化Kriging插值，所有離散數據點統一進行三維可視化成像，結果如圖8a）所示，圖8b）為沿z軸方向切片。

2.2 三維網格單元構建

理論上，三維網格構建就是以單類型數據作為建模數據來源構建三維有限元網格模型，并建立不同材料屬性與網格間的映射關系。三維網格構建的步驟如下：1）材料閾值分割;2）網格劃分及材料屬性映射。

2.2.1 材料閾值分割

材料閾值分割的過程：1）確定閾值Te;2）三維可視化模型[F（i，j，k）]中的離散數據點，其值對應不同的值域段;3）區分不同材料的分布區域[10]。可用如下公式表示，即

所以，對于材料閾值分割，最關鍵的一步是選擇合適的閾值T。峰值法[11]是一種全局閾值分割方法，假設不同材料對應不同的值域段，用多個正態分布函數表示，相鄰兩函數的曲線擬合會得到2個分離的峰值，那么兩峰之間的谷底對應的值就是分割閾值T。引入直方圖h（F），該直方圖通過統計不同值域段離散點的數目來確定，即

2.2.2 網格劃分及材料屬性映射

建立與數據點陣形成映射關系的三維網格模型，三維網格模型由8節點的六面體單元構成，每個數據點對應1個單元，數據點位于單元中心，單元沿三軸向的尺寸與數據點三軸向間距相同。對網格模型的節點編號，如圖9所示，模型第n個節點的坐標為

通過以上步驟，完成了三維網格模型重構，所得三維網格模型如圖10所示，網格步長設為0.25 m×0.25 m×0.25 m，共劃分121 545個節點，114 048個單元。該三維網格構建方法可以控制三維網格模型的單元規模，在保證分析精度的前提下提高計算效率。

3 數值分析

3.1 模型材料參數及荷載

采用商用有限元分析軟件Abaqus的有限元靜力分析模塊進行承載力數值模擬計算，研究巖溶地質構造穩定性。根據工程勘察報告，基巖為灰巖，且發育空溶洞和半填充溶洞，在數值分析中，將溶洞簡化為飽和黏土材料，各材料具體參數如表1所示，表中：ρ為材料重度，E為彈性模量，μ為泊松比。

在數值分析中對圖9所示模型的的位移邊界進行設置：約束模型前、后、左、右4面x軸和y軸方向的位移，約束模型底面x軸、y軸和z軸方向的位移，模型頂面為自由邊界。計算時，首先在模型中施加重力荷載，待初始應力平衡后，在模型頂面施加均布荷載，分別計算均布荷載為50 MPa、100 MPa、150 MPa、200 MPa這4種工況。

3.2 穩定性分析

圖11為均布荷載200 MPa工況下模型的應力等值云圖，從圖中可知，在模型頂面均布荷載作用下，在溶洞周邊巖體應力較大且存在應力集中，有塑性破壞風險，從而使溶洞頂板坍塌進而向上延伸，引起地基承載能力不足，影響上部基礎穩定性。

圖12為均布荷載200 MPa工況下模型的沉降等值云圖，從圖中可知模型頂面產生不均勻沉降，最大沉降量達110.73 mm，最小沉降量64.39 mm。

3.3 注漿加固效果評價

從以上分析可以得知，上部基礎穩定性受下伏溶洞影響很大，如不進行有效處治，存在安全隱患。綜合考慮現場施工條件和經濟因素，建議采用注漿填充對溶洞進行加固，提高基礎穩定性。根據現有的注漿材料的室內試驗報告，注漿填充材料的參數如表2所示。

通過數值模擬分析，計算求得注漿加固后巖溶地質構造的力學特性。圖13為注漿加固后均布荷載200 MPa工況下模型的應力等值云圖，對比圖10可知，溶洞周邊的應力集中明顯減小。

圖14為注漿加固后均布荷載200 MPa工況下模型的沉降等值云圖，模型頂面最大沉降量由110.73 mm減小至50.33 mm，下降了54.55%，最小沉降量由64.39 mm減小至43.57 mm，不均勻沉降現象明顯減弱。

沿模型頂面y軸方向中軸線，提取點A（4.5，4.5，-77）、B（4.5，9，-77）、C（4.5，13.5，-77）注漿前后的應力和沉降值，如圖15所示。

在A、B、C 3點，注漿后不同荷載作用下應力和沉降量都有明顯降低。分析其原因，當溶洞填充后，溶洞內部填充的注漿結石體對溶洞頂板有一定的支撐作用，不僅能有效減少不均勻沉降，而且還能使巖溶地質構造內應力分布更加均勻，減少應力集中，并能有效防止溶洞頂板破壞及進一步向上延伸破壞。

4 結論

1）本文提出了基于三維網格模型重構的巖溶地質構造穩定性數值分析方法，該方法根據彈性波CT探測結果，重構能精準反映巖溶地質分布狀態的三維網格模型，提高了數值分析結果的可靠性。

2）該方法中，對插值對象進行分區塊Kriging插值，避免了非相關切片的影響，提高建模精度的同時，實現巖溶地質彈性波速度構造三維可視化。

3）本文基于巖溶地質三維可視化模型數據結構的特點，設計三維網格構建方法，可控制三維網格模型的單元規模，在保證分析精度的前提下提高計算效率。

4）開展無錫市某超高層項目基礎落位區彈性波CT探測現場試驗，并重構該區域巖溶地質構造三維網格模型，通過數值分析獲取了該區域不同荷載作用下的應力和沉降狀況，對該項目基礎的設計施工有重要指導意義。同時對溶洞注漿加固處置后的應力和沉降狀況進行分析，結果表明：溶洞周邊巖體應力集中現象減弱，溶洞頂板穩定性提高，模型頂面最大沉降量下降了54.55%，不均勻沉降趨勢大大減小。

參考文獻：

[1]? ? 付敬，丁秀麗，張練. 巖溶回填改善地下洞室群圍巖穩定性的數值分析[J]. 長江科學院院報，2006，23（4）：47-50.

[2]? ? 唐前松. 巖溶區公路路基穩定性及處治效果三維數值分析[J]. 公路工程，2017，42（1）：143-146.

[3]? ? 蘇謙，黃俊杰，白皓，等. 鐵路站場巖溶地基加固補強技術與效果分析[J]. 巖土力學，2013，34（3）：776-782.

[4]? ? 李堯，李術才，劉斌，等. 鉆孔雷達探測地下不良地質體的正演模擬及其復信號分析[J]. 巖土力學，2017，38（1）：300-308.

[5]? ? 付代光，周黎明，肖國強，等. TSP預報隧道不良地質體應用研究[J]. 地球物理學進展，2016，31（1）：417-426.

[6]? ? 張永杰，曹文貴，趙明華，等. 巖溶區公路路基穩定性的區間模糊評判分析方法[J]. 巖土工程學報，2011，33（1）：38-44.

[7]? ? 張聰，陽軍生，謝亦朋，等. 富水巖溶越江盾構隧道注漿材料試驗與應用研究[J]. 巖石力學與工程學報，2018，37（9）：2120-2130.

[8]? ? 潘東東，李術才，許振浩，等. 巖溶隧道承壓隱伏溶洞突水模型試驗與數值分析[J]. 巖土工程學報，2018，40（5）：828-836.

[9]? ? 趙明階，金鵬. 層間孔隙巖溶影響隧洞穩定性的數值模擬研究[J]. 地下空間與工程學報，2017，13（3）：827-832.

[10]? 鐘江城，周宏偉，任偉光，等. 基于CT圖像灰度分布的含雜質煤體三值化方法[J]. 力學與實踐，2018，40（2）：140-147.

[11]? 彭啟偉，羅旺，馮敏，等. 改進二維Otsu法和果蠅算法結合的圖像分割方法[J]. 計算機應用，2017，37（S2）：193-197.