Brinkman-双重孔隙率-孔隙弹性模型及其不连续Galerkin方法

《Computers & Mathematics with Applications》：A Brinkman–dual–porosity–poroelasticity model and its discontinuous Galerkin methods

　　贾良边|瑞丽|陈张新
陕西师范大学数学与统计学院，西安，710062，中国

摘要
本文提出了一种Brinkman-双孔隙-多孔弹性模型，用于模拟填充管道和双孔隙可变形多孔介质区域中的流体流动问题。目的是阐明两种多孔弹性介质（即基质和微裂缝）对流体流动行为的影响。该模型采用不连续Galerkin有限元方法求解。这两种系统需要通过一组界面条件进行相互作用，包括一个无交换条件。在理论分析中，通过引入总压力来重新表述双孔隙多孔弹性系统，以避免可能由系统引起的锁定现象。空间采用不连续Galerkin方法离散化，时间采用后向Euler方法离散化。建立了半离散化和全离散化的数值方案。证明了弱形式和半离散化方案的适定性，分析了全离散化方案解的存在性和唯一性，并推导了半离散化和全离散化方案的误差估计器。对于数值模拟，验证了数值解的理论收敛速率和界面条件的正确性，并模拟了实际水力压裂中的流体流动情况，以验证该方法的有效性和精度。

引言
近年来，关于自由流动与多孔介质流动耦合问题的研究[1]、[2]、[3]、[4]、[5]、[6]、[7]、[8]、[9]、[10]、[11]、[12]、[13]引起了众多研究人员的关注。随着石油和天然气勘探开发技术的进步，非常规油气资源开采方法[14]、[15]的发展得到了极大推动，特别是在裂缝孔隙储层和致密储层的开采方面。裂缝孔隙储层的储层结构复杂，包括基质、微裂缝、宏观裂缝和喀斯特洞穴等，其尺度变化很大，具有显著的多尺度特性和非均匀性。研究人员发现，为了确保开采过程中管道的不可压缩性，需要向管道中注入支撑剂以防止管道闭合。此外，由于后期地质活动沉积的泥浆、沙子和方解石等物质可能进一步填充裂缝和喀斯特洞穴，这些储层的非均质性进一步加剧。因此，采用Brinkman方程[16]、[17]、[18]、[19]、[20]、[21]、[22]、[23]、[24]、[25]、[26]、[27]来描述这种填充流体流动情况。

在致密储层中，水力压裂过程会在水平井周围形成裂缝和改良区。这些区域通常由三种不同的介质组成：基质、微裂缝和水力裂缝[28]、[29]、[30]、[31]。在压裂过程中，储层内的应力场分布会发生显著变化，可能导致基质变形。因此，在非常规储层开采过程中，基质变形与流体流动之间的相互作用是一个关键考虑因素。对于天然裂缝储层，尤其是那些密集分布且无法单独表征的微裂缝，采用双孔隙方法[32]、[33]来模拟流体流动、流体传输和基质变形过程。该方法采用双孔隙多孔弹性介质，包括基质系统和微裂缝系统。这种双孔隙结构固有的地质复杂性给碳氢化合物开采带来了重大挑战。

应用于Stoke方程耦合和Biot系统[39]、[40]、[41]、[42]、[43]、[44]的流动模型[34]、[35]、[36]、[37]通常只关注多孔弹性介质中基质对流体流动行为的影响，而忽略了存在的微裂缝。此外，在实际非常规油气资源开采中，可能存在填充的自由流动区域。这种限制常常影响模型对复杂多孔弹性介质中流体流动的准确描述。为了解决这些问题，本研究采用Brinkman方程来描述填充管道内的流体动力学。这些方程基于Stokes方程[45]、[46]、[47]、[48]、[49]、[50]、[51]，考虑了复杂的渗透现象。此外，为了考虑基质的可变形性，并评估基质和微裂缝双孔隙多孔弹性介质对流体流动行为的影响，将双孔隙多孔弹性系统（DPP）[32]、[52]、[53]、[54]、[55]、[56]、[57]、[58]引入模型中，以更准确地模拟流体流动行为。

