多孔介质曲折度对膨润土衬垫渗透性能的影响

侯娟 滕宇阳 李昊 刘磊

摘要：膨润土的颗粒膨胀所引起的孔隙率及流径曲折度的改变是影响膨润土衬垫（Geo? synthetic Clay Liner，GCL）防渗性能的关键.通过 COMSOL 建立模拟膨润土中渗流的数值模型，研究了膨润土颗粒膨胀过程对 GCL 防渗性能的影响. COMSOL模型研究结果表明，流体通过膨润土粒间孔隙时，边界通道的最大速度明显小于非边界通道速度，这与泊肃叶公式的基本结论一致，初步验证了模型的准确性.同时，微观粒子速度分布规律与宏观流体场的速度分布规律基本一致，因此可以从微观的角度研究GCL 的宏观渗透规律.利用 COMSOL 中的粒子轨迹监测方法，分析了颗粒膨胀对膨润土孔隙率和流径曲折度的影响.研究结果表明，流径的曲折度随孔隙率的增加而减小，且呈指数函数关系.基于 GCL 试样的多孔介质毛管模型和修正泊肃叶公式，考虑 GCL 孔隙率与流径曲折度之间的指数函数关系，提出了能综合反映孔隙率与流径曲折度影响的 GCL渗透系数理论预测模型，并将理论预测结果与已有文献试验数据进行了对比，两者的比值介于1/5～5，验证了模型的准确性.

关键词：膨润土衬垫;渗透系数;COMSOL;曲折度;理论模型

中图分类号：TU411.4  文献标志码：A

Influence of Porosity and Tortuosity on Hydraulic Conductive of Geosynthetic Clay Liner

HOU Juan1，2，3?，TENG Yuyang1，LI Hao1，LIU Lei3

（1. School of Mechanics and Engineering Sciences，Shanghai University，Shanghai 200444，China;

2. School of Engineering，University of Virginia，Charlottesville，VA 22904，USA;

3. State Key Laboratory of Geomechanics and Geotechnical Engineering，Institute ofRock and Soil Mechanics，Chinese Academy of Sciences，Wuhan 430071，China）

Abstract：The change of porosity and the tortuosity induced by the swelling progress of Geosynthetic Clay Liner （GCL） influence the hydraulic conductive obviously. The model simulating the seepage in bentonite is established by COMSOL to investigate the effect of the expansion process of bentonite particles on the hydraulic conductive of GCL.The COMSOL simulation results showed that the maximum velocity of the boundary channel is less than that of the non-boundary channel， when the liquid passes through the pores between the bentonite particles， which is consistent with the Poiseuille equation and initially verifies the accuracy of the model. The microscopic particle velocity distri ? bution law is consistent with the velocity distribution law of macroscopic fluid field. Therefore， the macroscopic hy? draulic law can be investigated in the microscopic aspect. Meanwhile， the particle expansion on the porosity and the tortuosity of bentonite is analyzed by using the particle trajectory monitoring method in COMSOL. The results showed that the tortuosity decreases with the increase of porosity， and it is exponentially functional. Based on the capillary model of porous media and the Poiseuille Formula， the theoretical prediction model of the hydraulic conductive of GCL， which can reflect the influence of porosity and the tortuosity， is proposed， considering the exponential function relationship between the porosity and the tortuosity. The theoretical calculation prediction model for predicting the hydraulic conductive of GCL is proposed. The ratio of the theoretical and experiments is between 1/5 and 5， which confirms the accuracy of the theoretical model of GCL.

Key words：Geosynthetic Clay Liner（ GCL）;hydraulic conductivity;COMSOL;tortuosity;theoretical model

膨润土衬垫（Geosynthetic Clay Liner，GCL）是由两层土工织物间夹裹一层薄的膨润土构成的，近几十年来广泛应用于垃圾填埋场的顶部封场及底部衬垫系统，具有防渗性能好、厚度小以及自愈力强等优点.目前对于 GCL渗透系数的研究，主要集中在宏、微观试验研究方面[1-2]，Shackelford 等[3]对 GCL 进行了大量渗透试验，发现影响 GCL 渗透性的因素有孔隙比、膨润土级配、厚度等.刘志彬等[4]对膨润土内部微孔隙进行了定量分析，并与扫描电镜图像进行了对比.王宝等[5]通过试验发现孔隙比与渗透系数存在良好的线性关系，该结论与 Kang 等[6]用去离子水渗透 GCL结果相同. Mesri和 Olson[7]发现，渗透系数和孔隙比在对数坐标系下存在线性关系. Frataloc? chi[8]研究了应力对 GCL 渗透系数的影响，发现渗透系数的对数与孔隙比之间存在较好的线性关系，得到与 Petrov 等[9]相似的结论.何俊等[10]通过微观结构分析，推导了 GCL 的膨润土饱和渗透系数与孔隙比之间的计算公式.

