鄒 歐，薛 健，李 娜，王 鑫，楊 榮，李玉瓊*

（1.中國科學(xué)院大學(xué) 工程科學(xué)學(xué)院，北京 100049；2.中國科學(xué)院 力學(xué)研究所 流固耦合系統(tǒng)力學(xué)重點實驗室，北京 100190；3.中國科學(xué)院 力學(xué)研究所 非線性力學(xué)國家重點實驗室，北京 100190）

1 引言

土壤強度參數(shù)解譯在定量評估軍事裝備機動性和指導(dǎo)裝備輕量化結(jié)構(gòu)設(shè)計等方面具有重要作用，在地質(zhì)勘測、裝備全域機動、地外行星科考等重大需求領(lǐng)域中有廣泛的應(yīng)用背景［1-5］。傳統(tǒng)“取樣-室內(nèi)實驗”是典型的非原位測試方法，取樣造成的土體擾動——特別是原位應(yīng)力的解除——使得所采集試樣的力學(xué)性質(zhì)不能真實反映土壤原本的力學(xué)特征，從而影響對土壤承載性能及裝備車輛通行性的可靠評估［6］。因此，基于現(xiàn)場原位貫入觸探測試并結(jié)合相關(guān)理論實現(xiàn)土壤原位力學(xué)參數(shù)的精準、快速反演的需求日益緊迫。目前，土壤原位貫入測試方法主要有靜力觸探（Cone Penetration Test，CPT）和自由落體貫入觸探。靜力觸探能獲得土壤強度參數(shù)隨深度的變化情況，一般需要外載荷持續(xù)作用以保持勻速貫入（2.5 cm/s），設(shè)備安裝耗時長，操作費用高［7-8］，難以滿足復(fù)雜地質(zhì)環(huán)境和大范圍區(qū)域內(nèi)土壤力學(xué)參數(shù)的測量，限制了其在裝備機動性評估方面的應(yīng)用。自由落體貫入觸探技術(shù)是在靜力觸探基礎(chǔ)上發(fā)展起來的一種新型原位測試技術(shù)，它依靠自由落體貫入儀（Free Fall Penetrometer，F(xiàn)FP）以一定初始速度貫入土壤，通過傳感器實時記錄貫入過程中貫入儀的加速度和受到的阻力，來反演土壤強度參數(shù)。FFP 結(jié)構(gòu)輕便、操作簡單且較為經(jīng)濟，能夠?qū)崿F(xiàn)復(fù)雜環(huán)境和較大區(qū)域內(nèi)土壤強度參數(shù)的快速測量［7，9］。

目前，在海洋工程中，已有相關(guān)研究利用FFP 獲取土壤強度參數(shù)。Aubeny 等人［10］提供了基于數(shù)值模擬建立的不排水抗剪強度計算方法，給出了承載力系數(shù)的經(jīng)驗公式；Stoll 等人［11］通過對貫入實驗所得數(shù)據(jù)進行最小二乘擬合，得到了貫入阻力和最大減速度之間的關(guān)系。研究表明貫入過程中的阻力可表示為一個與貫入速度有關(guān)的阻力函數(shù)［12-14］。通過確定土壤強度參數(shù)與阻力函數(shù)之間的關(guān)系，結(jié)合貫入儀上傳感器記錄的加速度數(shù)據(jù)可反演土壤強度參數(shù)。True［15］給出了一個描述FFP 貫入海床土過程的半經(jīng)驗阻力函數(shù)，其中考慮了端部阻力、側(cè)壁摩擦阻力、土壤拖曳阻力以及土體浮力的作用。Chow 等人［9］利用離心機實驗中傳感器測量的加速度數(shù)據(jù)，反演了飽和軟黏土的不排水抗剪強度，并考慮了率效應(yīng)對軟黏土強度的影響。劉君等人［16］采用數(shù)值仿真的方法模擬了貫入儀在軟黏土中的運動，研究其所受阻力與貫入速度、土壤強度和界面摩擦系數(shù)之間的關(guān)系，提出了一種反演軟黏土不排水抗剪強度的方法。以上反演方法均基于True 提出的半經(jīng)驗阻力函數(shù)，其中包括承載力系數(shù)，界面摩擦系數(shù)，拖曳阻力系數(shù)等許多經(jīng)驗參數(shù)，當缺少場地土壤條件前期資料時，經(jīng)驗參數(shù)取值會對反演方法的準確性產(chǎn)生重要影響。此外，以上方法僅能針對海洋工程中常見飽和軟黏土的不排水抗剪強度進行反演。這表明，由于實際地質(zhì)環(huán)境的復(fù)雜性，基于自由落體貫入技術(shù)的沖擊貫入機理還有待更深入的研究，參數(shù)反演方法還需完善，其適用性還需拓展。

在貫入阻力的理論研究方面，現(xiàn)有方法主要包括承載力理論［17］，應(yīng)變路徑理論［18］和空腔膨脹理論［19-26］。承載力理論僅考慮土體破壞區(qū)域的整體平衡，而破壞區(qū)域內(nèi)各點并不處于平衡狀態(tài)。Janbu 等人［27］通過黏性土的靜力觸探實驗表明，利用承載力理論求得的貫入阻力相比于實測值低約40%。應(yīng)變路徑法只能近似求解飽和黏性土的情況［18］，對于砂土，由于貫入時的流場難以構(gòu)造，該方法并不適用。

