鄒海峰，蔡國軍，劉松玉，林 軍

（東南大學(xué) 巖土工程研究所, 江蘇 南京 210096）

1 引 言

巖土工程可靠度設(shè)計(jì)建立在對(duì)土體固有空間變異性的客觀、準(zhǔn)確評(píng)價(jià)基礎(chǔ)之上[1－3]。地質(zhì)統(tǒng)計(jì)學(xué)模型是模擬土體固有空間變異性的一種有效方法，最初應(yīng)用于采礦工程中，其后擴(kuò)展至各類工程學(xué)中。

目前地質(zhì)統(tǒng)計(jì)學(xué)已經(jīng)成為巖土工程空間變異性分析中常用的方法之一，Jones等[4]將不同孔壓靜力觸探CPTU參數(shù)的地質(zhì)統(tǒng)計(jì)學(xué)理論變差函數(shù)模型參數(shù)進(jìn)行了匯總，包括錐尖阻力、側(cè)壁摩阻力和孔隙水壓力，給出了模型參數(shù)的經(jīng)驗(yàn)取值范圍。Dawson等[5]將地質(zhì)統(tǒng)計(jì)學(xué)應(yīng)用美國 Northridge兩個(gè)沖積相地區(qū)的液化評(píng)估中，根據(jù)標(biāo)準(zhǔn)貫入試驗(yàn)（SPT）結(jié)果，估算了液化概率的空間分布規(guī)律。Murakami等[6]將地質(zhì)統(tǒng)計(jì)學(xué)方法應(yīng)用于獲取日本北 Kanto平原的地面沉降大小分布圖。Lee等[7]采用地質(zhì)統(tǒng)計(jì)學(xué)方法評(píng)價(jià)了次壓縮系數(shù)的空間不確定性。Saftner等[8]將地質(zhì)統(tǒng)計(jì)學(xué)方法用于定量評(píng)價(jià) Indiana西南部某地區(qū)CPT錐尖阻力的空間變異性。Dasaka等[9]采用地質(zhì)統(tǒng)計(jì)學(xué)結(jié)合隨機(jī)場理論，分析了香港風(fēng)化土地區(qū)決定樁基礎(chǔ)深度取值的3項(xiàng)關(guān)鍵指標(biāo)空間變異性特征。

受到設(shè)計(jì)成本的限制，常規(guī)工程中難以獲得大量具有代表性的原狀試樣，使得設(shè)計(jì)參數(shù)的空間變異性分析受到樣本容量的限制?？讐红o力觸探（CPTU）是一種高精度原位測試技術(shù)[10－11]，具有快速、經(jīng)濟(jì)和高重復(fù)性等優(yōu)點(diǎn)，非常適合于土體固有空間變異性分析。本文以宿－新高速公路建設(shè)為工程背景，進(jìn)行了CPTU測試，采用地質(zhì)統(tǒng)計(jì)學(xué)方法模擬了錐尖阻力 qt的固有空間變異性。宿-新高速沿線廣泛分布粉土和粉質(zhì)砂土，受到歷史上廢黃河演變導(dǎo)致沉積環(huán)境變化的影響，其工程性質(zhì)復(fù)雜多變，工程設(shè)計(jì)中應(yīng)當(dāng)充分考慮土體強(qiáng)度參數(shù)的空間分布規(guī)律。本研究有助于提高后續(xù)設(shè)計(jì)的可靠度。

2 地質(zhì)統(tǒng)計(jì)學(xué)模型

2.1 變差函數(shù)

獲得對(duì)區(qū)域化變量X的觀測值x(z)時(shí)，若空間點(diǎn)z只在一維z軸上變化，則將區(qū)域化變量x(z)在z，z+h兩點(diǎn)處的值之差的方差之半定義為區(qū)域化變量x(z)在z方向上的實(shí)驗(yàn)變差函數(shù)，假設(shè)N(h)是間距為h的所有點(diǎn)對(duì)的總數(shù)，則實(shí)驗(yàn)變差函數(shù)的計(jì)算公式[12]為

式中：N(h)為步距為h數(shù)據(jù)對(duì)的數(shù)量；x(zi)、x(zi+h)為間距為h的2個(gè)采樣點(diǎn)值。

