随着我国城市的快速发展，城市固体废弃物产量以每年8%~15%的速度增长，包围了全国2/3城市^[1]. 目前，填埋处置仍是我国生活垃圾的主要处置方式之一^[2]，截至2014年，垃圾填埋场的数量为2万多座^[3]. 美国垃圾组分中厨余垃圾质量分数约为12.5%^[4]，垃圾中水的质量分数为15%~40%^[5]. 与欧美国家不同，我国垃圾组分中厨余垃圾质量分数高，约为60%^[6]，水的质量分数为50%~70%^[7]，导致我国垃圾初期水的质量分数高，初期降解作用强烈的现象. 垃圾降解产生很多有毒有害的渗滤液，且我国填埋场普遍存在渗滤液液位雍高现象，水头通常高于10 m，其渗漏速率随水头增高而明显增大^[8]，存在污染土壤和地下水的风险，严重威胁城市环境安全.

衬垫系统是阻止渗滤液扩散的重要屏障，是保障填埋场安全服役的关键. 目前国内外大多数使用的是土工膜（geomembrane，GM）+压实黏土衬垫（compacted clay liner，CCL）+ 天然衰减层（attenuation layer，AL）或GM + 钠基膨润土防水毯（geosynthetic clay liner，GCL）+ AL组成的复合衬垫. 根据复合衬垫的层数，又分为单层复合衬垫与双层复合衬垫. 美国于20世纪90年代开始大量推广双层复合衬垫，规定全部危废填埋场及多数生活垃圾填埋场必须使用双层复合衬垫. 我国仍以单层复合衬垫为主，相关规范^[9-11]给出的2类常用单层衬垫结构为GM + CCL + AL和GM + GCL + AL. 由于GCL价格经济、施工便捷及渗透系数低，被成功应用于填埋场以替代CCL^[12-14]，尤其是在场地条件不允许、CCL不便实施时，GCL具有较好应用前景.

已有不少学者开展了GCL衬垫与CCL衬垫之间的等效分析. Foose等^[15]得到在半无穷条件下，污染物在复合衬垫中一维扩散模型的解析解. Rowe等^[16]根据污染物对受体含水层的影响，评估GM + CCL + AL和GM + GCL + AL替代复合衬垫的等效性，发现与GM/CCL/AL衬垫相比，GM/GCL/AL衬垫对地下水具有相同甚至更好的保护作用. 陈云敏等^[17]建立在层状土中的一维扩散解析模型. 谢海建等^[18]建立溶质在黏土中的一维扩散模型，在考虑有限厚度边界及土中溶质背景浓度的条件下得到了相应的解析解. 张虎元等^[19-20]从衬垫的污染物累积量和由吸附产生的污染物衰减2个角度研究填埋场防渗衬垫等效替代，但只涉及一层土质衬垫. 为了研究GCL/AL防渗层防污性能与CCL之间的等效性，谢海建等^[12]建立污染物在GCL/AL防渗层的一维对流-弥散模型，在假设污染物在GCL中的运移为稳态过程的基础上得到了其解析解. Cleall等^[21]得到了在半无限与有限厚度条件下，有机污染物在复合衬垫中一维扩散模型的解析解. 冯世进等^[22]考虑垃圾土沉降导致的成层性，运用SEEP/W有限元软件模拟渗滤液喷洒回灌. 吴珣等^[23]建立了在考虑降解作用时有机污染物在GM + CCL和GM + GCL完好复合衬垫中的一维扩散模型，并得到其解析解. 张文杰等^[24]探讨了常用的零浓度边界、零浓度梯度边界和半无限边界对应的试验条件，求解得到各边界条件下统一的有限厚度土层中对流-扩散-吸附解析解.

多层衬垫的污染物运移的解析解较复杂，即使采用数值模拟软件求解，仍然须获得大量的现场数据作为率定支撑，同时要求设计人员熟练使用数值模拟软件，对于工程设计而言较为不便. 因此，有必要提出等效简化多层衬垫计算模型以期能替代复杂解析方法和数值方法成为快速、经济、有效的设计和评估手段.

在上述研究基础上，提出基于时间矩卷积的多层衬垫污染物运移参数的等效分析模型，以期获得渗流速度与水动力弥散系数的等效表达式. 采用填埋场污染物运移分析程序POLLUTE v7^[25]以及国内外土柱实验数据^[26-28]验证该等效模型的合理性、可靠性与实用性. 以衬垫底部相对浓度、瞬时通量以及累计通量等3个参数作为评价指标，对填埋场GM + CCL + AL及GM + GCL + AL这2类不同复合衬垫进行等效分析计算.

在建立污染物衬垫等效计算模型时，主要基于以下假设：1）污染物在土中的主要迁移方式为考虑对流和分子扩散作用；2）污染物在土中的扩散与对流是一维的；3）下边界半无穷；4）指示性污染物为Cl⁻（Cl⁻在填埋场渗滤液中浓度高、无吸附性、运移快，常被作为指示性污染物用来模拟或计算污染物在衬垫系统中的运移^[29-30]）.