然而，仅用孔隙率往往存在不能充分体现土体内部粒间孔隙特征等问题，也不能反映颗粒尺寸本身对渗透系数的影响[11].自Kozeny[12]提出流径曲折度概念以来，Carman[13]考虑孔隙通道的不规则性，在 Poiseuille 定律的基础上，结合Kozeny公式得到了曲折度与渗透系数之间的表达式. Nooruddin 和 Hos? sain[14]根据曲折度与土颗粒之间的关系，改进了Kozeny-Carman 模型. Tsang[15]研究了渗流路径的曲折度对流场的影响. Murata 等[16]基于孔隙介质的渗透系数Kozeny-Carman 方程，推导得到了修正的流量关系表达式.此外，还有众多学者采用试验研究[17]、理论分析[10，18]以及数值模拟[4]等手段研究了多孔介质孔隙流径的曲折度与孔隙率之间的关系.

GCL 中起主要防渗作用的膨润土，是一种颗粒遇水会膨胀的特殊多孔介质.然而，尽管目前已有大量研究表明 GCL 中膨润土孔隙率与渗透系数呈正相关的关系[19]，但是，如何从细观角度考虑膨润土颗粒尺寸及孔隙特征对 GCL 渗透性能的影响，成为进一步研究和解释膨润土渗透机理的难点[9，20]. GCL 中曲折度是影响渗透系数的关键因素之一[21]，但对其膨胀过程中量化分析的研究比较鲜见，目前尤其缺乏对孔隙渗流通道的曲折度与孔隙率之间关系的研究.随着计算机技术的高速发展，数值模拟成为目前研究多孔介质曲折度的主要手段之一[22].其中， COMSOL在分析多物理场耦合以及在多孔介质中流体流动等复杂问题方面具有明显的优势[23]，因此，本文采用 COMSOL 软件，构建了包含孔隙率与曲折度等细观变量的多孔介质 GCL结构模型，从细观角度，研究土体的颗粒尺寸及孔隙特征等对 GCL渗透性能的影响.首先，通过构建流体物理场和粒子追踪场，模拟分析蠕动流场的稳态速度图和压力矢量图，得到宏观流场分布规律;然后，通过对比宏观蠕动流场与微观流动（水）粒子分布规律，探究宏观流场分布与微观流动粒子之间的关系;同时，通过实时追踪液体粒子的运移过程，得到液体在多孔介质内流径并通过数学变换得到曲折度与流径曲折度的表达式;最后，采用毛管模型描述粒间孔隙组成的土骨架[24]，并基于此提出能够综合考虑孔隙率与流径曲折度的 GCL渗透系数理论预测模型，验证模型的准确性.

1模型建立和参数确定

COMSOL 以有限元为基础，通过求解偏微分方程（单场）和方程组（多场）实现多场耦合的仿真模拟[25] .建模的过程一般分为模型假设、确定物理场控制方程、几何模型建立、边界与初始条件确定、网格划分以及最终求解6个部分.

1.1模型假设

参考张志红等[26]的研究，本文基本假设为：多孔介质内部为均一几何孔隙;过程等温;渗流液均质各向同性;忽略毛细管吸力作用;流体不可压缩且为单相稳定渗流.

1.2物理场控制方程1）流体物理场

由于 GCL 的渗透系数较小，在膨润土粒间孔隙通道中，形成斯托克斯流（Stokes Flow），由 Navier- Stokes连续性方程得蠕动流方程：

式中：v 为动力黏度;ρ为液体密度;P 为流体压力;u 为流体的速度矢量;?为拉普拉斯算符.

2）粒子追踪场

采用了欧拉-拉格朗日法对运动颗粒进行追踪，同时采用牛顿第二定律描述运动过程.欧拉-拉格朗日法采用式（2）记录每个粒子的不同状态.

式中：mp为颗粒质量;vp为颗粒速度;F 为颗粒受力.

1.3模型验证

首先采用宏观试验对 COMSOL数值模型进行验证.依据 Petrov 等[9]试验结果，建立 COMSOL 数值模型，并采用稳态对蠕动流场进行模拟，得到渗透系数分布图如图1所示.计算得到渗流出口的渗透系数 k 为1.57×10-11 m/s，与试验渗透系数比值为1.121.同时，分别建立其他试验结果的 COMSOL数值模型，得到的渗透系数与试验对比汇总如表1所示.由表1可知，该 COMSOL 数值模型能够较好地模拟 GCL 在去离子水渗透时的工况.其次采用宏观流场对 COMSOL 数值模型进行验证.在本文的仿真模拟过程中，宏观流场（图1所示）的规律与金磊等[27]的研究结果类似，即在孔隙通道半径大处流速增加明顯.同时，本文流场速度随孔隙通道半径增大而增大，其结果与梁越等[28]的研究结论一致.