實(shí)驗(yàn)變差函數(shù)本質(zhì)上是跟據(jù)離散樣本獲得的非連續(xù)分布函數(shù)。為得到連續(xù)的空間變異性評(píng)價(jià)，需要采用理論變差函數(shù)對(duì)實(shí)驗(yàn)變差函數(shù)進(jìn)行擬合。巖土工程中常用的理論變差函數(shù)模型見表1[4]。具體的選擇取決于實(shí)驗(yàn)變差函數(shù)的特征。

表1 常用理論變差函數(shù)模型[4]Table 1 Typical theoretical semi-variogram models[4]

變差函數(shù)有變程、基臺(tái)值和塊金值3個(gè)模型參數(shù)。變程是用來度量空間自相關(guān)結(jié)構(gòu)的最大距離，是變差函數(shù)達(dá)到某一穩(wěn)定值時(shí)的空間距離[4]?；_(tái)值是變差函數(shù)在變程處達(dá)到的平穩(wěn)值，趨近于樣本的方差[4]。塊金值表示區(qū)域化變量在小于采樣間距（或稱為微尺度）時(shí)非連續(xù)變異，取決于區(qū)域化變量的微尺度變異性、采樣或統(tǒng)計(jì)誤差和測試誤差[3，13]。對(duì)多維條件下的地質(zhì)統(tǒng)計(jì)學(xué)模型，可建立各個(gè)方向上的變差函數(shù)，從而簡便地?cái)U(kuò)展模型維數(shù)。

2.2 Kriging差值

Kriging插值是根據(jù)待估計(jì)點(diǎn)周圍的已知資料，以變差函數(shù)為工具，在平穩(wěn)性假定的基礎(chǔ)上，將待估計(jì)點(diǎn)處的參數(shù)值表述為變程內(nèi)所有已知點(diǎn)處參數(shù)值的線性加權(quán)平均，對(duì)待估計(jì)點(diǎn)作出最優(yōu)（即估計(jì)方差最?。?、無偏（即估計(jì)方差的數(shù)學(xué)期望為0）的估計(jì)。根據(jù)權(quán)重系數(shù)計(jì)算方法的不同，Kriging插值更進(jìn)一步分為簡單 Kriging、普通 Kriging和泛Kriging等[14]，其中普通Kriging具有計(jì)算簡便和可靠性高而更為常用，本文應(yīng)用了這一方法。

普通Kriging插值假定區(qū)域化變量x(z)的數(shù)學(xué)期望為未知常量E[ x( z)]=m，從而推導(dǎo)出插值加權(quán)系數(shù)的表達(dá)式。假定任一點(diǎn)z0處的估計(jì)值[14]為

3 試驗(yàn)資料

3.1 工程概況

試驗(yàn)段位于江蘇省宿－新高速公路宿遷一標(biāo)某施工段，地質(zhì)構(gòu)造主要為古北東向構(gòu)造處于活動(dòng)狀態(tài)的郯廬斷裂帶，是一條地震活動(dòng)帶。試驗(yàn)段地貌類型以廢黃河堆積沖積平原、堆積波狀平原為主，零星分布剝蝕殘丘。該地段全為第四系全新統(tǒng)及上、中、下更新統(tǒng)松散沉積物。場區(qū)表層為近期人工雜填土及耕值土，其下均為粉土、粉砂，地下水位3.8 m，試驗(yàn)段土層主要物理指標(biāo)見表2。

表2 試驗(yàn)場地土層的主要物理性質(zhì)指標(biāo)Table 2 Main physical characteristics indexes of soils

3.2 CPTU試驗(yàn)

多功能CPTU試驗(yàn)采用美國原裝進(jìn)口多功能數(shù)字式車載 CPTU 系統(tǒng)， 配備了最新的多功能測試探頭。探頭規(guī)格符合國際標(biāo)準(zhǔn)：錐角為60°，錐底直徑為35.7 mm，錐底截面積為10 cm2，側(cè)壁摩擦筒表面積為150 cm2，孔壓測試元件厚5 mm，位于錐肩位置（u2位置）。貫入速率為2 cm/s，沿深度每5 cm測試1組讀數(shù)。貫入過程中連續(xù)測讀錐尖阻力qt、側(cè)壁摩擦力fs，孔隙水壓力u2。

