随着人类活动不断向地下深入，基岩裂隙水受到日益严重的污染，引起人们的广泛关注，基岩裂隙水内溶质运移规律成为当前地下水科学的研究热点^[1-4]。天然岩体裂隙往往由未固结的松散物质充填，使得充填裂隙水流及溶质运移特征既不同于非充填裂隙，也与多孔介质有差异^[5-9]。

以往研究者较多地关注于非充填裂隙水流及溶质运移特征^[10-17]，近年来，充填裂隙水流及溶质运移也引起了一些学者的兴趣。陶同康^[5]、赵恺^[18]、Liu^[19]、郭保华^[20]、王鹏飞^[21]等研究了充填裂隙水流特性，发现充填裂隙的渗透性仍可用立方定律予以量化描述，不仅受制于裂隙开度和充填物的孔隙性，也与充填介质的物质构成密切相关。Qian等^[22]研究发现，相较于ADE和ADE–R模型，采用MIM模型对充填裂隙溶质运移穿透曲线的拟合效果更好。Huang等^[23]研究发现充填单管裂隙模型的浓度峰值到达时间快于双管裂隙模型，峰值浓度也高于双管裂隙模型；双管裂隙模型存在两个浓度峰值，且分支裂隙越长，第一浓度峰值到达时间越快。Dou等^[24]通过模拟非均匀填充裂缝中的溶质运移，结果表明填充介质的非均质性不仅导致了非均质速度场的产生，还导致了BTCs的早到和拖尾现象。Zhou等^[25]针对两种不同流动条件下的充填裂隙−基质系统中的溶质运移开展研究，发现在径向分散性和注入速率保持不变时，该系统对分散性最敏感，对流动/固定比率不敏感，对径向流动模型中的一级传质速率最不敏感。

上述对充填裂隙溶质运移试验研究多采用取样测试溶质浓度、绘制穿透曲线等手段，由于取样点往往只有一处，不能精细刻画充填裂隙内部溶质运移过程。近年来，由于自然电位（SP）方法在研究地下水流、传热传质等方面的明显优势，该方法成功应用于地下含水介质渗透系数估算^[26-28]、热量运移^[29-30]及多孔介质溶质运移参数估算^[31-33]。鉴于此，本文通过设计加工充填裂隙溶质运移试验模型，开展不同流速条件下溶质运移试验，结合不同测点SP实时动态监测数据，以期详细刻画充填裂隙内部溶质运移过程。

自然条件下，即使不向地下进行人工供电，利用仪器在地面两测点之间一般也能观测到电位差，即地球表层存在天然电场，一般称自然电场。根据引起自然电位方式的不同，可分为过滤电位、氧化还原电位和扩散电位^[34-35]。

过滤电位，也叫渗流电位，因固相对水溶液中正、负离子有选择的吸附作用，使正、负离子分布不均衡而产生，并最终趋于稳定。氧化还原电位由处于地下矿体或微生物的氧化还原作用而产生；同时，由于氧化作用使围岩溶液的pH值随之不同，进而引起浓度电位差，从而产生氧化还原电位^[36]。

当具有浓度差的两种溶液接触时，浓度梯度使得溶质由浓度高的地方向浓度低的地方运移，便产生扩散电位^[35,37]，其电流密度（ ${J_{\rm{s}}}$ ）可表示为^[38]：

式中： ${k_{\rm{b}}}$ 为玻尔兹曼常数； $T$ 为绝对温度； $t_i^ \pm $ 为水中第 $i$ 中溶质的Hittorf数； ${\sigma _{\rm{w}}}$ 为溶液电导率，S/m； ${q_i}$ 为第 $i$ 中离子的电荷数； ${C_i}$ 为离子浓度，g/L。

如果地下水中注入NaCl示踪剂，由溶质扩散形成的任意一点自然电位 ${\varphi} $ （单位：V）与该点溶质浓度关系可表示为^[38-39]：

式中： $n$ 为介质孔隙度； $C$ 为溶液浓度； ${C_0} $ 为初始浓度，g/L； ${\alpha ^*}$ 为耦合系数，对于NaCl溶液， ${\alpha ^*}$ =5.15×10^–3 V。

试验模型（图1）主要由4部分组成：进水系统、裂隙主体、出水系统、监测系统。进水系统包含恒压进水槽、进水缓冲箱；出水系统包含恒压出水槽、出水缓冲箱。进水槽和出水槽内设有溢流挡板，以保证进水端和出水端水压恒定；进水缓冲箱和出水缓冲箱内设有多孔分流板，为裂隙模型中的水流进出起到缓冲作用。裂隙主体尺寸为长1 200 mm×宽200 mm×隙宽5 mm。裂隙主体进、出水两端50 mm处分别设有溶质注入口及取样口。

图1(Fig. 1)

