劉亞磊，梁 杏，2，朱常坤，李 靜

(1.中國地質(zhì)大學(xué)(武漢)環(huán)境學(xué)院，湖北武漢 430074;2.中國地質(zhì)大學(xué)(武漢)生物地質(zhì)與環(huán)境地質(zhì)國家重點實驗室，湖北武漢 430074)

分形理論是由 Mandelbrot［1］提出，用來研究無規(guī)則圖形以及復(fù)雜結(jié)構(gòu)特征的方法。分形是度量上述土壤物理因素的一種手段，分形維數(shù)則是其具體的表現(xiàn)形式［2］。Bartoli［3］和 Tyler等人［4］發(fā)現(xiàn):反映土壤質(zhì)地和結(jié)構(gòu)的物理因素，如容重、粒徑分布、孔隙度及孔隙的連通狀況等，都表現(xiàn)出分形特征，這些物理因素對土壤的水力參數(shù)存在直接或者間接的影響。因此，國內(nèi)外學(xué)者［4～7］開展了利用土壤相關(guān)物理因素的分形維數(shù)研究土壤水力參數(shù)的工作，并在前人的基礎(chǔ)上不斷提出新的看法與改進手段:Tyler等［4］首先將分形理論成功地應(yīng)用于預(yù)測土壤水分特征曲線，提出了相應(yīng)分形模型并將其擴展到三維空間。Kravchenko［8］對Tyler-Wheatcraft模型進行改進，提出了一種根據(jù)土壤顆粒分析數(shù)據(jù)分段計算孔隙表面分形維數(shù)的方法，并分段估計土壤水分特征曲線，取得了較好的結(jié)果;Toledo等［9］利用分形方法建立水力傳導(dǎo)率模型;在上述基礎(chǔ)上，劉建立等［10］由土壤顆粒質(zhì)量分布曲線計算得出孔隙表面分形維數(shù)，并利用Burdine模型和Mualem模型預(yù)測非飽和水力傳導(dǎo)度，預(yù)測精度較高。分形理論的引入使得間接計算土壤水力參數(shù)多了一種簡便而且準(zhǔn)確的方法。

但是，王國梁等［11］認(rèn)為由不同粒徑的土壤顆粒質(zhì)量計算分形維數(shù)存在不合理的假設(shè)，提出由土壤顆粒體積的大小和數(shù)量來計算體積分形維數(shù);楊金玲等［12］采用激光衍射法與吸管法實測了土壤的粒徑分布，計算并比較顆粒的質(zhì)量分形維數(shù)和體積分形維數(shù)，發(fā)現(xiàn)二者間存在一定的線性相關(guān)關(guān)系。隨著科技的進步，激光衍射可以快速獲取土壤顆粒體積累積曲線，獨立于顆粒質(zhì)量法，可以得到土壤任意兩粒徑之間的體積百分含量。利用顆粒體積分布曲線計算的孔隙表面分形維數(shù)，可以避免斯托克斯公式中不合理的假設(shè)，探討此法預(yù)測土壤水分特征曲線在國內(nèi)相關(guān)研究中并未見報道。

本文采用土壤顆粒體積分布計算土壤孔隙表面分形維數(shù)，結(jié)合de Gennes分形模型，預(yù)測土壤水分特征曲線，通過與實測數(shù)據(jù)的對比分析，探討采用土壤顆粒體積分布曲線計算表面分形維數(shù)預(yù)測土壤水分特征曲線的合理性。

1 確定表面分形維數(shù)的原理

土壤的孔隙大小分布狀況對土壤的水分特征曲線有重要影響，而孔隙表面分形維數(shù)則是描述三維空間內(nèi)土壤孔隙表面不規(guī)則性的一種量度。為了確定此分維，劉建立［10］等基于 Kravchenko 和 Zhang［8］提出的一種由土壤粒徑分布曲線計算孔隙表面分形維數(shù)的方法，在孔隙體積與孔隙半徑的對應(yīng)關(guān)系下，建立了土壤顆粒質(zhì)量累積曲線與孔隙表面分形維數(shù)之間的函數(shù)關(guān)系。本文參考劉建立［10］等的推導(dǎo)過程進行公式推演，具體步驟如下:

