非等温条件下土体固结变形及溶质运移规律分析

田改垒 张志红

(北京工业大学城市与工程安全减灾教育部重点实验室， 北京 100124)

随着城市化和工业化的发展，垃圾处理已成为发达国家和发展中国家最严重的环境问题之一，针对该类问题，卫生填埋已成为首选方法[1-2].垃圾堆体生化反应产生的热量导致垫层上部土体温度升高，一般性填埋场甚至可以达到60 ℃[3]，温度的改变使得土体内部孔隙水压力、土体变形及溶质运移进程发生变化，进而对黏土垫层设计及服役性能评估产生重要影响.

为探究温度对多物理场耦合进程的影响，国内外学者从不同角度对耦合机理及过程进行描述并建立了相应的耦合模型.Yan等[4]基于有限元离散方法推导出三维热液(TH)耦合模型，该模型可用于探究裂隙岩体中的流体流动及热量传输问题.Chen等[5]建立了可分析水-空气-蒸汽三相流动的热-水-力(THM)耦合模型，该模型考虑了热传输和孔隙率变化对耦合进程的影响.Song等[6]建立了考虑热渗透及机械-热效应的热-水-力耦合模型，研究发现温度梯度引起的热渗透效应会显著加快孔隙水的传输.郭志光等[7]以湿-热弹性理论为基础，建立了非饱和热-水-力耦合模型，并给出了与温度相关的状态方程.随着研究的不断深入，考虑温度影响的热-水-力-化(THMC)耦合模型也逐渐得到发展.张玉军[8-9]和Gens等[10]针对高放射核废料处置库周围缓冲层建立了热-水-力-化耦合模型，其中张玉军[9]采用自主研发的有限元程序模拟了热、水、应力和迁移耦合作用对温度、水压力及污染物浓度分布的影响.Cleall等[11]建立了可考虑多组分反应性运输的非饱和热-水-力-化耦合模型，该模型考虑了温度对溶质扩散及气相水迁移的影响.Majid[12]对热、水力和化学行为进行了耦合分析，重点探究了升温对压实膨润土缓冲垫层地球化学反应的影响.Zheng等[13]建立的热-水-力-化耦合模型可用于描述溶质交叉扩散和热化学渗透现象.冯忠等[14]建立了考虑厌氧降解作用的耦合模型, 探究了温度对城市固废生化降解过程的影响，并通过温度系数对生化降解反应速率一级动力方程进行了修正.Liu等[15]基于局部平衡方程，建立了可模拟甲烷水合物解离过程的热-水-力-化模型.Wu等[16]建立可考虑力学固结和热固结作用的热-水-力-化耦合模型.张志红等[17]基于孔隙流体质量守恒方程、能量守恒方程和污染物质量守恒方程，建立了可考虑黏土自身属性及环境特征(温度和浓度)对溶质运输机制影响的热-水-力-化耦合模型.

国内外学者对考虑温度影响的耦合模型研究极大丰富了对温度作用机制的认识，但现有模型较少考虑温度对力学固结及化学-渗透固结过程的影响，也未量化其对土体固结变形及溶质运移进程的影响程度.基于此，本文将可反映热固结、力学固结及化学-渗透固结联合作用的热-水-力-化耦合模型与水-力-化耦合模型进行对比分析，探究了温度对孔隙水压力、土体变形及溶质运移进程的影响.

1 THMC耦合模型

1.1 模型假设

假设荷载瞬时施加，土体为均质、饱和、各向同性的线弹性体，孔隙流体流动和污染物运输分别满足达西定律和菲克定律，土体变形、孔隙流体流动、热量传输和污染物运移仅发生在竖直方向.

1.2 土体变形控制方程

基于应力平衡方程，考虑力学固结、热固结及化学-渗透固结对土体变形的影响，推导得到了均质弹性饱和土体的变形控制方程为

(1)

式中，E为弹性模量；ν为泊松比；wz为土体竖向位移；c为污染物的浓度；ψc为浓度引起的体积变化系数；ψπ为渗透压引起体积变化系数；R为普适气体常数；α为土体热膨胀系数；T为绝对温度；u为孔隙流体压力；γw为孔隙水重度；i=h0/l为水力梯度，h0为渗滤液水头，l为垫层厚度；n为孔隙率；ρf和ρs分别为孔隙流体密度和固体密度；g为重力加速度.

力学荷载、热荷载及化学荷载共同作用下的孔隙率表示为

(2)

1.3 孔隙流体流动控制方程

