0.   引言

土-水特征曲线（SWCC）定义了在不同吸力桩体下，土体中的水分含量随土吸力变化的状态方程。SWCC可用于估算非饱和土的渗透系数方程和抗剪强度方程^[1-2]。SWCC在脱湿和浸润两个过程中会表现出明显的滞后现象，导致非饱和土的水力-力学特性也会表现一定的滞后性^[3-5]。国内外学者提出了不同的模型用以评估SWCC滞后特性^[6-22]。Pham等^[23]将这些模型分类为两类：域模型和经验模型。经验模型的计算结果在很大程度上依赖于模型参数的选取，而这些模型参数往往是通过拟合或经验的方法确定，其本身没有实际的物理意义。本文在域模型的基础上，提出物理-统计模型用以测算SWCC滞后曲线。因此，本文并未对经验模型展开深入讨论。

域模型可进一步划分为独立域模型和非独立域模型。域模型采用两种类型的元素定义分布函数，即：α类元素：定义当吸力从${\psi _{{\text{wet,p}}}}$随$\delta ({\psi _{{\text{wet}}}})$变化而减小时，填充到孔隙中的水量（即土体水分增量）；β类元素：定义当吸力从${\psi _{{\text{dry,p}}}}$随$\delta ({\psi _{{\text{dry}}}})$变化而增大时，从孔隙中排出的水量（即土体水分减量）。Poulovassilis^[9]根据试验数据，提出确定分布函数的方法，如图 1所示。土体分别在D点和H点开始浸润，得到不同的浸润路径（即DC和HC），最终得到不同的浸润曲线。当土体饱和时，HC与DC汇集于同一点C，D点的含水率比H点高0.71，表明相比于DC路径，沿HC路径单位土体额外吸收0.71的水分。Poulovassilis^[9]根据HC和DC曲线的差异性，计算了这额外0.71水分在不同尺寸孔隙中的分布，如图 1（b）。沿路径DC，当土吸力从16 cm水头降至12 cm水头时，土体水分变化量为δV_12-16DC；沿路径HC，在吸力变化相同的情况下土体水分变化量为δV_12-16HC。对比发现δV_12-16HC比δV_12-16DC多0.09，表明额外0.71的水分有0.09分布在12~16 cm水头对应的孔径孔隙中。以此类推，额外0.71的水分在对应0~4，4~8，8~12，12~16，16~20 cm水头孔径孔隙的分布分别为0.04，0.18，0.36，0.09，0.04。

图  1  Poulovassilis^[9]的分布函数测定方法示意图

Figure  1.  Illustration of distribution diagram from Poulovassilis^[9]

下载: 全尺寸图片 幻灯片

Poulovassilis等^[13]通过对砂土做SWCC干湿循环试验，发现孔隙水体积在浸润过程中不但包含独立的α类元素，而且还存在非独立的β类元素。依据该发现，Poulovassilis等^[13]提出了SWCC滞后曲线的非独立域模型。总体而言，域模型（包括独立域模型和非独立域模型）能够描述SWCC滞后现象，但无法预测SWCC的滞后曲线，对干湿循环过程土体中的水分迁移机理也缺乏较为深入的研究。

Bear^[24]和Klausner^[25]指出，SWCC的滞后现象主要由4个主要因素主导，包括：留截空气、干湿过程中接触角的差异性（即“雨滴”效应）、“墨水瓶”效应和触变性恢复。Klute^[26]和Zhai等^[22]采用图 2对SWCC滞后曲线参数进行相关定义。1-S₀定义截留空气量，表示浸润结束后残留在土体中空气的含量。大多数室内SWCC试验，干湿循环的次数是有限的，触变性恢复对SWCC滞后现象的影响通常不明显。因此，本文着重研究截留空气、接触角差异和“墨水瓶”效应对SWCC滞后特性的影响，并在本文物理-统计模型的基础上提出相应的计算方法。

