王 策，張展羽,，曹德君，陳 于，齊 偉，馬 靚

基于VG-PENG收縮特征曲線與收縮各向異性的裂隙率估算模型

王 策1，張展羽1,2※，曹德君3，陳 于4，齊 偉2，馬 靚1

（1. 河海大學(xué)農(nóng)業(yè)科學(xué)與工程學(xué)院，南京 211100； 2. 河海大學(xué)水利水電學(xué)院，南京 210098；3. 南京市長江河道管理處，南京 210011； 4. 江蘇省農(nóng)村水利科技發(fā)展中心，南京 210029）

為量化農(nóng)田裂隙發(fā)育程度，考慮脫濕過程中土壤孔隙在基質(zhì)域、沉降域和裂隙域間轉(zhuǎn)化，該研究提出基于土壤收縮特征和收縮各向異性的裂隙體積比率（裂隙率）關(guān)于含水率的預(yù)測模型。該模型包括3個子模型：改進(jìn)VG型式的基質(zhì)域收縮特征VG-PENG模型，描述收縮各向異性的幾何因子Logistic模型，基于上述VG-PENG收縮特征模型和幾何因子模型的裂隙率預(yù)測模型。通過土壤收縮試驗和裂隙演化監(jiān)測試驗，采用圖像處理技術(shù)提取裂隙數(shù)據(jù)，評價了該模型的優(yōu)度及適用性。結(jié)果表明，VG-PENG收縮模型具有較好的連續(xù)性和明確的物理意義，可精確描述土壤收縮特征（2>0.98）；該研究引入Logistic曲線描述土壤收縮幾何因子，揭示了收縮過程中土壤橫向開裂和縱向沉降的各向異性機理，提出了脫濕初期縱向沉降（幾何因子趨近1）、中期主沉降-副開裂（幾何因子處于1～3之間）、后期趨于穩(wěn)定3個階段，Logistic模型可精確描述收縮幾何因子隨含水率變化；基于VG-PENG收縮模型和Logistic幾何因子模型，構(gòu)建了裂隙率關(guān)于含水率的演化模型，該模型呈“S”型曲線，取決于土壤收縮屬性及其各向異性特征，裂隙率模擬值和實測值吻合較好，呈顯著水平（2>0.90，<0.001）。該研究裂隙率預(yù)測模型修正了土壤收縮各向異性在裂隙率估算中造成的誤差，并突破性地將VG-PENG收縮特征曲線進(jìn)一步推演并應(yīng)用于裂隙率模擬，可方便、快捷地通過土壤收縮數(shù)據(jù)預(yù)測農(nóng)田裂隙率隨含水率演化規(guī)律，為膨縮土裂隙流研究提供理論依據(jù)和參數(shù)基礎(chǔ)。

農(nóng)田；土壤；干縮裂隙；收縮特征曲線；收縮幾何因子；Logistic模型；裂隙率模型

0 引 言

農(nóng)田干縮裂隙破壞土壤黏聚結(jié)構(gòu)[1]，阻隔根系延展[2]，其網(wǎng)狀大孔隙為水分、溶質(zhì)遷移提供優(yōu)先通道，減少水分或營養(yǎng)物質(zhì)在土壤表層持留時間[3-5]，增加農(nóng)業(yè)污染物遷移的危害。裂隙發(fā)育具有較強隨機性，導(dǎo)致變固相骨架土壤持水、導(dǎo)水特征難以精確模擬[6]。干縮開裂土壤中水分、溶質(zhì)遷移依賴于裂隙動態(tài)變化和土壤水力特征演化，因此，裂隙率量化及預(yù)測是裂隙流理論的前提，可為裂隙優(yōu)先流模型構(gòu)建提供參數(shù)基礎(chǔ)[7]。

土壤干濕縮脹是裂隙發(fā)育的根本原因，土壤內(nèi)黏土礦物（包括高嶺土、蒙脫石、蛭石等）在干濕過程中產(chǎn)生縮脹效應(yīng)，導(dǎo)致土顆粒、團聚體重新排列和孔隙體積變化，從而引起土壤體積收縮膨脹現(xiàn)象。通常采用孔隙比和含水率關(guān)系描述土壤收縮現(xiàn)象，即土壤收縮特征曲線，該曲線包括結(jié)構(gòu)收縮、線性收縮、殘余收縮和零收縮4個階段。收縮特征曲線作為必不可少的水力參數(shù)，用于預(yù)測的模型主要有多項式模型、三直線模型、sigmoid模型[8-9]等。Peng等[10]研究發(fā)現(xiàn)，基于土壤結(jié)構(gòu)（例如容重、團聚體結(jié)構(gòu)）等可將收縮特征曲線劃分為6種型態(tài)，且均可采用VG型式的收縮模型模擬。脫濕過程中，土顆?；驁F聚體間吸應(yīng)力導(dǎo)致土壤收縮效應(yīng)，當(dāng)顆粒間吸應(yīng)力大于土壤抗拉強度時產(chǎn)生裂隙[1]。裂隙形態(tài)量化是土壤持水、收縮特征與裂隙流模型研究的基礎(chǔ)，裂隙在不同質(zhì)地、黏土礦物類型、有機質(zhì)含量等因素下呈現(xiàn)出隨機網(wǎng)絡(luò)特征。目前，尚未存在公認(rèn)的裂隙指標(biāo)體系，通常采用裂隙率、長度密度、連通性等指標(biāo)對不同吸力值和含水率狀態(tài)下的裂隙動態(tài)演化進(jìn)行研究[6,11-15]。