Pachepsky［13］認(rèn)為孔隙體積增量 dVp(≤r)和 dVp(＞r)與孔隙半徑r之間存在如下關(guān)系:

而Perrier［14］等根據(jù)孔隙體積與孔隙半徑之間的關(guān)系提出孔隙體積增量的一種變化形式:

式中:E——歐氏幾何中的拓撲維數(shù)，E=3;

β——常數(shù);

DS——孔隙表面分形維數(shù)。

假設(shè)土壤孔隙為圓柱狀，土壤孔隙體積可表示為:

此處認(rèn)為土壤孔隙最小半徑為0，結(jié)合式(2)、(3)可求出孔隙半徑小于等于r的孔隙長度L(≤r)為:

Tyler等［5］認(rèn)為N個等價半徑為R的固體顆粒構(gòu)成的毛細管長度為l(R):

式中:R——顆粒粒徑;

N——顆粒粒徑為R時的顆粒數(shù)量;

D——彎曲毛細管的孔隙分形維數(shù)，Kravchenko等［8］認(rèn)為 D=DS-1。

Kravchenko等［8］假設(shè)土粒密度 ρs不變，認(rèn)為個數(shù)N可以由半徑為R的顆粒的質(zhì)量W(R)計算出來，即:

本文采用土壤顆粒體積V(R)計算相應(yīng)粒徑級別下的顆粒數(shù)量，避免了上式中關(guān)于土粒密度不變的假設(shè)，得到個數(shù)N的計算式:

式中:S——土壤顆粒形狀的系數(shù);

e——孔隙比，即孔隙體積與顆粒體積之比:

r是顆粒半徑為R時組成的土壤孔隙半徑。假設(shè)R最小為0，對式(5)進行積分可得由半徑小于等于R的顆粒組成的孔隙半徑小于等于r的孔隙累計長度:

綜合式(4)、(8)和(9)，可得半徑為R的顆粒累積V(R)關(guān)于土壤孔隙表面分形維數(shù)之間的關(guān)系式:

對上式進行積分則可得半徑小于等于R的土壤顆粒的累積體積，假設(shè)最小半徑為0:

式中:c——常數(shù)項。

對式(11)兩邊取對數(shù)得:

式(12)與Kravchenko和劉建立等人建立的方程是相似的，只是將土壤累積質(zhì)量W(≤R)替換為累積體積V(≤R)，此式的優(yōu)點是不必假設(shè)各級土壤顆粒密度不變，累積體積則可由激光衍射法測試得出。本文中所用顆粒半徑為相應(yīng)粒徑分級的上限值與下限值的算術(shù)平均值，通過上式與土壤顆粒體積的累積曲線，可計算土壤的表面分形維數(shù)Ds。

2 土壤水分特征曲線的分形模型

利用土壤的分形維數(shù)預(yù)測土壤水分特征曲線，眾多學(xué)者按照不同的理論提出多種預(yù)測模型。

de Gennes［16］根據(jù)土壤孔隙表面是由自相似性的孔洞互相嵌套組成或者團聚體連接而成的兩種模式下分別導(dǎo)出的預(yù)測模型如下:

式中:Ψ——負壓絕對值;

Ψa——進氣值;

θs——飽和含水量。

Perrier等［14］根據(jù)孔隙體積增量的變化形式，即式(2)，提出一種水分特征曲線的分形模型:

式中:V0——孔隙半徑為0時的孔隙體積;

V——總體積。

劉建立等［10］人假設(shè) V0/V=θs- θr(θr為殘余含水量)，代入式(14)即與Brooks-Corey模型一樣:

