一种用于阻尼不可压缩流的并行罚函数有限元方法

《Communications in Nonlinear Science and Numerical Simulation》：A parallel penalty finite element method for damped incompressible flows

　　

摘要
本研究提出了一种并行罚函数有限元方法，用于模拟由Navier-Stokes方程控制的阻尼稳态不可压缩流动。为了建立数值方法的理论基础，我们首先推导了罚函数有限元解的$H^1$和$L^2$误差估计，然后获得了定义在子域上局部解的先验局部估计。所提出的方法在全局复合网格上进行局部细化后，在重叠子域上独立计算局部解，使得该方法可以基于现有的顺序求解器轻松实现。利用局部先验估计这一理论工具，我们推导了算法近似解的误差界。此外，进行了若干数值实验以验证理论分析的正确性和所提方法的有效性。

引言
不可压缩Navier-Stokes方程是描述不可压缩流动的基本控制方程，具有极其广泛的应用范围。它们被采用于流体动力学与热传递耦合问题的建模（见[1]）、海洋动力学系统（见[2]），以及人体动脉和静脉中血液流动的数学表示（见[3]）。由流动阻力引起的阻尼效应可以有效地表征多孔介质流动、摩擦阻力和各种能量耗散等物理机制（见[4]、[5]）。因此，将阻尼项引入经典Navier-Stokes方程对于研究地球物理学和海洋声学中的不可压缩流动阻力问题具有重要意义（见[6]、[7]、[8]）。

设$Q$为$\mathbb{R}^n$（$n=2$或$3$）的一个有界子集，其边界$\Gamma=\partial Q$满足Lipschitz连续条件。以下是稳态不可压缩阻尼Navier-Stokes方程：
$$-\nu\Delta u+(u\cdot\nabla)u+\alpha|u|^{l-2}u+\nabla p=f,\quad \text{in }Q,$$
$$\nabla\cdot u=0,\quad \text{in }Q,$$
$$u=0,\quad \text{on }\Gamma.$$

其中，$u:Q\to\mathbb{R}^n$表示速度矢量场，$p:Q\to\mathbb{R}$和$f:Q\to\mathbb{R}^n$分别表示压力变量和给定的体积力函数。正参数$\nu$指运动粘度，与Reynolds数$Re$成反比。此外，项$\alpha|u|^{l-2}u$表示阻尼效应，其中$\alpha>0$为阻尼系数，指数$l$满足条件$2\leq l<\infty$。

学者们对此类方程进行了广泛的理论研究并获得了坚实的数值结果。Zhao和Li [9]研究了非线性阻尼二维Navier-Stokes问题的长时间动力行为。Cai和Lei [10]建立了受非线性阻尼控制的不可压缩流体强解和弱解的存在性。已证明阻尼项可以引起线性双曲问题的有限时间爆破解（见[11]、[12]、[13]）。此外，人们付出了大量努力来开发和分析相关问题的混合有限元方法。包括稳定化混合有限元格式（见[14]）、双网格和并行算法（见[15]、[16]、[17]）以及求解阻尼Navier-Stokes型变分不等式的数值策略（见[18]）在内的多种数值方法已在理论上建立并经过数值验证。

正如[19]、[20]、[21]、[22]、[23]、[24]、[25]中所讨论的，罚函数方法被用于松弛散度无约束条件并形成伪可压缩方案。对于系统(1.1)–(1.3)，我们采用$(u_\varepsilon, p_\varepsilon)$来近似原解$(u, p)$。这些罚变量满足罚函数形式的稳态阻尼Navier-Stokes方程（见[26]、[27]、[28]、[29]）：
$$-\nu\Delta u_\varepsilon+B(u_\varepsilon,u_\varepsilon)+\alpha|u_\varepsilon|^{l-2}u_\varepsilon+\nabla p_\varepsilon=f,\quad \text{in }Q,$$
$$\nabla\cdot u_\varepsilon+\varepsilon\nu p_\varepsilon=0,\quad \text{in }Q,$$
$$u_\varepsilon=0,\quad \text{on }\Gamma.$$

在此表达式中，$\varepsilon\in(0,1)$为罚系数，而$B(u,\varsigma)=(u\cdot\nabla)\varsigma+\frac{1}{2}(\nabla\cdot u)\varsigma$表示修正的双线性项。该罚函数方法最初由Temam在[30]中提出，并由Codina在[31]、Shen在[21]、Carey和Krishnan在[32]以及He等在[22]、[23]中进一步推广到Navier-Stokes方程。

在本工作中，我们开发了一种耦合求解策略，将罚函数型混合有限元离散化方法与完全重叠域分解算法相结合，用于由Navier-Stokes方程控制的阻尼不可压缩流动的高效数值模拟。所开发的数值方法采用满足离散LBB（inf-sup）稳定性约束的混合有限元空间。通过罚函数逼近来缓解不可压缩性约束，并采用并行域分解技术进行计算加速。在此并行域分解框架下，每个计算节点使用在其子域周围进行局部细化的全局网格，在自己的子域上独立计算局部罚解。该局部罚函数有限元模块自适应地校准与网格尺寸对应的罚系数，并获得最优收敛阶。此外，并行罚函数算法将并行空间分解结构与局部罚函数有限元技术相结合，在计算和并行效率方面实现了显著改进。

