多重干湿循环下的土壤裂缝网络几何与分形特征

刘 阳，张展羽，范世敏，王 策

(1.河海大学南方地区高效灌排与农业水土环境教育部重点实验室，南京 210098；2.河海大学水利水电学院，南京 210098)

0 引 言

干缩开裂是土壤中会出现的一个普遍现象，在农业中，土壤的干缩开裂具有十分广泛的影响，一方面能改善排水，改变土壤结构，对土壤渗透作用的改进，但另一方面由于水分和溶质通过土壤裂缝的快速运移，导致作物生长过程中的容易遭受水分胁迫[1]和养分胁迫，也容易导致地下水污染[2]。因此研究土壤裂缝对提高灌溉水利用率和作物产量，防治地下水污染等方面具有重要的理论与实践意义。

近年来，国内外多名研究者围绕干湿循环对土壤性质的影响开展了大量的工作。土壤干缩开裂受土壤特性影响，Horgan and Young[3]提出了基于裂缝几何特征的二维模型模拟裂缝干缩开裂过程，利用模型描述裂缝开展的多个参数，沈珠江等[4]基于非饱和土简化固结理论对黏土干湿循环过程中表面裂缝的演变过程进行了数值模拟。张展羽等[5-7]通过室内试验研究了土壤干缩裂缝发育规律，分析了不同含水率下土壤裂缝网络几何形态特征以及干燥增湿过程中裂缝的开闭规律，同时发现了作物根系通过限制裂缝拓宽对土壤开裂有一定的抑制作用，而盐分对裂缝发育起到促进作用。熊东红等[8]通过对不同类型土壤的研究分析，发现裂缝发育程度与土壤黏粒含量、容重正相关，而与有机质含量及总孔隙度负相关。通过干湿循环的作用，土壤某些特征会发生改变。张俊然等[9]对多次干湿循环过程中土-水特征曲线进行了试验研究及模拟。土壤裂缝与干湿循环之间也存在很大关系。唐朝生等[10]研究发现随着湿干循环次数的增加，土壤开裂临界含水率、裂缝率及试样干缩稳定厚度明显增加并趋于稳定，且裂缝变得曲折粗糙而成锯齿状。

由于土壤裂缝形成的复杂性，没有一个单一的分类方案能够充分的描述所有的裂缝图案的性质，而分形几何[11]提供了一个数学框架，以量化的方式来模拟裂缝网络，裂缝分形维数的大小直接反映出裂缝网络的复杂程度和自相似程度，有适当的土壤性质与裂缝分形特征的描述，在实际中可能会提供更快速的水分与污染物的运移。赵龙山[12]利用体积-表面积法对微地形分形维数进行了初步计算。胡晓蕾[13]发现同一质地的颗粒体积分形维数与土壤水分特征曲线分形模型中的孔隙大小分形维数呈良好的线性关系。YU.gm等[14]利用分形维数定量描述土壤裂缝网络的演化过程。J.U. Baer[15]研究确定了土壤干缩开裂的裂缝边缘分形维数和质量分形维数，并进行了土壤裂缝特征分析。彭瑞东等[16]对岩石孔隙的CT图像进行了分形维数计算，发现了相同的孔隙率下，孔隙结构越复杂，其分析维数越大。

本文对多次干湿循环影响下土壤干缩裂缝的形成、发展与闭合过程进行室内模拟，对裂缝基本几何特征以及分形维数进行了定量分析。研究成果将有助于揭示干缩裂缝开闭规律及干湿循环对裂缝形态的作用规律，有助于进一步探求裂缝网络开裂闭合机理及土壤水分在固相骨架变化时的非饱和运移规律。

1 材料与方法

1.1 材料与处理

