考慮土性空間變異特征的智能地層劃分方法

摘" 要：同一土層中的巖土參數(shù)會隨著位置改變而變化，在一定距離范圍內(nèi)存在一定程度的相關(guān)關(guān)系，該相關(guān)性隨著距離的增大而衰減，直至不相關(guān)。地層劃分是巖土工程勘察設(shè)計的重要內(nèi)容，在地層劃分中應(yīng)考慮土性的空間相關(guān)性，且有必要提高地層劃分的智能化程度。該研究發(fā)展一套能“自然”地考慮土性空間變異特征的智能地層劃分方法。首先分析某場地巖土工程勘察數(shù)據(jù)的空間變異性特征，并據(jù)此生成一個隨機場實例，構(gòu)建“完備地層”模型；對構(gòu)建的地層模型進行6個孔靜力觸探模擬，基于靜力觸探結(jié)果展示基于貝葉斯壓縮感知的數(shù)據(jù)擴展和基于支持向量機、高斯模型和隱馬爾科夫隨機場的地層智能劃分方法，將劃分結(jié)果與“完備地層”模型對比，論證地層劃分方法的準確性。

關(guān)鍵詞：土性變異性；地層劃分；相關(guān)距離；隨機場理論；靜力觸探

中圖分類號：TU43" " " 文獻標志碼：A" " " " " 文章編號：2095-2945（2024）36-0139-06

Abstract： Geotechnical parameters in the same soil layer will change with the change of location. There is a certain degree of correlation within a certain distance. This correlation will attenuate with the increase of distance until it is irrelevant. Stratum division is an important part of geotechnical engineering survey and design. The spatial correlation of soil properties should be considered in stratum division， and it is necessary to improve the intelligence level of stratum division. This study developed a set of intelligent stratigraphic classification methods that can \"naturally\" consider the spatial variation characteristics of soil properties. First， the spatial variability characteristics of geotechnical engineering survey data of a certain site are analyzed， and a random field example is generated based on this to build a \"complete stratum\" model; Six hole static penetration simulations were carried out on the constructed stratum model. Based on the static penetration results， data expansion based on Bayesian compressive sensing and intelligent stratum classification based on support vector machines， Gaussian models and hidden Markov random fields were demonstrated. The classification results were compared with the \"complete stratum\" model to demonstrate the accuracy of the stratum classification method.

Keywords： soil variability; stratum division; correlation distance; random field theory; static sounding

在巖土工程活動中，準確獲取地層分布是優(yōu)化工程設(shè)計、保障工程安全、提高工程質(zhì)量的關(guān)鍵一環(huán)。在實際巖土工程中，地層劃分多依賴于鉆探取樣和靜力觸探試驗。鉆探取樣可以直觀地獲取該鉆孔的土層分布情況；靜力觸探試驗獲取的數(shù)據(jù)，雖不能直觀提供土樣分層，但可以提供在深度方向近乎連續(xù)的數(shù)據(jù)點，參數(shù)沿深度變化剖面也有利于描述地下土體分層情況，并且根據(jù)所獲取的參數(shù)還能反映不同土層的力學特性。以上2種方法都只針對單一探測鉆孔獲取沿深度方向的地層數(shù)據(jù)，而真實的地層結(jié)構(gòu)往往在水平方向同樣呈現(xiàn)出變異性。由于時間和經(jīng)費限制，不論是鉆探取樣還是靜力探測的數(shù)量都很有限且分布稀疏，很難獲取未知區(qū)域地層的完整分布情況，這是一些嚴重工程事故的原因之一。而巖土材料的空間相關(guān)性、土性分類的不確定性和土層邊界復雜性又是地層劃分中不可忽視的關(guān)鍵因素，它們導致了地層劃分剖面不可避免的不確定性[1]。有必要在充分考慮實際巖土參數(shù)的隨機性及空間相關(guān)性的基礎(chǔ)上對地層進行劃分。