对于所构造数值方案的理论探索，第一步是建立与定义在指定子域上的罚函数有限元解相关联的局部先验估计。然后对局部算法和并行罚函数算法进行系统的误差分析。结果表明，所提出的方法在粗网格和细网格尺度之间具有合理的比例关系时可以实现最优收敛率。此外，构建了各种数值测试并执行以检验分析预测，并表明所开发的方法可以应用于实际工程任务。

本工作的主要创新和贡献总结如下。首先，提供了定义在子域上的罚函数有限元近似的局部先验估计的严格数学推导，这具有独立的研究价值，并为设计和分析罚函数离散化方案奠定了重要的理论基础。其次，构建了局部和并行罚函数离散化方法来模拟阻尼稳态不可压缩流动，并进行了相应的收敛性研究。第三，进行了多种典型数值算例，包括解析基准测试、方腔驱动流和绕圆柱流动，以验证理论结论的正确性和所提方法的有效性。

本文其余部分按以下结构安排。在第2节中，我们构建了阻尼稳态Navier-Stokes系统的罚函数格式，并相应地获得了对近似解的误差估计结果。第3节建立了与罚函数有限元数值解相关联的局部先验估计。在第4节中，在有限元离散化的基础上，讨论了带局部罚处理并行罚函数算法的构造和验证，并构建了相应的收敛性分析。在第5节中，通过数值算例验证所提方法的有效性和可靠性。最后，第6节总结了本工作。

预备知识
按照[33]的标准记号，对任意整数$m\geq 0$，采用Sobolev空间$W^{m,p}(Q)^n$及其范数$\|\cdot\|_{m,p,Q}$和半范数$|\cdot|_{m,p,Q}$。对于特殊情况$p=2$，我们用$H^m(Q)^n$代替$W^{m,2}(Q)^n$，其标准范数简记为$\|\cdot\|_{m,Q}$。此外，我们定义$H_0^1(Q)^n=\{u\in H^1(Q)^n: u|_{\partial Q}=0\}$，其中边界条件$u|_{\partial Q}=0$按迹的意义理解。应指出的是，在$H_0^1(Q)^n$上定义的范数和半范数互等价。相应地，我们使用

局部先验误差估计
在本节中，给出了定义在子域上的罚函数有限元解的先验局部界，作为局部算法和并行算法设计与分析的理论支持。

引理3.1（见[37]）设$\varphi\in C_0^\infty(Q)$且$\text{supp}\,\varphi\Subset Q_0\Subset Q$。则
$$\nu\|\varphi\theta\|_{1,Q}^2\leq\nu(\nabla\theta,\nabla(\varphi^2\theta))+c\nu\|\theta\|_{0,Q_0}^2\quad \forall\theta\in Y.$$

引理3.2（见[40]）设$V\Subset Q_0\Subset Q$，$f\in H^{-1}(Q_0)^n$且$\varepsilon=O(h_{Q_0})$。若假设A1–A4成立，且$(w_{\varepsilon h}, r_{\varepsilon h})\in Y_h(Q)\times R_h(Q)$满足
$$\nu(\nabla w_{\varepsilon h},\nabla v_{\varepsilon h})-(\nabla\cdot v_{\varepsilon h},r_{\varepsilon h})+b(u_{\varepsilon h},w_{\varepsilon h},v_{\varepsilon h})+b(w_{\varepsilon h},u_{\varepsilon h},v_{\varepsilon h})=(f,v_{\varepsilon h}),$$

带局部和并行实现的罚函数有限元方法
本节首先构造局部化罚函数有限元数值方案。在此基础上，通过将局部算法与完全重叠域分解策略相结合，进一步设计了并行罚函数有限元求解方法。

数值实验
本部分通过数值测试来确认理论分析的有效性并验证所开发算法的实际性能。具体而言，展示了三个数值算例，包括构造的具有精确解的基准测试、传统的方腔驱动流和经典的绕圆柱流动问题。数值测试基于公共软件FreeFem++ [58]构建。采用满足inf-sup条件的三角形Taylor-Hood $P^2$-$P^1$单元。

结论
针对带有非线性阻尼项的稳态Navier-Stokes方程，本研究基于完全重叠域分解技术构建了局部化并行罚函数有限元算法。在该框架下，实现了不同子域上的独立局部求解。通过采用全局构造网格并在子域周围进行局部细化，该方法具有编程方便和易于并行的特点。相关的理论证明和数值

资助
本工作得到重庆市自然科学基金（编号：ctc2021jcyjmsxmX1044）的资助。

CRediT作者贡献声明
袁付容：撰写——初稿，可视化，软件，方法论，形式化分析，数据整理。
尚跃强：撰写——审阅与编辑，验证，监督，方法论，概念化。

利益冲突声明
作者声明他们不存在已知的可能影响本文所报告工作的竞争性经济利益或个人关系。

袁付融 | 尚跃强
西南大学数学与统计学院，重庆，400715，中国