图  2  刚性骨架土体的SWCC滞后示意图^{[22, 26]}

Figure  2.  Diagram of hysteresis in SWCC of soil with rigid structure

下载: 全尺寸图片 幻灯片

1.   理论

本节介绍了修正毛细模型的矩阵表达形式，该模型可以计算在脱湿和浸润过程中，不同吸力状态下，不同尺寸系中的水分分布状态。基于该矩阵模型，详细阐释留截空气、“雨滴”效应和“墨水瓶”效应对SWCC滞后特性的作用。最终，提出物理-统计模型用以预测非饱和土SWCC滞后曲线。

1.1   修正毛细模型的矩阵表达

为了将土体中的水分含量实现量化计算，土体中不同尺寸的孔隙被简化为一系列半径不等的毛细管。Childs等^[27-29]首先采用该模型将SWCC等效为孔径分布函数，用以解释水分在土体中的迁移规律。在传统“毛细管束”模型（bundle of cylindrical capillaries, BCC）中，管中毛细水高度代表土吸力的大小，毛细水高度和毛细管半径遵循毛细定律。Zhai等^[30]指出“毛细管束”模型无法准确描述非饱和土含水率随土吸力的变化规律。尽管毛细水高度可定义对应的土吸力，但管中水的体积不能代表存储在毛细管的水分量。例如，定义毛细管半径为r₁, r₂, …, r_N，毛细管水头高度为h_c1, h_c2, …, h_cN。完全饱和状态下，单位土体中水分含量（V_{w, sat}）可采用下式计算（即将不同尺寸毛细管中水柱体积的累加量）：

基于“毛细管束”模型，当吸力从h_ci增加至h_ci+1时，毛细管r_i中存储的水会全部排出，导致土体水分含量减少（ΔV_{w, i}），如下所示：

由式（2）可知，土体水分随土吸力的变化量只是孔隙半径的函数，与土的类别没有关系。很明显，这一结论与实际SWCC试验结果相悖。Zhai等^[30]还指出，传统“毛细管束”模型假定毛细管柱相互独立，不相互连接，导致毛细管柱中的水分无法流至其它毛细管柱中。因此，“毛细管束”模型无法解释水分在土体中的迁移规律，有必要对“毛细管束”模型提出必要的修正。

本文采用图 3所示的矩阵描述修正毛细管模型。图 3中的首列表示孔隙的尺寸r₁, r₂, r₃, …, r_N（其中，r₁ > r₂ > r₃ > … > r_N）以模拟土体中不同半径的孔隙；首行代表土吸力的大小，ψ₁, ψ₂, ψ₃, …, ψ_N（其中，ψ₁ < ψ₂ < ψ₃ < … < ψ_N），以模拟土体所处的吸力状态。在脱湿过程中，将每一列不同尺寸孔隙中的水分量累加起来，如末行所示，就能得到单位体积的相对水分含量（即饱和度）。

图  3  在干燥过程中，土体水分在不同尺寸孔隙中的分布状态示意图

Figure  3.  Illustration of distribution of water in pores with different radii in drying process

下载: 全尺寸图片 幻灯片

当土吸力足够小时（如吸力为图 3中ψ₁），所有孔径都被水分填充，第ψ₁列每个单元格中的数值代表水分量在该尺寸孔隙中的占比（全列水分占比累加量等于饱和度）。以此类推，第ψ₁列至ψ₂₀列、每个单元中的数值代表该孔径中水分含量占土体饱和水分含量的百分比。当土吸力值为ψ_i、依据毛细定律，半径大于r_i的孔隙中的水会被排出。随着吸力的增大，模型中的累积含水率会随之减少，这与脱湿SWCC是完全一致的。

因Fredlund等^[31]模型，如式（3），能够较好定义SWCC并已被广泛应用，本文采用该模型计算每个单元格中的水分含量。

