高粘性岩土化学渗透效应数值模型试验分析

刘志猛，胡晓维

(中化学土木工程有限公司,江苏 南京 210031)

粘性岩土中富含大量黏土矿物，渗透系数较低，以往研究多认为其不透水或弱透水，因此在构建地下水流、溶质运移模型中往往会忽视流向，溶质运移期间弥散作用发挥主导作用。但黏土并非完全隔水，其在水头梯度下会呈现出越流现象[1]。部分学者经过试验获得了粘性岩土渗透系数、弥散度等相关数据，并在对流弥散模型中应用。近年研究表明，粘性黏土作为低渗透黏土介质，在性质上与生物半透膜具有一定的相似性，基于不同梯度产生的化学渗透、超滤现象即“耦合流”现象，其会对地下水流运动产生不同程度的影响，传统Darcy定律对地下水流动的刻画已经不能够满足研究需求[2]。目前，在土木工程、沉积盆地、环境科学等领域，粘性岩土已经被充分考虑，对粘性岩土化学渗透效应的分析有着重要的社会意义与应用价值[3-5]。目前我国关于粘性岩土化学渗透性质的研究相对较少，数值模型能够反映出粘性岩土剪应变、体应变及耦合效应，帮助判断不同应力路径下黏土弹塑性本构关系，掌握粘性岩土化学渗透规律。

1 化学渗透、渗透及耦合流原理

1.1 化学渗透与反渗透

浓度不同的2种溶液隔半透膜，浓度较低溶液的水分子及溶剂会经过半透膜流向浓度高的溶液，这一过程被称为“渗透”。通常渗透现象主要指的是非水头对水流的驱动现象，化学渗透及反渗透、超滤则是在化学梯度驱动下产生的渗透现象[6]。半透膜2端化学势能是导致化学渗透压的主要因素，假设渗透压用π(kPa)表示，其计算公式：

Δπ=vRTΔC，其中，解离常数用v表示；半透膜两端浓度差采用△c表示，假设NaCl浓度为0.1 mol/L，温度为25 ℃ ，经过计算可以发现半透膜2端渗透压力之差为495.7 kPa，当然该值为其理论上的差值，接近于50 m水头压力(超过大气压或任意标准压力的气体压力)。

反渗透即超滤，主要指的是高浓度一端的溶液经过较高渗透压力水头压力的作用，向逆行自然渗透方向进行反作用的过程。基于半透膜选择性渗透的特点，高浓度端会出现溶剂分子选择性聚集，此时半透膜的低压端获得更多的溶剂，被称为渗透液；浓缩的溶液则集中在高压端，被称为浓缩液[7-8]。

1.2 黏土膜性能

粘性岩土或粘质页岩半透膜性质与传统生物半透膜呈现出明显的差异性，其性质主要取决于黏土颗粒表面的双电层。通常，含水层透水性越好，孔隙水溶液离子受到双电层影响越小，其能够自由通过。若沉积物富含黏土物质，受到各方压力的影响，会缩小黏土颗粒间距，此时双电层会呈现出重叠现象，带电离子在溶液中会受到电荷排斥作用，无法顺利通过黏土介质[9]。不带电分子——水则不会受到排斥，会顺利通过黏土介质。作为天然半透膜，黏土会对溶质迁移进行一定的约束，但容许其具有迁移能力，一般可采用w表示，即化学渗透率。有学者在研究中用σ表示化学渗透率系数，认为当半透膜2端处于平衡状态，化学渗透所致的水头压力/理想中渗透压可以表示：

式中：J表示半透膜渗流速度；△P表示水头压力差变化值；△π表示理想渗透压变化值。

一般情况下，σ取值范围为0～1。当无选择通过介质，即σ=0时提示不存在化学渗透作用，介质中溶液均有对流-弥散；假设σ=1，提示其对所有溶质完全限制通过。通常，自然状态下粘土半透膜均属于非理想半透膜，σ介于二者之间。理想半透膜，则是人工制备下的高分子产物。

2 粘性岩土化学渗透数值模型构建