1.4几何模型建立

1.3节验证了 COMSOL 中的蠕动流场可对 GCL 的层流运动进行宏观模拟.由于在孔隙率不变时，等效渗透系数与孔隙尺寸的平方成正比[30-31]，考虑到有限元计算的时效性，采用几何相似性理论，将验证模型尺寸进行适当扩大，进行下列计算.同时，也有研究发现颗粒的几何外形[32]对其渗透系数的影响并不显著，因此，本模型构建边长为10 mm 的大正方形表示 GCL试样，交错布置9×9边长为1 mm 的小正方形表示膨润土颗粒，在孔隙率 n=0～1范围内进行优化计算.

已有学者[33-35]研究发现，在膨润土膨胀过程中的晶体膨胀阶段，是膨润土固相颗粒的膨胀间距不断增大的膨胀过程，即颗粒几何尺寸扩大而孔隙率减小的过程.同时采用正方形的放大命令，确保模型膨胀前后几何外形一致，也可保证模型达到所需的最小孔隙率，实现模拟孔隙率 n=0～1的全区域变化过程.因此通过 COMSOL 内逐渐放大小正方形的命令模拟颗粒膨胀过程，从而控制几何模型的孔隙率逐渐变小，部分工况如图2所示.选用蠕动流模块和流动粒子追踪模块进行瞬时计算，通过实时监测流动粒子的位置坐标，最终确定液体在多孔介质内流通的渗流路径.

1.5 边界与初始条件确定

蠕动流模块边界条件设定为小正方形（颗粒）外围及大正方形（膨润土试样）上下边界均为封闭边界，流体受压力 P 作用下，从左边界流入，右边界单向流出.

流动粒子追踪模块边界条件为粒子均匀排布在左边界，由蠕动流中的速度场控制粒子运动至右侧壁，最后附着在右侧壁上后其速度为0.上下边界、小正方形外围边界均为保持粒子动量守恒的粒子反弹壁，使得粒子在接触到上述边界后立即反弹，速度大小不变，不发生动能损失.整体 COMSOL 模型边界条件如图3所示.

工程实际中，垃圾填埋场中 GCL 衬垫上部渗滤液收集系统及第一层垃圾体上覆荷载一般约为14 kPa[36-37].因此，模型取入口边界压力为 P=14 kPa，该条件同时满足蠕动流场雷诺数小于1的要求.模型参数详见表2.

1.6 网格划分

考虑时效性，参考 Zimmerman[25]的研究，采用三角形网格进行划分计算，并分别对 n=0～1之间不同孔隙率的多孔介质模型进行模拟.

值得一提的是，在 GCL 渗透实验中，测得 GCL 试样的初始孔隙率一般为0.6[5，37].但当垃圾填埋场开始服役时，GCL 衬垫在吸水膨胀与上覆荷载作用下，孔隙率的数值将明显减小[38-39].同时，Seiphoori等[40] 对 MX-80膨润土试样进行膨胀试验后表明，不同压实状态下膨润土的孔隙率介于0.34～0.58之间，因此，以下主要以孔隙率 n=0.41的工况为例进行分析计算，具体网格划分如图4所示.

2仿真结果与分析

2.1 蠕动流计算结果

基于 COMSOL 软件自身的求解器的特性，一般认为稳态求解可用以验证瞬时求解的结果是否趋于稳定，而瞬时求解有利于分析宏观试验现象随时间的变化趋势，因此，将蠕动流场和粒子追踪模块进行耦合瞬时计算.

首先，采用稳态对蠕动流场进行模拟.入口边界选取充分发展的流动命令以减少边界效应，即当流体到达左边界（入口）时，水流趋于平稳，这样可以有效地减少因为边界加压而导致的入口处水流速度不均匀的现象.图5所示为模拟试样中流体的速度图与压力矢量图.在速度矢量图（图5（a））中，可观测到 A、B、D 点处路径的最大速度小于 C 點处最大速度.分析原因是 C点的水流速度源自上、下两侧竖向通道，类似于通道半径在此处增加.因此，根据泊肃叶公式，通道半径增加会导致水流速度增加，即通道个数直接影响流速.同时，由压力矢量图（图5（b））可见，粒间孔隙压力从左到右（出口）出现均匀递减的趋势.此外，图5（b）左边界下部压力分布（虚线）比上部压力分布（实线）更紧密，主要原因是几何模型孔隙通道布置不均匀，导致压力在上、下边界上出现不均匀分布.