目前常用土壤物理方法（收縮特征曲線）、斷裂力學(xué)理論方法和裂隙網(wǎng)絡(luò)隨機法3種途徑進(jìn)行裂隙模擬，例如，Vogel 等[16]將土壤間黏聚力簡化為彈性彈簧系統(tǒng)斷裂問題，提出了基于胡克彈簧系統(tǒng)的裂隙延展模型；Chertkov[17]在細(xì)觀尺度下對黏性顆粒和團聚體的收縮特征進(jìn)行模擬，可預(yù)測細(xì)觀尺度有限尺寸條件下的裂隙率隨含水率變化；在數(shù)值模擬和計算機技術(shù)大背景下，研究者也嘗試采用LEFM、離散元和有限元方法，對裂隙發(fā)育演化問題進(jìn)行模擬研究[18-22]。基于土壤收縮特征曲線，Stewart等[23]提出了團聚體、沉降域和裂隙域的孔隙隨含水率變化物理控制方程，然而未引入土壤收縮各向異性特征。Neely等[24]采用線脹率（Coefficient of Linear Extensibility，COLE）對Bronswijk裂隙體積預(yù)測模型進(jìn)行了改進(jìn)，合理地預(yù)測了田間原位裂隙體積，然而同樣忽略了收縮各向異性問題。在土壤裂隙開閉及其優(yōu)先流研究中，通常過于關(guān)注裂隙本身而忽略變形各向異性問題，例如Coppola等[25]構(gòu)建的裂隙流雙滲透模型中，僅考慮裂隙優(yōu)先通道隨含水率的動態(tài)愈合過程而忽略了土壤縱向膨縮。土壤膨縮各向異性或裂隙發(fā)育程度不均勻性關(guān)系到裂隙率預(yù)測及其優(yōu)先流模擬精度，是值得關(guān)注的問題。

本文聚焦土壤收縮特征及其各向異性問題，通過揭示孔隙在土壤基質(zhì)域、沉降域和裂隙域中的形式轉(zhuǎn)變，探索土壤收縮特性與變形各向異性特征關(guān)于含水率的動態(tài)變化，據(jù)此提出考慮收縮特征及其各向異性的裂隙率預(yù)測模型，并采用試驗和已有文獻(xiàn)數(shù)據(jù)進(jìn)行驗證。

1 裂隙率預(yù)測模型理論

土壤收縮或開裂問題的本質(zhì)是土壤孔隙結(jié)構(gòu)在形態(tài)與位置上的轉(zhuǎn)變。本文將脫濕過程土壤孔隙轉(zhuǎn)化作為切入點，首先將VG-PENG收縮特征曲線[10]引入該裂隙率預(yù)測模型，其次提出了基于Logistic函數(shù)的收縮各向異性模型，最終基于VG-PENG模型和Logistic各向異性模型，構(gòu)建了裂隙率預(yù)測模型，采用試驗數(shù)據(jù)進(jìn)行驗證。

1.1 模型假設(shè)和物理描述

1.1.1 裂隙率預(yù)測模型假設(shè)

1）收縮特征曲線是土壤基本屬性，模型僅考慮土壤固、液、氣三相物理轉(zhuǎn)化效應(yīng)；2）根據(jù)Stewart等[26]對三維裂隙形態(tài)研究，裂隙截面形態(tài)為垂直土壤表面的三角形或楔形，因此土壤收縮分量（即橫向開裂和縱向沉降）可明確區(qū)分；3）土壤開裂過程無外荷載擾動，不存在顆粒運移現(xiàn)象，因此土壤固相質(zhì)量維持恒定，固、液、氣轉(zhuǎn)化服從質(zhì)量守恒；4）土壤變形或開裂緩慢，因而既定含水率下處于熱力學(xué)平衡狀態(tài)，即收縮變形（或裂隙率）與含水率單調(diào)對應(yīng)，可忽略收縮時間效應(yīng)。

1.1.2 收縮（或開裂）物理描述

模型研究對象為裂隙發(fā)育活動區(qū)域，其總體積total由固相顆粒（Solid，s）、水分（Water，w）和氣體（Air，a）組成。研究區(qū)域土壤內(nèi)體積和質(zhì)量表示為

式中total、s、w、a和v分別為研究區(qū)土壤總體積、土顆粒體積、水分體積、氣體體積和孔隙總體積，cm3，下標(biāo)v表示孔隙（Void），total、s和w分別表示研究區(qū)土壤總質(zhì)量、土顆粒質(zhì)量和水分質(zhì)量，g。黏壤土顆粒集聚形成團聚體，團聚體集聚形成土壤基質(zhì)，土壤孔隙結(jié)構(gòu)呈現(xiàn)出雙孔隙特征[23,27-28]。干燥過程中土壤收縮本質(zhì)為孔隙體積v轉(zhuǎn)化為沉降sub和裂隙體積crack的過程，式（1）可改寫為