由于 θr對上式的影響不敏感，θr≈0，代入式(15)，其形式與式(13)完全相同。因此，選取式(13)作為本文預(yù)測土壤水分特征曲線的分形模型，采用Brooks-Corey模型作為擬合實測曲線的參考模型。

3 應(yīng)用實例

3．1 土壤樣品與粒徑分級

本文中所用到的材料為華北平原辛集地區(qū)土樣(表1)。

采用張力計法進行土壤水分特征曲線的實測;通過激光粒度儀對測試土樣進行顆粒分析，獲取土壤顆粒體積累積曲線。

本文所用樣品主要為細粒土，土壤顆粒最大半徑小于 1mm(表 1)。結(jié)合 Kravchenko［8］與郭中領(lǐng)［16］等人的研究成果，基于美國土壤質(zhì)地分類系統(tǒng)，進行以下幾種粒徑分級(表2)。對不同的分級方法分別進行表面分形維數(shù)的計算。

3．2 表面分形維數(shù)計算

由所獲取的顆粒體積累積曲線，分別按照表2所列的四種不同分級方法，以lg(Rave)、lg［V(≤R)］為X、Y軸(圖1)，基于最小二乘法的線性回歸，得出二者關(guān)系曲線的斜率K。

表1 土壤樣品顆粒分級情況Table 1 Particle-size class of samples

表2 不同粒徑分級方法Table 2 Different methods for grading particle-size class

結(jié)合式(12)有K與Ds的關(guān)系:

對式(16)進行求解，由于為三維空間，取兩個計算結(jié)果中較大的數(shù)值為DS，結(jié)果見表3。

由表3可以發(fā)現(xiàn)，Ds處于2.37～2.91之間。不同分級方法中，顆粒分級界限的最大半徑越小，累積體積對數(shù)與粒徑半徑對數(shù)之間的相關(guān)系數(shù)越大，即:D1＞D2＞傳統(tǒng)7級，說明利用顆粒體積分布計算表面分形維數(shù)同樣存在無標(biāo)度區(qū)間。計算結(jié)果顯示隨著樣品粘粒含量的減少，質(zhì)地越粗，表面分形維數(shù)逐漸降低的規(guī)律，與其他學(xué)者計算的分形維數(shù)規(guī)律是一致的。四種分級計算結(jié)果相比較，D3計算得到的表面分形維數(shù)處于D2與傳統(tǒng)7級之間，且與D2計算的結(jié)果間差距較小，說明增大分級密度對計算表面分形維數(shù)影響不大;采用不同粒徑分級方法計算會影響表面分形維數(shù)的大小，但不影響其反映的物理意義。

圖1 不同分級方法下K的擬合圖(部分)Fig．1 Slope coefficients of different particle-size grading methods

表3 計算得DS值與相關(guān)系數(shù)Table 3 DSand related coefficient of samples

3．3 實驗法獲取SWRC

本節(jié)選取辛集1～4號樣品進行土壤水分特征曲線的參數(shù)擬合與預(yù)測分析。首先將1～4號原狀樣通過張力計與稱重方式獲取試樣的土壤水分特征曲線，通過RETC軟件中的Brooks-Corey模型進行擬合求參，獲取進氣值Ψa，并與預(yù)測擬合結(jié)果對比。表4為所選4個樣品的參數(shù)擬合結(jié)果。

表4 土壤水分特征曲線擬合參數(shù)Table 4 The fitted parameters of SWRC

3．4 利用分形維數(shù)進行SWRC的預(yù)測

分別將D1和D2兩種分級計算的表面分形維數(shù)(表3)與表4中所列參數(shù)代入de Gennes預(yù)測模型(式14)，繪制預(yù)測土壤水分特征曲線(圖2);并采用均方根誤差RMSE作為衡量預(yù)測模型與實測值之間準(zhǔn)確性的計算標(biāo)準(zhǔn)(以含水量作為評價值)，見式(17):

式中:Xobs——實測值;