2.2 流动粒子追踪

在瞬时状态下，将蠕动流场和粒子追踪模块进行耦合计算.模拟过程中发现粒子追踪场在20000μs 时基本达到稳定流的速度状态，因此，仅监测流动粒子在 t=0～20000μs 内的位置坐标.图6给出了部分时刻流动粒子轨迹图.

对比图5（a）和图6（a）可知，微观流动粒子分布呈现与宏观流体场分布类似的规律.图6（a）中，E、F 点路径内的最大速度小于 G 点的最大速度（实线所示），分析其原因是 G点等效通道半径增加，导致速度增加;在相同的原因下，导致 H点的最大速度小于 I 点最大速度（虚线所示）.因此，图6（a）中宏观蠕动流场是细观粒子流动场作用的结果.

但同时对比发现，图6（a）中蠕动流场平均速度为1.53 m/s，因此，如果不考虑曲折度影响，粒子按直线规律运移的话，应该在 t=6535μs 时到达右边界.然而，如图6（d）所示，实际粒子运移中，当 t=6535μs 时，其仅运移到距右边界一定位置处（图6（d）中虚线所示），而在 t=8000μs（约为整体运移时间的前三分之一）时，大部分粒子才到达右边界，少数粒子甚至需要比这更长的时间（如图6（e）所示）.这是由于流径具有一定的曲折度，使得粒子实际运移路径大于直线距离.此外，观察图6（b）～（e）（粒子运移过程）可发现，在流体速度控制下，流动粒子先从左边界均匀分布，通过粒间孔隙通道，通道较窄，速度增加，随后发生分流，在上、下侧通道内速度减缓;继续通过粒间孔隙通道，通道收窄，速度再次增加.依次反复后，流动粒子呈不均匀状黏附在出口处.

2.3 GCL 渗透系数理论预测模型

GCL 衬垫中起主要防渗作用的膨润土是一种多孔介质. Kozeny[12]和 Carman[13]提出了预测多孔介质渗透率的半经验公式（Kozeny-Carman （ K-C ）方程），该式适用于多孔介质渗透系数较低的情况.首先， Nooruddin 和 Hossain[14]在 K-C 方程基础上，将土体通道简化为毛管模型[41]（图7）.长度为 L、截面面积为 A 的多孔介质土体（图7（a）），此时的固相颗粒为不规则分布，由孔隙组成的渗流通道也不均匀分布，均质土体的固相颗粒一般为土颗粒;采用等流量条件下渗流通道的等效毛管化，将土颗粒规则排布后，形成土颗粒以外的孔隙通道，等效后为直径2R 的通道（图7（b））;但由于土体的实际渗流路径如图7（c）中实线所示，即在土颗粒与土颗粒之间形成曲折的路径;假设该路径中孔隙通道直径为2r （ r≤R），如图7（d）所示，以此形成等效毛管模型来表示液体在多孔介质中的移动规律.

同时，由于土体的渗流路径 L 具有一定的曲折性. Collins[41]采用曲折度τ来反映多孔介质中平均有效流动路径长度 Li （m）与流体流动方向的直线路径 L （m）的比值，如公式（3）所示.不少学者[42-47]对曲折度τ进行进一步研究后发现，孔隙率与曲折度的关系可通过数学表达式进行表述，如表3所示.

由于膨润土具有自相似性，因此可将土颗粒通过单位体积法简化为单元体.同时，将简化后的单元体表示为毛管模型，以便通过渗流路径获得渗流的等效理论模型.本文等效模型作如下假设：假设在等效渗透通道模型中，总共有 M 根半径为ri、长度分别为 Li 的毛管束.结合修正的哈根-泊肃叶[48-49]方程可知，单位时间内通过 M 根毛管的渗流速度 q （m3/s）为：

式中：μ为渗滤液的动力黏度，Pa·s;ρw 为液体的密度，kg/m3;g 为重力加速度，m/s2;π=3.14159;ri为在 M 根弯曲毛管中，第i根毛管的通道半径，m;Li 为第i根毛管的渗流长度，m.

将式（3）代入式（4），得通过土体的总流量 Q （m3/s）为：

又由达西定律可知，水流通量 J （m/s）与水力梯度ΔH/L 成正比：

式中：k 为渗透系数，m/s;ΔH 为水力梯度;L 为渗流路径，m.

由于水流通量 J 是总流量 Q 与土体截面面积 A （m2）的比值（J=Q/A），将式（6）代入式（5），换算可得渗透系数 k 为：

值得一提的是，尽管式（7）中的曲折度τ是基于多孔介质本身性质的参数，其难于用试验手段直接测得，但可以通过孔隙率 n 计算获得（见表3）.