圖1為所選擇的5個(gè)CPTU 的qt測試剖面及相應(yīng)的坐標(biāo)。當(dāng)不同的地質(zhì)單元未被分開時(shí)，空間變異性分析結(jié)果將會(huì)趨向于高估土體的變異性，同時(shí)產(chǎn)生不合理的趨勢項(xiàng)估計(jì)，導(dǎo)致不平穩(wěn)的觀測資料。第③層粉質(zhì)砂土土層單元是實(shí)踐中需要重點(diǎn)考慮的抗液化加固對(duì)象，本文中選擇該單元為研究對(duì)象，深度范圍為5～15 m。

圖1 CPTU測試剖面與坐標(biāo)Fig.1 Piezocone profiles and coordinates

4 地質(zhì)統(tǒng)計(jì)學(xué)分析

4.1 趨勢項(xiàng)移除

由于CPTU剖面在豎直方向上表現(xiàn)出隨深度增加而增大的趨勢，水平方向上趨勢不明顯，本文選擇移除線性趨勢項(xiàng)來獲得弱平穩(wěn)的數(shù)據(jù)。二維條件下測量值Z(s)表達(dá)式[9]為

式中：m(si)為趨勢項(xiàng)；R(si)為殘差或稱為波動(dòng)分量；si為坐標(biāo)值，包含常量、x、y軸3項(xiàng)；β為趨勢項(xiàng)系數(shù)，采用最小平方根法進(jìn)行線性回歸分析獲得。圖 2為宿-新高速公路粉質(zhì)砂土層qt沿空間的趨勢項(xiàng)回歸分析結(jié)果，其趨勢項(xiàng)方程為

式中：x為里程坐標(biāo)；y為深度坐標(biāo)。

圖2 CPTU測試剖面的線性趨勢項(xiàng)移除Fig.2 Linear trend removal on CPTU profiles

從圖2中可以看出，CPTU錐尖阻力qt具有隨深度增加而增大的趨勢，沿里程方向趨勢并不明顯，這一結(jié)果并未超出意料，因?yàn)橼厔蓓?xiàng)的存在具有一定的物理意義[3]。在深度方向，隨著圍壓的增加，土體強(qiáng)度參數(shù)趨向于增大，qt具有明顯的增大趨勢；在水平方向上，同一土層大致經(jīng)歷了類似的沉積作用，趨于具有類似的強(qiáng)度特征。

4.2 變差函數(shù)擬合

對(duì)移除趨勢項(xiàng)之后的殘差進(jìn)行分析，計(jì)算其實(shí)驗(yàn)變差函數(shù)。圖3為4種理論變差函數(shù)模型對(duì)豎向?qū)嶒?yàn)變差函數(shù)的擬合結(jié)果。指數(shù)型模型對(duì)本文研究數(shù)據(jù)給出了最優(yōu)擬合度，采用該模型建立地質(zhì)統(tǒng)計(jì)學(xué)的理論變差函數(shù)，公式為

圖3 CPTU殘差剖面的變差函數(shù)計(jì)算及擬合Fig.3 Curves fitting for semivariogram of CPTU residuals

從圖3中可以看出，間距h < 1.05 m時(shí)，實(shí)驗(yàn)變差函數(shù)逐漸增大；h > 1.05 m后，在γ=4.45 MPa2附近穩(wěn)定地上下波動(dòng)，可以認(rèn)為殘差中不再包含趨勢項(xiàng)，表現(xiàn)出弱平穩(wěn)性，故線性趨勢項(xiàng)的移除是可以接受的。根據(jù)曲線擬合結(jié)果，殘差的塊金值為0.01，趨近于 0，表明在豎直方向上采樣誤差和 qt的測試誤差非常小。豎直方向的變程為3a=1.05 m，表明宿遷粉土qt資料在1.05 m的豎向間距范圍內(nèi)，數(shù)據(jù)具有較強(qiáng)的自相關(guān)性，而超過該范圍則幾乎不存在自相關(guān)特征。同時(shí)這一數(shù)值落在Jones等[4]推薦的取值范圍之內(nèi)，符合文獻(xiàn)資料。