基于上述假设，如图1所示的各层土质衬垫中污染物的运移方程^[29]可以描述为

(1) $\begin{array}{l} \dfrac{{\partial {C_i}}}{{\partial t}} = {D_i}\dfrac{{{\partial ^2}{C_i}}}{{\partial {z^2}}} - {v_i}\dfrac{{\partial {C_i}}}{{\partial z}};\; i = 1,2,\cdots,n. \end{array}$

式中：C_i为第i层土质衬垫的污染物浓度，D_i为第i层土质衬垫的水动力弥散系数，v_i为第i层土质衬垫的渗流速度，t为时间，z为空间坐标，n为土质衬垫的总层数.

式(1)的初始条件为

(2) ${C_i}\left( {z,0} \right) = 0.$

相应的边界条件为

(3) ${C_0}\left( {0,t} \right) = {{M\delta \left( t \right)}}/{q},$

(4) ${C_1}\left( {{0^ + },t} \right) = {C_0}\left( {0,t} \right),$

(5) ${C_i}\left( {L_i^ - ,t} \right) = {C_{i + 1}}\left( {L_i^ + ,t} \right),$

(6) ${C_n}\left( {\infty ,t} \right) = 0.$

式中：M为单位面积污染物注入质量，δ(t)为狄拉克函数，q为达西流速，L_i为第i层土质衬垫底部与土工膜底部的距离.

土质衬垫的污染物运移过程可以近似等效为如下过程. 将单位面积的溶质M在土层表面（z=0）瞬间注入的过程近似为狄拉克条件，其分布函数^[31]为

(7) ${C^\delta }\left( {z,t} \right) = \frac{M}{{\varphi z}}{C^{ - \delta }}\left( {z,t} \right).$

式中：C^−δ为单位孔隙体积的注入质量传递函数， $\varphi $为孔隙率. 式(7)的卷积形式为

(8) $C\left( {z,t} \right) = \int_0^t {\frac{{M\left( \tau \right)}}{{\varphi z}}{C^{ - \delta }}\left( {z,t - \tau } \right)} {\rm{d}}\tau .$

易知

(9) $\frac{{M\left( \tau \right)}}{{\varphi z}} = \frac{v}{z}{C_0}\left( {0,\tau } \right).$

则式(8)化简为

(10) $C\left( {z,t} \right) = \frac{v}{z}\int_0^t {{C_0}\left( \tau \right){C^{ - \delta }}\left( {z,t - \tau } \right)} {\rm{d}}\tau .$

对于成层土质衬垫，可以一般化为

(11) $ \begin{array}{l} {C_i}\left( {{z_i},t} \right) = \dfrac{{{v_i}}}{{{z_i} - {L_{i - 1}}}}\displaystyle\int_0^t {{C_{i - 1}}\left( {{L_{i - 1}},\tau } \right)C_i^{ - \delta }\left( {{z_i} - {L_{i - 1}},t - \tau } \right)} {\rm{d}}\tau ;\\ \;\;\;\;\;i = 1,2,\cdots,n.\\[-12pt] \end{array} $

当C₀(0,t)=1时，在半无穷均质土的C(z,t)的解析解^[32]如下：

(12) $C\left( {z,t} \right) = \frac{1}{2}\left[ {{\rm{erfc}}\;\left( {\frac{{z - vt}}{{2\sqrt {Dt} }}} \right) + \exp\; \left( {\frac{{vz}}{D}} \right){\rm{erfc}}\;\left( {\frac{{z + vt}}{{2\sqrt {Dt} }}} \right)} \right].$

联立式(10)、(12)，可以解得

(13) ${C^{ - \delta }}\left( {z,t,v,D} \right) = \frac{{{z^2}}}{{2vt\sqrt {{\text{π}} Dt} }}\exp \;\left[ {\frac{{ - {{\left( {z - vt} \right)}^2}}}{{4Dt}}} \right].$

令

(14) $\overline C (z,p) = \int_0^\infty {C(z,t){{\rm{e}}^{ - pt}}{\rm{d}}t} .$

对式(13)进行拉普拉斯变换，得到

(15) $\overline {{C^{ - \delta }}} \left( {z,p,v,D} \right) = \frac{z}{v}\exp\; \left[ {\frac{{z\left( {v - \sqrt {4Dp + {v^2}} } \right)}}{{2D}}} \right].$

对式(11)进行拉普拉斯变换，并利用其卷积性质，得到

(16) $\begin{array}{l} \overline {{C_i}} \left( {{z_i},p} \right) = \dfrac{{{v_i}}}{{{z_i} - {L_{i - 1}}}}\overline {{C_{i - 1}}} \left( {{L_{i - 1}},p} \right)\overline {C_i^{ - \delta }} \left( {{z_i} - {L_{i - 1}},p,{v_i},{D_i}} \right);\\ \;\;\;\;i = 1,2,\cdots,n.\\[-12pt] \end{array}$