1.恒压进水槽；2.进水缓冲箱；3.溶质注入口；4.电极，5.裂隙模型；6.取样口；7.出水缓冲箱；8.恒压出水槽。 图1 充填裂隙溶质运移模型示意图 Fig. 1 Schematic diagram of fissure solute transport physical model

裂隙主体内充填河砂，其颗粒级配曲线如图2所示。采用饱水法测得充填砂的孔隙度为0.362，利用达西渗透仪测得渗透系数为1.79 mm/s。

图2(Fig. 2)

图2 充填砂级配曲线 Fig. 2 Gradation curve of filled sand

监测系统由监测电极（直径2 mm的铜棒）、DT85G数据采集仪、电脑构成。参比电极设置在进水缓冲箱内，裂隙主体上布设等距的6个测量电极（记1^#～6^#），1^#测量电极设置在距溶质注入口100 mm处，每两个电极间隔150 mm。通过信号线将电极与数据采集仪连接，数据采集仪获取数据后传输并储存到电脑上。

试验选用化学性质稳定的NaCl作为溶质示踪剂。试验前配置已知浓度NaCl溶液，测试其电导率，标定NaCl溶液浓度（c）与电导率（E_c）关系，如图3所示。

图3(Fig. 3)

图3 NaCl溶液浓度与电导率关系曲线 Fig. 3 Relationship between NaCl concentration and electrical conductivity

1）水流试验

利用备好的砂样填充裂隙模型，调节上下游水箱，通过改变上下游水位差调整水流速度。试验时，在出水口处采用体积法测得流量，计算出水流平均流速，开展了平均流速为0.202～2.652 mm/s范围内的充填裂隙水流试验。

2）溶质运移试验

选取水流平均流速为0.247、0.431 和0.661 mm/s，开展了充填裂隙溶质运移试验。试验过程中，溶质注入后，立即在取样口处取样，采用电导率仪（DDS–307）测试其电导率，计算NaCl溶液浓度，绘制穿透曲线（BTC）。

溶质运移试验开始前，要求模型系统内水流稳定，首先进行示踪剂运移预试验，以确定溶质初始和最后到达时间，作为取样与监测时间间隔依据。连续两次溶质运移试验之间，充填裂隙试验用水冲洗足够长时间，直至NaCl浓度低于检出限。

3）自然电位监测

利用监测系统动态采集溶质运移试验过程中自然电位响应数据，采集时间间隔为5 s，试验用水为实验室内自来水，试验过程中室温基本稳定在25 ℃。

试验过程中获得的渗流速度（v）与水力梯度（J）数值见表1，绘制v−J关系，如图4所示。由图4可见，在试验范围内，渗流速度与水力梯度之间满足线性关系，说明水流符合达西流。同时，计算的充填裂隙渗透系数1.85 mm/s与多孔充填介质渗透系数1.79 mm/s也较为接近，表明河砂充填裂隙渗透性主要取决于充填介质的渗透性。于龙^[6]和速宝玉^[40]等采用两平行玻璃板模拟裂隙，结果也表明充填裂隙渗流特性不仅与裂隙的宽度有关，还与充填材料的颗粒组成、孔隙率、颗粒直径等有关。

表1(Tab. 1)

表1 水力梯度与渗流速度值 Tab. 1 Correspondence between seepage velocity and hydraulic gradient

表1 水力梯度与渗流速度值 Tab. 1 Correspondence between seepage velocity and hydraulic gradient J v/ (mm·s–1) J v/ (mm·s–1) 0.0417 0.0733 0.2540 0.4800 0.0833 0.1500 0.3170 0.5840 0.1210 0.2030 0.4040 0.7330 0.1580 0.2910 0.5130 0.9630 0.2080 0.3820

图4(Fig. 4)

图4 渗流速度与水力梯度关系 Fig. 4 Relationship between seepage velocity and hydraulic gradient

图5为不同平均流速（u）条件下充填单裂隙的溶质运移穿透曲线。据图5可知，不同流速下，穿透曲线大致呈对称状，流速较小时表现出一定的拖尾现象。随流速升高，对应的溶质峰值浓度逐渐升高，穿透曲线开口越小，穿透曲线的波峰右侧陡峭程度逐渐增加。

图5(Fig. 5)

图5 填充裂隙溶质运移穿透曲线 Fig. 5 Fracture filled solute transport breakthrough curves

考虑试验检测的精确性，当检测浓度超出浓度峰值的1%，即认为溶质已到达取样口。从上述穿透曲线上提取溶质初始到达和峰值到达时间，并绘制与水流平均流速关系，如图6所示。随着平均流速的增大，溶质初始到达、峰值到达取样口所需时间显著缩短。研究表明，在地下水中溶质运移峰值浓度区域是以地下水宏观平均流速往前推进^[41]。根据溶质通过取样口峰值浓度对应时间，计算平均流速依次为0.201、0.327、0.539 mm/s，小于根据流量数据计算的平均流速。

