一种用于分层异质介质中多裂纹扩展的梯度聚合能量变分XFEM方法

《COMPUTER METHODS IN APPLIED MECHANICS AND ENGINEERING》：A gradient-aggregated energy variational XFEM for multi-fracture growth in layered heterogeneous formations

　　**摘要**
在层状异质地层中，多裂纹扩展是一个关键问题，其特点是分支裂纹之间的竞争性传播，这在各种工程应用中都非常重要，尤其是在水力压裂中。然而，传统的基于能量的变分有限元方法（XFEM）结合基于能量的变分理论，在这种地质结构中会失败，因为界面处的能量梯度不连续性会导致搜索步骤不可靠，偏离能量下降方向。本文提出了一种新的变分XFEM方法，该方法能够聚合历史能量梯度，用于层状异质地层中的多裂纹扩展。在迭代过程中，该方法利用历史迭代路径中的梯度信息构建分段二次模型，该模型在非光滑区域使用能量最优的梯度聚合权重来确定搜索步骤。一旦迭代离开非光滑区域，搜索步骤可以退化为标准的牛顿形式，从而实现类似牛顿的快速收敛率。此外，该方法独特的空步机制确保了每次迭代时总能量的稳定和充分下降，从而提高了整体求解效率。通过四个基准裂纹扩展案例的数值结果验证和与实验数据的对比，证明了该方法的计算精度和预测可靠性。最后，通过对涉及层状异质地层中多裂纹扩展的水力压裂案例的能量变化和迭代收敛性的分析，显示了该方法在能量梯度不连续情况下的总能量最小化方面的优越性能。所提出的方法为模拟实际地质结构中多裂纹系统中各分支裂纹之间的复杂几何和物理相互作用提供了一个通用的XFEM框架。

**引言**
裂纹扩展是多种工程应用中遇到的基本问题，包括结构完整性评估、材料设计和地下能源提取[[1], [2], [3]]。特别是在油/气开采的致密地层中，通过高压流体注入来启动和维持裂纹扩展，以克服原位应力和宿主岩石的固有机械阻力[[4,5]]。由于岩石裂纹、流体流动和长期矿物沉积引起的层状异质性之间的强耦合，这一过程非常复杂[[6], [7], [8], [9], [10], [11]]。在层状异质地层中，不同材料层之间的丰富界面引发了多样的裂纹扩展模式（包括穿透、偏转和分支），并产生了更复杂的裂纹扩展演化[[12], [13], [14]]。为了研究层状异质地层中的裂纹扩展机制，过去几十年里发展了大量的理论和数值工作。大多数理论模型基于简化的裂纹几何形状来推导解析解，这限制了它们在更一般裂纹扩展情况下的应用[[15], [16], [17], [18]]。因此，开发了数值方法来处理裂纹扩展，如有限元方法（FEM）[[19,20]]、相场方法（PFM）[[21,22]]和黏聚区有限元方法（CEM）[[23,24]]。这些数值方法需要在裂纹附近进行局部网格细化，从而导致模拟裂纹扩展的计算成本较高。扩展有限元方法（XFEM）将先验解析解作为丰富函数引入近似空间，以表征裂纹体沿线的强不连续性和裂纹尖端附近的奇异性。通过这种方式，XFEM消除了围绕扩展裂纹进行局部网格细化的需要，显著提高了计算效率，并在处理裂纹扩展方面显示出优势[[25,26]]。