对于1，2，···，n层，均有

(17) $\left. {\begin{array}{{c}} {\overline {{C_1}} \left( {{L_1},p} \right) = \dfrac{{{v_1}}}{{{L_1}}}\dfrac{M}{q}\overline {C_1^{ - \delta }} \left( {{L_1},p,{v_1},{D_1}} \right)},\\ {\overline {{C_2}} \left( {{L_2},p} \right) = \dfrac{{{v_2}}}{{{L_2} - {L_1}}}\overline {{C_1}} \left( {{L_1},p} \right)\overline {C_2^{ - \delta }} \left( {{L_2} - {L_1},p,{v_2},{D_2}} \right)},\\ {\vdots}\\ \overline {{C_{n - 1}}} \left( {{L_{n - 1}},p} \right) = \dfrac{{{v_{n - 1}}}}{{{L_{n - 1}} - {L_{n - 2}}}}\overline {{C_{n - 2}}} \left( {{L_{n - 2}},p} \right)\times\\ \overline {C_{n - 1}^{ - \delta }} \left( {{L_{n - 1}} - {L_{n - 2}},p,{v_{n - 1}},{D_{n - 1}}} \right),\\ {\overline {{C_n}} \left( {z,p} \right)\!\! =\!\! \dfrac{{{v_n}}}{{z \!-\! {L_{n - 1}}}}\overline {{C_{n - 1}}} \left(\! {{L_{n - 1}},p} \right)\overline {C_n^{ - \delta }} \left( \!{{z_n} \!-\! {L_{n - 1}},p,{v_n},{D_n}} \!\right)\!.\!\!} \end{array}} \right\}$

累乘式（17），得到

(18) $\begin{split} \overline {{C_n}} \left( {z,p} \right) =& \dfrac{M}{q}\exp\; \left[ {\dfrac{{{L_1}\left( {{v_1} - \sqrt {4{D_1}p + v_1^2} } \right)}}{{2{D_1}}}} \right.{\rm{ + }}\\ &\dfrac{{\left( {z - {L_{n - 1}}} \right)\left( {{v_n} - \sqrt {4{D_n}p + v_n^2} } \right)}}{{2{D_n}}} + \\ &\left. {\displaystyle\sum\limits_{i = 2}^{n - 1} {\dfrac{{\left( {{L_i} - {L_{i - 1}}} \right)\left( {{v_i} - \sqrt {4{D_i}p + v_i^2} } \right)}}{{2{D_i}}}} } \right]. \end{split}$

Aris^[33]定义时间矩为

(19) ${m_j} = \int_0^\infty {{t^j}C\left( {z,t} \right){\rm{d}}t;j = 0,1,2,\cdots} .$

定义均值 $\overline \mu$及高阶中心矩μ_j分别为

(20) $\overline \mu = {{{m_1}}}/{{{m_0}}},$

(21) ${\mu _j} = \frac{1}{{{m_0}}}{\int_0^\infty {\left[ {t - \overline \mu } \right]} ^j}C\left( {z,t} \right){\rm{d}}t;j = 1,2,\cdots.$

二阶中心矩μ₂为

(22) ${\mu _2} = {{{m_2}}}/{{{m_0}}} - {\overline \mu ^2}.$

实际上，不必通过C(z,t)解求得m_j，也可以通过拉普拉斯域浓度解直接求得^[33]：

(23) ${m_j} = {\left( { - 1} \right)^j}\mathop {\lim }\limits_{p \to 0}\; \left[ {\frac{{{\partial ^j}\overline C \left( {z,p} \right)}}{{\partial {p^j}}}} \right].$

其中，零阶时间距为

(24) ${m_0} = \mathop {\lim }\limits_{p \to 0}\; \overline {{C_n}} \left( {z,p} \right).$

将式(18)代入式(24)，可以得到

(25) ${m_0} = {M}/{q}.$

一阶时间距为