目前关于粘性岩土化学渗透效应的研究主要存在2个方面问题：一方面为粘性岩土介质自身的复杂性；另一方面为工程应用强调模型尽可能简单。因此细观角度对渗透效应的研究尚存在不足。随着宏观土力学的发展，室内试验与数值模型相结合成为研究粘性岩土化学渗透效应的可靠路径[10-11]。通过室内试验能够帮助确定土样特性参数，提升数值模型的准确性与可靠性，二者相互验证能够对粘性岩土化学渗透问题本质予以揭示。

2.1 数值模型试验准备

2.1.1基本物理指标测定

研究采用体育馆深基坑土样，已经有扰动。由于本试验主要研究粘性岩土渗透特性，这需要对土样进行重塑，使其成为试验所需土样；按照已有研究[12-13]的试验规程结果，本试验土样的含水量、风干含水量及密度等设计如表1所示。

表1 土样含水量试验设计结果Tab.1 Water content test

标准试验直径用D表示，其参数为39.2 mm，取2 kg代表性土样，经过风干及碾碎处理，实施过筛，土粒粒径控制在0.1D以下；风干含水率如表2所示；密度试验的数据如表3所示。

表2 风干含水量试验Tab.2 Air drying water content test

表3 密度试验记录Tab.3 Density test records

在横坐标、纵坐标上分别设置含水率与圆锥下沉深度，然后绘制关系曲线。下沉深度17 mm含水率为标记，作为液限；2 mm位置为塑限。塑性指数计算公式：

IP=(wL-wP)×100=19.65>17

式中：IP表示塑性指数；WL表示液限；WP表示塑限。

式中：ρW表示水密度;ρd表示土干密度;Gs表示土体密度。

根据含水量对饱和度予以计算，公式

式中：Sr表示饱和度；Gs表示土体比重；e表示材料孔隙比；w表示含水量。

2.1.2试样制备

选择2 kg代表性土样，根据上述测得的含水率，计算需要加的数量，饱和含水率计算方法：

所需加水量：

在土料上喷洒需要加的数量并使其均匀分布，静置后放入塑料袋，保存于密闭容器，以20 h为宜，确保含水率均匀。将土料取出后对含水率进行复测，确保差值小于1%。土质量应根据干密度进行确定，将试样密度差值控制在0.02 g/cm3，经过抽气饱、反压饱制备黏土。

2.2 数值模型试验过程

2.2.1接触连接模型

接触连接模型假想颗粒接触间存在弹簧，其切向与法向刚度往往保持不变，弹簧的抗拉强度、抗剪强度强，颗粒的接触连接会对滑动现象产生影响，若存在接触连接，颗粒重叠量为0，提示受到了张力的影响[5]。接触连接强度能够为颗粒流程序赋予接触连接强度。当对颗粒产生作用的法向应力大于等于接触连接强度，会破坏接触连接，导致连接强度的消失。若切应力大于等于切向接触连接强度，那么也会破坏接触连接，颗粒间接触力则不会出现鲜明的变化，其关系如图1所示。研究经过调试确定法向与切向粘接强度均为30 kN。

图1 接触粘接模型Fig.1 Contact bonding model

2.2.2粒子的生成

可采用命令在颗粒流程序中生成颗粒密度、半径及数量等相关数据，构建线性或弧形墙，粒子在no-shadow侧产生，经no-shadow侧穿越墙体，可达到shadow，可利用墙对颗粒生成范围进行限定。采用半径扩大法获得试样，其能够将颗粒集合体置于密室状态。粒子与墙体之间无粘接，因此在渗流模拟时容易出现粒子相对滑动现象[15-16]。针对这一问题，研究引入不排水边界这一概念，对2侧、底部墙体粒子予以固定处理，将其作为模拟的试样边界[17]。试验时将土样上部作为一个渗流出口，保证数值模型与土工试验一致。研究根据土颗粒半径范围，对所要达到的孔隙率进行预测，对土颗粒总数进行计算，获得试样模型，然后对孔隙率进行计算；其相关参数如表4所示，其中3e4表示3倍的e的4次方，为切向黏结强度；1e7表示1倍的e的7次方，表示法向黏结强度。