因此，本文孔隙率 n 从0到1的一系列 COMSOL 模型，采用流动粒子追踪方法获得孔隙率全范围内流动粒子的运动轨迹，并统计计算相应的渗流路径以及曲折度，模型选取的试验验证参数如表4所示.由此计算得到不同孔隙率与曲折度之间的关系（图8中散点标示），可以看出，曲折度随着孔隙率的增加而减小.

由拟合模拟的一系列工况所得结果（图8中实线标示），可得到孔隙率与曲折度之間的数学表达式为：

将式（8）代入式（7）中，即可得到能综合反映孔隙率与流径曲折度影响的 GCL 渗透系数理论预测模型：

为验证式（9）的准确性，本文采用大量已有试验结果（Chai 等[50]、Jo 等[51-52]、Petrov 和 Rowe[53]、Chen 等[54]、Setz等[55]、Bradshaw 和 Benson[56]以及Seiphoor[57]）进行计算比较，结果如图9所示.以 Chai 等[50]数据为例，具体的计算过程如下.

由式（9）可得此时 GCL 的预测渗透系数 k：

本文采用对数坐标系下渗透系数偏差分析的表示方法[58-61]（如图9所示）.图中，横坐标为预测渗透系数值，纵坐标为试验渗透系数值.从图中可看出，预测值与试验值的范围为0.23～3.86.按照行业标准《钠基膨润土复合防水衬垫》（FZ/T 64036—2103），渗透系数值的误差常在一个数量级内进行表征，如规定渗透系数≤5×10-9～1×10-10 cm/s.同时，对数坐标系下渗透系数的偏差一般介于1/10～10[37，58-61].本文 GCL 的理论预测模型计算结果与试验值的比在1/5～5，小于一个数量级，说明该 GCL 理论预测模型可以较为准确地用于预测 GCL 的渗透系数，可以指导工程实际.

3结论与展望

本文通过 COMSOL模拟了 GCL在膨胀变形过程中，曲折度随孔隙率的变化而变化的规律，得到多孔介质渗流流径的细观模型，研究了孔隙率和流径曲折度对 GCL渗透系数的影响，主要结论如下：

1）当入口边界条件为均匀压力时，在蠕动流场中，流体通道半径的增加会导致流体速度增加，同时，流体压力从入口到出口处呈现均匀递减的趋势.

2） GCL试样中，微观流动粒子分布与宏观流体场分布规律基本一致，且大部分流动粒子在整体时间的前1/3段完成运移.

3） GCL试样中，流径的曲折度随孔隙率的增加而减小，且呈指数函数关系.

4）提出的能综合反映孔隙率与流径曲折度影响的 GCL 渗透系数理论预测模型，可以较为准确地预测试验结果，两者比值介于1/5～5.

需要注意的是，本文膨润土曲折度与孔隙率之间的关系主要是基于模拟均一几何尺寸和均一几何孔隙表征膨润土试样得到的.因此，研究结果可以从概念上定性地分析曲折度随膨润土颗粒膨胀的变化规律，同时反映曲折度对 GCL 渗透系数的影响.但是，后续还需要对大量不规则试样展开进一步的深入和系统的研究，才有可能将研究结果推广到定量分析以及实际应用中.

参考文献

[1]  BOUAZZA A. Geosynthetic clay liners[J]. Geotextiles and Geomembranes，2002，20（1）：3-17.

[2]  ROWE R K，MUKUNOKI T，LI M H，et al. Effect of freeze-thawon the permeation of arctic diesel through a GCL[J]. Journal of ASTM International，2004，1（2）：11741.

[3]  SHACKELFORD C D，BENSON C H，KATSUMI T，et al. Evaluating the hydraulic conductivity of GCLs permeated with non- standard liquids[J]. Geotextiles and Geomembranes，2000，18（2/3/4）：133-161.

[4] 劉志彬，施斌，王宝军.改性膨胀土微观孔隙定量研究[J].岩土工程学报，2004，26（4）：526-530.

LIU Z B，SHI B，WANG B J. Quantitative research on micropores of modified expansive soils[J]. Chinese Journal of Geotechnical Engineering，2004，26（4）：526-530.（In Chinese）

[5] 王宝，董兴玲.不同有效应力下矿山渗滤液对土工合成黏土衬垫渗透特性影响的试验研究[J].岩土力学，2017，38（5）：1350-1358.

WANG B，DONG X L. Hydraulic conductivity of mine leachate through geosynthetic clay liners under different effective stresses [J]. Rock and Soil Mechanics，2017，38（5）：1350-1358.（ In Chinese）