Hegazy等[2]將聚類分析應(yīng)用在地層界面分層的問題中，該方法能夠利用聚類分析測試出數(shù)據(jù)點之間的細微差別，從而發(fā)現(xiàn)地層類型或者性質(zhì)的劇烈變化這一特性，客觀地劃分土壤剖面中相似的數(shù)據(jù)組。夏艷華等[3]引入地層層序序列和水平集理論來解決地層界面稀疏特性，實現(xiàn)了地質(zhì)解釋和空間分割等技術(shù)。Ching等[4]利用小波分析法可以精確地識別剖面中變化點的位置，可以生成可視化和直觀方法的圖形和圖表，并指出依據(jù)貝葉斯方法劃分土層更加嚴謹，但由于貝葉斯方法的局限性，其基礎(chǔ)理論和貝葉斯方程計算復雜，導致其很難在實際工程中應(yīng)用。Wang等[5-6]采用期望最大化算法來修正已經(jīng)建立的隱馬爾可夫模型，此方法提供了相關(guān)突出特性的合理估計，能夠更加精準地識別土層分布情況，更加精準地確定不同土層之間的邊界劃分情況。曹子君等[7]提出了一種基于靜力觸探試驗數(shù)據(jù)自動劃分土層的方法，此方法客觀反映了不同土層的統(tǒng)計特性，能夠合理地反映識別出的土層邊界的可靠性，不僅能劃分土層剖面，而且定量地表示了土層剖面的不確定性。許程[8]基于概率論統(tǒng)計分析同一土層力學參數(shù)之間的線性相關(guān)關(guān)系，并基于隨機場理論計算不同土層的相關(guān)距離，根據(jù)各類資料總結(jié)鎮(zhèn)江地區(qū)不同地層年代所形成的土層分布規(guī)律。Molina-Gómez等[9]提出了一種基于聚類分析的多元統(tǒng)計方法，能夠客觀地對土壤剖面中相似的數(shù)據(jù)進行分類。Suryasentana等[10]研究了一種概率方法，通過使用貝葉斯檢測方法處理動力觸探數(shù)據(jù)，檢測不同土層之間的土壤數(shù)據(jù)趨勢的突變來圈定不同土層，實現(xiàn)了更可靠的土壤地層識別方案。

上述研究雖提出了眾多地層劃分方法，但其普適性驗證不充分，且劃分結(jié)果無“標準答案”可對照，對地層本身存在的空間變異性特征也未充分考慮。為此，本研究首先基于某場地巖土工程勘察數(shù)據(jù)分析了各層土的土性空間變異性量化指標，并據(jù)此生成該場地的一個隨機場實例，作為“完備地層”模型。對該“完備地層”，模擬其靜力觸探過程，基于靜力觸探結(jié)果進行地層智能劃分方法探索，將劃分結(jié)果與“完備地層”模型對比，探討了地層劃分的準確性。

1" 考慮空間變異性的“完備地層”模型構(gòu)建

1.1" 土性參數(shù)空間變異性分析

對某場地進行的巖土工程勘查，共實施了14個鉆孔取樣和6個靜力觸探，識別到六層土（自上而下分別為填土1、黏土1、粉質(zhì)黏土1、粉土粉砂、粉質(zhì)黏土2和黏土2），其彈性模量、黏聚力、內(nèi)摩擦角均值和標準差見表1。

土性相關(guān)距離是土性剖面隨機場模型的重要度量。在地質(zhì)變化與土體沉積過程中，土層因地質(zhì)、環(huán)境、物理、化學等因素形成了土體天然的變異性。對于土體內(nèi)的任意兩點，在一定距離內(nèi)存在一定的相關(guān)性，大于此距離則可認為這兩點之間不存在相關(guān)性，這個距離就稱之為相關(guān)距離。因測試上述力學指標的試樣的取樣深度已知，此處采用相關(guān)函數(shù)法、基于線指數(shù)型相關(guān)函數(shù)擬合得到6種土體的垂直相關(guān)距離，各個參數(shù)的相關(guān)距離非常接近（后文取其均值），擬合結(jié)果見表1。由于各孔間距較大，無法直接獲得水平方向的相關(guān)距離，在不失一般性的情況下，依據(jù)前人工作[11]，取10倍的垂直相關(guān)距離作為水平相關(guān)距離。