空腔膨脹理論是研究侵徹問題的重要模型。由于FFP 的剛度遠大于土壤，在貫入過程中變形較小，因此FFP 貫入土壤的過程可視作剛性彈侵徹靶體的過程。若將貫入儀看作由一系列與其內(nèi)切的球形空腔組成，則貫入儀與土體之間的相互作用就可看作這些球形空腔在無限大介質(zhì)中均勻膨脹時與土體介質(zhì)間的相互擠壓［28］。通過求解介質(zhì)中空腔膨脹的徑向應(yīng)力分布，就可得到彈體貫入過程中受到的貫入阻力。Bishop［19］首先給出了準靜態(tài)柱/球形空腔膨脹理論的控制方程，并利用它求解了尖錐體對金屬靶的壓入阻力問題。Chadwick［20］將這個理論應(yīng)用到理想土壤材料，結(jié)合Mohr-Coulomb 屈服準則，得到了土壤徑向位移的表達式。Cao 等人［21］假定土的本構(gòu)為修正劍橋模型，通過求解準靜態(tài)空腔膨脹理論控制方程，得到了土體內(nèi)的超孔隙水壓力和總應(yīng)力。Hill［22］則首先給出動態(tài)空腔膨脹理論控制方程，得到了不可壓縮理想彈塑性靶體對彈體的阻力。基于動態(tài)空腔膨脹理論，F(xiàn)orrestal 等［23-25］開展了針對于不同靶體材料（如巖石、混凝土及土壤等）的侵徹研究，大大拓寬了這一理論的應(yīng)用范圍。Shi 等人［26］針對干砂的沖擊貫入，建立了考慮壓縮性的動態(tài)空腔膨脹理論。動態(tài)空腔膨脹模型準確地刻畫了穩(wěn)定貫入過程中貫入阻力和貫入速度平方間的線性關(guān)系，展現(xiàn)了清晰的物理圖像，而土壤強度參數(shù)則反映在阻力系數(shù)當中，這為土壤強度參數(shù)反演提供了理論基礎(chǔ)。

本文基于沖擊貫入觸探技術(shù)和動態(tài)空腔膨脹模型，開展了土壤Mohr-Coulomb 參數(shù)（內(nèi)聚力和內(nèi)摩擦角）反演方法的研究。利用土壤沖擊貫入實驗中測得的阻力-速度曲線、加速度時程曲線、速度時程曲線以及貫入深度-初速度關(guān)系，對貫入阻力系數(shù)進行了擬合，進而對Mohr-Coulomb 參數(shù)進行了反演。最后，對參數(shù)可解域和參數(shù)敏感性進行了討論。

2 動態(tài)球形空腔膨脹模型

一般來講，在對土壤的沖擊貫入過程中，貫入儀可視作剛體。此時，動態(tài)球形空腔膨脹模型為計算彈-靶相互作用提供了一個簡單的物理圖像。這個模型將貫入儀看作是由一系列與其內(nèi)切的球形空腔組成，空腔球心在貫入儀的對稱軸上［28］。貫入儀與土體的相互作用可看作是這些球形空腔在無限大介質(zhì)中均勻膨脹時與土體介質(zhì)間的相互擠壓。這樣一來，貫入儀側(cè)表面上每一點的法向應(yīng)力可近似為球形空腔膨脹時內(nèi)表面的徑向應(yīng)力，它與空腔膨脹的速度有關(guān)，而空腔膨脹速度為貫入儀側(cè)表面上該點的法向速度分量［25，28］。由于球?qū)ΨQ性的引入，經(jīng)過以上簡化后，貫入儀側(cè)表面上的壓力分布就可轉(zhuǎn)化為一個一維問題進行求解。貫入儀圖示及球形空腔膨脹模型如圖1。

球形空腔的膨脹將導(dǎo)致土體介質(zhì)出現(xiàn)相應(yīng)的動態(tài)響應(yīng)，形成彈塑性分區(qū)［25］（如圖1（b））。記空腔膨脹速度為V、彈塑性邊界擴張速度為c、初始時刻空腔的半徑為零，則在t時刻塑性區(qū)的范圍為（Vt，ct），彈性區(qū)的范圍為（ct，+∞）。在彈性區(qū)內(nèi)，土壤的本構(gòu)關(guān)系采用Hooke 定律來描述，并且假定彈性區(qū)內(nèi)土壤不可壓縮。在塑性區(qū)內(nèi)，采用體積鎖定的假設(shè)［25］來描述土壤的應(yīng)力應(yīng)變關(guān)系，即發(fā)生塑性變形后土壤的體積應(yīng)變?yōu)槌?shù)η*，不隨平均應(yīng)力的增大而變化。在這個假定下，靜水壓力將不影響土壤的后繼屈服。記ρ0和ρ*分別是土壤的初始密度和鎖定密度，則：

圖1 貫入儀及球形空腔膨脹模型Fig.1 Penetrometer and spherical cavity expansion model

土壤的屈服準則采用Mohr-Coulomb 準則，在球坐標系下可寫為［25］：

其中：σr，σθ分別為徑向應(yīng)力和環(huán)向應(yīng)力（以壓為正）；p=（σr+2σθ）/3 為平均應(yīng)力；τ0，λ是與土壤內(nèi)聚力C和內(nèi)摩擦角φ有關(guān)的常數(shù)，且有［29］：

