0.   引言

黏土干缩开裂广泛存在于自然界中，会对土体的物理力学性质产生负面影响^[1-2]，进而导致边坡失稳、核废料处置库防渗屏障失效等工程问题^[3-4]。研究黏土干缩开裂的成因和规律，对于有效预防和治理相应工程灾害具有重要意义。早期学者们开展了土体膨胀性、土层边界条件、温度等土体自身与环境参数对黏土干缩开裂影响的室内外试验^[5-10]，揭示了黏土在各类因素作用下的干缩开裂性状。但通过试验仅能获得最终的裂缝结果，难以对裂缝形成的微观机理进行深入研究。随着计算机科学的发展，各类数值模拟方法开始被应用于其中，以更加深入探讨干缩裂缝形成的微观机理。有限单元法及其拓展方法较早地应用到黏土裂隙的模拟中，取得了较好的计算结果^[11-14]。但以有限元为代表的连续介质方法，模拟黏土的干缩开裂需要预先假定开裂位置，这与天然土体的散粒体特性、干缩开裂的随机性不符。而离散单元法将研究对象视为颗粒的集合体，能够较好地解决上述问题，因此许多学者们采用离散元对黏土的干缩开裂过程进行了研究。

Peron等^[15]首次基于Youssoufi等^[16]提出的颗粒粒径收缩经验公式，通过离散元进行了长条形土的开裂模拟，取得了与室内试验较为吻合的结果。Sima等^[17]、司马军等^[18]在该经验公式的基础上使用三维胶结模型开展了圆盘状黏土的开裂模拟，发现黏土开裂受到土样内部拉应力和内部缺陷的共同作用。林朱元等^[9]、Lin等^[19]同样基于该经验公式重点研究了土体干缩开裂的边界效应问题，发现土体裂缝的发育过程主要与顶面蒸发失水和底面摩擦两种边界条件相关。Le等^[20]则将离散元与土样的水分场、应力场进行耦合，发现黏土干缩开裂的位置会加剧此处水分的蒸发使得裂缝不断发展。

值得注意的是，上述离散元研究对于土颗粒之间的作用均采用了胶结模型，这在物理意义上与实际黏土颗粒间的作用力不相符，因为胶结一旦破坏便无法恢复，而黏土颗粒间的作用力由物理化学力、毛细吸力等组成^[21]，是可恢复的。唐朝生等^[1]发现，土样的基质吸力与抗拉强度是控制裂缝发育的两个关键力学指标。施斌等^[22]则指出黏性土的干缩开裂过程不仅存在毛细水的张力或土的基质吸力，还应当存在其它颗粒间作用力如物理化学力、离子静电引力等。因此，如何在干缩开裂过程中充分考虑不同含水率时黏土颗粒之间的作用力，是当前研究的重点与难点。

为了描述非饱和土中土颗粒间作用力随含水率的变化，Lu等^[23]提出了吸应力曲线(suction stress characteristic curve，简称SSCC)。SSCC描述的是土颗粒之间的吸应力与含水率之间的关系，其中吸应力由颗粒间物理化学力、胶结力、表面张力等组成。针对黏土中存在的毛细作用与吸附作用，Zhang等^[24]提出了统一有效应力方程，将吸应力明确区分为毛细吸应力与吸附吸应力。当前SSCC理论已被广泛应用于土体强度^[25]、干缩变形^[26]等方面的研究。

本文结合SSCC理论和统一有效应力方程，将黏土中的毛细作用与吸附作用进行区分，并基于此建立考虑黏土颗粒间吸引力随含水率变化的离散元接触模型，使用该模型开展了黏土干缩开裂的离散元模拟，并将其与室内试验结果进行对比。此外，在明确区分毛细作用与吸附作用的情况下，研究黏土裂缝在不同含水率阶段的发育情况、土颗粒的位移情况、土颗粒总接触数量的变化以及土颗粒间接触力的变化。

1.   黏土干缩开裂相关的非饱和土理论

黏土干缩开裂过程与非饱和土特性密切相关。黏土自身水分减少而发生干缩开裂时，土体由饱和状态过渡为非饱和状态或正处于非饱和状态。为了研究黏土干缩开裂问题，需要引入以下相关非饱和土理论。

在非饱和土研究中，土水特征曲线(soil water characteristic curve，简称SWRC) 描述了土壤基质吸力与含水率（或饱和度）之间的关系，目前较为常用的是VG (Van Genuchten) 模型^[27]，表达式如下：