试验土壤取自河海大学南方地区高效灌排与农业水土环境教育部重点实验室试验大田(31°86′N，118°60′E)的旱地土壤表层，土壤以黄褐土为主，土壤质地为壤质黏土，田间持水率 31.4%，凋萎系数 15.0%，土壤容重 1.38 g/cm3。将土样过 2 mm 筛后风干，按原容重回填至边长22.0 cm的方形有机玻璃容器，回填厚度10 mm(容器外壁有刻度尺)，均匀加水使土样完全饱和。以土样饱和为初始状态，将其放置于温度和相对湿度相对恒定的实验室内(温度为20±2 ℃，相对湿度为45%±10%)进行多次干湿循环处理。

表1 土壤试样物理性质Tab.1 Physical properties of the tested soil

干湿循环试验： 干燥过程中每隔 8 h 称重一次，同时采用佳能 EOS 60D 数码相机(分辨率 5 184×3 456)拍照记录土壤试样表面形态；当表面出现裂缝后，间隔时间缩短至 4 h。当 12 h 内(连续 3 次测量)试样质量变化在 0.3%内时，认为裂缝形态达到稳定状态，本次干燥过程结束，开始进行增湿试验；增湿过程中使用气压喷壶对试样表面喷水，每次均匀喷水约 20 g，耗时 3 min，喷水强度约为 8.3 mm/h。喷灌过程中保证土壤表面结构不被水分冲击破坏，并且在土壤表面不产生洼水或径流。一次喷水 30 min 后，土壤因增湿而产生的膨胀变形达到稳定，称取质量，拍照记录，并再次喷水；重复上述步骤，直至土壤裂缝完全闭合，增湿试验结束。为了保证试验条件的一致性，每次拍摄时试样、相机(三脚架固定)、电子秤放置在同一位置，同时遮蔽所有外界光源，仅用日光灯斜射照明。

1.2 数字图像处理技术

为消除边界效应，采用 Photoshop CS6截取出试样中间 16 cm×16 cm 为研究区域，并将截取后的图片转换为灰度图像。结合Matlab图像处理功能，通过图片二值化，去除杂点，骨架化等过程[6]，将土壤裂缝图像进行处理。

1.3 裂缝几何参数

(1)重量含水率θ(%)：试样中水分质量与干土质量之比。

(2)利用Matlab的图像分析及及形态学测量功能统计裂缝面积Acr(cm2)，裂缝骨架总长度Lcr(cm)，裂缝平均宽度Wcr(cm)。

(3)裂缝面积率Rcr，计算公式如式(1)。

(1)

式中:Rcr为裂缝面积率，%；Ai为第i条裂缝的面积，cm2；A0为研究区域总面积，cm2。

(4)连通性指数K：裂缝交叉点数与裂缝节点数之比，计算公式如式(2)。

(2)

式中：K为裂缝连通性指数；Nb、Ne分别为裂缝交叉点数与端点数。

1.4 分形特征算法

几何分形方法主要有周长-面积维数、盒维数和信息维数等。盒维数能精确地计算不规则物体的分维，且计算简单、操作方便、适用范围广。因此，本文采用盒维数法[17]来研究土壤裂缝的分形特征，该方法的具体思路如下:取边长为ε1(又称标度)的正方形网格(即小盒子),把研究区域覆盖起来，计算非空盒子数，记为N(ε1)。然后缩小盒子的尺寸ε2=ε1/2,得到非空盒子数N(ε2)，…，依次类推，当εi=εi-1/2时，统计出非空盒子数N(εi)，以loge(εi)作为横坐标、logeN(εi)作为纵坐标，作loge(εi)～logeN(εi)图，在一定区域内直线logeN(εi)=-Dloge(εi)+C的斜率D，即为研究区域的盒维数。

当r→0时，得到盒维数：

(3)

利用Matlab软件计算出研究区域的分析维数D，相关系数R2[18]。

2 结果与分析

2.1 视觉对比裂缝形态变化