[6]  KANG J B，SHACKELFORD C D. Consolidation of a geosynthetic clay liner under isotropic states of stress [J]. Journal of Geotechnical and Geoenvironmental Engineering，2010，136（1）：253-259.

[7]  MESRI G，OLSON R E. Mechanisms controlling the permeabilityof clays[J]. Clays and Clay Minerals，1971，19（3）：151-158.

[8]  FRATALOCCHI E. Hydraulic perfomance of GCLS with differentstress  history [ C]//Intemational  Workshop “Hydra-physica- mechanies of Landfills”. Grenoble，France：University of Joseph Fourier，2005：1-4.

[9]  PETROV R J，ROWE R K，QUIGLEY R M. Selected factors influencing GCL hydraulic conductivity[J]. Journal of Geotechni? cal and Geoenvironmental Engineering，1997，123（8）：683-695.

[10]何俊，施建勇.膨润土中饱和渗透系数的计算[J].岩石力学与工程学报，2007，26（ S2）：3920-3925.

HE J，SHI J Y. Calculation of saturated permeability coefficient of bentonite[J]. Chinese Journal of Rock Mechanics and Engi? neering，2007，26（ S2）：3920-3925.（In Chinese）

[11] PITTMAN E D. Relationship of porosity and permeability to various parameters derived from mercury injection capillary pressure curves for sand-stones [J]. AAPG Bulletin ， 1992，76（2）：191-198.

[12] KOZENY J. Uber kapillareleitung des wassersimboden[J]. Sitzungsberichte der Wiener Akademie der Wissenschaften，1927，136（2a）：271-306.

[13] CARMAN P C. Fluid flow through granular beds[J]. ChemicalEngineering Research and Design，1997，75：S32-S48.

[14] NOORUDDIN H A，HOSSAIN M E. Modified Kozeny-Carmencorrelation for enhanced hydraulic flow unit characterization[J]. Journal of Petroleum Science and Engineering，2011，80（1）：107-115.

[15] TSANG Y W. The effect of tortuosity on fluid flow through a singlefracture [J]. Water Resources Research ，1984，20（9）：1209-1215.

[16] MURATA S，MITSUISHI H，SAITO T. Characterization of fracture permeability by using a fractal model[ C]// SPE Asia Pacific Oil and Gas Conference and Exhibition. Melbourne，Australia： SPE，2002：1-8.

[17] DIAS R，TEIXEIRA J A，MOTA M，et al. Tortuosity variation in alow density binary particulate bed[J]. Separation & PurificationTechnology，2006，51（2）：180-184.

[18] XU P，YU B M. Developing a new form of permeability andKozeny-Carman constant for homogeneous porous media by means of fractal geometry[J]. Advances in Water Resources，2008，31（1）：74-81.

[19] CARMAN P C. Flow of gases through porous media[M]. London：Butterworths Scientific Publications，1956：182.

[20] SHAN H Y，LAI Y J. Effect of hydrating liquid on the hydraulicproperties of geosynthetic clay liners [J]. Geotextiles and Geo? membranes，2002，20（1）：19-38.

[21] ROWE R K. Geosynthetic clay liners：Perceptions and misconceptions [J].  Geotextiles and Geomembranes，2020，48（2）：137-156.

[22]郭建春，庄园，刘超.考虑非达西效应的酸蚀裂缝流场数值模拟[J].岩土力学，2015，36（11）：3315-3321.

GUO J C，ZHUANG Y，LIU C. Numerical simulation of flow field of acid etched fractures considering non-Darcy effect[J]. Rock and Soil Mechanics，2015，36（11）：3315-3321.（In Chinese）

[23]侯会明，胡大伟，周辉，等.考虑细观等效热学参数的高放废物处置库围岩应力-渗流-温度耦合模拟方法[J].岩土力学，2019，40（9）：3625-3634.

HOU H M，HU D W，ZHOU H，et al. Thermo-hydro-mechanical coupling simulation method of surrounding rock in high-level ra? dioactive waste repository considering effective meso-thermal pa ? rameters [J]. Rock and Soil Mechanics，2019，40（9）：3625-3634.（In Chinese）

[24] KOVACS G. Seepage Hydraulics[ M]. Amsterdam：Elsevier Scientific Publishing Company，1981：71.

[25] ZIMMERMAN W B J. COMSOL Multiphysics 有限元法多物理建模与分析[M].北京：人民交通出版社，2007：89-113.

ZIMMERMAN W B J. COMSOL Multiphysics modeling with finite element methods [ M]. Beijing： China Communications Press，2007：89-113.（In Chinese）

[26]张志红，师玉敏，朱敏.黏土垫层水力–力学–化学耦合模型研究[J].岩土工程学报，2016，38（7）：1283-1290.