表4 模型相关参数表Tab.4 Parameters related to the model

2.2.3模型微观参数选取

粘性岩土渗透效应建模涉及到多项微观参数，确保微观与宏观参数的对应存在较高的工作量，工作难度大。因此，首先应计算流固耦合流体参数对结果的敏感性，研究参数主要包括残余孔径、渗透系数、流体体积模数与距离缩放因子。土样是否会产生裂纹很大程度由土样内部应力状态决定，若内部应力对微观参数敏感性高，那么可以将其确定为影响结果的参数[18]。研究在试样中心确定一位置作为测量圆，半径为15 mm，对土颗粒平均应力进行记录，并绘制变化曲线；考虑到粒子集合体生成需要一定的运算时步进行平衡，研究对模型底部水头压力施压后的应力变化予以分析，对单一参数进行控制，保持其他参数值不变；对另一微观参数进行变化，对土样应力变化予以分析。从研究结果可以发现，流体体积模数、颗粒表面间法向距离及法向应力的取值处于合理范围内，能够将对计算结果的影响控制到最小。根据上述分析，流体计算参数取值情况如表5所示。

表5 流体计算参数取值情况Tab.5 Fluid calculation parameter values

2.2.4土样渗透判断标准

(1)土样未发生渗透：未见土样内部裂纹，或裂纹减少，但未产生渗透通道。土样孔隙率在监测过程中先表现为上升，然后逐渐趋于稳定。土样位移变量经过一定时间的波动，会逐渐稳定；系统不平衡力监测显示不运行后出现波动，但无突然增大。

(2)土样渗透标准：可见土样内部有明显的裂纹，有完整的渗透通道；可见土样孔隙率先出现波动，然后稳定；系统平均不平衡力监测变量随着运算时步的运行而不停的波动，但未出现突然增大的情况。

2.2.5数值模型的运行

载入Fish语言编程文件，通过对试样底部水头压力的施加，使得试样顶端与底部产生水头差。在模拟过程中，压力传递在相邻颗粒间的接触与颗粒“管道”相近，其与颗粒间接触法向位移呈现出正相关，在建立起数值模型时，通道间隙不容易发生变化，往往是粘接破坏后或颗粒发生位移后出现变化[19-20]。采用逐级加压法施加水头压力，在施加压力时遵循由小到大的原则，运行程序后明确有无渗透，若出现渗透破坏应逐级递减，根据差值法施加压力，明确精准的渗透压。

3 数值模型结果及分析

土样渗透破坏示意图，如图2所示。

图2 土样渗透破坏示意图Fig.2 Schematic diagram of seepage damage of soil sample

从图2可以看出，当出现渗透破坏时，在试样底部最先出现裂纹；然后逐渐向边壁延伸，发展为贯通渗透通道，且裂缝数量与时步无线性关系。这是因为裂缝尖段只有积聚一定的能量才能够扩展裂缝，形成新的裂缝。当能量变大，将会出现大量新裂缝。

对经过实验后出现渗透破坏的粘性岩土土样进行分析，对比了其渗透破坏位置与数值模型的差异，发现二者渗透通道均为底部裂纹连接的周边，并有明显裂纹，由此可以判断室内试验与数值模型具有较高的一致性。经过数值模型构建，其获得的抗渗允许坡降值与试验对比结果如表6所示。

表6 室内试验与数值模型结果比较Tab.6 Comparison of results between laboratory test and numerical model

由表6可知，数值模型与试验在抗渗允许坡降方面存在一定的差异性，但误差在可接受范围内。数值模型对粘性岩土渗透效应的规律反映与试验结果具有较高一致性，说明应用该方法是切实可行的。

4 结语

综上所述，数值模型下分析粘性岩土化学渗透效应，与室内试验具有较好的一致性，二者渗透通道位置基本一致，能够反映出粘性岩土化学渗透规律，具有较高的可行性与实用价值。但关于渗透性能与颗粒级配间的规律联系尚存在一定的不足，需要充分考虑地下水流与溶质运移的影响，这些都是今后粘性岩土渗透性能研究的方向。