式中v,agg、sub和crack分別表示基質(zhì)域中孔隙體積、沉降域體積和裂隙體積，cm3，下標(biāo)agg表示基質(zhì)域。

圖1a顯示了脫濕過程中土壤孔隙從基質(zhì)域孔隙轉(zhuǎn)化為裂隙和沉降的過程，實線外輪廓體積即為基質(zhì)域體積，由土顆粒體積、基質(zhì)內(nèi)水分和空氣體積組成。

注：表示體積，cm3，其下標(biāo)total、crack、sub、v,agg、a、s和w分別表示體積總量、裂隙體積、沉降域體積、基質(zhì)域孔隙體積、氣相體積、土顆粒體積和水分體積。

Note:is volume, cm3, and the subscript ‘total, crack, sub, v,agg, a, s and w’ denote ‘total volume, crack volume, subsidence volume, porosity volume of aggregated matrix, air phase volume, soil particles volume and water volume’, respectively.

圖1 土壤收縮過程中基質(zhì)域、裂隙域和沉降域隨含水率變化

Fig.1 Variations in matrix, crack and subsidence domain with varying water content during soil shrinkage

根據(jù)脫濕收縮最大路徑對土壤含水率和各孔隙率變量去量綱化可得：

且有

1.2 裂隙率預(yù)測模型

1.2.1 基質(zhì)域孔隙率模型——土壤收縮VG-PENG曲線

Peng等[29]研究發(fā)現(xiàn)，土壤持水特征曲線與收縮特征曲線具有相似性，均呈現(xiàn)出反曲“S”型特征，據(jù)此提出了基于VG模型的收縮特征曲線（VG-PENG模型）：

或

1.2.2 沉降域孔隙率模型

引入收縮幾何因子s量化土壤橫向收縮和縱向沉降的各向異性[30]。在開裂試樣中，橫向收縮包括裂隙和邊界孔隙。在室內(nèi)尺度，邊界效應(yīng)較小，可將其識別成一種特殊裂隙；在田間尺度，由于農(nóng)田土壤尺寸較大，邊界效應(yīng)可忽略不計，橫向收縮基本等同于干縮裂隙。收縮幾何因子s定義為[30]：

或

式中Δ為體積收縮量，cm3，0即為total。0和Δ分別為土層初始厚度和沉降高度，cm。V和H為第次測定時土壤試樣體積（cm3）和土層厚度（cm），見圖2。

注：、Δ、0和Δ分別為收縮前土壤總體積（cm3）、體積收縮量（cm3）、土樣初始厚度（cm）和試樣沉降厚度（cm）。

Note:, Δ,0and Δare total soil volume before shrinkage (cm3), change in soil volume (cm3), initial height of soil layer (cm) and vertical subsidence of soil layer (cm), respectively.

圖2 干燥過程中土壤縱向沉降和橫向收縮（或開裂）示意圖

Fig.2 Sketch of soil vertical subsidence and horizontal shrinkage (cracking) during soil drying

對式（11）變形可得：

將式（12）代入式（11）可得：

聯(lián)立式（10）和式（13）可得：

（14）

式（14）即為沉降域孔隙率預(yù)測模型。

1.2.3 裂隙率預(yù)測模型

聯(lián)立式（6）、式（9）和式（15）：

可得：

1.3 收縮幾何因子模型

因此，s可根據(jù)特征值劃分為以下5個區(qū)間：

1）s=1.0，土壤基質(zhì)僅存在沉降變形，無橫向收縮；

2）1.0

3）s=3.0，土壤收縮呈現(xiàn)出各向同性特征；

4）s>3.0，土壤基質(zhì)以水平收縮或開裂為主導(dǎo)，存在微弱縱向沉降；

s()～曲線通常呈現(xiàn)出存在上、下閾值的“S”型曲線形狀，具有初始值、增長區(qū)間、極值和最大斜率等特征參數(shù)，因此本文選取資源環(huán)境有限條件下的增長規(guī)律函數(shù)——Logistic模型[31]模擬s()曲線，公式如下：

式中s,max和s,min分別為s()的極大值和極小值，為相對含水率，和分別為平移系數(shù)和曲線增長速率系數(shù)，對應(yīng)于斜率極值處的相對含水率。

1.4 模型驗證與評估

1）收縮幾何因子Logistic模型的可行性及擬合優(yōu)度；

2）收縮特征曲線VG-EPNG模型的適用性；

3）基于上述2個模型，驗證裂隙率模型的預(yù)測能力。

2 材料與方法

擬采用室內(nèi)收縮試驗和開裂試驗，以及已公開發(fā)表文獻(xiàn)資料[32-37]中的實測數(shù)據(jù)，對上述模型進(jìn)行驗證。

2.1 供試土樣及其裂隙演化試驗