忽略污染物浓度对孔隙流体密度的影响，根据孔隙流体质量守恒推导得到孔隙流体流动控制方程为

(3)

式中，βf为孔隙流体压缩系数；ω为化学渗透效率系数；kT为热渗透系数；kh=kh0(n/n0)3为水力传导系数[18]，kh0为初始水力传导系数；M为溶质摩尔质量.

1.4 热传导控制方程

根据能量守恒方程，考虑孔隙流体流动对热量传输的贡献，推导得到热传导控制方程为

(4)

式中，Cps为土颗粒比热容；Cpf为孔隙流体比热容；λT为土体固有热传导率.

1.5 污染物运移控制方程

黏土颗粒具有吸附污染物离子的能力，污染物在黏土垫层中运移时不仅会随孔隙流体流动，还会被吸附到黏土矿物表面.基于液相和固相污染物质量守恒，推导得到的污染物运移控制方程为

(5)

式中，Kd为线性吸附系数；Dc=τ(1-ω)D0为污染物的有效扩散系数[19]；τ=n2为曲率因子，D0为开放水体中污染物扩散系数；ST为索雷特系数.

2 模型验证

为保证所建模型的正确性和有效性，忽略本文所建模型的化学作用，并将模拟结果与文献[20]的热固结试验进行对比.

初始条件T(z,0)= 20 ℃ (0≤z≤l；l=76 mm)；边界条件T(l,t)= 100 ℃ (t≥0)，u(z,t)= 4 MPa (t≥0).

(α2+β2T)(ln(pu))2

(6)

式中，α0、α1、α2、β1、β2、p均为常系数，α0=4.505×10-4℃-1，α1=9.156×10-5℃-1，α2=6.381×10-6℃-1，β1=1.2×10-6℃-2，β2=-5.766 × 10-8℃-2，p=15 MPa-1.

图1显示了土体孔隙率随温度变化的对比结果.由对比结果可知，所建模型预测结果与试验结果吻合度较高，保证了所建模型的正确性.

图1 孔隙率变化对比结果

3 HMC耦合模型

将本文所建THMC耦合模型中的温度因素忽略，则可得到与文献[21]相一致的HMC耦合模型.控制方程如下：

土体变形控制方程

(7)

孔隙流体控制方程

(8)

污染物运移控制方程

(9)

4 数值模拟

4.1 模拟参数及边界条件

选取1 m垫层进行分析，堆体荷载设为50 kPa，模拟所选特征溶质为CuCl2，将上边界质量浓度设置为135 kg/m3.

初始条件u(z, 0)=0，c(z, 0)=0，T(z, 0)=20 ℃，wz(z, 0)=0(0≤z≤l).

边界条件u(0,t)=0，c(0,t)= 135 kg/m3，T(0,t)=25,30,40,60 ℃，u(l,t)=0，c(l,t)=0，T(l,t)=20 ℃，wz(l,t)=0 (t≥0).

模型参数[21]D0=5×10-10m2/s，k0=1×10-10m/s，E=1.49×10-10Pa，ν=0.3，R=8.314 J/(mol·K)，ρs=2 700 kg/m3，ρf=1 000 kg/m3，γw=10 kN/m3，ω=0.005，n0=0.5，kd=0.814×10-3m3/kg，ψv= 2.5×10-9Pa-1，ψπ=5×10-7Pa-1，ψc=1.05×10-8m3/kg，M=135 g/mol，h0=0.

4.2 模拟结果

图2给出了不同温差条件下，固结作用时间分别为0.1、1、10和50 a时，HMC及THMC耦合模型所对应的孔隙水压力变化规律对比结果.由图2(a)～(d)可知，初始阶段(0.1 a)，较深土层(大于0.1 m)孔隙水压力为正，且相比于HMC耦合模型，THMC耦合模型所对应的正孔隙水压力大，但随着时间的推移，孔隙水压力表现为负，此时THMC耦合模型所对应的负孔隙水压力小.这是因为力学荷载及热荷载会引起正孔隙水压力，而化学荷载引起负孔隙水压力.初始阶段，力学荷载及热荷载引起的正孔隙水压力未完全消散且大于化学渗透流动引起的负孔隙水压力，故较深土层孔隙水压力表现为正.但随着力学固结和热固结的完成，化学渗透流动引起的化学-渗透固结占据主导地位，此时负孔隙水压力遍布整个土体.由于THMC耦合模型考虑了温差引起的热渗透和热扩散效应，污染物运移速率加快，因此负孔隙水压力的消散速率加剧.

(a) 5 ℃