土壤在若干次干湿循环中形成的裂缝形态会产生差异，4次干湿循环稳定状态下的土壤裂缝网络如图1所示。观察发现，随着土壤干湿循环次数的增加，形成裂缝的面积与宽度相对减小，但裂缝的总体形态与网络结构保持相似。研究发现，在黏土土壤的第一个干燥周期的收缩会导致土壤内部结构产生不可逆的变化，在这个周期中颗粒之间的凝聚力可能被破坏，使得土壤颗粒间的吸应力变得弱化。增湿过程中，由于水分的浸入，进一步减弱了重新排列后的土壤结构。随后的干燥将再次引起收缩，裂缝将出现在土壤胶结力最薄弱的位置。在多次干湿循环试验中，增湿过程会导致部分裂缝愈合，但是这些区域相对土壤基质团聚体仍较为薄弱，并在随后的循环周期中，这些原先开裂区容易重新开裂，所以在首次循环之后，裂缝的形态不会发生显著的改变。

图1 不同干湿循环次数稳定裂缝网络形态Fig.1 Stable crack morphology graphs with different drying-wetting cycles

2.2 裂缝几何参数规律

(1)裂缝面积率随含水率的变化规律。干湿循环过程中裂缝面积率与含水率的关系曲线如图2所示。S1～S4代表分别进行过干湿循环的次数。

图2 多次干湿循环裂缝面积率与土壤含水率关系曲线Fig.2 Curves of crack area ratio and moisture content with multiple drying-wetting cycles

图2中可以看出，每次干湿循环过程中，干燥过程与增湿过程形成一条“8”字形曲线:在含水率15%左右时，干燥过程与增湿过程的裂缝面积率达到一致；当含水率<15%时，干燥过程中裂缝面积率随含水率变化更缓慢；当含水率>15%时，增湿过程中裂缝面积率随含水率变化更缓慢。研究表明，土壤脱水与吸水过程中，相同含水率下裂缝形态有较大出入，其主要原因可联系到土壤水分特征曲线的回滞作用[19]。相同含水率下，脱湿过程中的吸力要大于吸湿过程，造成土颗粒间的吸应力在湿润与干燥过程中有所差异，膨胀收缩程度在相同含水率下是不可逆的，因此产生裂缝的速率也不尽相同。此外，土壤结构发生变化亦会造成脱水与吸水过程中裂缝形态的不同：在干燥过程中，土壤收缩，微观状态下的土壤细颗粒发生移动，导致土壤固相骨架的收缩，而在增湿过程中，土壤膨胀并非完全可逆，因此裂缝未按照原路径闭合。

干燥过程中，随着含水率的减少，裂缝面积率逐渐降低，过程大体分为3个部分：当40%≤含水率≤48%时，土壤开始出现裂缝，裂缝面积率变化缓慢；当10%≤含水率≤40%时，土壤裂缝发育迅速，裂缝面积率的增幅明显；当含水率≤10%时，裂缝面积率增加变慢，在含水率为5%左右达到最大值并保持稳定；增湿过程中，当3%≤含水率≤13%时，裂缝面积率减幅明显，当13%≤含水率≤30%时，裂缝闭合速度变缓，在含水率超过30%后裂缝迅速闭合。在4次增湿过程中，初始的裂缝面积率不同，但4条曲线的变化规律基本相似。4次干燥过程中土壤面积率与含水率的变化规律相似，其差异性体现在每经过一次干湿循环之后，含水率所对应的裂缝面积率相对减少，4次干燥过程的裂缝面积率最终稳定在14.09%、10.58%、9.43%、8.69%，说明干湿循环对裂缝的开裂有一定的抑制作用，并且抑制作用随着干湿循环的次数逐渐减弱，最终趋于稳定。

(2)裂缝长度随含水率的变化规律。干湿循环过程中裂缝长度与含水率的关系曲线如图3所示。

图3 多次干湿循环裂缝总长度与土壤含水率关系曲线Fig.3 Curves of crack length and moisture content with multiple drying-wetting cycles