式中：${S_{\text{e}}}$为有效饱和度；$\psi$为基质吸力（kPa）；$\alpha$为进气值吸力的倒数（kPa^-1）；$n$为孔隙分布参数；$m$为反映SWRC整体几何形状的参数；$S$为饱和度；${S_{\text{r}}}$为残余饱和度。

Lu等^[23]提出的SSCC理论则描述了土颗粒之间的吸应力${\sigma _{\text{s}}}$(kPa)与含水率之间的关系。Lu等^[28]通过热力学推导${\sigma _{\text{s}}}$的表达式为

为了进一步区分由毛细作用、吸附作用产生的吸应力，Zhang等^[24]提出了统一有效应力方程，并将吸应力区分为毛细吸应力$\sigma _{\text{s}}^{{\text{cap}}}$(kPa)和吸附吸应力$\sigma _{\text{s}}^{{\text{ads}}}$(kPa)：

对于吸附吸应力$\sigma _{\text{s}}^{{\text{ads}}}$：

式中：$\sigma _{\text{s}}^{{\text{dry}}}$为宏观尺度的土颗粒之间的范德华吸引力（kPa）；${f_{\text{ads}}}(w)$为无量纲尺度函数，

式中：${\text{erf}}$为误差函数；$\beta$为反映吸附吸应力强度的无量纲参数；$w$为质量含水率（%）；$w_{{\text{tran}}}^{{\text{SS}}}$为过渡质量含水率（%）。

对于毛细吸应力$\sigma _{\text{s}}^{{\text{cap}}}$：

式中：${\alpha ^{{\text{SS}}}}$为平均毛细吸应力值的倒数（kPa^-1）；${w_{\text{s}}}$为饱和质量含水率（%）；${n^{{\text{SS}}}}$为毛细吸应力相关的孔隙分布参数；${f_{\text{cap}}}(w)$为无量纲尺度函数，

图 1为根据式（4）~（8）得到的典型曲线。根据图 1可知，在不同的质量含水率区间下，毛细吸应力与吸附吸应力的大小相应改变，总吸应力在质量含水率较低时趋向于$\sigma _{\text{s}}^{{\text{dry}}}$。

图  1  典型的SSCC

Figure  1.  Typical SSCC

下载: 全尺寸图片 幻灯片

2.   离散元接触模型与颗粒粒径收缩公式

2.1   离散元接触模型

需要合适的颗粒接触模型以描述黏土颗粒间的力-位移响应。将黏土颗粒间的作用力分为法向物理接触力$F_{\text{n}}^{\text{l}}$(N)、切向物理接触力$F_{\text{s}}^{\text{l}}$(N)、吸引力${F_{\text{a}}}$(N)。颗粒间的物理接触力采用线弹性接触模型，$F_n^l$与$F_{\text{s}}^{\text{l}}$的计算公式如下：

式中：${k_{\text{n}}}$，${k_{\text{s}}}$分别为法向与切向接触刚度（N/m）；$d$为两颗粒表面之间的距离（m）；${({F_{\text{s}}})_{\text{o}}}$为一个计算时步内初始时刻的切向力（N）；$\Delta {\delta _{\text{s}}}$为一个计算时步内的相对剪切位移（m）；$\mu$为摩擦系数。

根据第1节所述可知，非饱和土中土颗粒间存在各种作用力，在离散元中将其统一为吸引力${F_{\text{a}}}$。由于黏土颗粒间的作用力（如范德华力、双电层斥力、毛细吸力等）均存在一定的作用范围，当两个土颗粒未产生物理接触时吸引力${F_{\text{a}}}$仍会存在，因此需要考虑作用距离${d_{\text{a}}}$(m)。当两颗粒表面间的距离$d$大于${d_{\text{a}}}$时，${F_{\text{a}}}$为0，仅当$d$小于等于${d_{\text{a}}}$时考虑${F_{\text{a}}}$的存在。此处为了简便考虑，将${F_{\text{a}}}$与$d$的变化关系视为线性，此时法向接触示意图如图 2所示，颗粒间的法向作用力${F_{\text{n}}} = F_{\text{n}}^{\text{l}} + {F_{\text{a}}}$。${F_{\text{a}}}$与吸应力${\sigma _{\text{s}}}$直接对应，当吸应力${\sigma _{\text{s}}}$随含水率变化时，${F_{\text{a}}}$随之相应改变。

图  2  考虑吸引力的离散元法向接触模型

Figure  2.  DEM normal contact model considering attraction

下载: 全尺寸图片 幻灯片