沖擊貫入過程中，質(zhì)量、動量及能量守恒定律均需得到滿足。對于貫入阻力和貫入深度的解答，僅利用質(zhì)量守恒和動量守恒關(guān)系即可得到。在球坐標系下，拉格朗日坐標表示的質(zhì)量守恒方程和動量守恒方程可寫為：

其中：r是拉格朗日坐標，u為徑向位移（以向外為正）?？涨槐砻嫣幍倪吔鐥l件為：

彈塑性界面處的Hugoniot 條件可表示為［24］：

其中，v是介質(zhì)的徑向速度，下標1、2 分別代表彈塑性界面處的彈性區(qū)和塑性區(qū)。

利用以上本構(gòu)關(guān)系和邊界條件求解控制方程（4）和（5）可解析地得到貫入儀端部表面的法向應(yīng)力分布為［25］：

其中：Vz是貫入儀在豎直方向的速度；A，B是與τ0/E，λ和η*有關(guān)的參數(shù)。根據(jù)λ的不同取值，A，B的表達式存在如下三種不同的情況［25］：

式（12）中的第2 式表明，在貫入過程中空腔膨脹速度V與彈塑性邊界擴張速度c的比值始終保持為常數(shù)。隨著貫入深度的增加，貫入速度減小，空腔膨脹速度V也逐漸減小，對應(yīng)的彈塑性區(qū)域的擴張逐漸變慢。

對于貫入儀相對土壤垂直貫入的情況，將式（9）在貫入儀端部積分，則得到貫入儀所受到的豎直方向阻力Fz為［25］

其中：

αs和βs是與土壤力學(xué)特性有關(guān)的系數(shù)。其中ψ=s/（2r）是貫入儀頭部的曲徑比，并且有θ0=arcsin［（2ψ-1）/（2ψ）］（見圖1（a））；μ為貫入儀與土壤間的滑動摩擦系數(shù)，二者之間的切應(yīng)力可表示為στ=μσn。如果忽略貫入儀與土體間的摩擦（μ=0），并令：

式（18）和式（19）可進一步簡化為：

由式（9）和式（17）可以看到，在當前模型假設(shè)下，貫入儀端部的應(yīng)力分布及貫入阻力由兩部分貢獻，一項是靜阻力項，另一項是流動阻力項（即速度平方項），而黏性效應(yīng)和附加質(zhì)量效應(yīng)［14］未被考慮。

貫入儀在豎直方向的貫入速度Vz可由牛頓第二定律計算得到，注意到式（17），有：

其中，m為貫入儀的質(zhì)量。積分上式即得到：

其中，V0=Vz（t=0）是貫入儀的初始貫入速度。貫入過程中貫入儀加速度大小為：

最終貫入深度P為：

3 土壤Mohr-Coulomb 參數(shù)的反演

3.1 Mohr-Coulomb 參數(shù)C，φ的反演方法

土壤Mohr-Coulomb 參數(shù)即土壤的內(nèi)聚力C和內(nèi)摩擦角φ。由于動態(tài)空腔膨脹模型刻畫了土壤力學(xué)特性參數(shù)對彈靶相互作用的影響，因此對C和φ的反演可以建立在動態(tài)空腔膨脹模型的基礎(chǔ)上。由式（3）可知：

上式表明，要反演內(nèi)聚力C和內(nèi)摩擦角φ，實際上就是要求得參數(shù)τ0和λ的值。由式（18）、式（19）（或?qū)τ诓豢紤]摩擦?xí)r的方程（21）和（22））以及A、B的表達式可見，參數(shù)τ0、λ包含在阻力系數(shù)αs、βs中。一旦 求得αs和βs，聯(lián)立方 程（18）和（19）（或?qū)τ诓豢紤]摩擦的情況，聯(lián)立方程（21）和（22）），即可求得參數(shù)τ0和λ的值。

在聯(lián)立方程（18）和（19）（或方程（21）和（22））求解參數(shù)τ0和λ時，需針對A、B取值的三種情況進行討論。由于參數(shù)A、B的表達式過于復(fù)雜，一般不能得到τ0和λ的解析表達式。此時，可采用如下幾何化方法進行求解。為簡便起見，忽略掉貫入儀與土體間的摩擦（即令μ=0），不影響對方法的討論。

（1）如果λ≠0 且λ≠3/4：

此時，由式（10）和式（11）可見，參數(shù)A、B均為（τ0，λ）的函數(shù)。根據(jù)式（21）和式（22），αs、βs也均是（τ0，λ）的函數(shù)。因此，由式（21）和式（22）可定義如下曲面：

此種情況下若存在τ0和λ的解，那么參數(shù)τ0和λ必是曲面Π1和Π2的公共零點。

（2）如果λ=3/4 或λ=0：

此時，由式（13）～式（16）可見，參數(shù)A、B僅是τ0的函數(shù)。因此αs、βs也僅是τ0的函數(shù)。于是，根據(jù)式（21）和式（22）可定義如下曲線：

此種情況下若存在τ0的解，則τ0必是曲線L1和L2的公共零點。

至此，對C、φ的反演就歸結(jié)為對阻力系數(shù)αs和βs的計算，它們可通過貫入阻力、貫入速度、貫入加速度以及最終貫入深度的測量結(jié)果計算得到。下面就針對這些情況進行具體說明。

3.2 基于阻力-速度關(guān)系進行的反演