由图3中曲线可以发现，在干燥过程中，当含水率 40%时，裂缝形成初始时期，裂缝总长度发育缓慢，当20% 含水率 40%时，裂缝总长度增长迅速，而当含水率低于20%时，裂缝长度基本保持稳定，此时裂缝结构基本形成，土壤结构逐渐收缩，裂缝开始拓宽，平均宽度增加。4次干燥过程最终的裂缝总长度基本一致，说明最终形成的裂缝中心骨架基本一致。

在增湿过程中，随着含水率的增高，裂缝长度减短，其中分为3个阶段：第一个阶段，裂缝长度随含水率的减少基本呈水平线性发展，此时裂缝处于横向闭合状态，裂缝长度基本不变；第二个阶段，在含水率达到25%左右，裂缝长度随着含水率的增加而迅速减短，直至第3个阶段，整体裂缝95%闭合，曲线变化缓慢。对比4条增湿过程的曲线图，我们发现随着干湿循环次数的增加，第一阶段向第二阶段的过渡点含水率逐渐减小，其原因可能是随着干湿循环次数增加，裂缝面积率减小，在初始裂缝长度基本一致的情况下，裂缝面积率小的试样将更快的结束横向闭合时期，渡过第一阶段，来到第二阶段。

(3)裂缝平均宽度随含水率的变化规律。干湿循环过程中裂缝平均宽度与含水率的关系曲线如图4所示。

图4 多次干湿循环裂缝平均宽度与土壤含水率关系曲线Fig.4 Curves of crack mean width and moisture content with multiple drying-wetting cycles

通过观察图4中裂缝平均宽度与含水率关系曲线可以看出:在含水率<25%时，其变化规律与图2中曲线相似，而图3中，裂缝的长度在这个时期变化缓慢，所以裂缝处于横向开裂(闭合)状态，裂缝宽度的变化直接线性影响裂缝的面积率变化；在含水率<25%时，干燥过程中，裂缝平均宽度随含水率的减少稳定的增加；增湿过程中，裂缝大量闭合，裂缝长度迅速减小，而其中刚闭合的裂缝对裂缝平均宽度产生了较大影响，导致其不规则的波动下降，直至裂缝网络完全闭合。

将4次干燥过程结束后的裂缝平均宽度进行对比，发现平均宽度逐渐减短，并且减幅有逐渐降低的趋势，由于裂缝长度相差不大，不难发现干湿循环对裂缝平均宽度和对裂缝面积有相似的抑制作用，并且抑制作用强度随着循环次数逐渐减弱。

(4)裂缝连通性随含水率的变化规律。干湿循环中连通性指数与土壤含水率的关系曲线如图5所示。

图5 多次干湿循环裂缝连通性指数与土壤含水率关系曲线Fig.5 Curves of crack connectivity index and moisture content with multiple drying-wetting cycles

在干燥过程中[图5(a)]，连通性系数大体趋势随着含水率的减小而增大。其中，曲线在12%<含水率<25%的区域中出现了不规律的部分，是由于裂缝骨架形成稳定后，土壤中出现的部分独立的细小裂缝，在未连通到裂缝骨架之前，对数据产生的突然变化，使得曲线产生波动，4次干燥过程最终裂缝形态的土壤连通性指数为0.74、0.84、0.94、0.92，说明在经过干湿循环的作用后，裂缝网络形成的更加完整，并且干湿循环的促进作用下，连通性指数在第3次循环后逐渐稳定。在增湿过程中[图5(b)]，连通性指数随着含水率下降而减小，在增湿初期，曲线的降幅较小，趋势平缓，在含水率达到12%左右，曲线开始不规则下降，在含水率=25%时，连通性指数变化加快，在最终达到完全闭合状态。

2.3 裂缝网络分形规律

分形理论的出现使得人们对事物之间关系的认识思维方式有线性进展到了非线性阶段，其研究范围是具有自相似性，复杂性，不规则性的现象，而分形维数就是对于自相似性的现象在几何性质上的尺度描述。我们在对土壤裂缝二维图像的分形研究中，可以得到一个有效的空间维数来D表示，它不是整数，而是一个介于1到2之间可以连续变化的数。