2.2   颗粒粒径收缩公式

对于黏土干缩导致的体积变形，普遍的做法是采用Youssoufi等^[16]与Peron等^[15]提出的颗粒粒径收缩经验公式，该公式将离散元中的颗粒视为土颗粒与水分的集合体，水分蒸发导致颗粒收缩，进而引起土样的体积变形。粒径收缩经验公式如下：

式中：$R$为$t$时刻的颗粒半径（m）；${R_0}$为初始时刻$t = 0$的颗粒半径（m）；${\alpha _{\text{R}}}$为材料参数；$\tau$为干缩的总时间(s)。式（11）需要试验数据进行相应换算后再代入计算，过程较为繁琐且物理含义不够明确。

本文采用Guo等^[10]的方法，先通过试验数据推导得到土样体积变形，再与颗粒粒径进行换算，该方法更为直接，物理含义也更为明确。在三维情况下土样饱和时，体积变形全部由水分蒸发引起，此时含水率与体积变形的关系、体积变形与颗粒粒径的关系如下：

式中：${w_0}$为初始时刻含水率（%）；${e_0}$为土样初始时刻孔隙比；${G_{\text{s}}}$为土样相对质量密度，${\varepsilon _{\text{v}}}$为体积应变，${R_{\text{w}}}$为含水率在$w$时颗粒的半径（m）。

联立式（12），（13）可以得到在饱和状态下，失水后颗粒半径表达式：

室内试验结果发现土体的干缩裂隙大部分是在饱和阶段产生^[1]，同时已有离散元模拟研究^{[9, 15, 18]}表明在干缩过程中将土样的饱和度较高阶段近似视为饱和，认为土样体积变形均是由失水收缩引起，接着在该近似情况下开展数值模拟能够得到与室内试验较为相符的结果。本文同样基于该近似开展后续计算分析工作。

3.   离散元干缩开裂模拟与验证

3.1   室内试验结果

本文根据Peron等^{[15, 29]}的干缩开裂试验结果建立数值模型。该试验所用土样为Bioley粉土，液限为31.8 %，塑限为16.9%，相对质量密度${G_{\text{s}}}$为2.71。试验时先让土样与蒸馏水混合，令其充分饱和后制备300 mm（轴向）×50 mm（横向）×12 mm（竖向）的泥浆样，在底部设置若干平行的小凹槽以限制土样在长轴方向的位移，接着控制相对湿度为40%，温度为18 ℃，并令土样开始失水干缩。

试验初始时刻的含水率${w_0}$为49%，含水率在24%开始出现裂缝，在22%之后不再有新的裂缝出现，且含水率在49%~22%内随时间线性减小。通过游标卡尺记录土样应变。试验典型最终裂缝如图 3所示。

图  3  室内试验平行组土样干缩开裂后的典型裂缝图像^[29]

Figure  3.  Typical crack images after desiccation cracking of clay sample in parallel laboratory tests^[29]

下载: 全尺寸图片 幻灯片

3.2   离散元模拟参数取值、土样成样与干缩

本节首先通过第1节非饱和土理论与Bioley粉土的实测数据^{[15, 29]}得到Bioley粉土的SSCC。接着通过2.1节的接触模型对SSCC中的吸应力进行考虑。最后生成与实际试验尺寸相符的离散元土样，进行干缩开裂的模拟。

Bioley粉土的VG模型结果见图 4，VG模型中的$\alpha =$0.00115 kPa^-1，$n =$1.33，$m =$-0.754，${S_{\text{r}}} =$0.066。通过式（2）计算土样的SSCC，结果如图 5所示，可以发现随着基质吸力的增大，吸应力初始时呈现出线性增加的趋势，再之后趋向于一稳定值。

图  4  Bioley粉土的SWRC ^[15]与实测质量含水率^[29]

Figure  4.  SWRC^[15] of Bioley silt and measured mass water content^[29]

下载: 全尺寸图片 幻灯片

图  5  Bioley粉土的吸应力-基质吸力曲线

Figure  5.  Suction stress-matric suction curve of Bioley silt

下载: 全尺寸图片 幻灯片

结合图 4的质量含水率-基质吸力与图 5的吸应力-基质吸力曲线，得到土样在对应质量含水率下的吸应力数据，如图 6所示。

图  6  Bioley粉土的SSCC

Figure  6.  SSCC of Bioley silt

下载: 全尺寸图片 幻灯片