在沖擊貫入過程中，最重要的物理量是貫入儀受到的貫入阻力Fz，它綜合反映了靶體介質(zhì)的物理力學(xué)性質(zhì)以及貫入儀與靶體材料之間的相互作用。動態(tài)空腔膨脹模型的核心就是揭示了貫入過程中的穩(wěn)定階段貫入阻力Fz（t）與速度平方間的線性關(guān)系（式（17））。因此，我們首先利用這個關(guān)系來反演土壤的Mohr-Coulomb參數(shù)。

實際當中，貫入阻力Fz（t）可由加速度傳感器測量得到；而貫入速度Vz（t）則可由加速度時程曲線積分得出。因此，利用測量數(shù)據(jù)集｛Fz，i，Vz，i｝并基于式（17）進行線性擬合，即可得到參數(shù)αs和βs的數(shù)值，并進而求得內(nèi)聚力C和內(nèi)摩擦角φ。為此，構(gòu)造誤差函數(shù)：

其 中，N是數(shù)據(jù)采樣點數(shù)目。得 到αs，βs的數(shù)值后，即可按照上節(jié)所述步驟進行C，φ值的反演。下面采用Forrestal 等人［25］報道的實驗結(jié)果來對這個反演方法進行詳細說明。

文獻［25］中，通過對不同采樣地點取樣的土壤進行三軸壓縮及單軸壓縮標定實驗，測得的土壤參數(shù)如表1 所示。需要說明的是采樣地點和現(xiàn)場貫入實驗地點的土層按照深度可分為三層［25］，第一層（Ⅰ：0～1.8 m）的彈性模量為120 MPa，鎖定體積應(yīng)變?yōu)?.17；第二層（Ⅱ：1.8～3.0 m）的彈性模量為150 MPa，鎖定體積應(yīng)變?yōu)?.10；第三層（Ⅲ：3.0～6.1 m）的彈性模量為210 MPa，鎖定體積應(yīng)變?yōu)?.13。為了簡化模型計算，文獻［25］忽略了土壤的分層，將其視作均勻介質(zhì)，并取三層土壤的平均彈性模量（160 MPa）和平均鎖定體積應(yīng)變（0.13）作為土壤基本參數(shù)。本文按照同樣的假定進行計算。

表1 Forrestal 實驗土壤參數(shù)［25］Tab.1 Soil parameters of Forrestal′s experiments［25］

以文獻［25］中的貫入實驗3 為例，貫入儀質(zhì)量為23.1 kg，彈體半徑r=47.6 mm，彈頭曲徑比ψ=3。忽略貫入過程中貫入儀與土體間的摩擦。由實驗測得的加速度數(shù)據(jù)計算得到貫入阻力和貫入速度后，給出的阻力-速度平方關(guān)系如圖2 所示。由圖2 可見，除初始階段貫入儀與土體碰撞造成的震蕩以及貫入結(jié)束階段的土體擾動外，在中間的穩(wěn)定貫入階段，貫入阻力與貫入速度的平方之間呈線性關(guān)系，這表明動態(tài)空腔膨脹模型揭示了貫入過程的關(guān)鍵特征。根據(jù)式（17），對穩(wěn)定貫入階段的阻力-速度關(guān)系進行線性擬合，結(jié)果如圖2 中紅線所示，擬合得到的阻力系數(shù)為αs=85.39 kN，βs=2.52 kg/m。圖中A、B兩點分別是擬合范圍的起點和終點，黑色圓點是參與擬合的數(shù)據(jù)。關(guān)于擬合范圍的選取，將在第4.2 節(jié)中進行詳細討論。

圖2 貫入阻力Fz和貫入速度平方Vz2間的關(guān)系［25］Fig.2 Relationship between penetration resistanceFzand penetration speed squareVz2［25］

由擬合得到的阻力系數(shù)αs、βs，可按照式（28）和式（29）畫出曲面Π1和Π2的圖像，如圖3。由圖3（a）、圖3（b）可見，曲面Π1與坐標平面Π1=0 有交線l1（曲面Π1零點的集合）；曲面Π2與坐標平面Π2=0 有交線l2（曲面Π2零點的集合）。將交線l1和l2畫在同一坐標系下，如圖4 所示，可見二者有公共交點（τ0=5.52 MPa、λ=0.58），顯然，這個交點就是曲面Π1和Π2的公共零點。按照第3.1節(jié)的分析，此時的τ0和λ就是所要求的解。此外，對于λ=3/4 和λ=0 的情況也應(yīng)當進行校驗。按照式（30）和式（31）分別畫出曲線L1和L2的圖像（如圖5（a）和圖5（b）所示），可以看到，對于這兩種情況，L1和L2都沒有公共零點，因此不存在τ0和λ的解。綜合以上分析，將τ0=5.52 MPa 和λ=0.58 代入式（27），就得到C、φ的反演結(jié)果分別為：

圖3 曲面Π1和Π2的圖像Fig.3 Images of surfacesΠ1andΠ2

圖4 λ≠0 且λ≠3/4 時交線l1和l2的圖像及其交點Fig.4 Intersection linesl1andl2and their intersection point underλ≠0 andλ≠3/4

圖5 λ=3/4 或λ=0 時曲線L1、L2的圖像Fig.5 LinesL1andL2whenλ=3/4 orλ=0