图6为干湿循环中分形维数随着土壤含水率变化的曲线。

图6 多次干湿循环下分形维数与土壤含水率关系曲线Fig.6 Curves of fractal dimension and moisture content with multiple drying-wetting cycles

4次干湿循环中干燥与增湿过程曲线均不重合，说明土壤裂缝分形维数的变化是一个不可逆的过程，并且交点在含水率12%～15%区域之间，在这个交点之前，干燥过程的分型维数大于增湿过程中的分形维数，在这个交点之后，增湿过程的分形维数大于干燥过程中的分型维数。在增湿过程中，对比4条曲线我们可以发现，分形维数随着干湿循环次数的增加呈“阶梯状”下降，在含水率35%之前，曲线坡度较缓，在含水率35%之后，曲线迅速下落；在4次干燥过程中，S1的曲线相对平滑，分形维数随着含水率的减少而稳定的增加，与之形成对比的S2～S4曲线规律相似，在干燥初期，分形维数增长迅速，在含水率达到30%左右时，分形维数增长幅度急剧放缓，直至裂缝形态最终稳定。

由图2对比发现，在含水率 30%时，4次干燥过程的裂缝面积率随含水率变化规律是相似的，但在同一时期，首次干燥过程分形维数的增长速度是明显慢于第2～4次干燥过程的，说明在这期间4次干燥过程中裂缝面积的增长速率虽然是相近的，但首次干燥过程的裂缝网络的复杂性却低于在第2～4次干燥过程。这是由于土壤随着含水率的减少逐渐形成一级裂缝、二级裂缝、在二级裂缝上再长出裂缝分叉，最终形成裂缝网络，是一个随着含水率渐进的过程；而在经历干湿循环后的第2～4次干燥过程中，初始状态裂缝网络虽然完全闭合，但由于土壤颗粒之间黏聚力的不可修复性，使得整个裂缝骨架的开裂不再是一个渐进的过程，而是一个几乎同时进行的过程，这就是分形维数在裂缝中心骨架形成之前增长迅速，而后增长的速率放缓的原因。

图7为4次干燥过程土壤最终状态的分形维数。

通过盒子数(Number of Box)和盒子尺寸(Size of Box)的对数线性回归来确定图像的分析维数，拟合方程显示了相当高的准确性(R2均大于0.995)，所测得的分形维数分别为1.645 52、1.579 67、1.558 91、1.545 18，随着干湿循环次数的增加是一个下降的趋势。由此可知，干湿循环的过程对裂缝的分形性质有着一定的抑制作用。

图7 多次干湿循环稳定裂缝网络图像的分形维数Fig.7 Fractal dimension of stable crack morphology graphs with different drying-wetting cycles

3 结 论

本文通过数字图像处理技术及Matlab测量统计功能，对土壤在不同干湿循环级数和不同含水率下的几何与形态参数进行定量分析，得出以下结论：

(1)随着干湿循环次数的增加，干燥过程和增湿过程中的裂缝面积率随含水率变化的规律几乎不受影响，但干燥稳定似的裂缝面积率呈递减趋势，减幅逐渐减小。干湿循环过程对土壤干缩裂缝面积有一定的抑制作用。

(2)土壤在干燥过程中形成的裂缝稳定长度几乎不受干湿循环作用的影响，在增湿过程中，随着干干湿循环次数的增加，裂缝长度随含水率变化从缓慢到迅速的过渡点不断提前。

(3)裂缝平均宽度随着含水率变化分两个阶段，第一个阶段与裂缝面积率变化规律类似，第二个阶段增湿过程中裂缝平均宽度呈波动减小趋势。随着干湿循环次数的增加，裂缝稳定平均宽度也逐渐减小。