土壤試樣取自河海大學(xué)江寧節(jié)水園區(qū)試驗場稻麥輪作大田(31°86′N，118°60′E)，該地區(qū)位于長江中下游，年均降雨量1 106 mm，年均蒸發(fā)量900 mm，土壤試樣為黏棕壤，隸屬于棕色石灰土亞類灰泥土屬，總孔隙度為53.5%，年平均日照總時間為2 017.2 h，≥10 ℃積溫為4 838.2 ℃。在表層0～40 cm取樣過篩后自然風(fēng)干，將其按田間原表層容重1.30 g/cm3回填至尺寸為22.0 cm（長）× 22.0 cm（寬）×10 cm（高）的有機玻璃容器內(nèi)，回填土壤質(zhì)量為699.6 g，試樣回填高度精準(zhǔn)控制在10 mm，試樣厚度采用容器四側(cè)內(nèi)壁刻度測定，試樣初始含水率為4.8%。預(yù)濕潤過程采用模擬噴灌將400 g蒸餾水均勻噴灑至土壤表面，達(dá)到飽和且有1 mm左右微薄層水，放置于室內(nèi)自然干燥（溫度（25±2）℃，相對濕度45%±10%）。脫濕過程中，每8 h拍攝記錄裂隙形態(tài)以及土壤層厚度（分別測定四邊土層厚度求平均），保證相機固定（Canon 60D，18～135 mm鏡頭）和光源均衡，設(shè)置2組重復(fù)。同時，采用相同試樣在同等環(huán)境中測定土壤收縮特征，按原容重將風(fēng)干土回填至200 cm3的環(huán)刀中（直徑70 mm，高度52 mm），將試樣自下而上飽和處理。脫濕過程中稱質(zhì)量測定試樣質(zhì)量變化，采用游標(biāo)卡尺（精度0.02 mm）測定試樣的直徑和高度（各方向測5次求平均值），脫濕結(jié)束后置于105 ℃烘箱內(nèi)烘干，并測定土顆粒質(zhì)量和干容重，后期推算孔隙率和體積含水比，并繪制其關(guān)系曲線。

2.2 圖像處理與參數(shù)提取

采用圖像批量處理與形態(tài)學(xué)算法，基于MATLAB R2018a自編Script腳本提取裂隙面積率。截取裂隙圖像中心20 cm×20 cm作為研究區(qū)域（像素尺寸為2 400×2 400），消除容器邊界效應(yīng)。圖像中裂隙區(qū)域顏色深度高于基質(zhì)域，可根據(jù)圖像灰度直方圖分布，采用自適應(yīng)閾值分割方法優(yōu)選臨界灰度值（），將大于該值的基質(zhì)域和小于該值的裂隙域區(qū)分，分別賦值1和0，輸出裂隙二值圖像（見圖3），統(tǒng)計裂隙像素數(shù)量比例，計算裂隙域面積率（cr）[12,38]。試樣為薄層土壤，裂隙在發(fā)育過程中存在貫穿現(xiàn)象，可近似地將假設(shè)中“裂隙剖面呈現(xiàn)出倒三角形或楔形”進(jìn)一步簡化為裂隙矩形剖面，根據(jù)土壤層厚度和表面裂隙面積率估算裂隙率，即：

式中cr,i和H分別為第次測定下的裂隙面積（cm2）和

在收縮特征曲線測定中，土壤孔隙比和體積含水比分別采用以下公式進(jìn)行計算：

式中s和w分別為土顆粒密度和純水密度，g/cm3，b為土壤干容重，g/cm3，為土壤含水率，g/g。上述參數(shù)均可通過環(huán)刀內(nèi)土顆粒質(zhì)量s，水分質(zhì)量w和測得的土壤體積V進(jìn)行推算。

2.3 文獻(xiàn)數(shù)據(jù)提取

表1 土壤試樣物理參數(shù)

注：為直徑；為高度；RH為相對濕度，NA表示暫無數(shù)據(jù)。粒徑<0.002、0.002～<0.02、0.02～2.0 mm分別為黏粒、粉粒和砂粒。土壤質(zhì)地以ISSS標(biāo)準(zhǔn)劃分。SD，取樣深度; SVR，飽和孔隙比；PD，顆粒密度；RVD，殘余孔隙比。

Note:is diameter;is height; RH is relative humidity; NA means not available. Particle size <0.002, 0.002-0.02 and 0.02-2.0 mm are clay, silt and sand. Soil texture is classified by ISSS standard. SD, sampling depth; SVR, saturated void ratio; PD, particle density; RVD, residual void ratio.

2.4 統(tǒng)計分析及模型評價指標(biāo)

擬采用決定系數(shù)（Coefficient of determination，2）、均方根誤差（Root Mean Square Error，RMSE）和平均絕對誤差（Mean Absolute Error，MAE），評估模型預(yù)測優(yōu)度，檢驗?zāi)M值與實測值之間的差異[39]。

3 結(jié)果與分析

3.1 收縮幾何因子模型

表2 收縮幾何因子模型擬合參數(shù)及擬合優(yōu)度

注：s,max和s,min分別表示基于Logistic模型擬合的收縮幾何因子s的最大值和最小值，和分別為擬合的平移系數(shù)和曲線增長速率系數(shù)。

Note:s,maxands,minare maximum and minimum shrinkage geometry factor fitted from Logistic model, respectively;andare the fitted translation coefficient and growing velocity factor, respectively.