图6(Fig. 6)

图6 不同平均流速对溶质初始、峰值到达取样口耗时的影响 Fig. 6 Effect of different average flow velocity on the time of initial and peak solute reaching the sampling port

图7为不同平均流速条件下，充填裂隙溶质运移过程中的自然电位响应曲线。由图7所知，相同平均流速下，沿水流方向从1^#～6^#测量电极处峰值SP绝对值逐渐减小，且依次出峰时陡峭程度逐渐减小，反映了在充填裂隙内溶质运移过程中浓度峰值到达情况。这是由于溶质在充填裂隙内沿着水流方向运移过程中，存在水动力弥散作用，随着溶质运移时间增长，水动力弥散作用越显著，导致溶质浓度逐渐降低。相应地，沿水流方向上，1^#～6^#电极处峰值SP绝对值也逐渐减小，峰也越来越平缓。

图7(Fig. 7)

图7 试验过程中自然电位响应曲线 Fig. 7 Natural potential response curves during the test

不同水流平均流速下，2^#和3^#电极自然电位响应曲线如图8所示。由图8可见，随着水流平均流速增加，溶质流过相同测量电极所需时间减少，峰值SP自然电位绝对值不断增大，这与上述溶质运移穿透曲线规律相一致。在试验过程中，水流作为示踪剂的运输者，流速增大，示踪剂运移速度也不断增大，示踪剂到达监测点的所需时间减少，示踪剂峰值到达时间也不断缩短。

图8(Fig. 8)

图8 不同平均流速下2#和3#电极自然电位响应曲线 Fig. 8 Natural potential response curves of 2# and 3# electrodes under different average flow velocity

为判断溶质运移过程中的扩散作用，引入无量纲数 Peclet（ $\rm Pe$ ），其值计算公式如下^[41]：

式中： $u$ 为试验过程中水流宏观平均流速，mm/s； ${d_{50}}$ 为充填介质平均颗粒粒径，试验中 ${d_{50}}$ =0.428 mm； ${D_{\rm{d}}}$ 为溶液分子扩散系数，根据文献[42]，结合试验条件取 ${D_{\rm{d}}}$ =1.50×10^–7 m²/s。

根据式（3），在试验流速范围内计算得 $\rm P{\rm{e}}$ 为0.705～1.890，可知在该试验条件下，机械弥散与分子扩散均不能忽略，分子扩散所起作用仍大于机械弥散，但随着流速增加，机械弥散作用逐渐增强。溶质运移过程中对流与水动力弥散同时存在，当平均流速增大时，对流作用增强，示踪剂到达监测点的时间明显缩短，机械弥散与分子扩散的作用时间大大缩短，溶液净弥散作用减弱，示踪剂浓度降低的程度小，浓度峰值反而高，对应的峰值自然电位绝对值也增加。

从图7自然电位响应曲线上提取初始到达和峰值到达时间，分别绘制其对应电极距离关系，如图9所示。由图9（a）可知，相同平均流速条件下，对溶质初始到达时间与运移距离数据进行线性拟合，效果良好，拟合相关系数R²≥0.993，说明溶质弥散前锋也呈匀速状态随着水流运移，溶质弥散前锋平均速度(u_初)分别为0.214、0.396、0.685 mm/s。由图9（b）可知，相同平均流速条件下，对溶质浓度峰值到达时间与运移距离数据进行线性拟合，效果良好，拟合相关系数均达R²=0.999，溶质弥散峰值平均速度(u_峰)分别为0.189、0.347、0.571 mm/s。由图9（a）、（b）可知，溶质弥散前锋与溶质弥散峰值浓度随水流动均呈匀速运移状态，但u_峰<u_初，由于速度差值恒定，导致溶质弥散范围随着时间的增加而逐渐扩大。

图9(Fig. 9)

图9 不同平均流速对溶质运移至各电极的影响 Fig. 9 Effect of different average flow velocity on the migration of solute to each electrode

图10(Fig. 10)

图10 不同 u 下峰值SP−距离的关系 Fig. 10 Effect of different u on the migration of solute to each electrode

表2为试验过程中分别采用溶质峰值浓度、峰值SP计算水流平均流速与流量数据计算平均流速对比结果。由表2可知，利用溶质浓度和峰值SP两种方法计算流速值结果基本一致，而采用流量数据计算流速值大于上述两种方法计算所得流速值，平均误差分别为20.4%和18.9%。此外，u_峰比u_c更接近流量计算平均流速，且随着流速的增大，误差减小。

表2(Tab. 2)

表2 试验过程中流速实测值与计算值对比 Tab. 2 Comparison of measured values and calculated values of flow velocity during the test