通过式（4）~（8）对图 6的吸应力数据点进行拟合，得到参数$\sigma _{\text{s}}^{{\text{dry}}} =$850 kPa，$\beta =$2.00，${\alpha ^{{\text{SS}}}} =$0.120 kPa^-1，${n^{{\text{SS}}}} =$1.40，$w_{{\text{tran}}}^{{\text{SS}}} =$15.8 %，${w_{\text{s}}} =$49.7 %（拟合曲线见图 6）。

由图 6可知，当含水率处于不同的区间时，毛细吸应力与吸附吸应力的大小会随之改变。根据二者的大小可以划分4个质量含水率区间：①纯毛细作用阶段，对应区间49.7%~32.7%，该阶段仅存在毛细吸应力；②毛细作用占优阶段，对应区间32.7%~27.0%，该阶段同时存在毛细吸应力与吸附吸应力，且毛细吸应力的值大于吸附吸应力；③吸附作用占优阶段，对应区间27.0%~9.5%，该阶段同时存在毛细吸应力与吸附吸应力，且吸附吸应力的值大于毛细吸应力；④纯吸附作用阶段，对应区间9.5%~0.0%，该阶段仅存在吸附吸应力。

Yao等^[30-31]的研究表明，范德华力、双电层斥力等微观力的相对作用范围可以达到黏土片厚度的约1~10倍。Soulie等^[32]的研究表明毛细力作用范围可达颗粒半径的约一倍大小。综合考虑上述作用力的作用范围，本文结合多次试算结果将${d_{\text{a}}}$确定为$1.0 \times {R_{\text{w}}}$。对于吸引力${F_{\text{a}}}$则按照中值粒径${d_{50}}$进行计算：

离散元模拟所采用的各参数取值见表 1，其中土颗粒半径选择为1 mm，这样可避免采用实际黏土颗粒尺寸所带来的过大计算量，同时已有研究表明该做法能够较充分反映干缩过程的土样变形^{[9-10, 15, 18, 20]}。

表  1  离散元模拟参数取值

Table  1.  Values of parameters for DEM simulation

参数	变量	取值
初始颗粒半径/m	${R_0}$	$1.0 \times {10^{-3}}$
阻尼系数	$\beta$	0.7
颗粒法向刚度/(N·m-1)	$k_{\text{n}}^{{\text{ball-ball}}}$	$1.0 \times {10^6}$
颗粒切向刚度/(N·m-1)	$k_{\text{s}}^{{\text{ball-ball}}}$	$8.0 \times {10^5}$
颗粒-墙法向刚度/ (N·m-1)	$k_{\text{n}}^{{\text{ball-wall}}}$	$1.0 \times {10^7}$
颗粒-墙切向刚度/ (N·m-1)	$k_{\text{s}}^{{\text{ball-wall}}}$	$8.0 \times {10^6}$
颗粒摩擦系数	${\mu ^{{\text{ball-ball}}}}$	0.20
颗粒-墙摩擦系数	${\mu ^{{\text{ball-wall}}}}$	0.0
吸引力/N	${F_{\text{a}}}$	$d_{50}^2{\sigma _{\text{s}}}$
吸引力作用范围/m	${d_{\text{a}}}$	$1.0{R_{\text{w}}}$

下载: 导出CSV 

| 显示表格

对于离散元成样，林朱元等^[9]研究发现，立方堆积生成的初始土样可以较好地模拟实际试验中相对均匀的泥浆样。但单一的立方堆积会使得颗粒之间的作用力仅存在于竖向与水平向，不存在斜向的作用力。考虑到实际试验的裂缝会沿斜向发展，本文采用立方堆积与六方堆积组合的方式来制备土样。实际计算结果证明，二者组合使用能获得较为理想的裂缝效果。

首先让颗粒在重力作用下沉积，生成300 mm×50 mm×14 mm（12 mm厚度为土样颗粒，2 mm厚度为底部边界颗粒）的离散元试样，共计26433个颗粒。接着固定底部边界颗粒在轴向的位移，以模拟试验土样在轴向受到的约束。生成的初始试样如图 7所示。参照试验条件同样令含水率从49%~22%线性减小。由图 4可知，在49%~22%含水率区间时土样饱和度均接近100%，可认为此时土样体积变形均是由失水收缩引起^{[9, 15, 18]}。试验土样呈现出扁平状使得水分能够较均匀地蒸发，且试验实测发现土样沿深度的失水速率梯度非常小^{[15, 29]}，故模拟过程令土样均匀失水。实际试验的初始孔隙比${e_0}$为1.05^[15]，失水过程的颗粒体积变化根据式（14）进行计算，颗粒间吸引力根据式（4）~（8），（15）进行计算。同时对含水率在32%~24%时离散元试样的轴向、竖向、横向应变进行记录。