分析了18組試驗數(shù)據(jù)Logistic收縮幾何因子模型的擬合優(yōu)度，2、RMSE和MAE（見表2）。結(jié)果表明，Logistic模型可有效模擬收縮幾何因子s()隨含水率動態(tài)變化，檢驗結(jié)果為極顯著水平（<0.001），均方根誤差RMSE較小，除第3和4組外，決定系數(shù)均大于0.95，均表明Logistic模型具有較高預(yù)測精度。其中第3組和4組數(shù)據(jù)為離散破壞性試驗數(shù)據(jù)，擬合數(shù)據(jù)的連續(xù)性不高，擬合優(yōu)度相對較低，而其他組均為同組數(shù)據(jù)在相同干燥路徑的連續(xù)監(jiān)測結(jié)果，因此數(shù)據(jù)連續(xù)性、光滑性較高。值得注意的是，s極大值、斜率增長速率等參數(shù)不同，s()～曲線形態(tài)在各土壤間呈現(xiàn)出差異性，該現(xiàn)象反映了s()在不同土壤屬性或試樣尺寸下的預(yù)測難度。

3.2 收縮特征曲線VG-PENG模型

本文基于第1～2組CW1和CW2收縮和開裂試驗，以及已發(fā)表文獻(xiàn)中第5組和第18組實測數(shù)據(jù)，對式（9）VG-PENG收縮模型進(jìn)行了驗證，評估了模擬值和實測值之間擬合優(yōu)度，結(jié)果詳見表3及圖5。結(jié)果表明，VG-PENG收縮特征模型較好地描述基質(zhì)域孔隙率隨含水率變化過程，各組2>0.98, RMSE<0.04，擬合效果顯著，其結(jié)果與文獻(xiàn)[10]相一致。此外，基于和兩參數(shù)的Stewart收縮模型同樣具有高預(yù)測精度，與VG-PENG模型均為連續(xù)可導(dǎo)的“S”型函數(shù)，為后續(xù)裂隙率預(yù)測提供了可靠的曲線模型基礎(chǔ)。相比之下，VG-PENG收縮模型源起于VG持水曲線模型，具有一定的應(yīng)用基礎(chǔ)。

表3 收縮特征曲線VG-PENG模型和Stewart模型擬合參數(shù)及其優(yōu)度

注：max和min分別表示基質(zhì)域孔隙率的最大值和最小值，、和分別是基于VG-PENG收縮模型的擬合參數(shù)，和分別為Stewart模型的擬合參數(shù)；

Note:maxandminare the maximum and minimum porosity of matrix domain, respectively;,andare the fitted parameters in VG-PENG model, respectively;andare the fitted parameters in Stewart model.

3.3 裂隙域孔隙率預(yù)測模型驗證

4 討 論

4.1 土壤收縮各向異性及裂隙發(fā)育機理

土壤干燥收縮屬性是裂隙發(fā)育的根本原因。在細(xì)觀上，脫濕過程中土壤孔隙水分排出，顆?；驁F聚體間形成彎液面，內(nèi)部產(chǎn)生的基質(zhì)吸力或毛管力驅(qū)使顆粒間產(chǎn)生相互吸引的吸應(yīng)力，土壤骨架結(jié)構(gòu)（包括顆粒和團聚體結(jié)構(gòu)）趨于緊密排列，孔隙體積減小，進(jìn)而產(chǎn)生宏觀收縮效應(yīng)。土壤收縮特征曲線可分為4個階段：1）結(jié)構(gòu)收縮，2）線性收縮，3）殘余收縮和4）零收縮[10,33]。結(jié)構(gòu)收縮發(fā)生于存在團聚體大孔隙結(jié)構(gòu)的原狀土內(nèi)，此時大孔隙迅速排水，團聚體骨架未見明顯變化；在線性收縮階段，孔隙比和體積含水比斜率值趨近于1，排出水分體積完全轉(zhuǎn)化為孔隙變形，此時團聚體內(nèi)孔隙趨于近飽和狀態(tài)；當(dāng)氣體進(jìn)入團聚體內(nèi)部時（進(jìn)氣點），孔隙體積變化量開始小于水分排出量，土壤團聚體開始收縮變形。

4.2 裂隙率預(yù)測模型優(yōu)度分析

然而，該模型假設(shè)土壤含水率縱向分布均勻，未考慮大田尺度層間含水率異質(zhì)性，因此需推廣其在田間原狀土中的應(yīng)用。此外，該模型僅考慮脫濕收縮過程，土壤變形時間效應(yīng)可忽略，然而在干濕交替過程中，土壤水力特征存在脫濕-吸濕滯后現(xiàn)象，能否采用膨脹特征曲線模擬裂隙閉合演化過程，并將該模型在干濕循環(huán)中進(jìn)一步推廣，可在后期進(jìn)行探究。

5 結(jié) 論

1）根據(jù)孔隙在基質(zhì)域、沉降域和裂隙域之間轉(zhuǎn)化過程，構(gòu)建了裂隙率關(guān)于含水率的動態(tài)演化模型。該模型包括3個子模型：基于VG-PENG收縮特征曲線的基質(zhì)域孔隙率模型；描述變形各向異性的收縮幾何因子Logistic模型；綜合VG-PENG模型和Logistic模型的裂隙率預(yù)測模型。