层状异质地层中的裂纹扩展表现为一个复杂的多裂纹系统，涉及不同层之间界面产生的分支裂纹之间的几何交叉和物理相互作用。在应用XFEM处理多裂纹扩展时，准确预测裂纹路径依赖于裂纹扩展准则来确定裂纹尖端的生长长度和弯曲角度[[27,28]]。基于线性弹性裂纹力学（LEFM），裂纹扩展准则可以分为两类：基于能量的准则（如最大能量释放率（MERR）准则[[29], [30], [31], [32]]和基于应力/应变的准则（如最大切向应力（MTS）准则[[33], [34], [35], [36]]）。在基于能量的准则中，裂纹生长评估的能量相关参数是从孤立裂纹周围的局部应力/应变场的解析解中得出的。对于多裂纹扩展问题，每个裂纹的扩展必须考虑它们相互作用产生的额外能量耦合项，而这些在之前的基于能量的准则中并未包括[[37,38]]。基于应力/应变的准则能够在LEFM理论内满足局部裂纹尖端应力/应变场的叠加原理，从而描述裂纹相互作用。然而，通过叠加所有裂纹尖端的应力场来确定最大周向应力的裂纹尖端扩展方向仍然是一项非平凡的任务[[39,40]]。因此，由于缺乏裂纹相互作用信息，基于能量和基于应力/应变的准则都无法描述多裂纹扩展。遵循多裂纹系统的能量最小化基本原理，能量变分理论被用来控制多裂纹扩展[[41], [42], [43]]。在这种理论中，对每个裂纹尖端生长长度进行二阶能量导数运算会产生一个Hessian矩阵。该矩阵的对角线元素代表每个分支裂纹的自导数项，而非对角线元素对应于交叉导数项[[43]]。交叉导数项定量描述了传统基于能量的准则中缺失的分支裂纹之间的能量耦合，使其特别适合处理多裂纹扩展问题。

多裂纹扩展的总能量最小化的数值实现严重依赖于二阶能量导数矩阵的确定性，该矩阵表征了总能量函数的局部曲率。正定的Hessian矩阵对应于凸能量函数，而负定或不定矩阵则表示非凸性，这需要不同的策略来寻找全局能量最小值[[44]]。基于这种分类，通过在XFEM中引入能量变分理论，在均匀域中对多裂纹扩展的数值建模取得了一些进展[[45], [46], [47], [48], [49], [50]]。Budyn等人首次建立了用于多裂纹扩展的能量变分XFEM（EV-XFEM）[[45]]。在这种方法中，通过计算二阶能量导数矩阵的最大子行列式来确定裂纹尖端的生长长度。当总能量函数是凸的时，从能量Hessian矩阵的最大子行列式获得的最优解正好对应于总能量函数的全局最小值。局部MTS准则仍然用于裂纹尖端的弯曲角度。作为将能量变分理论引入XFEM的开创性研究，这种方法已被广泛用于研究各种情景下的多裂纹扩展，包括Xu等人研究的包含宏观和微观尺度裂纹的情况[[46]]，以及Wang等人对多裂纹扩展的数值分析，考虑了新裂纹的成核[[47]]。基于Budyn等人工作的局限性（仅限于凸能量最小化），Sutula等人[[48], [49], [50]]将EV-XFEM扩展到了非凸情况，这也表明了多裂纹系统的强烈竞争性生长。提出的基于能量梯度的策略使用能量梯度来确定每次迭代中的最优步骤，从而在能量最小化迭代收敛时可以获得裂纹尖端的生长长度和方向。这些工作在均匀介质中的多裂纹扩展预测方面取得了巨大进展。然而，对于层状异质介质中更一般的多裂纹扩展情况，上述针对均匀情况的EV-XFEM家族（包括Budyn等人的工作和Sutula等人修改的工作）面临以下挑战：一方面，界面之间的材料不匹配在裂纹尖端引入了不同的应力奇异性[[51]]，结合固体-流体耦合，使得计算一阶和二阶能量导数项更加困难；另一方面，双材料界面在裂纹尖端弯曲角度处引入了能量梯度的不连续性，导致在非光滑区域附近出现次梯度集。由于次梯度集中存在多个可接受的梯度，依赖于唯一定义的梯度和Hessian矩阵的传统迭代方法（如牛顿方法）无法识别表征总能量下降的可靠搜索步骤[[52]]。因此，迭代过程可能会陷入局部最小值甚至发散[[53]]。因此，需要进一步发展EV-XFEM以考虑层状异质情况下的多裂纹扩展。