水平方向，由于CPTU鉆孔數(shù)量有限，難以直接采用實(shí)驗(yàn)變差函數(shù)和理論變差函數(shù)擬合方法進(jìn)行計(jì)算，采用鄰域的方法進(jìn)行估算。設(shè)置水平方向的容差角為30°，則沿x軸方向±30°范圍內(nèi)的數(shù)據(jù)均可認(rèn)為是處于水平方向，從而可以獲得多組不同間距時(shí)的實(shí)驗(yàn)變差函數(shù)值。鄰域的概念在地質(zhì)統(tǒng)計(jì)學(xué)中已經(jīng)得到廣泛應(yīng)用[14]，以克服采樣間距不等和采樣不規(guī)則的問題。同樣利用曲線擬合方法，得到水平方向的變程為4.05 m。注意到水平方向的變程大約為豎直方向的4倍，表明試驗(yàn)段內(nèi)土體的空間自相關(guān)結(jié)構(gòu)具有顯著的各向異性，而豎直方向上變異性的變化速率要遠(yuǎn)遠(yuǎn)高于水平方向，這也與以往的研究成果相符合[4]。

4.3 Kriging插值

根據(jù)殘差剖面得到空間變異性特征之后，即可采用普通 Kriging插值預(yù)測空間中未采樣點(diǎn)處的土體工程性質(zhì)的殘差值，然后與趨勢項(xiàng)相加，得到土工參數(shù)空間分布的最優(yōu)線性無偏估計(jì)。

圖4為研究區(qū)域內(nèi)qt值的空間分布規(guī)律，圖5為相應(yīng)的二維投影圖。圖中，不連續(xù)的數(shù)據(jù)點(diǎn)為實(shí)測資料。從圖中可以看出，整體上qt仍然服從隨深度增加而增大的趨勢。Kriging插值方法能夠相對(duì)準(zhǔn)確地估計(jì)不同空間點(diǎn)處的qt值，大部分?jǐn)?shù)據(jù)點(diǎn)均位于預(yù)測曲面上或在預(yù)測曲面附近變化，因此，Kriging插值估計(jì)結(jié)果能夠與實(shí)測資料形成較好的對(duì)比。另一方面，在某些采樣點(diǎn)處 qt變化劇烈，qt具有極大值或極小值。在這些采樣點(diǎn)處，預(yù)測值往往相應(yīng)地偏小或偏大。這是由于連續(xù)的理論變差函數(shù)描述的是整個(gè)空間范圍內(nèi)自相關(guān)結(jié)構(gòu)的變化規(guī)律。從某種程度上講，基于該理論變差函數(shù)的Kriging插值是一種平滑估計(jì)，也即更加符合整個(gè)空間的平均趨勢，而不是某些極值附近的波動(dòng)規(guī)律，因此可能會(huì)降低極值點(diǎn)處估計(jì)值的準(zhǔn)確性。

圖4 研究區(qū)域qt值的Kriging預(yù)測結(jié)果Fig.4 Kriging estimates of qt values at target site

圖5 Kriging預(yù)測結(jié)果平面投影Fig.5 2D distribution of Kriging estimates of qt values

從圖5中還可以看出，當(dāng)未知點(diǎn)與已知點(diǎn)間距較大（里程為710～720 m處）時(shí)，Kriging預(yù)測值表現(xiàn)出平緩的波動(dòng)性，這是由于土體的自相關(guān)性受到距離限制，未知點(diǎn)與已知點(diǎn)間距過大時(shí)超出變程范圍之后，難以采用已知點(diǎn)觀測值進(jìn)行準(zhǔn)確估計(jì)。此時(shí)，Kriging插值算法往往給出接近于整個(gè)空間范圍內(nèi)平均值的估計(jì)結(jié)果。圖5中當(dāng)未知點(diǎn)遠(yuǎn)離測試剖面時(shí)預(yù)測值在趨勢項(xiàng)附近平緩變化。一種提高預(yù)測可靠性的方法為增加測試資料，因此，Kriging插值方法能夠?yàn)楣こ虒?shí)踐的補(bǔ)孔提供有用信息。

5 結(jié) 論