本研究开发了一种新的Brinkman-双孔隙-多孔弹性模型，通过五个界面条件将Brinkman流动动力学与可变形多孔介质耦合起来。该模型主要考虑管道（如水平井或水力裂缝）内的流体流动，以及流体在基质区域和微裂缝区域中的相互作用。与现有的Stokes-Biot模型[35]、[37]、[59]（关注自由流动区域和单孔隙可变形多孔介质中的流体流动问题）相比，当前模型考虑了填充自由流动区域和双孔隙可变形多孔介质中的流体流动。双孔隙可变形多孔介质包括两个相互渗透的区域：基质区域和微裂缝区域，每个区域都表现出不同的水力和传输特性。这些区域通过质量交换项耦合，流体仅从基质传递到微裂缝，绕过了与管道在界面处的直接相互作用。采用两种不同的压力：基质压力和微裂缝压力。在这种情况下，界面处的流动交换仅发生在微裂缝和管道之间。因此，Brinkman流动与双孔隙可变形多孔介质流动之间的相互作用通过一系列界面条件实现，包括基质与管道之间的无交换条件，以及微裂缝与管道之间的四个界面条件。

在自由流动与多孔介质流动耦合问题的数值研究中，采用了多种有限元方法[59]、[60]、[61]、[62]、[63]、[64]、[65]、[66]、[67]，其中在多边形网格上使用的不连续Galerkin有限元方法（DG）在处理多物理场、多区域离散化问题[68]、[69]、[70]、[71]、[72]、[73]、[74]、[75]、[76]、[77]、[78]、[79]、[80]、[81]方面具有许多优势。本文采用加权不连续Galerkin有限元方法[82]、[83]、[84]来数值求解Brinkman-双孔隙-多孔弹性模型。本研究的主要目的是扩展在多边形网格上的不连续Galerkin有限元方法（DGFEM），用于模拟Brinkman流动和双孔隙可变形多孔介质流动的耦合动力学。该公式采用经典的速度-压力变量来表示Brinkman方程。对于双孔隙多孔弹性系统，我们从[63]以及[85]、[86]、[87]、[88]、[89]、[90]的工作中获得灵感，通过引入四个场变量（位移、基质压力、微裂缝压力和总压力）来重新表述系统，从而有效避免了多孔弹性中的锁定问题。

在本研究中，基于参考文献[59]的理论框架，我们引入了总压力，同时考虑了基质和微裂缝压力。我们提出了一个半离散数值方案，用于Brinkman-双孔隙-多孔弹性模型，该方案通过WDG方法进行空间离散化。通过将Galerkin近似理论与微分代数方程（DAE）的相关性原理相结合，证明了弱形式和半离散方案的有效性，并推导了与半离散方案相关的误差估计器。然后基于半离散数值方案，使用后向Euler方法进行时间离散化，建立了全离散数值方案，并推导了全离散方案的误差估计器。最后，基于不同类型的网格剖面，对Brinkman-双孔隙-多孔弹性模型进行了五次数值实验，验证了模型数值解的理论收敛速率，证明了所提方法的有效性。此外，通过模拟实际水力压裂问题（包括开孔完井问题、压裂后水平井中的流体流动问题以及带有套管的多级水力压裂问题），进一步证明了所提方法的有效性。

本文的其余部分分为以下几节：第2节介绍了Brinkman方程，并通过引入总压力重新表述双孔隙多孔弹性系统，发展了Brinkman-双孔隙-多孔弹性模型。随后，对模型施加了初始条件、边界条件和界面条件。同时推导了Sobolev空间中的弱形式，并定义了双线性形式。第3节根据WDG方法开发了半离散数值方案，并验证了模型的适用性。第4节提出了全离散方案，并在后向Euler方法下进行了分析，提供了证明该方案解唯一性的具体步骤。第5节提出了半离散和全离散方案的误差估计器。第6节基于不同类型的网格剖面进行了五次数值实验，展示了模型的特性和有效性。第7节对整篇文章进行了简要总结。