對于實驗中的土壤試樣，F(xiàn)orrestal 等人［25］給出了三軸實驗下剪切強度τ與平均應(yīng)力p之間的關(guān)系，如圖6，此處，剪切強度τ定義為試樣屈服時第一主應(yīng)力σ1與第三主應(yīng)力σ3的差［25］（τ=σ1-σ3）。注意到空腔膨脹模型的球?qū)ΨQ性，第一主應(yīng)力σ1就是空腔膨脹時的徑向應(yīng)力σr，第三主應(yīng)力σ3就是空腔膨脹時的環(huán)向應(yīng)力σθ（此時σ2=σ3=σθ），于是Mohr-Coulomb 準則（式（2））可以寫為：

根據(jù)上式，對圖6 中的數(shù)據(jù)進行線性擬合（如圖6 中紅線所示），可解得相應(yīng)的C、φ值分別為C=3.73 MPa、φ=3.48°。

圖6 剪切強度與平均應(yīng)力的關(guān)系［25］Fig.6 Relationship between shear strength and average stress［25］

對比式（36）所示的反演結(jié)果可見，對于內(nèi)聚力C的反演結(jié)果與三軸實驗的結(jié)果較為接近，而對內(nèi)摩擦角φ的反演結(jié)果卻與實驗值差別較大。這一定程度上與采用的擬合方法有關(guān)。接下來，我們分別基于貫入加速度數(shù)據(jù)和貫入速度數(shù)據(jù)對阻力系數(shù)進行擬合，并進一步討論擬合結(jié)果對方法的依賴。

3.3 基于加速度時程曲線和速度時程曲線進行的反演

貫入加速度與貫入阻力間僅相差m倍，原則上講，由貫入加速度進行阻力系數(shù)αs、βs的擬合與采用貫入阻力進行擬合沒有本質(zhì)差別。但是，在利用阻力-速度關(guān)系進行擬合時，貫入速度是由加速度的數(shù)值積分獲得，由于進行了數(shù)值積分操作，數(shù)據(jù)處理過程與直接利用加速度時程曲線進行擬合是不同的，因此擬合結(jié)果也將產(chǎn)生差別。此外，由式（17）可見，貫入阻力Fz與貫入速度的平方Vz2間呈線性關(guān)系，阻力系數(shù)αs和βs通過線性最小二乘擬合即可得到；而由式（25）可見，貫入加速度az與時間t呈現(xiàn)非線性關(guān)系，αs和βs只能通過非線性最小二乘擬合得到。

基于加速度時程曲線［25］和速度時程曲線對阻力系數(shù)的擬合如圖7。根據(jù)貫入加速度表達式（25），采 用Levenberg-Marquardt 方法［30］對文 獻［25］中貫入實驗3 給出的加速度數(shù)據(jù)進行非線性最小二乘擬合，如圖7（a）中紅線所示。擬合的初值設(shè)置為根據(jù)阻力-速度關(guān)系擬合得到的數(shù)值，即αs=85.39 kN，βs=2.52 kg/m。由加速度曲線擬合得到的阻力系數(shù)分別為αs=89.22 kN、βs=2.69 kg/m。由此得到τ0=5.61 MPa、λ=0.72，Mohr-Coulomb 參數(shù)為：

圖7 基于加速度時程曲線［25］和速度時程曲線對阻力系數(shù)的擬合Fig.7 Fitting of the resistance coefficients based on acceleration and velocity curves

與三軸實驗測量值相比，C值的相對誤差為29.2%，與上一節(jié)中基于阻力-速度關(guān)系反演得到的結(jié)果接近，φ值的反演誤差仍然較大。此外，圖7（a）給出了按照文獻［25］中取τ0=10 MPa、λ=0 時計算得到的加速度時程曲線（圖7（a）中藍線），將其與紅線相比較可以看到，τ0=10 MPa 和λ=0 并不是測量數(shù)據(jù)的最小二乘解。

按照相同的非線性最小二乘擬合方法，利用貫入速度表達式（24）亦可對阻力系數(shù)αs和βs進行擬合，如圖7（b）所示。此時，同樣將擬合的初值設(shè)置為根據(jù)阻力-速度關(guān)系擬合得到的數(shù)值，即αs=85.39 kN，βs=2.52 kg/m。由速度曲線擬合得到的阻力系數(shù)為αs=101.48 kN、βs=1.87 kg/m。相應(yīng)的τ0=7.76 MPa、λ=0.11，Mohr-Coulomb參數(shù)為：

將上式與三軸實驗測量結(jié)果對比，C值的相對誤差為2.14%，較前兩種方法的反演結(jié)果誤差大大減小，φ值的相對誤差也減小為9.77%。

總結(jié)起來，對比式（36）、式（38）和式（39）的反演結(jié)果可見，利用阻力-速度關(guān)系進行擬合的結(jié)果與利用加速度直接進行擬合的結(jié)果相近，而與利用貫入速度進行擬合的結(jié)果差別較大。更重要的是，相較于前兩者的擬合結(jié)果，利用貫入速度擬合得到的內(nèi)聚力C和內(nèi)摩擦角φ與三軸實驗給出的結(jié)果更加接近。以上結(jié)論并不是偶然的，它反映了貫入阻力系數(shù)對加速度、速度等物理量有不同的影響，這將在第4.2 節(jié)中進行更進一步的討論。此外，以上結(jié)論為土壤Mohr-Coulomb 參數(shù)反演方法的實際應(yīng)用提供了重要參照。

3.4 基于深度-初速度關(guān)系進行的反演