(26) ${m_1} = \dfrac{M}{q}\left[ {\frac{z}{{{v_n}}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {{L_i}\left( {\frac{1}{{{v_i}}} - \frac{1}{{{v_{i + 1}}}}} \right)} } \right].$

因此，

(27) $\overline \mu = \frac{z}{{{v_n}}} + \sum\limits_{i = 1}^{n - 1} {\left[ {{L_i}\left( {\frac{1}{{{v_i}}} - \frac{1}{{{v_{i + 1}}}}} \right)} \right]} .$

二阶时间距为

(28) $\begin{split} {m_2} =& \frac{M}{q}\left\{ {\left[ {\frac{{2{D_1}{L_1}}}{{v_1^3}} + \frac{{2{D_n}\left( {z - {L_{n - 1}}} \right)}}{{v_n^3}} + \sum\limits_{i = 2}^{n - 1} {\frac{{2{D_i}\left( {{L_i} - {L_{i - 1}}} \right)}}{{v_i^3}}} } \right]} \right.{\rm{ - }}\\ &\left. {{{\left[ {\frac{z}{{{v_n}}} + \sum\limits_{i = 1}^{n - 1} {{L_i}\left( {\frac{1}{{{v_i}}} - \frac{1}{{{v_{i + 1}}}}} \right)} } \right]}^2}} \right\}.\\[-20pt] \end{split}$

因此，

(29) ${\mu _2} = \frac{{2{D_1}{L_1}}}{{v_1^3}} + \frac{{2{D_n}\left( {z - {L_{n - 1}}} \right)}}{{v_n^3}} + \sum\limits_{i = 2}^{n - 1} {\frac{{2{D_i}\left( {{L_i} - {L_{i - 1}}} \right)}}{{v_i^3}}}. $

对于单层土质衬垫情况，也可以求得相应的 $\overline \mu $和μ₂：

(30) $\overline \mu = {z}/{{{v_{\rm{e}}}}},$

(31) ${\mu _2} = {{2z{D_{\rm{e}}}}}/{{v_{\rm{e}}^3}}.$

式中：v_e、D_e分别为等效流速及等效水动力弥散系数. 联立式(27)、(29)~(31)，可以得到v_e、D_e，取z=L_n，可以得到

(32) ${v_{\rm{e}}} = \dfrac{{{L_n}}}{{\dfrac{{{L_1}}}{{{v_1}}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {\frac{{{L_{i + 1}} - {L_i}}}{{{v_{i + 1}}}}} }},$

(33) ${D_{\rm{e}}} = \frac{{L_n^2\left[ {\dfrac{{{L_1}{D_1}}}{{v_1^3}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {\dfrac{{\left( {{L_{i + 1}} - {L_i}} \right){D_{i + 1}}}}{{v_{i + 1}^3}}} } \right]}}{{{{\left( {\dfrac{{{L_1}}}{{{v_1}}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {\dfrac{{{L_{i + 1}} - {L_i}}}{{{v_{i + 1}}}}} } \right)}^3}}}.$

式（1）在单层土质介质的情况下，已有学者给出了其相对浓度的解析解^[32]（式（12）），代入v_e与D_e，第n层土质衬垫的浓度解计算公式为

(34) $\frac{{{C_{\rm{b}}}}}{{{C_0}}} = \dfrac{1}{2}\left[ {{\rm{erfc}}\;\left( {\frac{{1 - {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{{T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right) + \exp \;\left( {{P_{\rm{L}}}} \right){\rm{erfc}}\;\left( {\dfrac{{1 + {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{{T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right)} \right],$

(35) ${T_{\rm{R}}} = {{{v_{\rm{e}}}{t_{\rm{b}}}}}/{{{L_n}}},$

(36) ${P_{\rm{L}}} = \dfrac{{{v_{\rm{e}}}{L_n}}}{{{D_{\rm{e}}}}} = \dfrac{{{{\left( {\dfrac{{{L_1}}}{{{v_1}}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {\dfrac{{{L_{i + 1}} - {L_i}}}{{{v_{i + 1}}}}} } \right)}^2}}}{{\dfrac{{{L_1}{D_1}}}{{v_1^3}} + \displaystyle\sum\limits_{i = 1}^{n - 1} {\dfrac{{\left( {{L_{i + 1}} - {L_i}} \right){D_{i + 1}}}}{{v_{i + 1}^3}}} }}.$

式中：C_b为防污屏障底部污染物击穿浓度，erfc为互补误差函数，T_R为时间因子，P_L为Pelect数，t_b为指示性污染物击穿时间.

瞬时通量计算式如下^[34]：

(37) ${J_{\rm{I}}} = \frac{1}{2}\frac{D}{L}n{C_0}\left[ {\frac{{{v_{\rm{e}}}L}}{D}\left( {{Q_1} + {Q_2}} \right) + {Q_3} - {P_{\rm{L}}}{Q_2}} \right].$

式中：

(38) ${Q_1} = {\rm{erfc}}\;\left( {\frac{{1 - {T_{\rm{R}}}}}{{2\sqrt {{{{T_{\rm{R}}}}}/{{{P_{\rm{L}}}}}} }}} \right),$

(39) ${Q_2} = \exp\; \left( {{P_{\rm{L}}}} \right){\rm{erfc}}\;\left( {\frac{{1 + {T_{\rm{R}}}}}{{2\sqrt {{{{T_{\rm{R}}}}}/{{{P_{\rm{L}}}}}} }}} \right),$

(40) ${Q_3} = \dfrac{{2\exp \;\left[ { - {{\left( {\dfrac{{1 - {T_{\rm{R}}}}}{{2\sqrt {{{{T_{\rm{R}}}}}/{{{P_{\rm{L}}}}}} }}} \right)}^2}} \right]}}{{\sqrt {{{{\text{π}} {T_{\rm{R}}}}}/{{{P_{\rm{L}}}}}} }}.$

可以得到累计通量的表达式为

(41) ${J_{\rm{A}}} = \int_0^t {J_{{\rm{I}}}{\rm{d}}t} .$

化简得到

(42) ${J_{\rm{A}}} = \frac{1}{2}\frac{D}{L}n{C_0}t\left[ {{P_{\rm{L}}}{Q_1} + {Q_3} + {Q_4} + \frac{{2\left( {1 - {P_{\rm{L}}}} \right) - {Q_2}}}{{{T_{\rm{R}}}}}} \right].$

式中：

(43) ${Q_4} = \frac{{{P_{\rm{L}}} - 1}}{{{T_{\rm{R}}}}}{\rm{erfc}}\;\left( {\frac{{{T_{\rm{R}}} - 1}}{{2\sqrt {{{{T_{\rm{R}}}}}/{{{P_{\rm{L}}}}}} }}} \right).$

复合衬垫达西流速表达式^[35]为

(44) ${v_{\rm{a}}} = {{mQ}}/{A}.$

式中：m为土工膜中单位面积缺陷孔洞的个数，Q为渗漏率，A为流速的横截面面积.

当出现土工膜与褶皱水力连通的漏洞情况时，渗漏率表达式^{[16, 29, 36]}为

(45) $Q = {{2L\left( {{h_{\rm{w}}} + {L_{\rm{s}}}} \right)\left[ {kb + {{\left( {k{L_{\rm{s}}}\theta } \right)}^{0.5}}} \right]}}/{{{L_{\rm{s}}}}}.$

式中：L为与漏洞相连的褶皱长度，b为褶皱宽度，θ为土工膜与黏土衬垫界面的导水系数，h_w为土工膜上方水头，L_s为衬垫厚度，k为CCL或GCL的渗透系数.

式(6)为零浓度Dirichlet条件，工程中也有衬垫下方含水层流速不断减小或远处有不透水基岩的情况，则边界条件应采用零通量Neumann条件：

(46) ${{\partial {C_n}\left( {\infty ,t} \right)}}/{{\partial z}} = 0.$

将式(46)替换为边界条件式(6)，也可以得到与式(32)~(34)一样的结果，推导过程不赘述.

基于上述等效模型，采用POLLUTE v7软件对GM + GCL + AL型的衬垫结构进行模拟验证. 我国填埋场渗滤液水头较高^[8]，如上海老港填埋场主水位曾达到8 m，滞水位达25 m^[37]，远超《生活垃圾填埋场污染控制标准》（GB 16889—2008）^[9]中规定的运行期内渗滤液液位限值0.3 m，因此水头h_w分别选取0.3、5.0、15.0 m. 文献[38]指出与褶皱水力连通漏洞产生的渗漏率是平坦土工膜上漏洞的10~5 000倍，当褶皱长度为10 m时，前者为后者的9.37倍. 本研究对与褶皱水力连通漏洞频率取1个/hm²，褶皱宽度为0.2 m，褶皱长度为10 m，不考虑平坦土工膜上的漏洞.

《生活垃圾卫生填埋处理技术规范》（GB 50869—2013）^[10]指出GCL防渗层渗透系数不应大于5.0×10⁻¹¹ m/s. 但不少学者发现GCL易发生阳离子交换而削弱其吸水性及膨胀性，并导致其渗透系数增大近一个量级^[39-40]. Mery等^[41]指出在土工膜下伏的GCL服役一段时间后，由于离子交换作用，其渗透系数一般约为10⁻¹⁰ m/s. 此外，Meer等^[40]发现干湿循环作用也会使其渗透系数增大，在与离子交换的共同作用下，其渗透系数量级为10⁻⁸~10⁻⁶ m/s. 鉴于上述原因以及GCL目前在国内施工质量良莠不齐的现状，从设计安全性考虑，考虑实际运行中更恶劣的情况，GCL渗透系数参考Rowe等^[16]的取值，为2×10⁻¹⁰ m/s. 其他衬垫材料参数取值如表1^[16]所示.

在不同水头条件下，服役时刻t=100 a时，GCL复合衬垫(GM + 7 mm GCL + 2 m AL)的不同深度的相对浓度等效解与POLLUTE v7解较拟合，如图2所示. 在不同水头条件下，GCL复合衬垫（GM + 7 mm GCL + 1 m AL）不同时刻的相对浓度的等效解与POLLUTE v7结果较拟合，如图3所示. 在不同水头条件下，GCL复合衬垫（GM + 7 mm GCL + 1 m AL）不同时刻的累计通量等效解与POLLUTE v7解较吻合，如图4所示.

综上，本研究所提出的等效模型解与两层衬垫POLLUTE v7解在不同水头工况下的各时刻及各深度的相对浓度均一致，且在不同时间段的累计通量也一致，较拟合. 这充分说明了该等效模型的合理性和可靠性.

Sharma等^[28]开展多层土溶质运移的土柱实验，实验试剂采用NaCl溶液，实验采用脉冲型边界条件：

(47) ${C_1}\left( {0,\;t} \right) = \left\{ {\begin{array}{*{20}{l}} {{C_{\rm{0}}}{\rm{, \;0}} \leqslant t \leqslant {t_{\rm{0}}}}; \\ {{\rm{0, \;\; \; \; }}t{\rm{ > }}{t_{\rm{0}}}}. \end{array}} \right.$

式中：t₀=30 min. 多层土自上而下依次为中粒砂、粉土与细砂，三者渗流速度均为0.02 cm/s，厚度均为20 cm，水动力弥散系数依次为0.017、0.021、0.042 cm²/s. 通过本研究的等效模型，可以求得等效厚度L_n=60 cm，等效渗流速度v_e=0.02 cm/s以及等效水动力弥散系数D_e=0.026 cm/s. 将式(47)替换边界条件式(4)，运用拉式变化，可以求得第n层土质衬垫浓度解^[42]. 当0≤t ≤ t₀时，浓度解同式(34)；当t > t₀时，浓度解如下：

(48) $\dfrac{{{C_{\rm{b}}}}}{{{C_0}}} = \dfrac{1}{2}\left[ \begin{gathered} {\rm{erfc}}\;\left( {\dfrac{{1 - {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{{T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right) + \exp \;\left( {{P_{\rm{L}}}} \right){\rm{erfc}}\;\left( {\dfrac{{1 + {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{{T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right) \\ {\rm{ - erfc}}\;\left( {\dfrac{{1 - \Delta {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{\Delta {T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right){\rm{ - }}\exp \;\left( {{P_{\rm{L}}}} \right){\rm{erfc}}\;\left( {\dfrac{{1 + \Delta {T_{\rm{R}}}}}{{2\sqrt {\dfrac{{\Delta {T_{\rm{R}}}}}{{{P_{\rm{L}}}}}} }}} \right) \\ \end{gathered} \right],$

(49) $\Delta {T_{\rm{R}}} = {{{v_{\rm{e}}}\left( {{t_{\rm{b}}}{\rm{ - }}{t_0}} \right)}}/{{{L_n}}}.$

与Cl⁻作为示踪剂的三层土运移土柱实验数据相比，本研究等效模型和实验结果较拟合，平均相对偏差约为12%，如图5所示. 由此可见，本研究等效模型可以应用至不同的成层结构. 采用本研究等效模型拟合，过程较为简便，无须进行复杂的数值模拟计算. 这也从另一方面说明了本研究等效模型假设的合理和可靠性.

杨艳等^[26-27]通过室内土柱实验研究Cl⁻在层状土中的运移特性，该实验土柱高为20 cm并采用0.11 mol/L KCL溶液恒定注入. 土柱分3层层状土，自上而下依次为杨陵塿土、渭河砂土与杨陵塿土，厚度分别为7.2、5.0、7.8 cm. 杨陵塿土和渭河砂土的渗流速度依次为1.58×10⁻⁴、8.33×10⁻⁴ cm/s，水动力弥散系数依次为2.67×10⁻⁵、8.00×10⁻⁵ cm²/s. 通过本研究等效模型，可以求得等效厚度为20.0 cm，等效渗流速度为2.00×10⁻⁴ cm/s以及等效水动力弥散系数为1.77×10⁻⁴ cm/s. 如图6所示，本研究等效解与该土柱实验数据较拟合.

随着GCL的广泛应用，GCL与CCL之间的等效性研究愈发重要.《生活垃圾卫生填埋场防渗系统工程技术规范》^[11]允许填埋场采用GCL来替代CCL，且GCL下方必须设置一定厚度的AL，并规定AL渗透系数须小于10⁻⁷ m/s，但未对AL厚度L_AL提出要求. HDPE膜的厚度不应小于1.5 mm，CCL的厚度不得小于75 cm^[11]. GCL的厚度通常为7 mm^[36].

在工程设计上可以通过等效计算模型对上述2类不同的复合衬垫设计结构进行等效计算，从而对GCL下伏AL厚度L_AL提出要求，以便工程设计人员在不使用复杂数值模拟软件的前提下知晓在何种情况下，这2类衬垫结构可以等效替换.

对服役时间为100 a的上述2类衬垫结构进行等效计算，等效计算从3个等效指标进行计算，分别是使得两者在服役时刻的相对浓度（C_b/C₀）、瞬时通量（J_I）、累计通量（J_A）相等.

如图7~9所示分别为不同水头下CCL、GCL衬垫的AL层厚度与相对浓度、瞬时通量、累计通量等3个等效指标的关系曲线. 如图7所示，当服役时刻t=100 a时，对任意CCL衬垫结构(GM + 0.75 m CCL + AL)中AL层的厚度，都有唯一的GCL衬垫结构(GM + 7 mm GCL + AL)中AL层的与之对应，并使两者相对浓度相等. 此外，水头由0.3 m增加至15.0 m，相对浓度也增大，GCL衬垫结构中的AL等效厚度随之减少. 当CCL衬垫AL层厚度为1.0 m，且水头为15.0 m时，GCL衬垫相应的土质衬垫厚度减少了0.066 m. 如图8所示，当服役时刻t=100 a时，CCL衬垫与GCL衬垫的瞬时通量随AL层厚度的增加而减少. CCL衬垫的底部瞬时通量小于GCL的底部瞬时通量，两者的差距随AL层厚度的增大而增大. 如图9所示，当服役时刻t=100 a时，2类衬垫的累计通量随AL层的增大而减少. 当AL层厚度相同时，CCL衬垫的累计通量小于GCL的累计通量，两者的差距随AL层厚度的增大而减小.

综上所述，可以通过图7~9计算得到相应的与CCL衬垫等效的GCL衬垫的AL层厚度. 但在实际工程中，仍希望得到一个最保守的等效指标，即在不同水头工况下，等效后的GCL衬垫结构在服役时刻的相对浓度、瞬时通量以及累计通量均小于或等于CCL衬垫结构的瞬时通量和累计通量. 因此，考察不同水头工况下3个等效指标所计算的等效AL层厚度L_e的差别.

如图10~12所示分别为在服役时刻为100 a时，GCL衬垫(GM + 7 mm GCL + L_e m AL)与CCL衬垫（GM + 0.75 m CCL + 1.0、2.0、3.0 m AL）相对浓度、瞬时通量、累计通量的等效AL层厚度. 近年来，国内在役填埋场监测发现渗滤液水位雍高甚至个别出现超高水位的现象. Zhang等^[43]发现最大压力水头高至16 m，水位监测标高最高达50 m. 张海华等^[44]监测到某填埋场主水位标高为55.85~59.30 m. 因此将等效分析的水头范围考虑为0~60 m.

如图10所示，针对1 m AL的CCL衬垫的情况，相对浓度、累计通量及瞬时通量计算所得的等效AL层厚度均随水头的增大而减小. 意味着GCL衬垫的防污能力随水头的增大而显得更重要，因此GCL的环境保护能力在高水头工况下具有优越性. 在相同水头下，基于相对浓度计算所得的等效厚度大于基于累计通量计算所得的等效厚度，大于基于瞬时通量计算所得的等效厚度，因此，相对浓度是更保守的等效指标. 采用相对浓度计算所得的等效AL层厚度也可以保证GCL衬垫的瞬时通量以及累计通量均小于CCL衬垫的瞬时通量及累计通量.

此外，通过线性拟合，由相对浓度计算所得的等效AL层厚度与水头线性负相关：

(50) ${L_{\rm{e}}} = - 0.002\;8{h_{\rm{w}}} + 1.721\;5,{\rm{ }}({R^2} = 0.993\;0).$

因此，GCL在复合衬垫中的作用随水头增大而增强. 通过比较图中曲线斜率可知，采用GCL复合衬垫对污染物瞬时通量的减少效果最好.

如图11所示，在2.0 m AL的CCL衬垫的情况下，等效AL层厚度与水头仍线性负相关，其中相对浓度等效AL层厚度与水头的线性拟合解析式为

(51) ${L_{\rm{e}}} = - 0.010\;4{h_{\rm{w}}} + 2.649\;4,{\rm{ }}({R^2} = 0.991\;4).$

对比式（50）、（51）以及图10、11可知，当CCL衬垫AL层由1.0 m增至2.0 m后，GCL衬垫对污染物的削减效果更强.

如图12所示，在3.0 m AL的CCL衬垫情况下，对比式（50）、（51），GCL衬垫对污染物的阻隔效果会随水头增加以及土质衬垫厚度的增加而更好. 等效AL层厚度与水头仍线性负相关，其中相对浓度等效AL层厚度与水头的线性拟合解析式如下：

(52) ${L_{\rm{e}}} = - 0.016\;9{h_{\rm{w}}} + 3.575\;9,{\rm{ }}({R^2} = 0.991\;9).$

综上，从工程设计安全角度考虑，当服役时间为100 a时，采用相对浓度的等效方法能保证等效后的GCL衬垫的瞬时通量与累计通量均小于CCL衬垫的瞬时与累计通量. 因此，下面给出基于100 a击穿时间的2类衬垫的等效对比.

如表2所示，CCL复合衬垫的相对浓度、瞬时通量及累计通量是通过式（32）、（33）等效计算所得. 如表3所示，GCL复合衬垫的相对浓度、瞬时通量及累计通量也是通过式（32）、（33）求得，并通过试算AL层厚度使其与表2中底部相对浓度相等.

当CCL + AL总厚度为1.75 m时，等效的GCL + AL总厚度为1.56~1.73 m，采用GCL复合衬垫时，总厚度可以减少0.02~0.19 m. 当CCL + AL总厚度为2.75 m时，采用GCL复合衬垫，总厚度可以减少0.08~0.70 m. 当CCL + AL总厚度为3.75 m时，采用GCL复合衬垫，总厚度可以减少0.14~1.15 m. 因此，随设计水头与AL厚度增加，GCL复合衬垫具有更好的防污性能.

基于Cl⁻击穿时间100 a这一标准，当水头为15.0 m时，采用GCL复合衬垫可以使总厚度减少0.45 m. 实际垃圾填埋场的水头是变化的，从安全设计的角度，采用GCL衬垫结构应选用AL层厚度最大的情况进行等效设计. 此外，当CCL施工难度大或无法就地取材时，可以参考采用上述等效厚度的GCL + AL替代CCL + AL.

如图13所示，在相同的水头下，通过等效计算使得在服役时刻的GM + 0.75 m CCL + 1.0 m AL与(L_e−0.75) m CCL + 1.0 m AL的相对浓度相等，从而获得GM的等效CCL厚度L_e. 通过拟合，GM的等效CCL厚度与服役时刻呈对数线性关系：

(53) $\lg {L_{\rm{e}}}{\rm{ = }}\left\{ {\begin{array}{*{20}{c}} {0.619\;0\lg \;t - 1.249\;8,{R^2} = 0.994\;5,{h_{\rm{w}}} = 0.3\;{\rm{ }}{\rm{m}}} ;\\ {0.511\;8\lg\; t - 0.295\;8,{\rm{ }}{R^2} = 0.999\;5,{h_{\rm{w}}} = 5.0\;{\rm{ }}{\rm{m}}}; \\ {0.496\;1\lg \;t - 0.011\;8,{R^2} = 0.999\;8,{h_{\rm{w}}} = 15.0\;{\rm{ }}{\rm{m}}.} \end{array}} \right.$

随着服役时刻的增加，GM的等效CCL厚度增大，意味着服役时间越长，土工膜的防污能力更为重要. 当服役时间为100 a、水头为0.3 m时，GM的等效CCL厚度为0.94 m；当水头为5.0 m时，等效CCL厚度为5.29 m；当水头为15.0 m时，等效CCL厚度为9.49 m. 可见在高水头时土工膜在填埋场地下水土环境保护方面的重要性.

POLLUTE v7是基于污染物在复合衬垫中的半解析模型进行数值模拟求解. POLLUTE v7中多层土质衬垫之间的边界条件采用浓度连续及通量连续条件. 本研究等效模型考虑了浓度连续，将每层土的污染物运移过程近似以式(11)表达，利用卷积性质通过累乘式(17)，求得拉式域的浓度解. 这样处理是为了得到较为简化的等效渗流速度和扩散系数表达式. 将本研究解和POLLUTE v7的结果及国内外土柱实验结果进行对比，发现本研究解与数值分析软件和实验结果的均较拟合，也说明本研究模型假设的合理性及计算结果的可靠性.

本研究底部衬垫浓度解析解（即式(34)）基于半无穷的Dirichlet或Neumann底部边界条件求得. 在填埋场衬垫初步设计和评价时，在现场水文地质条件掌握不充分的情况下，采用本研究解能进行更方便和高效率的计算.

（1）给出基于时间矩卷积的等效分析模型，并在此基础上给出多层土质衬垫的等效渗流速度和水动力弥散系数的表达式. 采用填埋场污染物运移分析程序POLLUTE v7及国内外土柱实验验证该等效模型的合理性与可靠性.

（2）从相对浓度、瞬时通量以及累计通量3个角度对GM + CCL + AL及GM + GCL + AL这2类不同复合衬垫进行等效计算，发现采用相对浓度的等效方法最为合适，也能保证等效后的GCL衬垫的瞬时通量与累计通量均小于CCL衬垫的瞬时与累计通量. 通过拟合，发现等效AL层厚度与水头呈负相关关系，GCL在复合衬垫的中的作用随水头增大而增长. 采用GCL复合衬垫对污染物瞬时通量的减少效果最好.

（3）基于Cl⁻击穿时间100 a这一标准，当水头为15.0 m时，采用GCL复合衬垫可以使土质衬垫的总厚度减少0.45 m.

（4）通过相对浓度等效，得到GM的等效CCL厚度，发现土工膜防污能力随水头增加而增加，并与服役时间呈对数线性增长关系. 当水头为15.0 m时，GM的等效CCL厚度为9.49 m.

（5）该等效模型可以对成层土质衬垫进行简化计算，为工程设计提供参考，便于一线工程设计人员应用. 本研究成果可以为填埋场衬垫结构初步设计、等效性评价和拟合试验数据等提供参考. 由于本研究等效模型中未考虑污染物吸附及生物降解的作用，应进一步推导包含污染物吸附与生物降解作用的等效模型.