ZHANG Z H，SHI Y M，ZHU M. Coupled hydro-mechanical- chemical model for clay liner[J]. Chinese Journal of Geotechni? cal Engineering，2016，38（7）：1283-1290.（In Chinese）

[27]金磊，曾亚武，程涛，等.基于 LBM-IMB-DEM 分析岩土多孔介质的渗透淤堵特性[J].岩土工程学报，2021，43（5）：909-917.

JIN L，ZENG Y W，CHENG T，et al. Seepage clogging character? istics of rock and soil porous media using LBM-IMB-DEM simula? tion method[J]. Chinese Journal of Geotechnical Engineering，2021，43（5）：909-917.（In Chinese）

[28]梁越，陈鹏飞，林加定，等.基于透明土技术的多孔介质孔隙流动特性研究[J].岩土工程学报，2019，41（7）：1361-1366.

LIANG Y，CHEN P F，LIN J D，et al. Pore flow characteristics of porous media based on transparent soil technology[J]. Chinese Journal of Geotechnical Engineering，2019，41（7）：1361-1366.（In Chinese）

[29] AHN H S，JO H Y. Influence of exchangeable cations on hydraulic conductivity of compacted bentonite [J]. Applied Clay Sci? ence，2009，44（1/2）：144-150.

[30]古成中，吴新跃.有限元网格划分及发展趋势[J].计算机科学与探索，2008，2（3）：248-259.

GU C Z，WU X Y. A review of FEM and trend of development [J]. Journal of Frontiers of Computer Science & Technology，2008，2（3）：248-259.（In Chinese）

[31]李守巨，劉迎曦，孙伟，等.多孔材料孔隙尺寸对渗透系数影响的数值模拟[J].力学与实践，2010，32（4）：12-17.

LI S J，LIU Y X，SUN W，et al. Numerical simulation of relation ? ship between hydraulic conductivity and porosity size for porous geotechnical materials[J]. Mechanics in Engineering，2010，32（4）：12-17.（In Chinese）

[32]周力沛，唐晓武，程冠初，等.双尺度渐进展开计算黏土渗透率影响因素研究[J].岩土工程学报，2018，40（7）：1205-1211.

ZHOU L P，TANG X W，CHENG G C，et al. Influencing factors for calculating clay permeability using asymptotic expansion method[J]. Chinese Journal of Geotechnical Engineering，2018，40（7）：1205-1211.（In Chinese）

[33] LAIRD D A. Influence of layer charge on swelling of smectites[J]. Applied Clay Science，2006，34（1/2/3/4）：74-87.

[34] LIU L C. Prediction of swelling pressures of different types ofbentonite in dilute solutions[J]. Colloids and Surfaces A：Physi? cochemical and Engineering Aspects，2013，434：303-318.

[35]项国圣.含盐环境中膨润土的膨胀理论及膨胀衰减机理研究[ D].上海：上海交通大学，2007：15-19.

XIANG G S. Theory of swelling propertied and mechanism on swell attenuation of bentonite in solution[ D]. Shanghai：Shanghai Jiaotong University，2007：15-19.（In Chinese）

[36] SARABIAN T，RAYHANI M T. Hydration of geosynthetic clayliners from clay subsoil under simulated field conditions [J]. Waste Management，2013，33（1）：67-73.

[37]王宝，陈彬，窦桐桐，等.水泥固化体淋滤液对 GCL 防渗性能的影响[J].岩土工程学报，2019，41（2）：390-396.

WANG B，CHEN B，DOU T T，et al. Influences of stabilization/ solidification product leachates on hydraulic performance of geo ? synthetic clay liners[J]. Chinese Journal of Geotechnical Engi? neering，2019，41（2）：390-396.（In Chinese）

[38] SUMMER M E. Interdisciplinary aspects of soil science. Handbook of Soil Science. Section G[ M]. Boca Raton USA：CRC Press，2000：8.

[39]何俊.基于微观分析的填埋场 GCL 中运移规律的理论研究[ D].南京：河海大学，2006：45-76.

HE J. Theoretical study of transport through GCL in landfill based on microscopic analysis [ D]. Nanjing：Hohai University，2006：45-76.（In Chinese）

[40] SEIPHOORI A，LALOUI L，FERRARI A，et al. Water retentionand swelling behaviour of granular bentonites for application in Geosynthetic Clay Liner （ GCL） systems [J]. Soils and Founda? tions，2016，56（3）：449-459.

[41] COLLINS R E. Flow of fluids through porous materials [ M]. Reinhold，New York：Reinhold Publishing Corp，1961：250-256.