阻力系數(shù)αs、βs的數(shù)值還可以通過對同一狀態(tài)的土壤進行多次貫入實驗后測得的貫入深度P和初始貫入速度V0來擬合。為此，根據(jù)式（26）構(gòu)造誤差函數(shù)：

其中，N是貫入實驗次數(shù)。

可解得：

將式（44）代入式（43）并令：

實驗中，彈體質(zhì)量為286 g，彈頭長度為140 mm，曲徑比ψ=8.54［31-32］。土壤的基本參數(shù)如表2 所示。

表2 何翔實驗的土壤參數(shù)［31］Tab.2 Soil parameters of He Xiang’s experiments［31］

表3 中列出了同一剛性彈以不同初速度多次貫入同一狀態(tài)土壤后測得的貫入深度結(jié)果。圖8為的圖像及貫入深度-貫入初速度關(guān)系。將表3 中的9 個實驗數(shù)據(jù)代入式（45）可畫出函數(shù)f的圖像（圖8（a）），并得到的零點為1.09×10-3s2/m2，聯(lián)立式（41）和式（44）可解得αs=276.92 N，βs=0.30 kg/m。取鎖定體積應(yīng)變η*=0.13，可得到τ0=52.51 kPa，λ=1.16，Mohr-Coulomb 參數(shù)為：

表3 何翔貫入實驗結(jié)果［31］Tab.3 Results of He Xiang’s penetration experiments［31］

圖8 的圖像及貫入深度-貫入初速度關(guān)系Fig.8 Image of the functionand the relationship between penetration depth and initial speed

反演得到的內(nèi)聚力C的相對誤差為58.74%，內(nèi)摩擦角φ的相對誤差為3.07%?？梢姡? 個實驗數(shù)據(jù)整體擬合反演得到的內(nèi)聚力C與三軸實驗測量結(jié)果相差較大，而得到的內(nèi)摩擦角φ與三軸實驗測量結(jié)果一致。

此外，對表3中所列數(shù)據(jù)進行分組擬合。根據(jù)文獻［31］，表3 中所列實驗結(jié)果對應(yīng)三種不同的貫入初速度V0，分別為150 m/s（序號1～3），200 m/s（序號4～6）和250 m/s（序號7～9），在每個速度等級下取1 個數(shù)據(jù)點，總計3 個數(shù)據(jù)點作為一組數(shù)據(jù)（總共27 組），來擬合αs、βs的值，結(jié)果見表4。分組擬合及反演得到的內(nèi)聚力C的平均值及標準差為29.45±20.69 kPa，內(nèi)摩擦角φ的平均值及標準差為25.00±17.88°。此時，分組擬合并反演得到的C的平均值較整體擬合結(jié)果有所提高，但相比61 kPa 仍差別較大；而對于φ的反演，整體擬合結(jié)果與實驗測得的30°更加接近。

表4 27 組數(shù)據(jù)的計算結(jié)果Tab.4 Calculation results of the 27 data sets

以上分析表明，利用貫入深度-初速度關(guān)系進行反演時，內(nèi)摩擦角φ的反演結(jié)果更加可靠。

4 討論

4.1 參數(shù)的可解域

土 壤Mohr-Coulomb 參 數(shù)C、φ并 不是任意的，具有一定范圍，相應(yīng)的阻力系數(shù)αs和βs也對應(yīng)了一定范圍，稱為參數(shù)的可解域。只有擬合得到的αs和βs落在可解域范圍內(nèi)，才能在動態(tài)空腔膨脹模型框架下進行C、φ值的反演。

圖9（a）和 圖9（b）分別給出了Forrestal［25］和何翔［31］等人報道的實驗結(jié)果對應(yīng)的參數(shù)可解域。圖9（a）中給出了利用阻力-速度關(guān)系、加速度時程曲線和速度時程曲線擬合得到的αs、βs值，其均落在參數(shù)的可解域內(nèi)。圖9（b）中的分組擬合值為第3.4 節(jié)中利用貫入深度-初速度關(guān)系對27 組數(shù)據(jù)擬合得到的αs、βs值；整體擬合值為利用9 個實驗數(shù)據(jù)擬合得到的αs、βs值。這些擬合結(jié)果也都在參數(shù)的可解域內(nèi)，因此可以利用動態(tài)空腔膨脹理論來進行C、φ值的反演。

圖9 動態(tài)空腔膨脹模型下參數(shù)的可解域Fig.9 Solution domain determined by dynamic cavity expansion model

實際的地質(zhì)環(huán)境往往十分復(fù)雜，數(shù)據(jù)噪聲也難以避免，根據(jù)測量數(shù)據(jù)擬合得到的αs和βs并不總能反演得到Mohr-Coulomb 參數(shù)C、φ的數(shù)值。此時，就需要對偏離模型的過程進行更細致的分析

參數(shù)的可解域是土壤Mohr-Coulomb 參數(shù)對阻力系數(shù)作出的限制，它反映了土壤物理力學(xué)特性對貫入阻力的影響。對于給定土壤，Mohr-Coulomb 參數(shù)總是確定的，因此在理想情況下，阻力系數(shù)的實驗值必定落在可解域內(nèi)。但在實際當中，由于噪聲的影響，或簡化的模型未能涵蓋貫入過程的復(fù)雜細節(jié)，阻力系數(shù)的擬合值可能超出參數(shù)的可解域。此時，就需要對沖擊貫入的物理過程進行更加詳細的分析，并對相應(yīng)的數(shù)據(jù)處理手段進行優(yōu)化。