式中：S_r为饱和度；C_r为对残余吸力的预估值，对于大部分砂土可取1500 kPa；a_f，n_f，m_f为模型参数，可以通过将式（3）和试验数据进行拟合得到。由此可知，参数a_f，n_f，m_f本身没有实际的物理意义。

图 3中第r_i行的水分含量如下式计算得到

当土体处于完全饱和状态时，水分完全填充土体中的所有孔隙，在土吸力增大的过程中，水分依次先从大的孔隙排出，再到小的孔隙并逐步排出土体。因此，图 3不仅描述了在脱湿过程水分在土体中的残留状态，也反映了土体不同尺寸孔隙的分布状态（即孔径分布函数）。这些不同尺寸的孔隙在土体中是随机相连的，连接概率取决于相邻孔隙的概率分布密度。当土体在浸润过程中，孔隙的空间连接状态（取决于不同孔径孔隙的连接概率）会直接影响水分对干燥孔隙的填充。

1.2   截留空气的量化计算

如图 3所示，当土吸力从ψ₁增加至ψ₁₁时，孔径r₁至r₁₀的孔隙中的水被完全排出。当土体开始浸润时，由于留截空气，当土吸力从ψ₁₁减小至ψ₁时，土体饱和度无法达到100%。在脱湿过程中，只有孔径r₁至r₁₀中的水分被排出，孔径r₁₁至r₂₀始终保持被水填充的状态。因此，留截空气只可能残留在孔径r₁至r₁₀的孔隙中。留截空气在r₁至r₁₀的孔隙中的分布形态可能会有多种，其中均匀分布最为简单、直观，也易于计算。作者采用这一基本假定，发现计算结果与文献中的试验数据吻合度较高，因此该假定较为合理。因此，本文采用均匀分布这一基本假定。如果浸润过程结束测定土体饱和度为S₀，则截留空气的单位体积可表达为(1-S₀)，单元体中孔径r₁至r₁₀的孔隙总体积为$1 - {S_{\text{r}}}({\psi _{10}})$。据此，留截空气在孔径r₁至r₁₀的孔隙残留比例（P_af）可用下式表示：

换言之，在浸润过程中，孔径为r₁至r₁₀的孔隙体积被水分填充的比例（P_wf）可用下式表示：

式中：${S_{\text{r}}}({\psi _{10}})$为土吸力为ψ₁₀时对应的饱和度，$1 - {S_{\text{r}}}({\psi _{10}})$为当土吸力从ψ₁增加至ψ₁₀时从单位土体中排出的总水量（即孔径r₁至r₁₀的孔隙体积总和）；${S_0} - {S_{\text{r}}}({\psi _{10}})$为当土吸力从ψ₁₀增加至ψ₁时填充至r₁至r₁₀的孔隙中单位水的体积。

在浸润过程中，考虑留截空气在干燥孔隙中均匀分布，留截空气量在不同尺寸孔隙中残留量可用该尺寸的孔隙体积乘以P_af进行计算。换言之，在浸润过程中，水分能够填充干燥孔隙的体积可用该尺寸的孔隙体积乘以P_wf进行计算。式（5）表明残留在土体中的留截空气总量与土样中干燥孔隙量及其分布有关。换言之，浸润过程初始吸力越大（即干燥孔隙孔径分布越广），留截空气总量越大，反之越小。

1.3   “雨滴”效应对SWCC滞后影响的量化计算

浸润过程水压驱使水气分界面移动，弯液面与土颗粒表面形成的接触角为推进接触角（θ_advancing）。脱湿过程气体驱使水气分界面移动，弯液面与土颗粒表面形成的接触角为后退接触角（θ_receding）。一般情况下，推进接触角会大于后退接触角，该现象通常被描述为“雨滴”效应。Zhai等^[22]定义浸润过程相对应接触角为α_wet，脱湿过程相对应接触角为α_dry。因接触角的直接测量难度较大（尤其对颗粒粒径分布较广的土体），Zhai等^[22]提出采用常数k来描述接触角α_dry和α_wet的差异性，如下式所示：