Xi——預(yù)測值，i=1，2，3，…，N。

由圖2可以很直觀地發(fā)現(xiàn)，采用D1分級方法計算出的表面分形維數(shù)預(yù)測的土壤水分特征曲線與實測值擬合得好，而通過D2進行的計算結(jié)果相對較差(表5)。

圖2 兩種分形維數(shù)下土壤水分特征曲線的預(yù)測值與實測值Fig．2 Prediction and observed points of SWRC by using two kinds of surface fractal dimension

表5 預(yù)測與實測的RMSE值Table 5 The RMSE of prediction and observed points

對D1分級方法的預(yù)測結(jié)果進行分析:最大的均方根誤差出現(xiàn)在4號土樣，為0.0105cm3/cm3，而其他三個土樣的RMSE均小于0.01cm3/cm3，也就是說，采用此種分級方法計算得出的表面分形維數(shù)對預(yù)測土壤水分特征曲線具有較高的適用性，盡管3號土樣在高吸力段的預(yù)測值與實測值略有偏差，但其他試樣的預(yù)測結(jié)果在整個吸力范圍內(nèi)與實測值較匹配。D2的預(yù)測結(jié)果與實測值只有在低吸力段匹配較好，隨著吸力的增加，預(yù)測值與實測值出現(xiàn)較大偏差，土樣預(yù)測值與實測值之間的RMSE≥2.11E-02cm3/cm3。

4 結(jié)論

(1)利用土壤顆粒體積累積曲線計算表面分形維數(shù)，避免了土壤顆粒密度不變的假設(shè)，改善了計算過程中的缺陷。

(2)相同試樣由不同顆粒分級方法計算得出的表面分形維數(shù)大小不同;相同分級方法下，試樣的表面分形維數(shù)隨著土壤粘粒的減少而變小，符合一般規(guī)律;即顆粒分級方法的不同不影響其反映的物理意義。

(3)按照D1計算的表面分形維數(shù)預(yù)測土樣吸力大于進氣值的曲線趨勢精度高，誤差小，這是由于土壤負壓超過進氣值后，土壤水分特征曲線主要受由半徑小的顆粒組成的細小孔隙控制，從側(cè)面說明土壤孔隙對水力參數(shù)的影響。因此，此法可以用于預(yù)測更高吸力值情況下的曲線變化趨勢;而D2的情況，其預(yù)測的效果較差，誤差偏大，但可以滿足田間粗略估計的需要;應(yīng)用土壤顆粒體積累積曲線在無標(biāo)度區(qū)間內(nèi)計算的表面分形維數(shù)更適合以de Gennes模型來預(yù)測細粒土的土壤水分特征曲線。

(4)筆者認(rèn)為，利用顆粒體積累積曲線計算表面分形維數(shù)的方法比較適用于相對均質(zhì)的土樣，可以避免由于參與顆粒分析的樣品量少而不具有代表性的不利影響。在無標(biāo)度區(qū)間計算的表面分形維數(shù)預(yù)測結(jié)果較好，這可能是由于大于0.1mm的土壤顆粒組成的孔隙并不具有很強的自相似性從而影響分形維數(shù)的計算。因此，利用土壤顆粒體積累積曲線按照無標(biāo)度區(qū)間計算表面分形維數(shù)，從而預(yù)測土壤水分特征曲線是可行的，此法具有較明確的物理意義。

［1］ Mandelbrot B B.The fractal geometry of nature［M］.San Francisco:W H Freeman，1983.

［2］ 任權(quán)，王家鼎，袁中夏，等.高速鐵路地基黃土微結(jié)構(gòu)的分形研究［J］.水文地質(zhì)工程地質(zhì)，2007，34(6):76-82.［Ren Q，Wang J D，Yuan Z X，et al.Fractal study on loess microstructure of one high-speed railway foundation［J］.Hydrogeology ＆ Engineering Geology，2007，34(6):76-82.(in Chinese)］