(4)随着干湿循环次数的增加，裂缝网络形成不断完善，连通性增强，但在多次干湿循环后，作用效果减弱，连通性指数在第3次循环之后趋于稳定。

(5)采用盒维数法对裂缝网络分形特征分析可知，土壤经历过第一次开裂后，在原裂缝处开裂的速率加快，开裂现象较为集中，分形维数增长迅速，形成裂缝骨架网络时所对应的含水率提前；裂缝稳定形态的分形维数随着干湿循环次数的增加是减小趋势，说明干湿循环对裂缝网络分形特征起抑制作用。

[1] 田洪艳，周道玮，李质馨，等. 土壤胀缩运动对草原土壤的干扰作用[J]. 草地学报, 2003,(3):261-268.

[2] Greve A, Andersen M S, Acworth R I. Investigations of soil cracking and preferential flow in a weighing lysimeter filled with cracking clay soil[J]. Journal of Hydrology, 2010,393(1-2SI):105-113.

[3] Horgan G W, Young I M. An empirical stochastic model for the geometry of two-dimensional crack growth in soil (with Discussion) [J]. Geoderma, 2000,96(4):263-276

[4] 沈珠江，邓 刚. 黏土干湿循环中裂缝演变过程的数值模拟[J]. 岩土力学, 2004,(S2):1-6.

[5] 张展羽，王 策，朱成立，等. 土壤干缩裂缝网络几何特征[J]. 地球科学(中国地质大学学报), 2014,(10):1 565-1 572.

[6] 张展羽，朱文渊，朱 磊，等. 根系及盐分含量对农田土壤干缩裂缝发育规律的影响[J]. 农业工程学报, 2014,(20):83-89.

[7] 李文杰，张展羽，王 策，等. 干湿循环过程中壤质黏土干缩裂缝的开闭规律[J]. 农业工程学报, 2015,(8):126-132.

[8] 熊东红，杨 丹，李佳佳，等. 元谋干热河谷区退化坡地土壤裂缝形态发育的影响因子[J]. 农业工程学报, 2013,(1):102-108.

[9] 张俊然，许 强，孙德安. 多次干湿循环后土-水特征曲线的模拟[J]. 岩土力学, 2014,(3):689-695.

[10] Tang C, Cui Y, Shi B, et al. Desiccation and cracking behaviour of clay layer from slurry state under wetting-drying cycles[J]. Geoderma, 2011,166(1):111-118.

[11] Mandelbrot B. The Fractal Geometry of Nature. W. H. Freeman and Company[M]. 1982.

[12] 赵龙山，梁心蓝，张青峰，等. 基于微DEM的黄土坡耕地微地形分形维数计算[J]. 节水灌溉, 2010,(10):5-7.

[13] 胡晓蕾，张富仓，赵 旭. 杨凌包气带深层黄土颗粒分形特性的实验研究[J]. 节水灌溉, 2009,(3):1-4.

[14] Yu G M, Xie H P, Zhao J F, et al. Fractal evolution of a crack network in overburden rock stratum[J]. Discrete Dynamics in Nature and Society, 2000,5(1):47-52.

[15] Baer J U, Kent T F, Anderson S H. Image analysis and fractal geometry to characterize soil desiccation cracks[J]. Geoderma, 2009,154(1-2):153-163.

[16] 彭瑞东，杨彦从，鞠 杨，等. 基于灰度CT图像的岩石孔隙分形维数计算[J]. 科学通报, 2011,(26):2 256-2 266.

[17] 孙博玲. 分形维数(Fractal dimension)及其测量方法[J]. 东北林业大学学报, 2004,(3):116-119.

[18] 杨书申，邵龙义. MATLAB环境下图像分形维数的计算[J]. 中国矿业大学学报, 2006,(4):478-482.

[19] Mualem Y, Beriozkin A. General scaling rules of the hysteretic water retention function based on Mualem's domain theory[J]. European Journal of Soil Science, 2009,60(4):652-661.