1.2" “完備地層”模型構(gòu)建

隨機變量的不可數(shù)無窮集構(gòu)成了隨機場，可用一系列離散的隨機數(shù)近似描述。隨機場處理離散問題的方法有：局部平均法、中心點法、Karhunen-Loeve（K-L）展開法等。其中，K-L展開法由于優(yōu)越的解析性而在隨機場分析中廣泛應(yīng)用，其實質(zhì)是將隨機場分解為一系列不相關(guān)的隨機變量和確定系數(shù)（特征函數(shù)、特征值）；與其他隨機場離散方法相比，其對于任意類型的隨機場均收斂，展開階數(shù)相同時具有最小的均方誤差。既有文獻表明，K-L展開對規(guī)則的矩形區(qū)域有較好的收斂性能，對不規(guī)則區(qū)域的收斂性不佳。參考前人經(jīng)驗，對典型地質(zhì)剖面內(nèi)的各層土體，先單獨在整個剖面計算其K-L展開結(jié)果，再根據(jù)指定的地層分界面對各矩形區(qū)域進行剖切合并，由此構(gòu)成二維的地層剖面。二維剖面數(shù)據(jù)點維度為50×110，即由5 500個像素點組成矩陣，結(jié)果如圖1所示。圖中橫向黑虛線表示地層分界線，縱向黑虛線表示后續(xù)模擬靜力觸探孔位置。

圖1是基于實際地層有限的地勘數(shù)據(jù)、指定的地層界面，利用隨機場理論“擴展”出的具有完備信息（包括彈性模量、黏聚力和內(nèi)摩擦角）的“完備地層”模型。這里的“完備”是指地層模型具有后續(xù)數(shù)值模擬所需的所有幾何與力學參數(shù)。其理論價值在于地層界面位置完全已知，力學參數(shù)完全已知；可基于該地層模型模擬地勘過程并獲取地勘參數(shù)（如靜力觸探曲線），然后基于地勘數(shù)據(jù)進行地層劃分，將該劃分結(jié)果與完全已知的“正確”地層界面位置對比，實現(xiàn)從理論上探討基于靜力觸探曲線的智能地層劃分方法的可靠性。

2" 靜力觸探過程的數(shù)值模擬

本節(jié)基于“完備地層”模型，提取各個靜力觸探孔附近的地層剖面數(shù)據(jù)，鏈接到ABAQUS軟件中，將隨機場理論與有限元分析相結(jié)合，實現(xiàn)靜力觸探的有限元模擬。模擬結(jié)果將用作后續(xù)地層智能劃分的依據(jù)。

本節(jié)共模擬6個孔（孔7—孔12）的靜力觸探，孔位如圖1所示。此處采用耦合歐拉-拉格朗日法模擬，將土體視作歐拉體，網(wǎng)格保持不變，材料可在網(wǎng)格內(nèi)發(fā)生任意流動；將觸探桿、錐設(shè)為拉格朗日體，材料變形與網(wǎng)格變形保持一致。該方法的優(yōu)點是網(wǎng)格不發(fā)生扭曲，并能準確地計算材料邊界，有效解決了有關(guān)大變形、材料破壞和流體材料等諸多問題。

土體模擬范圍是以探頭貫入路徑為中心、邊長為0.5 m和深度為25 m的長方體，模型底部為固定端，四周為法向約束，頂端為自由端。探頭采用變形體模擬，錐角為60°、錐底直徑為3.57 cm、探桿長度30 m，貫入速度為20 mm/s。貫入前，探頭錐尖位于土體頂面處。土體采用摩爾-庫倫本構(gòu)模型描述參數(shù)按照圖1賦值。采用ABAQUS內(nèi)置的通用接觸模擬土體與觸探錐和桿之間的相互作用，接觸算法為歐拉-拉格朗日接觸。切向接觸摩擦系數(shù)為0.2，接觸面法向行為定義為“硬”接觸。另對探頭變形體施加剛體約束，限制探頭本身在貫入的變形。

