具有变开口的非线性混合维断裂模型的混合有限元方法

《Applied Mathematics Letters》：Mixed finite element method for nonlinear hybrid-dimensional fracture model with varying aperture

　　费孙|洪星睿 中国山东省济南市山东大学数学学院，邮编250100 摘要 本文提出了一种新的非线性混合维裂缝多孔介质流动模型，该模型中裂缝的开口大小会随裂缝内的流体压力而变化。这一模型更为贴近实际场景。文中采用混合有限元方法来求解这一非线性耦合系统，最后通过一些数值实验来验证该近似方案的正确性以及所提模型的有效性。引言 由于在环境与能源问题中的广泛应用，如地下水管理、碳捕获与储存、污染物处理等，研究裂缝多孔介质中的流体流动模式越来越受到关注。裂缝与多孔介质在渗透率及几何尺寸上的巨大差异，给准确建模和模拟裂缝介质中的流体流动带来了很大挑战。由于裂缝与多孔介质的渗透率不同，裂缝既可以作为促进流体流动的快速通道，也可能成为阻碍流体流动的地理屏障。由于裂缝的开口尺寸远小于周围基质的尺寸，若要精确模拟流体流动，则需要过度细化局部网格，从而带来较高的计算成本。因此，2005年Martin等人将裂缝视为n维多孔介质中的（n?1）维界面，进而提出了结合达西流的混合维裂缝模型[1]，并且该模型在解决不同类型裂缝问题方面的优势也得到了验证。此后，这一混合维模型得到了进一步发展，出现了一系列用于求解此类问题的高效方法[2],[3],[4],[5],[6],[7]。值得注意的是，由于混合有限元方法能够同时以相当高的精度近似标量量（如压力）及其通量（如流速），该方法在工程应用中被广泛使用。在[1]中，混合有限元方法是在符合要求的网格上应用的；而在[8]中，则通过在裂缝穿过的三角形上引入不连续基函数并进行弱接口条件处理，将最低阶Raviart–Thomas混合有限元方法应用于非匹配网格上求解该模型。需要指出的是，在上述所有文献中，裂缝的开口尺寸都是恒定的。但实际上，裂缝开口会受到岩石变形及裂缝内流体压力所施加的应力影响，而开口尺寸的变化又会影响裂缝的孔隙度与渗透率特性，为这一问题增添了更多复杂性。目前，只有部分工程和地质领域考虑了裂缝开口尺寸的变化[9],[10]。本文在经典的混合维裂缝模型中考虑了裂缝开口尺寸的变化。为简化问题，我们做了一些假设：首先，作用在裂缝上的应力会导致裂缝发生正应变和剪切变形，本文仅考虑正应变；其次，仍将周围基质视为多孔介质而非多孔弹性介质，即外部力为恒定值，这样一来，裂缝开口尺寸仅取决于裂缝内的压力。基于经典裂缝模型[1]并考虑上述因素，本文提出了一种新的非线性混合维裂缝模型，该模型中裂缝开口尺寸取决于裂缝内的压力。与一些仅考虑压力连续接口条件的工程和地质模拟不同，我们采用了适用于描述不同类型裂缝的通用Robin型接口条件。据作者所知，这是首篇考虑具有开口尺寸变化的稳态非线性混合维裂缝模型的论文。我们需要强调的是，在实际生产中，若简化处理，即假设裂缝开口尺寸恒定，可能会导致模拟结果不准确，尤其是对于深部地质构造中的裂缝。因此，对于实际应用而言，与恒定开口尺寸模型相比，我们所提出的模型具有更大优势。为获得具有相当精度的压力和流速数据，我们采用混合有限元方法对这一非线性耦合系统进行离散化，具体使用了Raviart–Thomas–Nedelec元素。需要注意的是，在我们之前的工作[11]中，已研究了非稳态情况，并采用了线性化隐式有限元方法来求解压力问题，实际上只需求解线性方程组即可。而本文研究的是稳态耦合问题，需要求解非线性方程组，同时还需考虑流速这一重要参数，这给计算带来了较大挑战。最后，本文给出了两个数值实验，其中一个用于验证近似方案的正确性，另一个用于展示所提模型的性能。本文的其余结构如下：第2节介绍模型问题，第3节阐述弱形式表达及混合有限元近似方案，第4节通过数值实验来验证所提模型及近似方案的准确性与性能。章节节选 非线性混合维裂缝模型 设Ω?Rd（d=2,3）为一个有界且凸的多边形/多面体区域，其利普希茨连续边界?Ω代表多孔基质，该区域被裂缝区域Ωγ分为Ω1和Ω2两部分。由于Ωγ的宽度远小于整个区域的尺寸，因此裂缝Ωγ通常被视为(d?1)维界面γ=?Ω1∩?Ω2[1]。设ni表示Ωi的外法向量，在界面γ上，n1=?n2。我们给出一些符号以便更好地描述这一混合维模型。混合有限元近似方案 本节除了介绍混合有限元近似方案外，还给出了模型问题((2.1))–((2.5))的弱形式表达。首先，我们定义一些空间：Wi={vi∈H(div;Ωi):vi?ni∈L2(γ)}，Wγ=H(divτ;γ)，W=W1×W2×Wγ，Mi=L2(Ωi)，Mγ=L2(γ)，M=M1×M2×Mγ，同时还定义了双线性形式aξ(q;u,v)=∑i(Ki?1ui,vi)Ωi+(α(q)uγ,vγ)γ+∑i(β(q)[ξui?ni?(1?ξ)ui+1?ni+1],vi?ni)γ，b(u,w)=∑i(??ui,wi)Ωi+(?τ?uγ,wγ)γ?∑i(ui?ni,wγ)γ，显然，对于固定的q，aξ(q;u,v)关于u和v是线性的。那么弱形式表达为：寻找(u,p)数值试验 本节通过两个数值实验来验证混合有限元方案的正确性以及所提模型的性能。在这两个示例中，计算域均为矩形，即Ω1=[?1,0)×[0,1]，Ω2=(0,1]×[0,1]，界面γ={x=0}×[0,1]。Ω的三角剖分是通过将形状相同的矩形一分为二来形成的，且采用符合要求的网格。我们考虑裂缝中的达西定律，各耦合条件均为非线性形式，进而应用混合有限元方案。结论 本文的主要贡献在于提出了开口尺寸随压力变化的非线性混合维裂缝模型。该模型更贴近现实场景，尤其是对于压力极高的深部地质构造中的裂缝。与原始模型相比，考虑开口尺寸的变化会因其非线性特性而给求解带来很大困难。我们采用混合有限元方法来求解这一非线性耦合系统，从而获得压力和流速数据。致谢 本研究得到国家自然科学基金项目编号12131014的支持。