本研究由 Pan Yuhua（第一作者，重庆大学）、Liu Xianshan（通讯作者，重庆大学）、Song Yulin（重庆大学）、Cao Yiting（重庆大学）和 Wei Ning（中国科学院武汉岩土力学研究所）合作完成。研究论文题为“Impact of shear displacement on non-fickian solute transport within diverse 3D interconnected cross-fracture”，发表于 International Communications in Heat and Mass Transfer 期刊，发表时间为2025年（卷169，文章号109489），已于2025年8月22日在线发布。
本研究隶属于地质工程与岩石力学领域，具体聚焦于裂隙岩体中的流体流动与溶质运移这一交叉科学问题。该问题对于放射性废物处置、非常规油气开采、二氧化碳地质封存以及地热储层优化等关键地下工程至关重要。在这些工程活动中，外部力学的扰动会改变岩体中断裂网络的结构完整性，可能使原本独立的裂隙网络相互连通。
天然岩石裂隙网络由大量单个裂隙及其交叉点构成。裂隙交叉点（cross-fractures）是实现不同流体和溶质在相连裂隙间混合与再分配的关键位置。传统研究常将单一裂隙简化为光滑平行板模型，忽略了天然裂隙表面的粗糙度与网络的复杂性。实际中，粗糙的裂隙表面会导致曲折的流动路径和溶质运移异常（即非费克运移，Non-Fickian transport）。同时，增强型地热系统等工程通过水压剪切（hydroshearing）来刺激裂隙网络，此过程伴随着裂隙的剪切滑移（shear-induced slip），会显著改变裂隙开度分布，增强流动通道化，进而深刻影响流场和溶质运移行为。
尽管已有研究探讨了剪切对单一裂隙或交叉裂隙的影响，但对于在外部扰动下，经历剪切滑移的三维（3D）互联交叉裂隙网络中非费克溶质运移的机制仍缺乏深入探索。特别是，剪切如何通过改变空隙结构和交叉点几何形态来控制异常运移，尚不明确。因此，本研究的核心目标是：系统研究在外部扰动引起的剪切滑移过程中，三维互联交叉裂隙内的非费克溶质运移行为，揭示剪切位移、交叉几何结构对流动涡旋、溶质侵入、混合与分布的影响机制。
本研究流程严谨，包含几何模型构建、数值模拟和结果分析三大主要部分。
第一，构建三维粗糙交叉裂隙几何模型。 1. 生成粗糙裂隙表面：研究团队采用湖北芝麻白花岗岩，制备了边长为50毫米的标准立方体试样。使用配备定制夹具的WADJ-600伺服岩石剪切流变试验机，通过巴西劈裂法产生天然粗糙裂隙面。随后，利用Wiiboox 3D Reeyee物体扫描仪（单次测量精度≤0.1毫米）对裂隙面形态进行数字化扫描，获取高精度点云数据（图1，图2）。 2. 量化裂隙表面粗糙度：为客观选择代表性裂隙面，研究采用了多种统计参数进行量化。除了常用的节理粗糙度系数（JRC）的统计关联参数 Z2 （剖面一阶偏差的均方根）外，重点采用了基于分形理论的参数组合。通过计算裂隙剖面分形维数（D） 和分形参数（Kv），并以其乘积 D × Kv 作为综合粗糙度评价指标。该指标能同时反映空间变异性和幅度，比单一统计参数更优。研究沿裂隙面的x和y两个正交方向提取多个剖面进行计算，最终从所有候选表面中筛选出三个具有不同D × Kv值的代表性粗糙表面（编号：No.2-lower, No.4-lower, No.6-lower）用于构建模型（表1）。 3. 构建互联交叉裂隙（多交叉裂隙）模型：本研究创新性地构建了“多交叉裂隙”模型，以模拟两个交叉裂隙在外部力作用下相互连接形成的更复杂网络。具体方法是：一个水平裂隙分支（记为F1，来源于No.2-lower表面）同时作为两个交叉裂隙（Cross-fracture No.1 和 No.2）的共享水平分支。两个垂直裂隙分支分别来自No.4-lower（记为F2，属于交叉裂隙No.1）和No.6-lower表面（记为F3，属于交叉裂隙No.2）。通过改变第二个交叉裂隙（No.2）的二面角，构建了三种不同的几何组合：45° + 45°、45° + 90° 和 45° + 135° （图6）。 4. 施加剪切位移：为了模拟剪切滑移，研究对水平裂隙（F1）的下表面以及两个垂直裂隙分支施加沿x方向的水平剪切位移（δx）。设置了三个剪切水平：0毫米（未剪切）、1毫米和1.7毫米。将此剪切方法应用于上述三种几何组合，最终得到了 9个 不同的三维几何模型用于后续分析（图8）。 5. 分析模型几何特征：研究详细统计了各裂隙的平均开度、开度标准差、相对粗糙度和接触比（表2）。发现随着剪切位移增大，水平裂隙F1的平均开度减小，标准差增大，相对粗糙度和接触比显著增加，开度分布呈现双峰特征，表明空隙结构复杂性和各向异性增强（图9）。此外，重点量化了两个交叉点（Intersection No.1 和 No.2）的横截面积（反映空隙体积，图10）以及交叉点与相邻分支裂隙之间的界面接触面积（图11，图12）。这些几何参数是后续分析流场和溶质分布的基础。
第二，计算流体动力学与溶质运移数值模拟。 1. 控制方程与求解器：流体流动遵循纳维-斯托克斯方程（Navier-Stokes equations） 和连续性方程，假设为等温、稳态、不可压缩层流。溶质运移由对流-扩散方程（advection-diffusion equation） 控制。所有控制方程均使用商业软件 COMSOL Multiphysics 6.1 进行直接数值求解。 2. 边界与初始条件：在入口（垂直裂隙F2和F3的入口）施加层流入口边界条件，采用与压力无关的质量流量，以确保不同入口面积下的流量一致。出口设置为恒定压力条件。所有裂隙壁面为无滑移、不透水壁。初始时刻，整个计算域溶质浓度为零。瞬态溶质运移模拟通过在F2和F3入口瞬间施加恒定浓度源（C0 = 1 mol/m³）来启动（图13）。 3. 网格划分与求解配置：计算域采用自由四面体单元进行离散，网格单元数量超过158万，保证了计算精度（图14）。采用多物理场耦合框架：首先进行稳态计算，求解流场；然后将稳定的速度场和压力场导入瞬态求解器，模拟溶质运移过程，实现流体动力学与质量输运的双向耦合。模拟中使用的流体物性、扩散系数、流速（对应Peclet数 Pe = 13,074， Reynolds数 Re = 10）等关键参数详见表3。
第三，数据处理与分析方法。 研究对数值模拟结果进行了多维度、定量化的分析： 1. 流场分析：通过截面速度云图、速度流线图、特定截面（y=15, 25, 35 mm）的速度统计分布，分析流速大小、通道化以及流体在交叉点的再分配。 2. 涡量场分析：绘制涡量等值线图，分析涡旋的强度、空间分布及其与几何结构的关系。同样在特定截面对涡量大小进行区间统计，量化不同剪切位移下高强度涡旋区域的比例变化（见附录B）。 3. 溶质浓度场分析：可视化不同时刻的溶质浓度空间分布，观察溶质的前缘推进、在不同分支中的分配。绘制各出口（F1出口、F2出口、F3出口）的溶质穿透曲线（BTCs），分析早期到达、拖尾效应等非费克运移特征，并计算稳态时的归一化浓度。 4. 混合比例分析：在交叉点上下游定义分析截面（图33），计算混合比（Mixing Ratio， MR），用于量化溶质在交叉点的混合程度。MR=0.5对应于完全混合模型，偏离0.5则趋向于流线分配模型。通过对比不同模型的MR值及其差值（|δMR|），评估剪切和几何组合对混合模式的影响。
1. 剪切位移加剧流场不均匀性与通道化。 结果显示，所有模型的流速场均表现出高度不均匀性，高速流主要集中于开度较大的优势通道中（图15）。随着水平裂隙剪切位移（δx）的增加，流场的不均匀性进一步加剧。统计表明，高流速区域（如u > 0.05 m/s）的面积比例显著增大（图18）。剪切导致水平裂隙开度分布双峰化，产生了更多接触区域和大开度区域共存的复杂结构，这是流场异质性增强的根本原因。速度流线图（图19-21）清晰显示，流体倾向于绕过曲折的接触区，沿优势通道流动。剪切位移增大会改变交叉点的几何结构，进而影响流体从一个交叉点到另一个交叉点的再分配量，最终导致各出口分支流量比例发生变化。
2. 涡旋场的形成与演化受空隙结构复杂性控制。 涡量场分析（图22）表明，涡旋强度随剪切位移增大而增强。在开度均匀的区域，涡旋核心离散但相互连通（图23①）。而在开度受限、接触曲折的区域，流动出现分隔，形成离散的涡旋域（图23②）。二维截面观察常显示为封闭流线的孤立涡旋，但三维分析揭示这些涡旋与主流动通道之间仍存在流体动力学连通（图23①a,b）。定量统计（图26及附录B）证实，高强度涡旋（如涡量 > 350 s⁻¹）的区域比例随剪切位移显著增加。这表明，剪切引起的空隙结构复杂化是促进流体旋转运动、形成更强更广泛涡旋的关键因素。
3. 溶质侵入呈现两阶段模式，剪切影响溶质分配。 溶质浓度场的时空演化（图28-30）揭示了在粗糙裂隙中溶质侵入的两个截然不同的阶段：（i）对流主导阶段：溶质首先侵入优势流动通道，实现快速突破。（ii）扩散控制阶段：当优势通道相对饱和后，由浓度梯度驱动，溶质逐渐侵入周围的小空隙结构。溶质的空间分布始终遵循流线轨迹。
剪切位移对溶质在出口分支间的分配有系统性影响。随着δx增大，F2出口的稳态溶质浓度单调增加，而F3出口的浓度则下降（图31）。这种不对称的再分配并非简单地将溶质从水平分支转移到垂直分支，因为水平分支（F1）出口的归一化浓度本身也因剪切而大幅下降（降幅达39.7%-44.3%）。这意味着更多的溶质被滞留在裂隙网络内部。
4. 交叉点几何界面控制溶质混合与分配。 溶质在交叉点的分配与其界面接触面积强相关（图12）。对于流入交叉点（Intersection No.2），剪切位移增大了水平流入方向（2-l）的界面面积，减小了垂直流入方向（2-r）的面积。这导致水平流入速度相对垂直流入速度增大。在对流主导机制下，更高流速的流入分支会将更多流体（及溶质）分配至其对向的流出分支。因此，水平流入速度优势促使更多溶质流向水平分支（F1），进而减少了流向F3出口的溶质。对于流出交叉点（Intersection No.1），剪切的影响相反，最终导致更多溶质流向F2出口。溶质混合比（MR）的计算结果（图34，表5）进一步证实，所有交叉点的混合模式均显著偏离完全混合模型（MR≠0.5），且剪切位移增大使Intersection No.1的混合模式更趋向于流线分配模型（|δMR_1|增大）。研究强调，准确预测交叉点混合模式需要综合分析出口方向的界面面积和入口方向的流速大小。
5. 涡旋区域实现动态捕集-释放平衡，是溶质滞留的关键机制。 溶质穿透曲线（BTCs）显示出明显的非费克特征：早期到达和拖尾现象（图31）。垂直裂隙（F2， F3）BTCs拖尾不明显，而水平裂隙（F1）的拖尾随剪切显著增强，且在δx=1.7 mm时出现了一个“伪稳定”阶段（初始稳定后出现浓度短暂爬升）。结合涡旋场和开度分布的分析，本研究提出了一个创新性假说：高强度涡旋核心是溶质滞留的主要场所。然而，这种滞留并非单向的永久捕获。在稳态流动条件下，每个涡旋能够建立一个溶质储存容量的动态平衡。溶质可以“艰难地”侵入涡旋核心并被旋转流体携带，但当涡旋内的溶质浓度超过其旋转输运能力时，溶质又会被排出。因此，涡旋区通过动态捕集-释放平衡（dynamic entrapment-release equilibration） 机制维持着溶质储存，导致出口溶质总量始终低于进口总量，且该差值随剪切位移增大（表4）。这从机理上解释了BTCs的拖尾现象和溶质滞留。
本研究通过构建真实粗糙的三维多交叉裂隙模型，并结合数值模拟，系统揭示了剪切位移对裂隙网络中非费克溶质运移的影响机制，主要结论如下： 1. 裂隙空隙结构的复杂性控制着涡旋场的组织。均匀开度促进连通涡旋，而受限开度和曲折接触则导致流动分隔和离散涡旋域。 2. 粗糙裂隙中的溶质侵入遵循两个阶段：先是对流主导的优势通道快速突破，随后是浓度梯度驱动的小空隙扩散侵入。 3. 涡旋区域对溶质的滞留是通过动态捕集-释放平衡实现的，而非单向捕获。每个涡旋在稳态流下具有平衡的溶质储存容量。 4. 裂隙交叉点的溶质分布和混合模式主要由相交分支间的界面接触面积控制，需要结合出口方向界面面积和入口方向流速进行综合分析。
科学价值：本研究深化了对复杂裂隙网络中非费克运移机理的理解，特别是明确了剪切作用通过改变几何结构来影响涡旋发育、溶质分配和滞留的详细路径。提出的“涡旋动态平衡”假说为解释裂隙介质中的溶质延迟和储存提供了新的理论视角。
应用价值：研究成果可直接应用于需要精确预测地下污染物迁移、地热流体循环、油气驱替效率或二氧化碳封存安全性的工程领域。考虑剪切历史和多交叉点几何效应的模型，能够更准确地评估实际地质条件下流体和溶质的运移行为，为相关工程的设计、风险评估和优化提供更可靠的科学依据。
研究附录提供了丰富的补充信息，包括不同采样间隔（dx）对粗糙度参数计算影响的敏感性分析（附录A），以及所有模型在不同截面下涡量值分级统计的完整数据表（附录B）。这些内容为其他研究者复现或深化本研究提供了详细的数据基础和方法参考。同时，论文中详尽的符号命名表（附录C）也增强了文章的可读性和规范性。