图  7  离散元模拟生成的初始土样

Figure  7.  Initial soil sample generated by DEM simulation

下载: 全尺寸图片 幻灯片

3.3   离散元模拟结果与室内试验的对比

离散元的干缩开裂模拟结果如图 8所示，图 8中展示了不同含水率时颗粒沿轴向的位移。各含水率意义分别为49.0%（初始状态）、32.7%（纯毛细作用阶段结束）、27%（毛细作用与吸附作用相等）、24.0%（开始出现裂缝阶段）、22.0%（不再出现新裂缝阶段）。

图  8  离散元模拟得到的不同含水率情况下土颗粒的轴向位移

Figure  8.  Axial displacements of soil particles under different water contents obtained by DEM simulation

下载: 全尺寸图片 幻灯片

在图 8中，若相邻颗粒间的颜色差越大，说明颗粒的位移差越大，颗粒间接触会集中断开形成裂缝，如图 8（e）中的位置Ⅰ、Ⅱ所示。从图 8（a）~（e）可以发现，从初始状态的49.0%含水率至27.0%含水率，离散元试样颗粒间色差非常小，均未出现裂缝。在含水率为24.0%时，试样的左端出现一处色差较大的位置，说明此时出现裂缝，与实际试验出现裂缝时的含水率相符。随着含水率降低至22.0%，颗粒位移不断增大，颗粒间色差加剧，最终沿轴向形成显著的7条裂缝，裂缝的数量和形态与图 3的实际试验结果比较相符。

实际试验与离散元模拟得到的土样在32%~24%含水率区间时的竖向、轴向、横向应变如图 9所示。离散元模拟的竖向、轴向、横向应变趋势与试验吻合，最终值与试验接近。

图  9  室内试验^[29]与离散元模拟的土样应变-含水率关系

Figure  9.  Strain-water content relationship of soil samples in laboratory tests and DEM simulation

下载: 全尺寸图片 幻灯片

离散元模拟得到的裂缝发育历程、最终裂缝形态和土样应变都与室内试验较为相符，证明了前述理论建立的离散元模型的合理性。

3.4   干缩开裂过程的细观分析

（1）颗粒位移分析

为了进一步分析土颗粒在不同含水率阶段的具体位移情况，对离散元模拟中土颗粒轴向（即裂缝分布方向）位移绝对值的最大值、平均值进行统计，并分别除以各自的最终值进行归一化，得到归一化百分比如图 10所示。根据3.2节土样含水率区间，将图 10位移情况进行分段，离散元模拟的含水率区间为49.0%~22.0%，共包含3个阶段：纯毛细作用阶段、毛细作用占优阶段、吸附作用占优阶段。从图 10可以发现，最大值随含水率的降低迅速增大，不同阶段之间的变化趋势差异较小，而平均值在纯毛细作用阶段变化较小，在毛细作用占优阶段有所增加，到最后的吸附作用占优阶段急剧增加。在纯毛细作用阶段结束与毛细作用占优阶段结束时，最大值对应的归一化百分比分别为32.7%，56.7%，平均值对应的归一化百分比分别为13.3%，27.4%。从图 8可知最大值集中在土样的轴向两端，所以其变化仅能反映土样在干缩开裂过程整体变化比较平稳。

图  10  土颗粒轴向位移绝对值的最大值、平均值的归一化百分比

Figure  10.  Normalized percentages of maximum value and average value of absolute axial displacement of soil particles

下载: 全尺寸图片 幻灯片

对轴向位移绝对值的平均值进行分析发现，在吸附作用占优阶段之前，平均值对应的归一化百分比仅有27.4%，余下约70%的平均位移均在吸附作用占优阶段完成。根据前述裂缝图像结果，在吸附作用占优阶段，即27%含水率以前，土样没有形成裂缝，而在含水率区间为27%~22%的吸附作用占优阶段，土样会形成并发展裂缝。结合二者可知，土样在吸附作用占优阶段，土颗粒的轴向位移平均值会迅速增加，土样迅速生成并发展裂缝；在此之前，虽然含水率变化范围较大（含水率区间为49%~27%），但此时土颗粒的轴向位移变化较小，土样没有形成裂缝。这也说明，颗粒轴向位移平均值的变化能较好地反映土样裂缝的发展历程。