探頭、探桿采用八結(jié)點線性六面體的單元（C3D8R），減縮積分，沙漏控制。土體采用八節(jié)點線性歐拉六面體單元（EC3D8R）劃分網(wǎng)格。為減少計算時間和提高計算精度，探桿貫入范圍附近采用較密的網(wǎng)格劃分方式，外周采用相對稀疏的劃分方式。網(wǎng)格劃分如圖2所示。試算表明，圖中的網(wǎng)格劃分密度能得出合理的模擬結(jié)果。模型采用顯示動力學求解器進行求解計算，時間步長0.000 1 s，并開啟大變形選項。

從模擬結(jié)果中提取比貫入阻力沿深度的變化，并與地勘報告靜力觸探曲線對比，如圖3所示。由圖3可知，基于隨機場理論考慮彈性模量、黏聚力和內(nèi)摩擦角3種巖土參數(shù)的空間變異性，靜力觸探貫入模型能夠較好地模擬錐尖與土體之間的相互作用，比貫入阻力隨深度的變化趨勢與原位該處靜力貫入試驗結(jié)果基本一致，表明此處的數(shù)值模擬方法合理、可靠。

3" 地層智能劃分方法

本研究提出的智能地層劃分方法包括2個主要步驟：①基于貝葉斯壓縮感知（Bayesian Compressive Sensing，BCS）方法，對有限的地勘數(shù)據(jù)進行空間擴充，使得場地剖面各處均有擴充的勘察數(shù)據(jù)；②基于數(shù)據(jù)聚類識別各土類的判別邊界，據(jù)此對地層剖面各點判別其土類，形成地層劃分結(jié)果。以下對上述步驟逐一展開介紹。

3.1" 基于貝葉斯壓縮感知的地層數(shù)據(jù)擴充

BCS是一種廣泛應(yīng)用于雷達信號成像、人臉識別、無線通信、機器學習等領(lǐng)域的信號處理技術(shù)，其最大的優(yōu)點是以采集到的稀疏信號進行插值和展開，實現(xiàn)對原始數(shù)據(jù)的重構(gòu)。在巖土工程勘察中，靜力觸探數(shù)據(jù)可認為是對地層的一次稀疏采樣，本研究采用Fang等[12]發(fā)展的BCS算法，基于少量觸探數(shù)據(jù)（此處為第2節(jié)模擬的靜力觸探結(jié)果）完成對地層剖面數(shù)據(jù)擴充。

勘察剖面長度為110 m、深度為25，靜力觸探孔共6個，從左到右各孔水平間隔分別為23、21、24、23、19 m。垂直分辨率為0.5 m，水平分辨率為1.0 m，則待擴充的場地數(shù)據(jù)維度為110×50，而已知稀疏數(shù)據(jù)為6個鉆孔模擬得到的比貫入阻力，維度為6×50，經(jīng)BCS算法擴充后的比貫入阻力分布如圖4所示。由圖4可知，擴充的結(jié)果水平向呈現(xiàn)帶狀分布特征，比貫入阻力總體沿深度方向的變化趨勢和對應(yīng)的模擬數(shù)據(jù)變化趨勢基本一致。

3.2" 智能地層劃分方法