4.2 擬合區(qū)間的選取及對反演結(jié)果的進一步討論

由第3 節(jié)的分析可以知道，要準確地反演C、φ的數(shù)值，關(guān)鍵是要準確地擬合阻力系數(shù)αs和βs。第3.1 節(jié)至3.3 節(jié)利用貫入時程曲線進行擬合時，擬合數(shù)據(jù)是穩(wěn)定貫入階段貫入阻力與貫入速度平方呈線性關(guān)系的部分（圖2 中A、B之間的數(shù)據(jù)）。數(shù)據(jù)擬合范圍對阻力系數(shù)的擬合及Mohr-Coulomb 參數(shù)的反演有重要影響，在這一節(jié)中就來對數(shù)據(jù)擬合范圍的選取進行討論。

為此選取擬合窗口起點為圖2 中的點A（對應(yīng)時刻為t=4.84 ms），改變擬合窗口終點tend進行阻力系數(shù)αs和βs的擬合以及Mohr-Coulomb 參數(shù)的反演，結(jié)果如圖10 所示。由圖可見，不論擬合窗口的大小如何改變，基于阻力-速度關(guān)系反演的結(jié)果與基于加速度曲線反演的結(jié)果接近，而與基于速度曲線的反演有較大的差別；并且，對于內(nèi)摩擦角φ，前兩者的反演結(jié)果始終大于后者。

圖10（a）顯示了兩個關(guān)鍵擬合范圍。第一個范圍對應(yīng)三條擬合曲線的交點（tend=t2）。對于內(nèi)聚力C的反演，當tend＜t2時，基于阻力-速度關(guān)系和加速度曲線反演得到的結(jié)果小于利用速度曲線反演的結(jié)果；當tend＞t2時，基于阻力-速度關(guān)系和加速度曲線反演得到的結(jié)果大于利用速度曲線反演的結(jié)果。在tend=t2的擬合范圍下，三種擬合方法得到的C的反演結(jié)果一致。此時反演得到的數(shù)值（C=4 MPa）比三軸實驗結(jié)果（C=3.73 MPa）大7.24%。另一個關(guān)鍵擬合范圍對應(yīng)于圖10（a）中速度曲線擬合結(jié)果的最小值點tend=t1（這一點也是圖10（b）中速度曲線擬合結(jié)果的最大值點），此時，由速度曲線擬合并反演得到的內(nèi)聚力C=3.81 MPa，內(nèi)摩擦角為φ=3.14°，與三軸試驗結(jié)果一致。這一范圍就是在第3.1 節(jié)至3.3 節(jié)中擬合時所采用的范圍，圖10 中時刻tend=t1就對應(yīng)圖2 和圖7 中的點B。

圖10 擬合區(qū)間大小對C、φ反演結(jié)果的影響［25］Fig.10 Influence of the fitting interval on the inversion results ofCandφ［25］

以上分析為反演方法的實際應(yīng)用提供了一個可行方案：對于沖擊貫入實驗測得的數(shù)據(jù)，用速度曲線進行參數(shù)C（φ）的擬合，改變區(qū)間的大小進行搜索，使得C（φ）的擬合值最小（大）時對應(yīng)的區(qū)間即為最終的擬合區(qū)間，相應(yīng)的C、φ值即為最終的反演結(jié)果。

由圖2、7、10 可以看到，擬合區(qū)間終點tend一般并不能覆蓋整個貫入過程的數(shù)據(jù)段，圖10 中能夠反演得到Mohr-Coulomb 參數(shù)的最大tend為27.80 ms，這是由于參數(shù)可解域的限制。當tend大于27.80 ms 后，擬合得到的αs和βs已經(jīng)超出了參數(shù)的可解域，這表明超過27.80 ms 后的貫入過程可能存在模型忽略掉的細節(jié)。在利用動態(tài)空腔膨脹模型進行參數(shù)反演時，假定了土壤介質(zhì)的初始狀態(tài)是均勻的，這一假設(shè)保證了模型的球?qū)ΨQ性，而在現(xiàn)場貫入實驗中，土壤實際上是分層的。圖11 給出了文獻［25］中貫入實驗3 的加速度時程曲線并根據(jù)貫入儀進入不同土層（Ⅰ、Ⅱ、Ⅲ）的時刻進行了劃分。由第3.2 節(jié)的說明可知，這三層土壤的彈性模量和鎖定體積應(yīng)變是不同的，特別是第三層土壤的彈性模量明顯高于前兩層；此外，第三層土層的厚度也幾乎是前兩層的兩倍，貫入儀大部分的運動時間在第三層土壤中。這表明，第三層土壤對貫入過程的貢獻應(yīng)當比前兩層土壤的大，然而平均化的假定抹平了這些差別。當擬合區(qū)間為圖11 中的AB范圍時，貫入儀在三個土層中的運動時間相差不大，此時，平均化的假定基本適用；而當擬合區(qū)間范圍逐漸擴大，比如對于AC段或AD段，此時第三個土層的貢獻逐漸增大（但仍與前兩層可比擬），平均化的假定將產(chǎn)生誤差，此時可能需要考慮土壤分層的影響?？傊瑴蚀_地反演土壤的力學(xué)參數(shù)，需要模型能夠準確地把握貫入過程的關(guān)鍵細節(jié)。

