Analysis Model of Geotechnical Materials Based on Relative Humidity Field and Variable Modulus Constitutive Model
-
摘要: 随着含水率（湿度）的增加，红层黏土（岩）会表现出体积膨胀和吸水软化的双重特征，这是诱发红层滑坡和路基隆起等地质灾害的主要原因。为研究含水率对受荷岩土体变形的影响，本文在湿度应力理论的基础上，对含水率归一化处理，提出基于相对湿度的变模量本构，建立岩土膨胀分析模型；以固结黏土和红层黏土岩为研究对象，开展膨胀和吸水软化试验。结果表明：膨胀试验中，黏土岩和固结土的侧限膨胀率分别为0.345和0.205，膨胀应力分别为225.0和112.5 kPa，黏土岩具有更强的膨胀能力；固结土吸水效率更高，能够更快地达到膨胀平衡状态。软化试验中，干燥状态下黏土岩的压缩模量较大，但吸水过程中模量的衰减也较快，饱和状态下固结土与黏土岩的压缩模量基本相同。基于试验数据，对本文模型的稳态参数进行反演，并基于已标定的分析模型对侧限膨胀试验进行模拟，一方面，验证本文方法在稳态分析时的适用性；另一方面，给出了相对湿度场平均扩散率K的标定方法。以圆形洞室遇水问题为算例，进行了稳态理论计算和瞬态数值计算，以说明本文模型在应用于岩土膨胀大变形分析时的优越性。本文分析模型参数易获得，数值模拟结果准确性高，可为红层膨胀岩和膨胀土地区工程建设提供理论指导。Abstract: With the increase of water content (humidity), red bed clay (rock) shows the characteristics of volume expansion and water-softening simultaneously, which is the main reason for geological hazards, such as landslides and subgrade uplifts in red bed regions. To understand the influence of water content on the deformation regulation of loaded geo-materials, a water-softening model based on relative humidity (normalized water content) was proposed and a corresponding expansive analysis model was established. Firstly, taking consolidated clay and red bed clay rock as research objects, expansion, and water-softening tests were carried out to calibrate the steady parameters of the model. The expansive tests show that the confining expansion rate of clay rock and consolidated soil are 0.345 and 0.205, and the expansion stress is 225.0 kPa and 112.5 kPa respectively, which reflected that the clay rock has a stronger expansion capacity than the consolidated clay. Meanwhile, the test also shows that the consolidated soil absorbs water more quickly and can reach the expansion equilibrium state faster. The water-softening test shows that the compression modulus of clay rock is larger than consolidated clay in the dry state, and the decrease of compression modulus with water content is also faster. The compression modulus of consolidated clay and clay rock is basically the same in the saturated state. Then, the confining expansion test was numerically simulated based on the calibrated analysis model, and the applicability of the method in terms of steady analysis was verified. Meanwhile, the average diffusivity K (which is a transient parameter) was calibrated by the comparison of tests and numerical results. Finally, taking the problem of circular cavern encountering water as an example, the steady theoretical analysis and the transient numerical calculation were carried out to illustrate the superiority of the model in expansive large deformation conditions. In summary, the parameters of the novel model are easy to obtain, and the numerical result is accord with tests, which can guide the engineering practice in the red bed expansive clay or swelling rock areas.
-
中国红层地区膨胀土、膨胀岩的分布极为广泛，吸水后产生的膨胀变形对工程结构的安全十分不利，特别是对于边坡、路基和隧道等交通工程，由岩体膨胀导致的地质灾害较多[1–3]。
岩土材料的吸水膨胀是一个十分复杂的流固耦合过程，与矿物成分、应力状态及浸水条件等因素均有关[4–6]。目前，从宏观角度直接建立的膨胀模型都是基于一些数学上的假设[7]，如：经典的1维Huder–Amberg膨胀模型假设膨胀率与膨胀压力的对数呈线性关系[8]；Gysel[9]提出的膨胀模型则假设膨胀力与吸水率服从指数关系。这一类模型属于半经验模型，模型参数的确定一般需要数据拟合，缺乏理论基础。缪协兴等[10–11]借鉴了热弹性力学的分析方法，提出湿度应力理论，并给出了控制方程组。白冰等[12–13]给出了湿度应力场的理论证明，认为该理论在力学模型上是完备的；同时，也指出湿度应力场虽然是一种合理的理论模型，但需要根据实际情况进行必要的修正和完善。为此，在湿度应力理论框架下，研究人员开展了大量拓展。如：陈有亮等[14]推导了考虑膨胀应力和剪胀效应的深埋隧道弹塑性解；陈钒等[15]建立了考虑时间效应的膨胀模型；任松等[16]得到了考虑时间效应的隧道弹–膨胀解析模型。然而，为方便理论推导，现有研究多假设固体相的变形是弹性或者弹塑性的，这一假设忽略了土水或水岩作用下岩土材料的软化效应，是一个明显的不足。且从矿物成分的角度来看，膨胀岩、膨胀土多含有蒙脱石和伊利石等亲水性黏土矿物，一般也有较高的软化系数，遇水后变形模量会持续下降[17–19]。对于有侧限或围压的情况，完全软化状态的岩土体变形模量一般会趋于定值[20]。吴建勋等[21]研究了围岩软化对隧道围岩自稳时间的影响，认为围岩软化作用不能忽略。综上所述，在湿度应力的理论框架内，如果能够充分考虑变形模量随湿度的变化，就能够更好地为膨胀岩地区的地下工程提供指导。
本文以压力固结的红层黏土（后文简称固结土）和红层黏土岩（后文简称黏土岩）为研究对象，首先，开展受荷膨胀试验，通过侧限膨胀率，确定变形模量随含水率的变化规律；然后，基于归一化的相对湿度建立了变模量本构，给出模型参数的确定方法，并通过室内试验和相应数值模拟的对比，验证了本文膨胀分析模型的可靠性；最后，给出了湿度平均扩散系数的标定方法，以圆形洞室遇水为算例，进行理论计算和数值模拟，讨论了本文模拟方法的优越性。
1. 湿度应力场的改进
1.1 湿度应力场的控制方程
湿度应力理论是对复杂水–岩土体物理膨胀问题的一种简化研究方法，借鉴了热弹性问题的经典理论。湿度应力场的控制方程包含湿度场的扩散方程，应力场的平衡微分方程、几何方程，以及材料的本构方程[10–12]。
岩土体中的水分迁移过程十分复杂，与土或岩石的渗透性、饱和度、孔隙结构乃至应力状态等因素都密切相关，湿度场的引入是对水分迁移过程的一种简化。湿度是含水率（或含水量）的表征。根据质量守恒，单位体积含水量的增加应与水分的净流入量相等；而根据达西定律，渗流速度与湿度梯度成正比，那么湿度场的控制方程实际是一个无源扩散方程，各向同性条件下可以表示为：
$$ {\qquad \frac{{\partial W}}{{\partial t}} - \kappa \left( {\frac{{{\partial ^2}W}}{{\partial {x^2}}} + \frac{{{\partial ^2}W}}{{\partial {y^2}}} + \frac{{{\partial ^2}W}}{{\partial {{\textit{z}}^2}}}} \right) = 0 }$$ (1) 式中：W为湿度场，对应的场变量为含水率w，是坐标x、y、z和时间t的函数；κ为湿度扩散率，在各向同性条件下为常数。
应力场的平衡方程、几何方程为：
$$ \frac{{\partial {\sigma _{ij}}}}{{\partial {x_i}}} + {f_j} = 0 $$ (2) $$ {\varepsilon _{ij}} = \frac{1}{2}\left( {\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}}} \right) $$ (3) 式中，εij为应变分量，σij为应力分量，ui和uj为位移分量，xi和xj为坐标分量，fi为体力分量。
湿度应力场与固体力学场的区别主要表现在本构模型方面，除荷载外，湿度的变化也会引起材料的体变，无约束状态时[10]：
$$ {\varepsilon _{ij}} = \alpha {\delta _{ij}}W $$ (4) 式中：α为膨胀系数；δij为克罗内克记号，仅当i=j时取1，其余取0。
有约束状态时，湿度场的变化会产生附加应力，附加应力与应变间满足广义虎克定律，总的本构模型就可以描述为[10]：
$$ {\qquad {\varepsilon _{ij}} = \frac{1}{E}\left[ {\left( {1 + \nu } \right){\sigma _{ij}} - \nu {\delta _{ij}}\Theta } \right] + \alpha {\delta _{ij}}\Delta W }$$ (5) 式中， Θ为应力张量的第1不变量，E、
$\nu $ 为弹性模量和泊松比。湿度应力场理论的数学实质是在初边界条件下求解由式（1）～（3）和式（5）组成的偏微分方程组。
1.2 基于相对湿度的扩散方程
湿度场与温度场具有相似性，但也存在一些实质性的区别，如：湿度必然有一个上限指标，即饱和含水率。当岩土体达到饱和后，不会自发地继续吸水，湿度也不会继续增加。同时，岩土体会固结排水，湿度也不仅仅是空间位置和时间的函数，还是应力状态的函数，严格地说，在湿度应力理论中，湿度场和应力场是双向耦合的。为简化耦合模式，本文对含水率行归一化处理，引入相对湿度的概念，将湿度场转化为相对湿度场来描述水分在岩土体的运移规律，相对湿度ϖ 的定义为：
$$ \varpi = \frac{w}{{{w_{{\text{sat}}}}}} $$ (6) 式中，w为含水率，wsat为饱和含水率。对于相对湿度，采用体积含水率和质量含水率是等效的。
相对湿度场的控制方程仍可简化为一个无源的扩散方程，即：
$$ {\qquad \frac{{\partial \varpi }}{{\partial t}} - K\left( {\frac{{{\partial ^2}\varpi }}{{\partial {x^2}}} + \frac{{{\partial ^2}\varpi }}{{\partial {y^2}}} + \frac{{{\partial ^2}\varpi }}{{\partial {{\textit{z}}^2}}}} \right) = 0} $$ (7) 式中，K为相对湿度的扩散率，是求解瞬态问题的基础，理论上可通过非饱和土的渗流试验进行标定[22–23]。
1.3 考虑软化效应的本构方程
文献[8]在提出湿度应力场时就强调了岩石软化问题的重要性，即在考虑岩土体吸水膨胀的同时，不应忽略其弹性模量的降低。因此，式（5）中的弹性模量应当是含水率的函数，进而也是相对湿度的函数，而非一个固定的材料参数，即：
$$ E = E\left( w \right) = E\left( \varpi \right) $$ (8) 那么，式（5）也不再是线性的偏微分方程，应力的求解难度也将显著增加。
综上所述，基于相对湿度并考虑软化效应的湿度应力场瞬态分析模型控制方程为：
$$ \left\{ {\begin{array}{l} {\dfrac{{\partial \varpi }}{{\partial t}} - K\left( {\dfrac{{{\partial ^2}\varpi }}{{\partial {x^2}}} + \dfrac{{{\partial ^2}\varpi }}{{\partial {y^2}}} + \dfrac{{{\partial ^2}\varpi }}{{\partial {{\textit{z}}^2}}}} \right) = 0,} \\ {{\varepsilon _{ij}} = \dfrac{1}{{E\left( \varpi \right)}}\left[ {\left( {1 + \nu } \right){\sigma _{ij}} - \nu {\delta _{ij}}\Theta } \right] + \alpha {\delta _{ij}}\Delta \varpi ,} \\ {{\varepsilon _{ij}} = \dfrac{1}{2}\left( {\dfrac{{\partial {u_i}}}{{\partial {x_j}}} + \dfrac{{\partial {u_j}}}{{\partial {x_i}}}} \right),} \\ {\dfrac{{\partial {\sigma _{ij}}}}{{\partial {x_i}}} + {f_j} = 0} \end{array}} \right. $$ (9) 由式（9）可知，相对湿度扩散率K、膨胀系数α的取值和弹性模量E与相对湿度ϖ间的函数关系E(ϖ)都需要通过试验确定。对于稳态问题，不考虑具体的湿度场变化过程，不涉及湿度扩散方程，那么可以不对相对湿度扩散率K进行准确标定[10–11,15]。
2. 湿度应力场稳态参数的确定
以采自湖南株洲的红层黏土（固结土，压实系数大于0.95，编号A）和采自甘肃天水的红层黏土岩（编号B）为研究对象，开展物理特性试验，结果见表1。
表 1 土样的基础物理参数Table 1 Basic physical parameters of samples土样
编号液限 塑限 缩限 含水率
/%孔隙率
/%干密度
/(g·cm–3)A 39.3 13.5 11.2 12.6 16.8 1.75 B 32.5 12.7 9.2 10.7 15.1 1.62 由表1可知，固结土和黏土岩的基础物理参数较为接近，且均处于硬塑状态，具有较强的可比性。
2.1 膨胀力及侧限膨胀试验
膨胀试验通过WG型三联高压固结仪完成，试验装置测得土样A、B的膨胀压力时程曲线如图1所示。
图 1 室内试验得到的土样膨胀压力时程曲线Fig. 1 History curves of expansive pressure in laboratory test
下载:
全尺寸图片
由图1可知：固结土的吸水速度较快，在30 min时膨胀压力已达到稳定；而黏土岩的膨胀压力在45 min时才达到稳定。分别选取试验中位值112.5 kPa和225.0 kPa作为固结土和黏土岩的膨胀应力进行后续分析。
在侧限膨胀试验中，以膨胀应力为参考值，固结土组的法向压力分别设定为0、12.5、25.0、50.0、75.0和100.0 kP；黏土岩组的法向压力分别设定为0、25、50、100、150和200 kPa，各进行3次平行试验，结果见表2和3。
表 2 固结土侧限膨胀率Table 2 Confined expansive rates of consolidated clay法向压力
σn/kPa侧限膨胀率 第1组 第2组 第3组 中位值 0 0.205 0.200 0.220 0.205 12.5 0.166 0.175 0.160 0.166 25.0 0.135 0.130 0.130 0.130 50.0 0.085 0.085 0.095 0.085 75.0 0.055 0.045 0.045 0.045 100.0 0.010 0.010 0.005 0.010 表 3 黏土岩侧限膨胀率Table 3 Confined expansive rates of clay rock法向压力
σn/kPa侧限膨胀率 第1组 第2组 第3组 中位值 0 0.330 0.355 0.345 0.345 25 0.290 0.275 0.280 0.280 50 0.210 0.225 0.235 0.225 100 0.130 0.140 0.130 0.130 150 0.055 0.055 0.050 0.055 200 0.010 0.015 0.015 0.015 2.2 稳态参数的标定方法
假设膨胀变形结束后试样已处于饱和状态，引入相对湿度后，膨胀系数α与无荷载作用时的侧限膨胀率δep0为同一参数的不同表述。膨胀系数可以通过侧限膨胀试验直接得到，固结土和黏土岩的膨胀系数分别为0.205和0.345。
基于文献[24]的修正Huder–Amberg模型，对表2中的数据进行拟合，得到侧限膨胀率ε与法向压力间的函数关系如式（10）所示，曲线如图2所示。
图 2 侧限膨胀率−法向压力曲线Fig. 2 Confined expansive ratio−normal stress curves下载:
全尺寸图片
$$ \varepsilon = {\varepsilon _{\text{b}}}\left[ {1 - \frac{{\lg \left( {{\sigma _{\text{n}}} + {\sigma _{\text{e}}}} \right)}}{{\lg \left( {{\sigma _0} + {\sigma _{\text{e}}}} \right)}}} \right] $$ (10) 式中：εb为膨胀参数；σ0为膨胀应力；σn为施加的法向压力；σe为侧限膨胀试验的预加荷载，可取为σe=10 kPa。
由图2可知，试验值与理论曲线的一致性较好，可以对不同法向压力对应的侧限膨胀率进行预测。
饱和含水率可以通过初始孔隙率和侧限膨胀率按式（11）计算：
$$ {w_{{\text{sat}}}} = \frac{{{e_0} + {\delta _{{\text{ep}}}}}}{{1 + {\delta _{{\text{ep}}}}}} $$ (11) 式中，e0为初始孔隙率，δep为侧限膨胀率。
根据式（11）可以在法向压力–饱和含水率坐标系内得到一系列散点，通过线性拟合建立两者的函数关系，如图3所示。由于在相同的法向压力条件下黏土岩具有更大的膨胀变形，因此也对应更大的饱和含水率。
图 3 法向压力–饱和含水率曲线Fig. 3 Normal stress−saturated water content curves下载:
全尺寸图片
根据表2的数据，将不同法向压力状态视为饱和试样侧限排水压缩时的中间状态，按式（12）反演得到不同含水率阶段时的弹性模量：
$$ E = - \frac{{\Delta {\sigma _{\text{n}}}}}{{\Delta {\delta _{{\text{ep}}}}}} $$ (12) 式中，Δσn为侧限膨胀试验中法向压力的增量，Δδep为对应的侧限膨胀率增量。
对图3中的饱和含水率指标进行归一化处理，可以计算得到不同相对湿度下的压缩模量，如图4所示。侧限膨胀试验的相对湿度变化范围较小，因此本文还补充了塑限、缩限和干燥状态下的侧限压缩试验。
图 4 弹性模量–相对含水率曲线Fig. 4 Elastic modulus-relative humidity curves下载:
全尺寸图片
假设弹性模量与相对湿度呈线性关系，那么基于图4中的数据即可建立变模量本构如式（13）所示：
$$ E = {E_0}\left( {1 - \varpi } \right) + {E_1} $$ (13) 式中，E1为饱和弹性模量，E0为饱和过程弹性模量的减小量，E0+E1为干燥弹性模量。
由图4可知：相对湿度相同时，黏土岩的压缩模量大于固结土；但随着含水率的增加，其压缩模量衰减速度也更快，在饱和时的压缩模量基本相同。
3. 模型可靠性验证
由于相对湿度场与温度场的控制方程具有相同的形式，在数值模拟过程中可以通过在0～1范围内变化的温度场来等效模拟，并定义弹性模量随温度变化实现变模量本构。
3.1 验证性模型
通过膨胀应力和侧限膨胀率对本文模型的可靠性进行验证，建立的几何模型如图5所示。
图 5 数值模拟中的几何模型Fig. 5 Geometric models in numerical simulation下载:
全尺寸图片
如图5所示，侧限膨胀模型包括试样和试验装置（环刀、垫片）两部分，其中，试样采用温变弹塑性材料进行模拟，试验装置采用刚体材料进行模拟。模型的环刀部分节点均为固定节点，试样底部施加法向位移约束。试样和环刀间不共节点，采用无摩擦接触；垫片底面与试样顶面也不共节点，采用点–面双向接触。在模拟侧限膨胀试验时，在垫片上表面施加均布压力，并对垫片法向位移数据进行监测和输出；在模拟膨胀应力试验时，对垫片顶部节点进行固定约束，监测试样与垫片间的接触力并输出。
初始时刻，将模型的相对湿度设置为0；定义上、下表面为恒定湿度边界，相对湿度设置为1；模型采用显式方法按瞬态问题进行求解，求解器为LS-solver（R11.0）。由于膨胀应力和侧限膨胀率均为稳态指标，因此并不需要准确的相对湿度扩散率Κ。为加速求解，计算时将相对湿度扩散率设置为10 m2/s。数值模拟得到无荷载的侧限膨胀率时程曲线如图6所示，膨胀压力时程曲线如图7所示。由于采用了放大的相对湿度扩散率，图6和7中横轴T对应数值计算时的虚拟时长，而不是实际时间t。
图 6 无荷载侧限膨胀率时程曲线Fig. 6 History curves of confined expansive rates下载:
全尺寸图片
图 7 膨胀压力时程曲线Fig. 7 History curves of expansive pressure下载:
全尺寸图片
如图6和7可知：当虚拟时长达到0.14 s时，膨胀应力模型的相对湿度达到稳定；而在膨胀率模型中，黏土岩和固结土相对湿度达到稳定的虚拟历时分别为0.16和0.19 s，略长于膨胀应力模型。这是由于模型采用了相同的相对扩散率，而膨胀率模型中体积膨胀增加了湿度扩散的距离。数值模拟得到的膨胀应力和侧限膨胀率与室内试验结果在数值上基本一致，也验证了本文方法在求解稳态问题时的可靠性。
3.2 相对湿度扩散率的标定
根据湿度场控制方程（7），增大湿度扩散率后，虚拟时长T和实际的渗透时间t之间的变换是线性的，且缩小的比例因子与湿度扩散率的放大比例相同。将图6的数值曲线与室内试验结果对比，通过相对湿度场达到稳定的历时来间接地标定相对湿度扩散率Κ。以固结土为例，膨胀压力稳定历时为30 min，而模拟中对应的虚拟历时为0.14 s，相对湿度扩散率放大比例为12 857，那么实际湿度扩散率应当为预设值的1/12 857，即7.8×10–4 m2/s；同理，可以计算得到黏土岩的相对湿度扩散率为6.0×10–4 m2/s。严格地说，本方法标定的相对湿度场湿度扩散率Κ是湿度场的平均扩散率，而非瞬态值。
3.3 模型的验证
除无荷载侧限膨胀试验外，本文进行了有法向压力的侧限膨胀数值模拟，得到的受荷侧限膨胀率时程曲线如图8所示，其中横坐标为已转换为实际时间t。
图 8 有荷载的侧限膨胀率时程曲线Fig. 8 History curves of loaded confined expansive rates下载:
全尺寸图片
由图8可知，在有法向荷载作用时，数值模拟得到的稳定侧限膨胀率在数值上与室内试验的一致性仍较高，但模拟的侧限膨胀率时程曲线与试验结果存在一定差异，主要表现为膨胀变形稳定前实测的膨胀率大于模拟值。这是由于试验对象为非饱和的黏土（岩），其湿度扩散率与含水率密切相关，随含水率的增加，湿度扩散率减小[25]，而模拟过程中并未考虑这一问题。然而，由平均扩散率引起的差异较小，且会随时间的增加而持续减小，对模拟结果的影响并不明显。
4. 算 例
4.1 稳态理论解
文献[26]给出了用湿度应力理论求解圆形洞室遇水的算例，简化后的力学模型如图9（a）所示。图9（a）中：ρ为岩土体的密度；H为埋深；a为隧道半径；b为影响区范围，取b=5a；σr=a、σr=b分别为隧道壁和影响区边缘的应力边界条件，取σr=a=0，σr=b=ρgH。
图 9 圆形洞室模型及稳态理论解Fig. 9 Model of circular chambers with continuous internal water erosion and its theoretical solution下载:
全尺寸图片
根据弹性力学小孔问题，考虑5倍半径区域的影响，并忽略重力梯度，采用均布的压应力ρgH代替，孔内侧为自由边界无荷载，并假设湿度场为线性分布，即：
$$ W = {W_0}\left( {a/r} \right) $$ (14) 式中，W0为洞室内边界湿度，a为洞室半径。
若采用相对湿度来描述，则为：
$$ \varpi = a/r $$ (15) 对应的线弹性理论解为[26]：
$$ \left\{ {\begin{array}{l} {{\sigma _r} = - \rho gH\left( {1 - \dfrac{{{a^2}}}{{{r^2}}}} \right) - \dfrac{{E\alpha \cdot a}}{{1 - \nu }} \cdot \dfrac{{r - a}}{{{r^2}}},} \\ {{\sigma _\theta } = - \rho gH\left( {1 + \dfrac{{{a^2}}}{{{r^2}}}} \right) - \dfrac{{E\alpha \cdot {a^2}}}{{\left( {1 - \nu } \right){r^2}}},} \\ {{\tau _{r\theta }} = 0} \end{array}} \right. $$ (16) 式中，σr、σθ为极坐标下的应力分量，约定以受拉为正，受压为负。
作为对比，本文采用相同的简化条件，并考虑由湿度场引起的软化效应，即，将式（13）代入式（16），得到对应的理论解为：
$$ \left\{ \begin{array}{l} {\sigma }_{r}=-\rho gH\left(1-\dfrac{{a}^{2}}{{r}^{2}}\right)-\left({E}_{1}+\dfrac{a}{r}{E}_{0}\right)\dfrac{\alpha \cdot a\left(r-a\right)}{\left(1-\nu \right){r}^{2}}, \\ {\sigma }_{\theta }=-\rho gH\left(1+\dfrac{{a}^{2}}{{r}^{2}}\right)-\left({E}_{1}+\dfrac{a}{r}{E}_{0}\right)\dfrac{\alpha \cdot {a}^{2}}{\left(1-\nu \right){r}^{2}}, \\ {\tau }_{r\theta }=0 \end{array} \right.$$ (17) 以固结土中半径a=2 m的圆形洞室为例，本文及文献[22]模型在5倍洞室半径的稳态理论解如图9（b）所示。对于式（16）的模型，为增加可比性，取弹性模量E=0.5[E1+(1/5)E0]，即本文模型变形模量变化范围内的中间值。
由图9可知，本文模型计算得到的洞室周边环向压应力分量较文献[22]的模型小，而在较远处环向应力分量则较大，这是由于本文模型考虑了围岩遇水软化引起的应力重分布。
4.2 瞬态分析
本文模型也可适用于瞬态问题的数值分析，瞬态分析模型可采用显–隐式方法求解，即：首先将相对湿度设置为0，施加应力边界条件σr=b=ρgH，通过隐式迭代法进行应力初始化；然后，定义内边界相对湿度为1，外边界相对湿度为0，并对应力初始化的模型进行显式的重启动求解。
围岩的材料参数在第2、3节中已经进行了标定，模拟中将渗透系数放大864 000倍进行求解，那么虚拟时长0.1 s即相当于真实时间1 d。仍取洞室半径a=2 m，埋深H=10 m。求解得到洞室在120 d内相对湿度的变化过程如图10所示。为便于比较，图10中给出了ϖ=0.167的参考线。
图 10 数值模拟得到的相对湿度场变化规律Fig. 10 Evolution of relative humidity field in numerical simulation下载:
全尺寸图片
由图10可知：随着时间的增加，固结土和黏土岩内相对湿度场的变化均趋缓，在第60天至第120天阶段，黏土岩内相对湿度ϖ =0.167参考线移动距离仅为0.2 m；相同阶段，固结土内参考线的移动距离已小于0.1 m，但均未达到理论上的稳态。可见，岩土体吸水膨胀问题实际是一个持续的过程，相对湿度场在短时间内无法达到稳态，本算例作为瞬态问题分析更具有工程意义。
第3、15、30、60和120天时的σr和σθ随r的变化规律如图11所示。
图 11 极坐标下应力分量的变化规律Fig. 11 Evolution of stress components in the polar coordinate system下载:
全尺寸图片
由于本文模型可以考虑围岩软化引起的应力重分布，在图11中，随着时间的增加，洞室周边的环向应力有所下降。对比固结土与黏土岩可以发现，黏土岩虽然具有更大的膨胀性，但由于其遇水软化明显，洞室周边岩体实际应力水平低于相同条件下的固结土洞室。且黏土岩内应力下降比例较高，在第3天至第60天阶段，环向最大压应力由385下降为331 kPa，下降比例为14.0%；而该阶段固结土内的环向最大压应力下降比例仅为9%，这是因为黏土岩软化引起的应力重分布更加明显。因此，进行膨胀岩（土）洞室稳定分析时，除围岩膨胀性外也应关注其软化性。由图11还可知，在围岩软化尚未稳定前，本文数值模拟得到的最大压应力已小于理论模型，这是由于式（17）的理论模型是基于小变形假设推导的，未考虑围岩大变形产生的卸载，计算应力偏大。
5. 结 论
本文模型通过引入归一化的相对湿度场和变模量本构，同时考虑了膨胀岩（土）的体积膨胀和吸水软化，相较于传统的湿度应力理论，能较好地体现围岩软化产生的应力重分布，计算结果更加合理。就模型的使用而言，本文模型参数标定过程简单，物理意义明确，数值模拟易于实现，具有较强的可推广性。
以典型的红层膨胀岩和膨胀土为研究对象，就膨胀岩遇水问题具体分析得到以下结论：
1）固结土和黏土岩的膨胀应力分别为112.5和225 kPa，无荷载侧限膨胀率分别为0.205和0.345，膨胀应力和变形均较大；就红层地区的地下结构遇水的问题而言，采用大变形分析更为合理。
2）相对湿度相同时，黏土岩变形模量大于固结土；但随着含水率的增加，黏土岩的变形模量衰减速度更快，软化现象更明显。因此，黏土岩洞室遇水后，由软化效应引起的应力重分布更加明显。
3）黏土岩和固结土的相对湿度场平均扩散率分别为6.0×10–4和7.8×10–4 m2/s；对于工程尺度的膨胀岩（土）遇水的问题，相对湿度场达到稳态的历时较长，采用瞬态分析更具工程价值。
-
图 1 室内试验得到的土样膨胀压力时程曲线
Fig. 1 History curves of expansive pressure in laboratory test
下载:
全尺寸图片
图 2 侧限膨胀率−法向压力曲线
Fig. 2 Confined expansive ratio−normal stress curves
下载:
全尺寸图片
图 3 法向压力–饱和含水率曲线
Fig. 3 Normal stress−saturated water content curves
下载:
全尺寸图片
图 4 弹性模量–相对含水率曲线
Fig. 4 Elastic modulus-relative humidity curves
下载:
全尺寸图片
图 5 数值模拟中的几何模型
Fig. 5 Geometric models in numerical simulation
下载:
全尺寸图片
图 6 无荷载侧限膨胀率时程曲线
Fig. 6 History curves of confined expansive rates
下载:
全尺寸图片
图 7 膨胀压力时程曲线
Fig. 7 History curves of expansive pressure
下载:
全尺寸图片
图 8 有荷载的侧限膨胀率时程曲线
Fig. 8 History curves of loaded confined expansive rates
下载:
全尺寸图片
图 9 圆形洞室模型及稳态理论解
Fig. 9 Model of circular chambers with continuous internal water erosion and its theoretical solution
下载:
全尺寸图片
图 10 数值模拟得到的相对湿度场变化规律
Fig. 10 Evolution of relative humidity field in numerical simulation
下载:
全尺寸图片
图 11 极坐标下应力分量的变化规律
Fig. 11 Evolution of stress components in the polar coordinate system
下载:
全尺寸图片
表 1 土样的基础物理参数
Table 1 Basic physical parameters of samples
土样
编号液限 塑限 缩限 含水率
/%孔隙率
/%干密度
/(g·cm–3)A 39.3 13.5 11.2 12.6 16.8 1.75 B 32.5 12.7 9.2 10.7 15.1 1.62 表 2 固结土侧限膨胀率
Table 2 Confined expansive rates of consolidated clay
法向压力
σn/kPa侧限膨胀率 第1组 第2组 第3组 中位值 0 0.205 0.200 0.220 0.205 12.5 0.166 0.175 0.160 0.166 25.0 0.135 0.130 0.130 0.130 50.0 0.085 0.085 0.095 0.085 75.0 0.055 0.045 0.045 0.045 100.0 0.010 0.010 0.005 0.010 表 3 黏土岩侧限膨胀率
Table 3 Confined expansive rates of clay rock
法向压力
σn/kPa侧限膨胀率 第1组 第2组 第3组 中位值 0 0.330 0.355 0.345 0.345 25 0.290 0.275 0.280 0.280 50 0.210 0.225 0.235 0.225 100 0.130 0.140 0.130 0.130 150 0.055 0.055 0.050 0.055 200 0.010 0.015 0.015 0.015 -
[1] 谭罗荣,张梅英,邵梧敏,等.灾害性膨胀土的微结构特征及其工程性质[J].岩土工程学报,1994,16(2):48–57. Tan Luorong,Zhang Meiying,Shao Wumin,et al.Nicrostructure characteristics and engineering properties of expansive soil[J].Chinese Journal of Geotechnical Engineering,1994,16(2):48–57 [2] Mao Dawei,Nilsen B,Lu Ming.Analysis of loading effects on reinforced shotcrete ribs caused by weakness zone containing swelling clay[J].Tunnelling and Underground Space Technology,2011,26(3):472–480. doi: 10.1016/j.tust.2011.01.004 [3] Puppala A J.Performance evaluation of infrastructure on problematic expansive soils:Characterization challenges,innovative stabilization designs,and monitoring methods[J].Journal of Geotechnical and Geoenvironmental Engineering,2021,147(8):04021053. doi: 10.1061/%28ASCE%29GT.1943-5606.0002518 [4] 戴福初,董文萍,黄志全,等.南水北调中线段原状膨胀土抗剪强度试验研究[J].工程科学与技术,2018,50(6):123–131. doi: 10.15961/j.jsuese.201800242 Dai Fuchu,Dong Wenping,Huang Zhiquan,et al.Study on shear strength of undisturbed expansive soil of middle route of south-to-north water diversion project[J].Advanced Engineering Sciences,2018,50(6):123–131 doi: 10.15961/j.jsuese.201800242 [5] Rao S M,Thyagaraj T,Thomas H R.Swelling of compacted clay under osmotic gradients[J].Géotechnique,2006,56(10):707–713. doi: 10.1680/GEOT.2006.56.10.707 [6] Wang Jinshu,Wu Guang,Zhao Zhiming,et al.Experimental study on factors influencing the expansive properties of disturbed swelling rock[J].Soil Mechanics and Foundation Engineering,2020,57(1):84–91. doi: 10.1007/s11204-020-09641-8 [7] Jamsawang P,Nuansrithong N,Voottipruex P,et al.Laboratory investigations on the swelling behavior of composite expansive clays stabilized with shallow and deep clay-cement mixing methods[J].Applied Clay Science,2017,148:83–94. doi: 10.1016/j.clay.2017.08.013 [8] 秦幼林,许崇帮,高晓静,等.膨胀岩本构力学模型研究进展及展望[J].公路交通科技,2021,38(6):23–31. Qin Youlin,Xu Chongbang,Gao Xiaojing,et al.Study progress and prospect of constitutive mechanical model of swelling rock[J].Journal of Highway and Transportation Research and Development,2021,38(6):23–31 [9] Gysel M.Anhydrite dissolution phenomena:Three case histories of anhydrite Karst caused by water tunnel operation[J].Rock Mechanics and Rock Engineering,2002,35(1):1–21. doi: 10.1007/s006030200006 [10] 缪协兴,杨成永,陈至达.膨胀岩体中的湿度应力场理论[J].岩土力学,1993(4):49–55. doi: 10.16285/j.rsm.1993.04.009 Miao Xiexing,Yang Chengyong,Chen Zhida.Humidity stress field theory in swelling rock mass[J].Rock and Soil Mechanics,1993(4):49–55 doi: 10.16285/j.rsm.1993.04.009 [11] 缪协兴.湿度应力场理论的耦合方程[J].力学与实践,1995,17(6):22–24. [12] 白冰,李小春.湿度应力场理论的证明[J].岩土力学,2007,28(1):89–92. doi: 10.16285/j.rsm.2007.01.017 Bai Bing,Li Xiaochun.Proof of humidity stress field theory[J].Rock and Soil Mechanics,2007,28(1):89–92 doi: 10.16285/j.rsm.2007.01.017 [13] 白冰.CO2煤层封存流动–力学理论及场地力学稳定性数值模拟方法[D].武汉:中国科学院武汉岩土力学研究所,2008. Bai Bing.CO2 Fluid-mechanical theory of coal seam sealing and numerical simulation method of site mechanical stability[D].Wuhan:Institute of Rock and Soil Mechanics,Chinese Academy of Sciences,2008. [14] 陈有亮,刘耕云,杜曦,等.考虑膨胀应力和剪胀的深埋隧道弹塑性解[J].岩土力学,2020,41(8):2525–2535. doi: 10.16285/j.rsm.2019.1799 Chen Youliang,Liu Gengyun,Du Xi,et al.Elastoplastic solution for a deep-buried tunnel considering swelling stress and dilatancy[J].Rock and Soil Mechanics,2020,41(8):2525–2535 doi: 10.16285/j.rsm.2019.1799 [15] 陈钒,吴建勋,任松,等.基于湿度应力场理论的硬石膏岩膨胀试验研究[J].岩土力学,2018,39(8):2723–2731. doi: 10.16285/j.rsm.2017.0544 Chen Fan,Wu Jianxun,Ren Song,et al.The swelling tests of anhydrite based on the humidity stress field theory[J].Rock and Soil Mechanics,2018,39(8):2723–2731 doi: 10.16285/j.rsm.2017.0544 [16] 任松,欧阳汛,吴建勋,等.含时间效应的硬石膏隧道弹-膨胀解析模型[J].浙江大学学报(工学版),2018,52(5):896–905. Ren Song,Ouyang Xun,Wu Jianxun,et al.Elastic-swelling analytical model of anhydrite rock considering time-dependent effect[J].Journal of Zhejiang University(Engineering Science),2018,52(5):896–905 [17] 张善凯,冷先伦,盛谦.卢氏膨胀岩湿胀软化特性研究[J].岩土力学,2020,41(2):561–570. doi: 10.16285/j.rsm.2019.0358 Zhang Shankai,Leng Xianlun,Sheng Qian.Study of water swelling and softening characteristics of expansive rock[J].Rock and Soil Mechanics,2020,41(2):561–570 doi: 10.16285/j.rsm.2019.0358 [18] 张嘉凡,程树范,王焕,等.西部弱胶结软岩细观结构及水理特性试验[J].煤田地质与勘探,2020,48(3):116–121. Zhang Jiafan,Cheng Shufan,Wang Huan,et al.Experiment of mesostructure and hydrologic characteristics of weakly cemented soft rocks in Western China[J].Coal Geology & Exploration,2020,48(3):116–121 [19] Wild K M,Wymann L P,Zimmer S,et al.Water retention characteristics and state-dependent mechanical and petro-physical properties of a clay shale[J].Rock Mechanics and Rock Engineering,2015,48(2):427–439. doi: 10.1007/s00603-014-0565-1 [20] 张坤勇,徐娜,陈恕,等.膨胀土完全软化强度指标试验研究[J].岩土工程学报,2020,42(11):1988–1995. Zhang Kunyong,Xu Na,Chen Shu,et al.Experimental study on fully softened shear strength of expansive soil[J].Chinese Journal of Geotechnical Engineering,2020,42(11):1988–1995 [21] 吴建勋,任松,范金洋,等.考虑软化效应的隧道围岩自稳时间解析模型[J].工程科学与技术,2017,49(3):68–75. doi: 10.15961/j.jsuese.201600338 Wu Jianxun,Ren Song,Fan Jinyang,et al.Analytical model for stand-up time of rock mass surrounding tunnel with consideration of softening effect[J].Advanced Engineering Sciences,2017,49(3):68–75 doi: 10.15961/j.jsuese.201600338 [22] 甘永德,刘欢,贾仰文,等.膨胀性土壤饱和水分运动参数计算模型[J].工程科学与技术,2018,50(2):77–83. doi: 10.15961/j.jsuese.201700680 Gan Yongde,Liu Huan,Jia Yangwen,et al.Swelling soil saturated water movement parameters calculating models[J].Advanced Engineering Sciences,2018,50(2):77–83 doi: 10.15961/j.jsuese.201700680 [23] 曾铃,李光裕,史振宁,等.降雨入渗条件下非饱和土渗流特征试验[J].中国公路学报,2018,31(2):191–199. doi: 10.19721/j.cnki.1001-7372.2018.02.020 Zeng Ling,Li Guangyu,Shi Zhenning,et al.Experiment on seepage characteristics of unsaturated soil under rainfall infiltration[J].China Journal of Highway and Transport,2018,31(2):191–199 doi: 10.19721/j.cnki.1001-7372.2018.02.020 [24] 邓建华,黄醒春.考虑时间因素的膏溶角砾岩膨胀本构关系[C]//第六届全国土木工程研究生学术论坛论文集.北京:中国土木工程学会,2008:128–133. Deng Jianhua,Huang Xingchun.Study on time-dependent swelling constitute relation of gypsum breccia[C]//The 6th National Civil Engineering Graduate Academic Forum.Beijing:China Civil Engineering Society,2008:128–133. [25] 翟张辉,张亚国,李同录,等.考虑非饱和土吸力效应的圆柱孔扩张有效应力解析[J].工程科学与技术,2020,52(3):100–106. doi: 10.15961/j.jsuese.201900402 Zhai Zhanghui,Zhang Yaguo,Li Tonglu,et al.Effective stress solution for the cylindrical cavity expansion in unsaturated soil considering suction effect[J].Advanced Engineering Sciences,2020,52(3):100–106 doi: 10.15961/j.jsuese.201900402 [26] 缪协兴.用湿度应力场理论解圆形硐室遇水作用问题[J].岩土工程学报,1995,17(5):86–90. doi: 10.3321/j.issn:1000-4548.1995.05.014