当k=1时，脱湿和浸润过程的接触角相等。Bechmann等^[32]、Eral等^[33]回顾了当前有关推进接触角和后退接触角的相关理论和试验结果，发现当固有接触角相同时，推进接触角大于后退接触角，因此在大部分工况下k值是大于1的。另外，颗粒的矿质成分、颗粒大小以及表面形态都可能影响接触角的测量结果。因此，接触角的试验测量难度较大，采用数学模型描述脱湿和浸润过程的接触角差异性可能是一个较为合理的选择。当土吸力从${\psi _i}$增加到${\psi _{i + 1}}$时，水分从孔径${r_i}$的孔隙排出，而当土吸力从${\psi _{i + 1}}$减少到${\psi _i}$时，因${\alpha _{{\text{dry}}}}$和${\alpha _{{\text{wet}}}}$的差异性，水分无法填充孔径为r_i的孔隙。脱湿和浸润过程，弯液面在毛细管的极限平衡状态可以用下式表示：

式中：T_s为水的表面张力。

当土吸力大于${\psi _{{\text{dry}}}}$时，平衡态被打破，水分排出毛细管；当土吸力小于${\psi _{{\text{wet}}}}$时，平衡态才能被打破，水分填充毛细管。从式（8），（9）可知，对于相同孔径的毛细管，因${\alpha _{{\text{wet}}}}$大于${\alpha _{{\text{dry}}}}$，排水的驱动土吸力会大于浸水的驱动土吸力。换言之，如果排水的驱动土吸力为${\psi _i}$，只有当土吸力降至${\psi _i}$/k时，浸润过程的平衡状态才能被打破，水分才可以填充孔径为${r_i}$的孔隙。

由此可见，“雨滴”效应对SWCC滞后特性的影响主要体现在水分对孔隙填充的延缓作用上，可用常数k进行计算。Zhai等^[22]建议采用浸润过程的部分试验数据校验k值，校验的k值可用以预测其他浸润曲线。本文沿用了Zhai等^[22]这一方法确定k值。

1.4   “墨水瓶效应”对SWCC滞后影响的量化计算

Taylor^[34]用毛细管模型（如图 4），形象地描述了“墨水瓶”效应。依据开尔文毛细定律，图 4（a）中毛细高度h_c1可通过毛细管半径r₁计算求得。当土吸力超过毛细高度h_c1时，水从半径r₁的管A中排出。如果土吸力减小，水将再次注入管A。可能出现如图 4（b）所示的场景，较大孔径的毛细管B连接在A1和A2之间。尽管土吸力减少至可以使水分填充半径r₁的毛细管，但因管B阻碍了水分迁移路径，水分无法填充管A2。如果土吸力持续减小，直至水完全填充管B，管B中的水分随即填充管A2（如图 4（c）），则“墨水瓶”效应将消失。由此可见，“墨水瓶”效应是动态的，不同于留截空气。另外，如图 4（b）所示，“墨水瓶”效应对SWCC滞后特性的影响取决于管B和管A2连接的概率，其连接概率又依赖于管B和管A的孔径分布密度f(r)。

图  4  “墨水瓶”效应示意图^[34]

Figure  4.  Illustration of phenomenon of "ink-bottle" effect^[34]

下载: 全尺寸图片 幻灯片

当吸力从ψ₁₁减小至ψ₁₀/k时，水分可以填充孔径为r₁₀的孔隙。孔径r₁至r₉的孔隙处于干燥状态，这次干燥孔隙可能阻断部分水分的迁移路径，导致水分无法100%填充半径为r₁₀的孔隙。阻塞的概率（P_blocking）取决于干燥孔隙的与目标孔隙连接的概率，可采用下式计算：

因此，考虑“墨水瓶”效应，孔径为r₁₀的孔隙体积被水分填充的可能性可由下式确定：

式中：${S_{\text{r}}}({\psi _i})$为土吸力为${\psi _i}$时对应的饱和度。