綜合鉆孔取樣目測、物理力學特性測試等結(jié)果，可將土類與靜力觸探比貫入阻力進行關(guān)聯(lián)，如圖5所示，該“深度-比貫入阻力”二維平面內(nèi)以點的顏色區(qū)分不同土類。為更加準確地在圖5中劃分不同土類的分界面，實現(xiàn)數(shù)據(jù)聚類，首先采用基于線性核函數(shù)的支持向量機（Support Vector Machine，SVM）算法，獲得圖5中各土類的線性邊界。SVM能尋找出邊界線使各土類數(shù)據(jù)點盡可能地分布在邊界的兩側(cè)，結(jié)果如圖5所示（各聚類分區(qū)以色塊顏色區(qū)分）。在上述擴充的地層靜力觸探比貫入阻力（圖4）基礎(chǔ)上，利用SVM聚類結(jié)果可判斷勘察斷面上各點所屬的土類，結(jié)果如圖6所示。可以看出，本方法推測的二維剖面的土層分布基本符合原始的土層分布（黑虛線為原始邊界）。第2層土（黏土1）和第3層土（粉質(zhì)黏土1）相對于其他土層的估計準確度略低。是由于該土層相對其他土層面積偏小，對插值結(jié)果最為敏感。且在進行數(shù)據(jù)統(tǒng)計時，該土層數(shù)量也相對較少，準確度也受到統(tǒng)計不確定的影響。綜上所述，總體地層的分布較為清晰，分類較為準確，但土層突變處劃分粗糙，與原始地層分布還存在一定的誤差。經(jīng)計算整個研究場地的二維剖面的5 500個像素點分類準確的個數(shù)為5 203個，準確率為94.6%。

在SVM聚類基礎(chǔ)上，進一步采用高斯模型（Gaussian Model， GM）算法對上述聚類進行概率運算，其實質(zhì)是給出圖5中各數(shù)據(jù)點屬于某類土的概率。樣本數(shù)據(jù)x（深度，比貫入阻力）是多維數(shù)據(jù)，高斯分布用下列概率密度函數(shù)描述

式中：u為期望；？撞為協(xié)方差；D為數(shù)據(jù)維度。高斯模型參數(shù)？茲=（？滋，？撞）由期望和協(xié)方差組成。對于高斯模型中參數(shù)？茲的值可以通過極大似然函數(shù)估計

假設(shè)每個數(shù)據(jù)點是獨立的，N表示樣本點的總數(shù)。似然函數(shù)L（？茲）可以由概率密度函數(shù)給出

最終結(jié)果如圖5中各等概率橢圓所示：在每個聚類分區(qū)中，距離同族橢圓中心越近則概率越大。將上述擴充的地層中各點靜力觸探數(shù)據(jù)帶入即得各點屬于各土類的概率，取概率最大值對應(yīng)的土類為改點所屬土類，結(jié)果如圖7所示。可以看到，第4層土的上、下界面的識別準確性比圖6有明顯改善，土層劃分準確性提高到95.6%。

SVM+GM算法對于地層分界處的細節(jié)劃分仍有所欠缺，這是由于各個點之間僅根據(jù)圖5的聚類進行土類確定，相鄰點之間的聯(lián)系并未考慮。為了考慮相鄰點在土類屬性上的潛在聯(lián)系以更加符合實際情況，進一步采用隱馬爾可夫隨機場（Hidden Markov Random Field，HMRF）對分類模型進行進一步細化，具體算法詳見Zhao等[13]，細化結(jié)果如圖8所示。結(jié)果表明，經(jīng)過HMRF優(yōu)化后，識別準確率提高為97.35%，對土層與土層之間的邊界處像素點劃分更加精細準確。第1層土（填土1）和第6層土（黏土2）與原始地層分布基本一致，第2層土（黏土1）上邊界處誤判為填土1的點更少。第4層土（粉土粉砂）的上下分界處有少量位置不準確。

綜上，SVM+GM+HMRF算法具有最好的地層劃分精度。其中，SVM主要實現(xiàn)地層初步劃分，而GM、HMRF則逐步提高劃分的精細程度。BCS數(shù)據(jù)擴充和地層劃分方法中涉及的參數(shù)都是基于地勘數(shù)據(jù)自動識別獲取，無需人為干涉，因而可實現(xiàn)整個地層劃分的智能化。此外，土性參數(shù)的空間變異性特征“自然地”體現(xiàn)在BCS數(shù)據(jù)擴充和HMRF中。

4" 結(jié)論