此外，从图2中还可以观察到，随着垫层上下边界温度差的增加，2种耦合模型所对应的曲线间距逐渐加大，可见随着温差的增加，温度对孔隙水压力的影响逐渐增强.当模拟时间为50 a、温差为40 ℃时，与HMC耦合模型负孔隙水压力峰值(-6.2 kPa)相比，考虑THMC耦合时的负孔隙水压力峰值(-3.6 kPa)减小了42%.

图3给出了不同温差条件，HMC及THMC耦合模型所对应的沉降量变化规律对比结果.由图3(a)～(d)可知，初始阶段THMC耦合模型所对应的土体沉降量小于HMC耦合模型所对应的土体沉降量，但沉降峰值却比HMC耦合模型所对应的沉降峰值大.这是因为THMC耦合模型考虑了土体的热膨胀变形即负向沉降，负向沉降抵消了力学荷载引起的正向沉降，从而减小初始阶段土体变形.但由于THMC耦合模型考虑了土体热固结效应，因此其所对应的沉降峰值较HMC耦合模型大.此外，从图3中还可以看到，随着时间的延长，土体沉降开始出现回弹，且相比于HMC耦合模型，THMC耦合模型所对应的回弹量大.这是因为随着污染物的运移，垫层上下边界污染物浓度梯度减小，土体所受化学渗透压力降低，此过程可视为卸载过程，因此土体表现出一定程度的回弹.由于THMC耦合模型考虑了温度引起的热渗透及热扩散效应，污染物传输速率加快，故THMC耦合模型所对应的沉降回弹量增大.

(a) 5 ℃

还需注意的是，当温差为40 ℃时， THMC耦合模型所对应的土体沉降可分为3个阶段：负向增长阶段(Ⅰ)、正向增长阶段(Ⅱ)、回弹阶段(Ⅲ).计算发现模拟时间为50 a时，THMC耦合模型所对应的沉降回弹量(4.3 mm)约为HMC耦合模型所对应的沉降回弹量(0.8 mm)的5.4倍.

图4给出了不同温差条件，固结作用时间分别为0.1、1、10和50 a时，HMC及THMC耦合模型所对应的污染物浓度变化规律对比结果.由图4(a)～(d)可知，THMC耦合作用下的污染物运移速率比HMC耦合作用下的污染物运移速率快，这主要是因为垫层上下边界温差会引起热渗透及热扩散现象的发生，由于热渗透、热扩散方向与对流、弥散方向一致，因此污染物的运输加快.

(a) 5 ℃

此外，从图4中还可以看到，随着温差的增加，2种耦合机制所对应的污染物运移速率差异逐渐增大，在图中表现为2条曲线之间的间距增加.这主要是因为随着温差的增加，热渗透和热扩散作用逐渐增强.

击穿时间是评估垃圾填埋场防渗垫层服役性能的重要指标，具体指垫层底部污染物积累浓度达到规定极限浓度所需要的时间.根据《中华人民共和国地下水质量标准》(GB 14848—2017)所规定的Ⅲ类水标准来确定CuCl2的极限质量浓度为2.1 mg/L，即当垫层下部CuCl2累积质量浓度达到2.1 mg/L时就认为污染物已击穿黏土垫层.表1给出了2种耦合机制影响下污染物的击穿时间差值.根据表1计算结果可知，当垫层上下边界温差为40 ℃时，相比于HMC耦合模型(击穿时间为33 a)，THMC耦合模型所对应的击穿时间缩短了11.4 a.可见，在评估垃圾填埋场服役性能和服役年限时必须考虑温度影响.

表1 污染物击穿时间差

5 结论

1) 垫层上下边界温差会引起正孔隙水压力，同时加快负孔隙水压力消散，此外随着温差的增加，温度对孔隙水压力的影响逐渐增强.当温差达到40 ℃、模拟时间为50 a时，考虑温度影响时的负孔隙水压力峰值减小42%.

2) 垫层上下边界温差会使土体产生负向沉降，同时加大沉降峰值及后续沉降量的回弹，且随着温差的增加，其对土体变形的影响逐渐增强.当温差为40 ℃时， THMC耦合模型所对应的土体变形分为负向增长阶段、正向增长阶段、回弹阶段3个阶段.模拟时间为50 a时，考虑温度影响的沉降回弹量扩大了5.4倍.

3) 垫层上下边界温差可显著加快污染物运输，缩短污染物的击穿时间，且随着温差的增加，温度对污染物运移进程的影响逐渐增强，当温差为40 ℃时，考虑温度影响的污染物击穿时间缩短了11.4 a.