“墨水瓶”效应对SWCC滞后特性的影响主要体现在干燥的大孔隙可能阻挡水分填充目标孔隙。当这些大孔隙也被水填充时，“墨水瓶”效应也随之消失。“墨水瓶”效应对SWCC滞后特性的影响可采用式（11）进行量化计算。式（11）表明，“墨水瓶”效应对SWCC滞后特性的影响依赖浸润过程的土吸力大小。土吸力越小，“墨水瓶”效应越低。如果土吸力足够小（即$S({\psi _i})$接近1），“墨水瓶”效应可以完全消除。

1.5   考虑留截空气、“雨滴效应”、“墨水瓶效应”计算浸润过程不同尺寸孔隙中的水分量

为了精确计算不同孔径的孔隙中水分分布状态，需综合考虑留截空气、“雨滴”效应和“墨水瓶”效应，如图 5所示。白色区域代表留截空气所占的体积，其体积在后续脱湿和浸润过程中都不会产生改变，可用式（5）进行计算。灰色区域表示“墨水瓶”效应，因干燥的大尺寸孔隙阻挡部分水分的迁移路径，导致水分无法完全填充目标孔隙。当土吸力减小足够小时（水分能够填充干燥的大尺寸孔隙），水分便可完全填充目标孔隙，最终“墨水瓶”效应得以消除。“墨水瓶”对水分在孔隙中分布的影响可用式（11）进行计算。另外，因“雨滴”效应导致水分对孔隙填充的延缓，无法在图 5中具体表达。在浸润过程，需要一个更低的驱动土吸力才能将水分填充目标孔隙。浸润过程的驱动土吸力可用式（7）进行计算。

图  5  润湿过程中特定孔隙不同比例示意图

Figure  5.  Illustration of different proportions of specific pore in wetting process

下载: 全尺寸图片 幻灯片

最终，浸润过程中，不同尺寸孔隙中的水分残余量可进行图 5所示的分类计算，并可得到不同吸力状态的水分分布，如图 6所示。当土体从吸力ψ₁开始脱湿，增至ψ₁₁时，孔隙r₁至r₁₀中水分被排出，统计不同吸力状态对应土体的残余含水率，可以得到初始脱湿曲线。当土体从吸力ψ₁₁开始浸润时，降至ψ₃时，孔径为r₃的孔隙因“雨滴”效应，无法被水分填充，因留截空气和“墨水瓶”效应，被水分填充的孔隙也并非百分百被填充，具体填充量可采用式（5），（7）进行计算。当土体从吸力ψ₃开始再次脱湿，此时，孔径r₃至r₁₀的孔隙并未完全被水填充，因此排出的水分量会小于初始脱湿过程的排水量。本文提出的物理-统计模型可以准确描述初始脱湿曲线，并依据初始脱湿曲线测算出浸润主线。模型中的k值可以通过浸润主线的试验数据进行校验，校验后的k值可用以预测脱湿扫描曲线和浸润扫描曲线。

图  6  湿润过程中孔隙中的储水条件示意图

Figure  6.  Illustration of water storage condition in pores in wetting process

下载: 全尺寸图片 幻灯片

2.   试验数据验证

为了验证本文提出的统计模型，笔者查阅了国内外相关文献，并搜集了12组试验数据，包含不同类别土质的脱湿和浸润SWCC试验数据。其中国内文献数据6组包括：西安黄土（1）^[35]、西安黄土（2）^[36]、青藏黏土^[36]、全风化凝灰岩^[37]、泾阳黄土^[38]、河南南阳原状黄土^[39]；国外文献数据6组包括：砂-高岭土混合物SK-5^[5]（砂-高岭土比例15∶85）^[5]、SK-10（砂-高岭土比例35∶65）^[5]、SK-17（砂-高岭土比例55∶45）^[5]、砂土Ⅳ^[5]、砂土Ⅴ^[5]、砂^[10]。调研发现国内文献对二次脱湿的SWCC试验数据相对较少。因此，选取的国外6组试验数据包括初次脱湿、首次浸润和二次脱湿的相关试验数据。