［3］ Bartoli F.Structure and self-similarity in silty and sandy:the fractal approach［J］.Soil Sci，1991，42:167-185.

［4］ Tyler S W，Wheatcraft S W.Fractal scaling of soil particle-size distribution:Analysisand limitations［J］.Soil Sci.Soc.Am.1992，56:362-369.

［5］ Tyler S W，Wheatcraft S W.Application of fractal mathematics to soil water retention estimation［J］.Soil Sci Soc Am，1989，53:987-99.

［6］ 楊培嶺，羅遠培，石元春.用粒徑的重量分布表征的土壤分形特征［J］.科學(xué)通報，1993，38(20):1896-1899.［YANG P L，LUO Y P，SHI Y C.Fractal feature ofsoilon expression by weight distribution of particle size［J］.Chinese Science Bulletin，1993， 38(20):1896-1899.(in Chinese)］

［7］ 楊建，陳家軍，楊周喜，等.松散砂?？紫督Y(jié)構(gòu)、孔隙分形特征及滲透率研究［J］.水文地質(zhì)工程地質(zhì)，2008，35(3):93-98.［YANG J，CHEN J J，YANG Z X，et al.A study of pore structure，pore fractal feature and permeability of unconsolidated sand［J］.Hydrogeology＆ Engineering geology，2008，35(3):93-98.(in Chinese)］

［8］ Kravchenko A，Zhang R D.Estimating the soil water retention from particle 2 size distributions:A fractal approach［J］.Soil Science，1998，163:171-179.

［9］ Toledo P G，Novy R A，Davis H T，et al.Hydraulic conductivity of porous media at low water content［J］.Soil Science Society of America Journal，1990，54:673-679.

［10］ 劉建立，徐紹輝，劉慧，等.確定田間土壤水力傳導(dǎo)率的分形方法［J］.水科學(xué)進展，2003，14(4):464-469.［LIU J L，XU S H，LIU H，et al.Fractal method for estimating field hydraulic conductivity function from soilparticle-size distribution［J］.Advances In Water Science，2003，14(4):464-469.(in Chinese)］

［11］ 王國梁，周生路，趙其國.土壤顆粒的體積分形維數(shù)及其在土地利用中的應(yīng)用［J］.土壤學(xué)報，2005，42(4):545-550.［WANG G L，ZHOU S L，ZHAO Q G.Volume fractal dimension of soil particles and its applications to land use［J］.Acta Pedologica Sinica，2005，42(4):545-550.(in Chinese)］

［12］ 楊金玲，李德成，張甘霖，等.土壤顆粒粒徑分布質(zhì)量分形維數(shù)和體積分形維數(shù)的對比［J］.土壤學(xué)報，2008，45(3):415-419.［YANG J L，LI D C，ZHANG G L，et al.Comparison of mass and volume fractal dimensions of soil particle-size distributions［J］.Acta Pedologica Sinica，2008，45(3):415-419.(in Chinese)］

［13］ Pachepsky Y A，Shcherbakov R A，Korsunskaya L P.Scaling of soil water retention using a fractal model［J］.Soil Science，1995，159:99-104.

［14］ Perrier E，Rieu M，Sposito G，et al.Models of water retention curve for soils with a fractal pore size distribution［J］.Water Resources Research，1996，32:3025-3031.

［15］ de Gennes P G.Partial filling of a fractal structure by a wetting fluid［C］//Adler D.Physics of Disordered Materials.New York:Plenum Press，1985:227-241.

［16］ 郭中領(lǐng)，符素華，張學(xué)會，等.土壤粒徑重量分布分形特征的無標(biāo)度區(qū)間［J］.土壤通報，2010，41(3):537-541.［GUO Z L，F(xiàn)U S H，ZHANG X H，et al.Scale-ree Domain of Fractal Characteristic of the Soil Particle-size Distribution［J］.Chinese Journal of Soil Science，2010，41(3):537-541.(in Chinese)］