[42]孙军杰，田文通，刘琨，等.基于泊肃叶定律的土体渗透系数估算模型[J].岩石力学与工程学报，2016，35（1）：150-161.   SUN J J，TIAN W T，LIU K，et al. Estimation model of soil perme ? ability coefficient based on Poiseuille′s law[J]. Chinese Journal of Rock Mechanics and Engineering，2016，35（1）：150-161.（ In Chinese）

[43] IVERSEN N，J?RGENSEN B B. Diffusion coefficients of sulfateand methane in marine sediments：Influence of porosity[J].Geo? chimica et Cosmochimica Acta，1993，57（3）：571-578.

[44] PISANI L. Simple expression for the tortuosity of porous media[J]. Transport in Porous Media，2011，88（2）：193-203.

[45] KOPONEN A，KATAJA M，TIMONEN J. Tortuous flow in porousmedia[J]. Physical Review E，1996，54（1）：406-410.

[46] ATTIA A M. Effects of petrophysical rock properties on tortuosityfactor[J]. Journal of Petroleum Science and Engineering，2005，48（3/4）：185-198.

[47] KHABBAZI A E，HINEBAUGH J，BAZYLAK A. Determining theimpact of rectangular grain aspect ratio on tortuosity-porosity cor? relations of two-dimensional stochastically generated porous media [J]. Science Bulletin，2016，61（8）：601-611.

[48] YU B M，CHENG P.  A fractal permeability model for bidispersed porous media [J]. International Journal of Heat and Mass Transfer，2002，45（14）：2983-2993.

[49] XIAO B Q，TU X，REN W，et al. Modeling for hydraulic permeability and Kozeny-Carman constant of porous nanofibers using a fractal approach[J]. Fractals，2015，23（3）：1550029.

[50] CHAI J C，PRONGMANEE N. Barrier properties of a geosynthetic clay liner using polymerized sodium bentonite[J].Geotex? tiles and Geomembranes，2020，48（3）：392-399.

[51] JO H Y，KATSUMI T，BENSON C H，et al. Hydraulic conductivity and swelling of nonprehydrated GCLs permeated with single- species salt solutions[J]. Journal of Geotechnical and Geoenvi? ronmental Engineering，2001，127（7）：557-567.

[52] JO H Y，BENSON C H，SHACKELFORD C D，et al. Long-termhydraulic conductivity of a geosynthetic clay liner permeated with inorganic salt solutions[J]. Journal of Geotechnical and Geoenvironmental Engineering，2005，131（4）：405-417.

[53] PETROV R J，ROWE R K. Geosynthetic clay liner （ GCL）-chemical compatibility by hydraulic conductivity testing and fac ? tors impacting its performance[J]. Canadian Geotechnical Jour? nal，1997，34（6）：863-885.

[54] CHEN J，FANG Y G，GU R G，et al. Study on pore size effect oflow permeability clay seepage [J]. Arabian Journal of Geosci? ences，2019，12（7）：1-10.

[55] SETZ M C，TIAN K，BENSON C H，et al. Effect of ammonium onthe hydraulic conductivity of geosynthetic clay liners[J].Geotex? tiles and Geomembranes，2017，45（6）：665-673.

[56] BRADSHAW S L，BENSON C H. Effect of municipal solid wasteleachate on hydraulic conductivity and exchange complex of geo ? synthetic clay liners[J]. Journal of Geotechnical and Geoenviron ? mental Engineering，2014，140（4）：04013038.

[57] SEIPHOORI A，FERRARI A，LALOUI L. Water retention behav?iour and microstructural evolution of MX-80 bentonite during wet? ting and drying cycles[J].Géotechnique，2014，64（9）：721-734.

[58] LI Q，CHEN J N，BENSON C H，et al. Hydraulic conductivity ofbentonite-polymer composite geosynthetic clay liners permeated with bauxite liquor[J]. Geotextiles and Geomembranes，2021，49（2）：420-429.

[59] BENSON C H，CHEN J N，EDIL T B，et al. Hydraulic conductiv?ity of compacted soil liners permeated with coal combustion prod ? uct leachates[J]. Journal of Geotechnical and Geoenvironmental Engineering，2018，144（4）：04018011.

[60] CHEN J N，SALIHOGLU H，BENSON C H，et al. Hydraulic con?ductivity of bentonite – polymer composite geosynthetic clay liners permeated with coal combustion product leachates[J]. Journal of Geotechnical and Geoenvironmental Engineering，2019，145（9）：04019038.

[61] SHACKELFORD C D，SEVICK G W，EYKHOLT G R. Hydraulicconductivity of geosynthetic clay liners to tailings impoundment so ? lutions [J].  Geotextiles and Geomembranes，2010，28（2）：149-162.