圖11 貫入過程的加速度時程曲線以及貫入儀進入不同土層對應(yīng)的時間區(qū)間［25］Fig.11 Acceleration curve and the time interval corresponding to the penetration instrument entering different soil layers［25］

4.3 參數(shù)的敏感性

在求解沖擊貫入問題時，對于貫入阻力和貫入深度的求解屬于正問題。此時，彈性模量和鎖定體積應(yīng)變等參數(shù)對于計算結(jié)果的影響不大［25］。然而，對于通過測量貫入阻力求解靶體介質(zhì)力學(xué)特性參數(shù)這樣的反問題，反演結(jié)果往往具有參數(shù)敏感性。這種敏感性不僅與材料參數(shù)本身有關(guān)還與阻力系數(shù)有關(guān)。

圖12 和圖13 分別給出了文獻［25］和文獻［31］中，土壤彈性模量E和鎖定體積應(yīng)變η*對C、φ值反演的影響。由圖可見，隨著彈性模量E的增大，反演得到的C值減小、φ值增大；隨著鎖定體積應(yīng)變η*的增大，反演得到的C、φ值均增大。

圖12 利用速度曲線［25］進行反演時，參數(shù)E、η*對C、φ值反演的影響，其中αs=101.48 kN、βs=1.87 kg/mFig.12 Influence of parametersEandη* on the inversion results ofCandφ，respectively，whereαs=101.48 kN andβs=1.87 kg/m（using the velocity curve［25］）

由圖12（a）和圖13（a）可見，彈性模量E對于C、φ值反演結(jié)果的影響較小。在圖12（a）中，當彈性模量E的變化范圍為135 MPa 時，C值的變化大約僅有1 MPa，而φ值的變化約3°。在圖13（a）中，當彈性模量E的變化范圍為45 MPa 時，C值的變化不到2 kPa，而φ值的變化小于1°。這就是說，C、φ值的反演結(jié)果對彈性模量E不敏感。相反地，由圖12（b）和圖13（b）可見，鎖定體積應(yīng)變η*對于C、φ值的反演結(jié)果卻有很大影響，即C、φ值的反演結(jié)果對參數(shù)η*是敏感的。

圖13 利用貫入深度-初速度數(shù)據(jù)［31］進行反演時，參數(shù)E、η*對C、φ值反演的影響，其中αs=276.92 N、βs=0.30 kg/mFig.13 Influence of parametersEandη* on the inversion results ofCandφ，respectively，whereαs=276.92 N andβs=0.30 kg/m（using the penetration depth-initial velocity data［31］）

此外，從圖12（b）和圖13（b）還可以看到，不同阻力系數(shù)導(dǎo)致根據(jù)相同范圍內(nèi)的η*反演得到的C、φ值也不相同。這說明反演結(jié)果的參數(shù)敏感性也與阻力系數(shù)αs、βs有關(guān)。

5 結(jié)論

本文基于動態(tài)空腔膨脹模型，利用土壤沖擊貫入實驗中測得的阻力-速度關(guān)系、加速度時程曲線、速度時程曲線，以及貫入深度-初速度關(guān)系，建立了對土壤Mohr-Coulomb 參數(shù)（內(nèi)聚力C和內(nèi)摩擦角φ）進行反演的方法。其中，基于前三者建立的反演方法是利用單次貫入過程中的時程曲線進行的反演；而基于貫入深度-初速度關(guān)系建立的反演方法則是在對同一狀態(tài)的土壤進行多次貫入時采用。

在貫入過程中的穩(wěn)定貫入階段，貫入阻力與貫入速度平方之間呈線性關(guān)系，準確地擬合線性系數(shù)（阻力系數(shù)）αs和βs是反演Mohr-Coulomb 參數(shù)的關(guān)鍵。本文揭示了利用動態(tài)空腔膨脹模型進行Mohr-Coulomb 參數(shù)反演時存在一個關(guān)于參數(shù)的可解域。對于測得的實驗數(shù)據(jù)，只有首先保證擬合得到的阻力系數(shù)落在可解域范圍內(nèi)，C、φ值的反演才能在動態(tài)空腔膨脹模型的框架下進行。

本文分析表明利用阻力-速度關(guān)系進行反演的結(jié)果與利用加速度時程曲線進行反演的結(jié)果相近，而這兩者與利用速度時程曲線進行反演的結(jié)果差別較大。這本質(zhì)上是由于阻力系數(shù)αs、βs分別對加速度和速度的貢獻不同。因此，這三種反演方法對C、φ值反演的準確性也是不同的。其中，利用速度時程曲線進行參數(shù)反演的結(jié)果最為準確。利用Forrestal 等人［25］報道的土壤原位貫入實驗結(jié)果進行的方法驗證表明，利用速度時程曲線對土壤內(nèi)聚力反演的相對誤差為2.14%，對土壤內(nèi)摩擦角反演的相對誤差為9.77%。

此外，Mohr-Coulomb 參數(shù)的反演具有參數(shù)敏感性，這不僅與材料參數(shù)本身有關(guān)還與阻力系數(shù)αs和βs有關(guān)。C、φ值的反演結(jié)果對鎖定體積應(yīng)變η*敏感，而對彈性模量E不敏感。

本文建立的反演方法可為復(fù)雜地質(zhì)環(huán)境下土壤Mohr-Coulomb 參數(shù)的快速確定提供基礎(chǔ)，并為大范圍內(nèi)的地質(zhì)力學(xué)信息勘測提供一個新的途徑。