（1）宿-新高速公路某試驗(yàn)段內(nèi)的 qt具有隨深度增加而增大的趨勢，然而隨里程的變化并不顯著，可以通過移除二維線性趨勢項(xiàng)獲得具有弱平穩(wěn)性的殘差剖面。

（2）指數(shù)型理論變差函數(shù)能夠精確描述宿新高速公路某試驗(yàn)段內(nèi)尖錐阻力 qt殘差的空間變異性特征，豎直向和水平向的變程分別約為1.20 m和4.05 m，表明土體空間變異性具有顯著的各向異性。

（3）普通Kriging插值能夠準(zhǔn)確地描述qt在整個(gè)試驗(yàn)段空間范圍內(nèi)的整體變化規(guī)律，提供未采樣點(diǎn)處的最優(yōu)線性無偏估計(jì)。然而，當(dāng)極值點(diǎn)變化劇烈或未知點(diǎn)遠(yuǎn)離采樣點(diǎn)時(shí)，估計(jì)值可靠度略微降低。

[1] CHRISTIAN T J, BAECHER G. Unresolved problems in geotechnical risk and reliability[C]//Geo-Risk 2011. [S. l.]:[s. n.], 2011: 50－63.

[2] PHOON K K, CHING J Y. Beyond coefficient of variation for statistical characterization of geotechnical parameters[C]//Keynote lecture of Geotechnical and Geophysical Site Characterization 4. Florida: CRC Press,2012: 113－130.

[3] JAKSA M B. Modeling the natural variability of over-consolidated clay in Adelaide, South Australia[C]//Proceedings of the Second International Workshop on Characterisation and Engineering Properties of Natural Soils. London: Taylor and Francis, 2007: 2721－2751.

[4] JONES A L, KRAMER S L, ARDUINO P. Estimation of uncertainty in geotechnical properties for performancebased earthquake engineering[R]. Berkeley: University of California, 2002.

[5] DAWSON K M, BAISE L G. Three-dimensional liquefaction potential analysis using Geostatistical interpolation[J]. Soil Dynamics and Earthquake Engineering, 2005, 25(5): 369－381.

[6] MURAKAMI S, YASUHARA K, SUZUKI K, KOMINE H. Reliable land subsidence mapping using a spatial interpolation procedure based on geostatistics[J]. Soils and Foundations, 2006, 46(2): 123－134.

[7] LEE W J, KIM D H, CHAE Y G, RYU D W. Probabilistic evaluation of spatial distribution of secondary compression by using Kriging estimates of geo-layers[J].Engineering Geology, 2011, 122(3-4): 239－248.

[8] SAFTNER D A, GREEN R A, HRYCIW R D, FADDEN M F, et al. Characterizing spatial variability of cone penetration testing through geostatistical evaluation[C]//GeoRisk 2011, Geotechnical Risk Assessment &Management. Atlanta, GA: [s.n.], 2011: 428－435.

[9] DASAKA S M, ZHANG L M. Spatial variability of in situ weathered soil[J]. Géotechnique, 2012, 62(5): 375－384.

[10] CAI G J, LIU S Y, PUPPALA A J. Comparison of CPT charts for soil classification using PCPT data: Example from clay deposits in Jiangsu province, China[J].Engineering Geology, 2011, 121(1-2): 89－96.

[11] CAI G J, LIU S Y, PUPPALA A J. Liquefaction assessments using seismic piezocone penetration (SCPTU)test investigations in Tangshan region in China[J]. Soil Dynamics and Earthquake Engineering, 2012, 41: 141－150.

[12] 李黎, 王永剛. 地質(zhì)統(tǒng)計(jì)學(xué)應(yīng)用綜述[J]. 勘探地球物理進(jìn)展, 2006, 29(3): 163－169.LI Li, WANG Yong-gang. Review of application of geostatistics[J]. Progress in Exploration Geophysics,2006, 29(3): 163－169.

[13] JAKSA M B, BROOKER P I, KAGGWA W S.Inaccuracies associated with estimating random measurement errors[J]. Journal of Geotechnical and Geoenvironmental Engineering, 1997, 123(5): 393－401.

[14] CHILES J P, DELFINER P. Geostatistics: modeling spatial uncertainty[M]. New York: John Wiley & Sons,1999.