为了解决上述挑战，本文提出了一种基于梯度聚合的能量变分XFEM（GEV-XFEM），用于层状异质地层中的多裂纹扩展，通过以下关键改进：首先，从Zhong等人[[54]]那里继承了全面的裂纹尖端丰富框架，以描述XFEM中双材料界面处的裂纹尖端应力奇异性，为准确计算总能量最小化提供了基础；其次，开发了一种新的基于束牛顿的耦合迭代方案，以考虑界面处的梯度不连续性，提供裂纹尖端生长长度和弯曲角度的能量最优搜索步骤。具体来说，通过束牛顿方法解决具有梯度不连续性的总能量最小化问题。该方法收集历史迭代路径中的梯度信息，并建立带有历史能量梯度聚合的分段二次束模型，以捕捉总能量函数的非光滑特性。通过将束模型最小化问题重新表述为凸二次规划（QP）问题，可以通过求解凸QP问题的Karush-Kuhn-Tucker（KKT）方程来确定裂纹尖端生长长度和弯曲角度的最优搜索步骤，这对应于历史梯度的能量最优聚合权重。此外，该方法的空步机制确保了迭代过程中的稳定和可靠的搜索步骤以及足够的能量下降。通过这些改进，GEV-XFEM在处理一般裂纹情况（包括均匀和异质情况）的多裂纹扩展时实现了高精度和稳健的收敛性。

**本文结构**
第2节介绍了层状异质地层中多裂纹扩展的控制方程及其XFEM离散化。第3节介绍了基于梯度聚合的能量变分XFEM（GEV-XFEM），包括用于解决总能量最小化问题的基于束牛顿的耦合迭代方案。第4节通过两个基准裂纹示例验证了所提出方法的精度。第5节通过涉及层状异质地层中多裂纹扩展的水力压裂案例分析了总能量变化和迭代过程，以说明所提出方法的优越性。第6节对所提出的GEV-XFEM得出了一些有价值的结论。在总能量导数计算的第三部分中，对于分层异质地层中多裂缝扩展的应用，梯度聚合能量变分XFEM（GEV-XFEM）依赖于对总能量函数的一阶和二阶导数的准确高效计算，本文通过刚度导数方法实现了这一点。这里建立了一个双层模型，其中一条裂缝从左中间边缘开始扩展，以验证刚度导数方法在总能量计算中的有效性。

**结果与讨论**
在本节中，研究了一系列在具有不同弹性模量、断裂韧性和层取向的分层异质地层中发生多裂缝扩展的水力压裂案例。对于所有测试案例，都实施了由Sutula等人改进的用于均匀介质的能量变分XFEM（EV-XFEM）[[48], [49], [50]]以及梯度聚合能量变分XFEM（GEV-XFEM）。通过对总能量变化和迭代过程的比较分析，...

**结论**
基于传统牛顿方法且依赖于多裂缝扩展平滑性假设的能量变分XFEM（EV-XFEM）在分层异质地层中面临两个主要挑战：复杂的裂缝尖端应力奇异性以及双材料界面处的能量梯度不连续性。本文提出了一种用于此类地质结构中多裂缝扩展的梯度聚合能量变分XFEM（GEV-XFEM）。

**致谢**
本工作得到了中国国家自然科学基金（12572234, U24B6004, 12388101, 124B2053 和 12232016）、中国教育部基础与交叉学科创新计划（JYB2025XDXM301）、青年教师科研创新能力支持项目（SRICSPYF-ZY2025050）以及中央高校基本科研业务费（WK2090000087）的共同支持。数值计算是在CRediT平台上完成的。

**作者贡献声明**
钟一伦：撰写——原始草案、数据收集与整理、概念构思。
余浩：撰写——审稿与编辑、监督、项目管理、资金筹集。
王权：资源与方法论支持。
陈秀媛：验证与软件实现。
孙继云：撰写——审稿与编辑、可视化处理。
Sundararajan Natarajan：撰写——审稿与编辑、形式化分析。
吴恒安：撰写——审稿与编辑、监督、项目管理工作。