在整篇文章中，符号?和?用于表示不等式在乘法常数范围内成立，这些常数与网格大小和时间步长无关。然而，它们可能取决于基本问题所固有的物理参数和与基本网格相关的形状规则常数。

**模型片段**
**耦合Brinkman-双孔隙-多孔弹性模型**
在本研究中，我们研究了在有界域Ω?Rd（d=2表示空间维度）上的Brinkman流动与双孔隙可变形多孔介质流动耦合问题，该域包括一个填充流体域ΩF和一个双孔隙可变形多孔介质域Ωp。这两个域都具有Lipschitz连续边界。此外，Ωp包含基质区域和微裂缝区域。设Γ=?ΩF∩?Ωp，其中?Ωi(i=F,p)是Ωi的边界，并具有外部单位法向量nΩi。

**弱形式的适定性**
首先，将计算域Ωi(i=F,p)划分为若干不相交的开多边形元素Ki，Ωi￣=?Ki∈ThiKi￣。设hi,K表示Ki∈Thi中元素的直径，he表示Ki中边e的长度。我们定义网格参数hi=maxKi∈Thihi,K。在此背景下，设EFo为ThF中所有内部边e的集合，EFD为ThF中属于Dirichlet边界的所有边e的集合，EΓ为...

**全离散方案**
在本节中，使用后向Euler方法对方程（3.7）-（3.12）进行时间离散化，以推导出全离散方案。时间步长固定为Δt=TN，tn=nΔt，n=0,1,?,N，为离散时间。f是一个通用标量或向量场，我们设置fn=f(tn)和?tnf=fn?fn?1Δt。在定理2的证明中构建的离散初始数据用作时间依赖问题的起点。在每个时间迭代n=1,2,?,N中，我们寻求解(uhFn,phF...

**误差分析**
在本节中，我们提出了半离散和全离散数值方案的误差估计器。这一过程的第一步是引入两个引理[69]、[100]。

**引理3（插值不等式）**
对于任何非负整数kuF,kpF,kηp,kpf,kpm,k?p，对于任何vF∈XF, qF∈QF, ξp∈Xp, qf∈Qf, qm∈Qm, ψp∈Wp，存在插值多项式ΠuFvF∈PkuF(K),ΠpFqF∈PkpF(K),Πηpξp∈Pkηp(K),Πpfqf∈Pkpf(K),Πpmqm∈Pkpm(K),Π?pψp∈Pk?p(K)，使得∥vF?ΠuFvF∥L2(K)+hF∥?(v...

**数值示例**
本节展示了使用所提方法对Brinkman-双孔隙-多孔弹性模型进行的五次数值测试结果。首先给出初始值ηp0(x), pm0(x), pf0(x)，并使用L2-正交投影找到ηhp0,phm0,phf0。然后，对于n≥1，使用（4.7）计算序列(uhFn,phFn,ηhpn,?hpn,phfn,phmn)，使其满足给定方案。在以下示例中，我们取kuF=kηp=kpm=kpf=l，kpF=k?p=l?1。首次测试验证了理论计算结果的准确性。在本研究中，我们采用了Brinkman-双孔隙率-孔隙弹性模型，该模型旨在模拟填充管道和基体区域以及微裂缝双孔隙率可变形多孔介质中流体的复杂流动现象。首先，引入了总压力这一参数，以重新构建双孔隙率孔隙弹性系统，从而克服系统可能出现的锁定现象。同时，还确定了Brinkman方程耦合过程中的初始条件、边界条件以及五个界面条件。