2）VG型式的VG-PENG收縮特征曲線模型包括2個估算參數(shù)和3個擬合參數(shù)，可精確模擬多種類型土壤收縮特征，擬合效果顯著（2>0.98, RMSE<0.04）。

3）引入Logistic模型描述土壤收縮幾何因子，多種土壤數(shù)據(jù)表明，土壤收縮呈現(xiàn)出各向異性：脫濕初期僅縱向沉降，脫濕中期主要表現(xiàn)為縱向沉降，同時伴隨著水平收縮或開裂，脫濕后期趨于穩(wěn)定。Logistic模型可精確模擬收縮幾何因子隨含水率變化過程，實測值和模擬值呈顯著相關(guān)性（2=0.964，<0.001）。

4）基于VG-PENG收縮模型和幾何因子Logistic模型，構(gòu)建了裂隙率預(yù)測模型。實測值和模型模擬值之間呈極顯著水平，該模型補足了收縮各向異性導(dǎo)致的裂隙率預(yù)測誤差，可為后續(xù)裂隙優(yōu)先流模擬提供參數(shù)基礎(chǔ)。

[1] Yesiller N, Miller C J, Inci G, et al. Desiccation and cracking behavior of three compacted landfill liner soils[J]. Engineering Geology, 2000, 57(1/2): 105-121.

[2] Mawodza T, Burca G, Casson S, et al. Wheat root system architecture and soil moisture distribution in an aggregated soil using neutron computed tomography [J]. Geoderma, 2020, 359(2): 113988.

[3] Allaire S E, Roulier S, Cessna A J. Quantifying preferential flow in soils: A review of different techniques[J]. Journal of Hydrology, 2009, 378(1/2): 179-204.

[4] Zhang Z B, Peng X, Zhou H, et al. Characterizing preferential flow in cracked paddy soils using computed tomography and breakthrough curve[J]. Soil & Tillage Research, 2015, 146(2): 53-65.

[5] Beven K, Germann P. Macropores and water flow in soils revisited[J]. Water Resources Research, 2013, 49(6): 3071-3092.

[6] Wang C, Zhang Z Y, Liu Y, et al. Geometric and fractal analysis of dynamic cracking patterns subjected to wetting-drying cycles[J]. Soil & Tillage Research, 2017, 170(6): 1-13.

[7] Vogel H J, Hoffmann H, Roth K. Studies of crack dynamics in clay soil-I. Experimental methods, results, and morphological quantification[J]. Geoderma, 2005, 125(3/4): 203-211.

[8] Groenevelt P H, Grant C D. Curvature of shrinkage lines in relation to the consistency and structure of a Norwegian clay soil[J]. Geoderma, 2002, 106(3/4): 235-245.

[9] 呂殿青，邵明安. 土壤干濕收縮特征研究進(jìn)展[J]. 土壤通報，2003，34(3)：225-228.

Lv Dianqing, Shao Ming’an. A review of soil shrinkage characteristics[J]. Journal of Soil Science, 2003, 34(3): 225-228. (in Chinese with English abstract)

[10] Peng X, Horn R. Identifying six types of soil shrinkage curves from a large set of experimental data[J]. Soil Science Society of America Journal, 2013, 77(2): 372-381.

[11] DeCarlo K F, Shokri N. Effects of substrate on cracking patterns and dynamics in desiccating clay layers[J]. Water Resources Research, 2014, 50(4): 3039-3051.

[12] Wang C, Zhang Z Y, Qi W, et al. Morphological approach to quantifying soil cracks: Application to dynamic crack patterns during wetting-drying cycles[J]. Soil Science Society of America Journal, 2018, 82(4): 757-771.

[13] Tang C S, Shi B, Liu C, et al. Experimental characterization of shrinkage and desiccation cracking in thin clay layer[J]. Applied Clay Science, 2011, 52(1/2): 69-77.

[14] 張中彬，彭新華. 土壤裂隙及其優(yōu)先流研究進(jìn)展[J]. 土壤學(xué)報，2015，52(3)：477-488.

Zhang Zhongbin, Peng Xinhua. A review of researches on soil cracks and their impacts on preferential flow[J]. Acta Pedologica Sinica, 2015, 52(3): 477-488. (in Chinese with English abstract)

[15] 李文杰，張展羽，王策，等. 干濕循環(huán)過程中壤質(zhì)黏土干縮裂縫的開閉規(guī)律[J]. 農(nóng)業(yè)工程學(xué)報，2015，31(8)：126-132.