表2 试验过程中流速实测值与计算值对比 Tab. 2 Comparison of measured values and calculated values of flow velocity during the test 流量计算平均流速u/(mm·s–1) 溶质峰值浓度计算流速 峰值SP计算流速 初始SP计算流速 uc/(mm·s–1) 误差/% u峰/(mm·s–1) 误差/% u初/(mm·s–1) 误差/% 0.247 0.201 18.6 0.189 23.5 0.214 13.4 0.431 0.327 24.1 0.347 19.5 0.396 8.1 0.661 0.539 18.5 0.571 13.6 0.685 –3.6 平均相对误差/% — 20.4 — 18.9 — 6.0

表2中同时列出了由初始SP数据计算流速（u_初）的数据，发现u_初更加接近由流量数据计算的流速数据。这与以往认为溶质运移峰值浓度区域是以地下水宏观平均流速往前推进的认识^[41]存在一定差异。

从自然电位响应曲线上提取各电极处自然电位峰值，绘制不同平均流速下，自然电位峰值与溶质运移距离的关系，对峰值SP与运移距离进行对数函数拟合（R²≥0.977），如图10所示。根据对数方程特性可知，斜率k>0且不断减小，因此，峰值SP变化率在不断减小，呈现峰值SP绝对值前期降低快而后期降低缓慢的变化趋势。这也符合水动力弥散特点。水动力弥散分为机械弥散和分子扩散，机械弥散主要是由于水流微观运动速度在方向和大小上有别于水流宏观运动速度，因而产生的扩散作用；分子扩散主要是由于示踪剂浓度分布不均匀而引起的弥散，体系内某一组分的浓度梯度驱使示踪剂具有向四周辐射的特点，这与体系的微观运动速度的运动方向不一致，造成组分偏离体系的扩散现象。当平均流速相同时，机械弥散相同，水动力弥散变化主要受分子扩散作用影响，进样口注入示踪剂，浓度梯度达到最大，分子扩散作用最明显，此时自然电位峰值绝对值最大；随着示踪剂沿水流方向运移，分子扩散作用导致示踪剂浓度显著降低，浓度梯度迅速减小，从而导致扩散作用逐渐下降，监测自然电位变化较为明显；此后分子扩散速度逐渐降低，沿着水流方向示踪剂浓度越来越趋于稳定，当运移距离足够远时，示踪剂浓度无限趋近于0。此阶段自然电位的变化逐渐降低，且数值逐渐趋近于试验用水的自然电位值。

图11(Fig. 11)

图11 不同u下SP峰值-距离的关系图 Fig. 11 Peak-distance relationship diagram under different u

对比上述两种溶质运移监测方法，传统取样监测时需取出一定量的溶液，这将对充填裂隙内部溶质运移产生一定扰动，且取样点少，取样监测值不能反映瞬时溶液浓度。采用SP方法进行溶质运移研究时，由于实时监测无需取样，弥补了传统取样监测的不足，同时能够实时地检测到溶质运移过程中自然电位响应数值，将溶质运移实时状态可视化，能更便捷地掌握充填裂隙内溶质运移的情况，提高试验过程中监测的准确性。

采用自然电位法，结合传统取样测试，通过改变充填裂隙上下游水头差，调整水流速度，开展了充填裂隙溶质运移试验，根据自然电位响应特征探讨了充填裂隙内部溶质运移过程。

1）在试验流速0.0733～0.963 0 mm/s范围内，充填裂隙渗透系数为1.85 mm/s，与充填多孔介质渗透系数1.79 mm/s较为接近，水流符合达西定律。

2）试验过程中水流实际平均流速为0.247、0.431、0.661 mm/s时，取样口处溶质峰值浓度分别为5.14、5.50、5.78 g/L，溶质初始到达和峰值到达时间均减少，其中初始到达时间分别为4 640、2 640、1 560 s，而峰值到达时间分别为5 480、3 360、2 040 s。

3）3种水流实际流速条件下，溶质运移过程中自然电位响应呈现出相似特征，沿水流方向从1^#到6^#测量电极处SP峰值绝对值逐渐减小，且依次出峰时陡峭程度逐渐减小；此外，均反映出溶质弥散前锋与溶质弥散峰值浓度点随水流动均呈匀速运移状态，其中溶质弥散前锋相应的平均流速分别为0.214、0.396、0.685 mm/s，溶质弥散峰值相应的平均流速为0.189、0.347、0.571 mm/s。

4）基于峰值SP计算的水流实际平均流速与采用溶质峰值浓度计算结果基本一致，但基于峰值SP计算结果更接近实际平均流速（流量计算结果），准确度提高1.5%；随着流速不断增加，采用峰值SP计算的宏观平均流速误差降低9.9%。

5）试验发现，与峰值SP值相比，采用初始SP值计算的充填裂隙水流宏观平均流速更接近由流量数据计算的水流宏观平均流速。