本文發(fā)展了基于靜力觸探結(jié)果對地層進行智能劃分的方法，主要結(jié)論如下。

1）基于線指數(shù)型相關(guān)函數(shù)可擬合得到土性參數(shù)空間變異性的相關(guān)距離，運用K-L展開法能得到合理的隨機場實例，考慮參數(shù)空間變異性的數(shù)值模擬能夠較好地模擬錐尖與地基土體之間的相互作用，模擬得到的比貫入阻力分布與原位靜力貫入試驗的結(jié)果基本一致。

2）基于貝葉斯壓縮感知理論對模擬數(shù)據(jù)進行擴充，可可靠地獲取整個研究剖面的比貫入阻力分布，其數(shù)值在地勘剖面上呈現(xiàn)帶狀分布，與原模擬數(shù)據(jù)趨勢吻合，數(shù)據(jù)擴充方法可行。

3）采用支持向量機（SVM）+高斯模型（GM）+隱馬爾可夫隨機場（HMRF）算法相結(jié)合對擴充后的地層數(shù)據(jù)進行地層劃分，模型準確率可達97.35%，各地層位置整體準確性高，部分層位附近少量點位劃分準確性稍低。驗證了本文地層劃分方法的可行性。該方法只需輸入擴充的地層數(shù)據(jù)，其余參數(shù)可自動從輸入數(shù)據(jù)中學習得到，從而支持智能化地層劃分的實現(xiàn)。

參考文獻：

[1] 鄭碩.基于靜力觸探的土層自動劃分方法及不確定性表征[D].武漢：武漢大學，2017.

[2] HEGAZY Y A， MAYNE P W. Objective site characterization using clustering of piezocone data[J]. Journal of Geotechnical and Geoenvironmental Engineering， 2002，128（12）：986-996.

[3] 夏艷華，白世偉.復雜三維地層模型地層劃分研究[J].安徽理工大學學報（自然科學版），2011，31（3）：35-39.

[4] CHING J， WANG J S， JUANG C H， et al. Cone penetration test （CPT）-based stratigraphic profiling using the wavelet transform modulus maxima method[J]. Canadian Geotechnical Journal， 2015，52（12）：1993-2007.

[5] WANG H， WELLMANN J F， LI Z， et al. A segmentation approach for stochastic geological modeling using hidden Markov random fields[J]. Mathematical Geosciences， 2017，49（2）：145-177.

[6] WANG X R， WANG H， LIANG R Y， et al. A hidden Markov random field model based approach for probabilistic site characterization using multiple cone penetration test data[J]. Structural Safety， 2018，70：128-138.

[7] 曹子君，鄭碩，李典慶，等.基于靜力觸探的土層自動劃分方法與不確定性表征[J].巖土工程學報，2018，40（2）：336-345.

[8] 許程.基于不確定性理論的鎮(zhèn)江第四紀土層層序劃分研究[D].南京：東南大學，2021.

[9] MOLINA-G？魷MEZ F， VIANA D F A， FERREIRA C， et al. Defining the soil stratigraphy from seismic piezocone data： A clustering approach[J]. Engineering Geology， 2021，287：106-111.

[10] SURYASENTANA S K， LAWLER M， SHEIL B B， et al. Probabilistic soil strata delineation using DPT data and Bayesian changepoint detection[J]. Journal of Geotechnical and Geoenvironmental Engineering， 2023，149（4）：6-23.

[11] 吳振君.土體參數(shù)空間變異性模擬和土坡可靠度分析方法應(yīng)用研究[D].武漢：中國科學院武漢巖土力學研究所，2009.

[12] FANG Y， WU J J， HUANG B M. 2D sparse signal recovery via 2D orthogonal matching pursuit[J]. Science China Information Sciences， 2012，55（4）：889-897.

[13] ZHAO T Y， WANG Y. Non-parametric simulation of non-stationary non-gaussian 3d random field samples directly from sparse measurements using signal decomposition and markov chain monte carlo[J]. Reliability Engineering and System Safety，2020，203（1）：107-108.