Li Wenjie, Zhang Zhanyu, Wang Ce, et al. Propagation and closure law of desiccation cracks of loamy clay during cyclic drying-wetting process[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2015, 31(8): 126-132. (in Chinese with English abstract)

[16] Vogel H J, Hoffmann H, Leopold A, et al. Studies of crack dynamics in clay soil - II. A physically based model for crack formation[J]. Geoderma, 2005, 125(3/4): 213-223.

[17] Chertkov V Y. Physical modeling of the soil swelling curve vs. the shrinkage curve[J]. Advances in Water Resources, 2012, 44(10): 66-84.

[18] Hirobe S, Oguni K. Coupling analysis of pattern formation in desiccation cracks[J]. Computer Methods in Applied Mechanics and Engineering, 2016, 307(10): 470-488.

[19] Thi Dong V, Pouya A, Hemmati S, et al. Numerical modelling of desiccation cracking of clayey soil using a cohesive fracture method[J]. Computers and Geotechnics, 2017, 85(5): 15-27.

[20] Sima J, Jiang M, Zhou C. Numerical simulation of desiccation cracking in a thin clay layer using 3D discrete element modeling[J]. Computers and Geotechnics, 2014, 56(3): 168-180.

[21] Sanchez M, Manzoli O L, Guimaraes L J N. Modeling 3-D desiccation soil crack networks using a mesh fragmentation technique[J]. Computers and Geotechnics, 2014, 62(8): 27-39.

[22] 沈珠江，鄧剛. 黏土干濕循環(huán)中裂縫演變過程的數(shù)值模擬[J]. 巖土力學(xué)，2004，25(增刊2)：1-6，12.

Shen Zhujiang, Deng Gang. Numerical simulation of crack evolution in clay during drying and wetting cycle[J]. Rock and Soil Mechanics, 2004, 25(Z2): 1-6, 12. (in Chinese with English abstract)

[23] Stewart R D, Rupp D E, Najm M R A, et al. A unified model for soil shrinkage, subsidence, and cracking[J]. Vadose Zone Journal, 2016, 15(3): 1-15.

[24] Neely H L, Morgan C L S, McInnes K J, et al. Modeling soil crack volume at the pedon scale using available soil data[J]. Soil Science Society of America Journal, 2018, 82(4): 734-743.

[25] Coppola A, Comegna A, Dragonetti G, et al. Simulated preferential water flow and solute transport in shrinking soils[J]. Vadose Zone Journal, 2015, 14(9): 1-22.

[26] Stewart R D, Najm M R A. Field measurements of soil cracks[J]. Soil Science Society of America Journal, 2020, 84(5): 1462-1476.

[27] Li X, Zhang L M. Characterization of dual-structure pore-size distribution of soil[J]. Canadian Geotechnical Journal, 2009, 46(2): 129-141.

[28] Wang G, Wei X. Modeling swelling-shrinkage behavior of compacted expansive soils during wetting-drying cycles[J]. Canadian Geotechnical Journal, 2015, 52(6): 783-794.

[29] Peng X, Horn R. Modeling soil shrinkage curve across a wide range of soil types[J]. Soil Science Society of America Journal, 2005, 69(3): 584-592.

[30] Bronswijk J B. Shrinkage geometry of a heavy clay soil at various stresses[J]. Soil Science Society of America Journal, 1990, 54(5): 1500-1502.

[31] 王海軍，劉藝明，張彬，等. 基于Logistic-GTWR模型的武漢城市圈城鎮(zhèn)用地擴展驅(qū)動力分析[J]. 農(nóng)業(yè)工程學(xué)報，2018，34(19)：248-257.

Wang Haijun, Liu Yiming, Zhang Bin, et al. Analysis of driving forces of urban land expansion in Wuhan metropolitan area based on Logistic-GTWR model[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2018, 34(19): 248-257. (in Chinese with English abstract)

[32] Sanchez M, Atique A, Kim S, et al. Exploring desiccation cracks in soils using a 2D profile laser device[J]. Acta Geotechnica, 2013, 8(6): 583-596.

[33] Peng X, Horn R. Anisotropic shrinkage and swelling of some organic and inorganic soils[J]. European Journal of Soil Science, 2007, 58(1): 98-107.

[34] Chertkov V Y, Ravina I, Zadoenko V. An approach for estimating the shrinkage geometry factor at a moisture content[J]. Soil Science Society of America Journal, 2004, 68(6): 1807-1817.

[35] Chertkov V Y. The shrinkage geometry factor of a soil layer[J]. Soil Science Society of America Journal, 2005, 69(6): 1671-1683.

[36] Huang C, Shao M, Tan W. Soil shrinkage and hydrostructural characteristics of three swelling soils in Shaanxi, China[J]. Journal of Soils and Sediments, 2011, 11(3): 474-481.

[37] 劉平，張虎元，嚴(yán)耿升，等. 土建筑遺址表部土體收縮特征曲線測定[J]. 巖石力學(xué)與工程學(xué)報，2010，29(4)：842-849.

Liu Ping, Zhang Huyuan, Yan Gengsheng, et al. Determination of soil shrinkage characteristic curve of surface soil on ancient earthen architectures[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(4): 842-849. (in Chinese with English abstract).

[38] 張展羽，朱文淵，朱磊，等. 根系及鹽分含量對農(nóng)田土壤干縮裂縫發(fā)育規(guī)律的影響[J]. 農(nóng)業(yè)工程學(xué)報，2014，30(20)：83-89.

Zhang Zhanyu, Zhu Wenyuan, Zhu Lei, et al. Effects of roots and salinity on law of development for farmland soil desiccation crack[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2014, 30(20): 83-89. (in Chinese with English abstract)

[39] 朱磊，馬榮，范東峻，等. 考慮相關(guān)長度的土壤開裂模型改進(jìn)及其參數(shù)對裂隙形態(tài)的影響[J]. 農(nóng)業(yè)工程學(xué)報，2018，34(24)：123-131.

Zhu Lei, Ma Rong, Fan Dongjun, et al. Soil cracking improved model considering with correlation length and effect of its parameters on crack morphology[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2018, 34(24): 123-131. (in Chinese with English abstract)

Crack porosity estimation model based on VG-PENG shrinkage characteristic curve and soil shrinkage anisotropy

Wang Ce1, Zhang Zhanyu1,2※, Cao Dejun3, Chen Yu4, Qi Wei2, Ma Liang1

(1.,,211100,; 2.,,210098,; 3.210011,; 4.,210029,)

Soil shrinkage and cracking are essential behavior during dehydration in the environmental geotechnical engineering of farmland. Their quantification can be necessary to determine the soil physical-hydraulic parameters and crack preferential flow in soils. In the present study, a new predicting model was proposed for the crack ratio with respect to soil water content using a VG-based soil shrinkage characteristic curve and soil shrinkage anisotropy factor. The migration and transition of pores were also considered in the solid-liquid-gas phase system in soil matrix-subsidence-crack domains. This model included a VG-based soil shrinkage characteristic model, a shrinkage anisotropy model using a Logistic curve, and the soil cracking ratio model. An indoor experiment was conducted to investigate the soil shrinkage characteristics and cracking behavior. The crack ratio was determined using image processing techniques and morphological features. The experimental data was used to evaluate the fitting of a model for cracking ratio evolution. The results showed that the VG-based shrinkage model (VG-PENG model),three fitting parameters, and two estimated parameters, well predicted the soil shrinkage characteristics in various types of soils (2>0.97, RMSE<0.04). The Logistic model was first introduced into the expression of soil shrinkage geometry, which was previously used to describe the growth principles under limited resources. Quantification of soil shrinkage anisotropy showed that the soil shrinkage was highly anisotropic. The soil shrinkage exhibited only subsidence with shrinkage geometric factor approximately equivalent to 1 in the early phase of the soil dehydration. The shrinkage showed mainly vertical subsidence with a light horizontal shrinkage (or cracking) in the middle phase, with shrinkage geometric factor varying between 1 and 3. The shrinkage geometry factor tended to be stabilized in the late phase, indicating a residual state of shrinkage. The anisotropic shrinkage with rapid change occurred in the relative water content of 0.3-0.7. The logistic shrinkage anisotropy model well predicted the shrinkage geometric factor with respect to the water content. A new model was also proposed to predict the evolution of crack ratio with respect to water content. The curve of the model showed sigmoid characteristics, depending highly on shrinkage properties and anisotropy. The simulated data showed better agreement with the experimental one, indicating an extremely significant level (2=0.974,<0.001). Since the water content within the soil layer was assumed evenly distributed, this model was considered to be appropriate in a relatively limited height of the soil layer. The evolution of crack ratio was predicted from a perspective view of soil physics rather than a mechanical view. Consequently, the soil shrinkage anisotropy was fully integrated into the modelling of the cracking ratio. A significant innovation was also made to apply the VG-type shrinkage characteristic curve to crack ratio prediction. The crack porosity prediction belonged to the field of soil physics to describe the evolution of cracking ratio using the shrinkage curve and geometry factor. The proposed model well predicted the cracking ratio in surface soils with high accuracy and convenience. A further investigation was also needed to explore the efficacy of the crack ratio model on undisturbed soils, and the effects of substrate properties on cracking behaviors. The research can provide a promising theoretical basis and parameter prediction for soil water-solute movement in soils with variable-solid phase and preferential flow in soil physics and hydrology.

farmlands; soils; cracks; shrinkage characteristic curve; shrinkage geometry factor; Logistic model; crack porosity model

2020-12-10

2021-01-31

國家自然科學(xué)基金面上項目（51879071）；江蘇省自然科學(xué)基金青年基金項目（BK20200523）；中央高?；究蒲袠I(yè)務(wù)費專項資金資助（B200201017）；江西省重點研發(fā)一般項目（20203BBGL73226）

王策，博士，講師，研究方向為農(nóng)田裂隙及其優(yōu)先流理論、高效灌排理論。Email：wangce@hhu.edu.cn

張展羽，博士，教授，研究方向為節(jié)水灌排及農(nóng)業(yè)水資源高效利用。Email：zhanyu@hhu.edu.cn

10.11975/j.issn.1002-6819.2021.07.014

S278

A

1002-6819(2021)-07-0112-10

王策，張展羽，曹德君，等. 基于VG-PENG收縮特征曲線與收縮各向異性的裂隙率估算模型[J]. 農(nóng)業(yè)工程學(xué)報，2021，37(7)：112-121. doi：10.11975/j.issn.1002-6819.2021.07.014 http://www.tcsae.org

Wang Ce, Zhang Zhanyu, Cao Dejun, et al. Crack porosity estimation model based on VG-PENG shrinkage characteristic curve and soil shrinkage anisotropy[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2021, 37(7): 112-121. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2021.07.014 http://www.tcsae.org