先采用Fredlund等^[31]的方程对初次脱湿的SWCC试验数据进行拟合，得到SWCC模型参数（a_f，n_f，m_f），列于表 1。将模型参数代入式（4），构建与之对应的矩阵模型，以便掌握土体的孔径分布。后续采用式（6），（7），（11）进行计算首次浸润的SWCC曲线。因接触角的相关参数k的取值有待进一步研究，本文将式（6），（7），（9）的计算结果与首次浸润试验数据进行校验，以确定k值。最后，利用校验的k值，预测二次脱湿的SWCC曲线，模型计算结果和试验数据如图 7所示。

表  1  土体土-水特征曲线的拟合参数和k值

Table  1.  Fitting parameters of SWCCs and values of k for soils

序号	土体名称	初始干密度${\rho _{\text{d}}}$/(g·cm-3)	围压/kPa	af /kPa	nf	mf	Cr/kPa	k
(a)	西安黄土（1）	1.51	0	9.02	1.33	0.41	1500	1.72
(b)	西安黄土（2）	1.39	0	32.69	1.02	0.45	1500	1.79
(c)	青藏黏土	1.46	0	2.91	1.23	0.36	1500	2.01
(d)	全风化凝灰岩	—	—	191.82	1.91	2.19	1500	1.74
(e)	泾阳黄土	1.23	0	8.32	1.19	0.56	1500	1.00
(f)	河南南阳原状黄土	1.70	0	3000	0.01	0.20	1500	2.25
(g)	SK-5	1.50	0	67.78	1.32	0.78	1500	1.47
(h)	SK-10	1.67	0	72.24	1.72	0.29	1500	2.30
(i)	SK-17	1.86	0	65.39	1.29	0.72	1500	2.30
(k)	Sand Ⅳ	1.77	0	3.47	5.34	1.46	1500	1.87
(l)	Sand Ⅴ	1.77	0	2.79	6.52	1.20	1500	1.55
(m)	Sand	—	0	2.83	5.97	1.73	1500	1.52

下载: 导出CSV 

| 显示表格

图  7  本文提出的统计方法对土体的预测湿润和干燥扫描曲线比较

Figure  7.  Comparison between predicted wetting and drying scanning curves for soils using proposed statistical method

下载: 全尺寸图片 幻灯片

图 7（a）~（f）表明，本文模型可以较好测算国内非饱和土的浸润曲线，对比模型结果与试验数据，发现决定系数R²基本大于95%。模型的测算结果精度在很大程度上依赖参数k的取值，可采用浸润试验部分数据校验k值。另外，需要指出，本文所提出的理论和计算式都是基于土骨架是不变形的（即孔径分布函数在脱湿和浸润的过程中是恒定不变的）这一基本假定。如果孔径分布在脱湿和浸润过程中产生变化，需要对本文理论做必要修正，该修正工作位子啊本文研究范畴之内。

图 7（g）~（l）表明校验的k值，可以准确预测后续的脱湿扫描曲线。对比模型结果与试验数据，对比的决定系数R²都大于99%。因此，本文提出的物理-统计模型可有有效评估非饱和土SWCC的滞后特性。

3.   结语

本文提出了一种采用矩阵形式的修正毛细模型，该模型可以准确描述水分在土体中不同尺寸孔隙中的分布状态。结合统计学理论，量化计算留截空气、“雨滴”效应和“墨水瓶”效应对SWCC滞后特性的作用。基于矩阵模型，计算在浸润过程，水分在不同尺寸孔隙中的残余量。最终提出测算SWCC滞后曲线的物理-统计模型。通过对比模型计算结果和国内外文献的试验数据，认定本文模型的理论较为完善，预测结果是可靠的。另外，考虑体积变形的SWCC滞后模型有待进一步的研究。