（2）颗粒总接触数量分析

裂缝的形成伴随着土颗粒间的脱开，对干缩过程中土颗粒总接触数量进行分析，更能反映裂缝的内在微观机理。离散元模拟的土样颗粒总接触数量与含水率的关系如图 11所示。在纯毛细作用阶段，总接触数量变化较小；在毛细作用占优阶段，总接触数量急剧减少；在吸附作用占优阶段，总接触数量略微减少后开始增加。结合图 8可知，在吸附作用占优阶段之前，土样并没有形成裂缝，但是总接触数量却大幅度减少，说明此时脱开的接触分布较为均匀，没有出现局部化情况（即形成裂缝）。随着含水率进一步降低，总接触数量进一步减少，接触集中脱开，土样开始形成裂缝，并裂开成若干小块。此后，各小块内部不再形成新的裂缝，已断开的接触会随着土样的进一步收缩而相应恢复。这与室内试验观察到的“收缩核”现象相符^[1]，即土样产生裂缝裂开成若干小块后，周围土颗粒会集中向小块区域进行收缩。

图  11  土颗粒间总接触数量与含水率的关系

Figure  11.  Relationship between total contact number of DEM soil samples and water content

下载: 全尺寸图片 幻灯片

（3）颗粒间接触力分析

裂缝的产生与颗粒间的接触力密切相关，干缩过程中土样应力集中的位置更容易导致颗粒接触断开形成裂缝^[9]。本文模拟的吸应力随含水率的变化如图 12所示，可以发现吸附吸应力在32.7%~22%含水率区间的毛细/吸附作用占优阶段呈指数型增长，使得总吸应力在该阶段急剧增加。吸应力的迅速增加直接使得颗粒间的吸引力迅速增加，进而导致颗粒间的物理接触力在上述阶段同样迅速增加，如图 13所示。

图  12  干缩过程吸应力与含水率的关系

Figure  12.  Relationship between suction stress and water content during desiccation

下载: 全尺寸图片 幻灯片

图  13  土颗粒间接触力与含水率的关系

Figure  13.  Relationship between contact force of soil particles and water content

下载: 全尺寸图片 幻灯片

为了进一步分析32.7%~22%含水率区间内颗粒间接触力与裂缝间的关系，截取图 8中位置Ⅰ、Ⅱ处颗粒物理接触力的力链进行分析，如图 14所示，可以发现：含水率为32.7%时，颗粒的力链分布表明此时颗粒间的受力较为均匀；随含水率降低，颗粒间物理接触力迅速增加，而土样底部沿轴向为固定边界，土颗粒收缩过程中沿轴向位移受到限制，土样因此会产生应力集中现象；当应力集中处土颗粒的受力无法平衡时，颗粒间的接触开始集中断开进而产生裂缝。

图  14  土颗粒间物理接触力力链与含水率的关系

Figure  14.  Relationship between physical contact force chain of soil particles and water content

下载: 全尺寸图片 幻灯片

4.   结论

本文通过SSCC理论和统一有效应力方程对黏土中的毛细作用与吸附作用进行区分，建立考虑黏土颗粒间吸引力随含水率变化的接触模型，同时考虑黏土失水过程中土样体积的变化，对黏土的干缩开裂过程开展了离散元模拟，并与室内试验结果进行对比，得到以下4点结论。

（1）结合SSCC理论、统一有效应力方程与颗粒粒径收缩公式的离散单元法能够较好地模拟黏土的干缩开裂过程。模拟得到的土样应变、裂缝发育情况和最终裂缝形态都与实际试验较为吻合。

（2）在黏土干缩开裂过程中，毛细与吸附作用共同发挥效果，且随着含水率的降低，吸附作用的影响逐渐超过毛细作用。对于模拟土样，土颗粒沿裂缝分布方向的平均位移集中在吸附作用占优阶段，占其最终值的73%，土样在此阶段会形成并发展裂缝。

（3）土颗粒的总接触数量会随着含水率的降低先减少，然后再急剧减少，但土颗粒间脱开的接触分布较为均匀，此时土样没有裂缝形成。此后土颗粒间接触脱开呈现局部化趋势，裂缝开始形成。最后土颗粒的总接触数量会随着土样的进一步收缩而略微恢复。

（4）土样进入毛细/吸附作用占优阶段后，土颗粒间的物理接触力会迅速增加导致土样内部产生应力集中现象。当应力集中处颗粒的受力无法平衡时，颗粒间接触会因此集中断开进而形成裂缝。