基于电荷守恒的多层混合有限元方法用于求解稳态无感应磁流体动力学方程

《Applied Numerical Mathematics》：Multi-level mixed finite element methods based on charge-conservative for solving the stationary inductionless MHD equations

　　夏彦德|杨云波
中国云南省昆明市云南师范大学数学系，邮编650500

摘要
本文基于斯托克斯、奥森和牛顿线性化方法，构建并分析了三种用于稳态无感应磁流体力学方程的多级混合有限元方法。这些算法在嵌套网格层次上运行：非线性问题仅在最粗的网格上求解，而通过求解线性化的无感应磁流体力学问题，在逐渐细化的网格上获得近似解。与传统的单级混合有限元方法相比，所提出的多级策略在合适的网格层次下能显著降低计算成本。本文严格证明了所提方法的稳定性及最优误差估计，数值实验也验证了理论结果，并展示了这些方法的效率。

引言
磁流体力学（MHD）用于描述存在电磁场时的导电流体的行为。其控制方程由流体动力学的纳维-斯托克斯方程和电磁学的麦克斯韦方程组构成[1]–[4]。特别是当流体中电流产生的磁场可忽略，或磁雷诺数较低时，可适用无感应磁流体力学（IMHD）模型。这一简化模型被广泛应用于工程领域，如磁流体力学泵、钢铁连续铸造以及晶体生长装置等[5]–[8]。

设Ω?Rd（d=2,3）为具有利普希茨连续边界Γ??Ω的有界域，我们研究的IMHD方程如下：
?Re?1Δu+(u·?)u+?p?κJ×B=f，inΩ；
divu=0，inΩ；
J+???u×B=g，inΩ；
divJ=0，inΩ。
边界条件为u=0，J·n=0 on Γ，其中u为流体速度，J为电流密度，p为流体压力，?为电势，n为?Ω上的单位外法向量，Re为流体雷诺数，κ为耦合数。f和g表示外部力项，B为给定的外加磁场。对于d=3的情况，u=(u1,u2,u3)，J=(J1,J2,J3)，B=(B1,B2,B3)，f=(f1,f2,f3)，g=(g1,g2,g3)；对于d=2的情况，u=(u1,u2,0)，J=(J1,J2,0)，B=(0,0,B3)，f=(f1,f2,0)，g=(g1,g2,0)。

近年来，人们越来越关注使用无散有限元方法求解MHD方程[9]–[13]。对于无感应磁流体力学模型而言，电流密度的无散约束?·J=0是一条重要的物理定律，可通过后处理或混合有限元方法来实现。例如，Ni等人[14]–[15]提出了适用于无感应磁流体力学方程的守恒电荷有限体积格式，并通过后处理计算Jh。许多研究还致力于构建能够保证无散性的混合有限元方法，参见[16]–[19]及其中的参考文献。针对IMHD系统(1.1)，张和丁[20]提出了一种守恒电荷的有限元近似方法，通过H(div, Ω)×L2(Ω)-一致对离散化电流密度和电势，从而得到严格无散的离散电流密度。

多级策略也受到了广泛关注[21]–[22]。实际上，多级算法是从两级方法发展而来的[23]–[24]。一般两级方法的核心思想是：首先在较粗的网格上求解较小的非线性问题，然后在较细的网格上求解较大的线性问题，最后在较粗的网格上进行线性修正。两级方法已广泛应用于纳维-斯托克斯方程的研究，参见Girault[25]–[26]、何等人[27]–[28]以及Layton[29]–[30]的工作。戴和陈[31]提出了另一种两级策略，省略了第三步中的粗网格修正，而是求解具有与第二步相同刚度矩阵但右侧不同的线性化系统。这种变体也得到了广泛研究，参见[32]–[34]。对于MHD方程，Layton等人[35]–[36]研究了两级牛顿有限元算法，即先在较粗的网格上求解非线性MHD系统，再在较细的网格上求解线性问题。张等人[37]进一步扩展了这一思路，提出了用于稳态不可压缩MHD方程的两级耦合修正和独立并行修正有限元方法。最近，周等人[19]将这一理念进一步应用于稳态IMHD方程的研究。尽管多级方法比标准单级方法高效得多，但其误差估计的推导通常比两级方法更为复杂。早期用于稳态纳维-斯托克斯方程的多级有限元方法由Layton[38]以及Layton和Leferink[21]提出。随后，Layton等人[39]提出了一种通用的J级牛顿基方法，该方法采用一系列网格：完整的非线性方程仅在最粗的网格上求解，而在更细的网格上则通过求解线性化的纳维-斯托克斯问题获得近似解。他们通过两个数值示例证明了该算法的高效率。后来，何和李[40]提出并分析了用于稳态不可压缩纳维-斯托克斯方程的多级间断有限元方法。最近，邱[41]将这种多级策略扩展到稳态不可压缩MHD方程，提出了三种基于斯托克斯、奥森和牛顿线性化的多级混合有限元算法。

本文遵循[41]的框架，提出了三种用于稳态IMHD系统(1.1)的多级算法。我们的核心思想是生成一系列离散迭代序列，用以逼近非线性问题(1.1)的解。具体而言，这些算法利用网格层次结构：完整的非线性IMHD系统仅在最粗的网格上求解，而在更细的网格上，则通过基于斯托克斯、奥森和牛顿线性化方法，围绕前一步网格的解来求解线性化的IMHD问题，从而获得近似解。在有限元近似方面，我们采用H(div, Ω)×L2(Ω)-一致对离散化电流密度和电势，确保离散电流密度严格无散。在精确解具有合理正则性的假设下，我们为速度建立了最优H1误差估计，为其他变量建立了最优L2误差估计。数值模拟表明，所提出的多级混合有限元方法相比单级方法能大幅节省计算时间。

本文的其余部分结构如下：第2节介绍索伯列夫空间的符号及连续问题的初步结果；第3节阐述离散问题的有限元离散化及初步估计；第4节提出三种多级混合有限元方法，并分析其稳定性和收敛性；第5节给出数值实验，以展示性能并验证理论收敛速率；最后，第6节总结主要研究成果。

在本文中，C（以及C1、C2等）表示可能在不同出现处发生变化的通用正常数，这些常数仅取决于方程的系数和域Ω。

数学预备知识
首先回顾索伯列夫空间的符号。记Wk,p(Ω)（1≤p≤∞）为具有通常范数‖·‖k,p和半范数|·|k,p的（整数k∈N或分数k>0）索伯列夫空间。空间Hk(Ω)为具有范数‖·‖k和半范数|·|k的k阶标准希尔伯特空间。通常，L2(Ω)配备内积和范数uv=∫Ωuvdx，∥u∥0=∥u∥L2(Ω):=(u,u)1/2，对所有u,v∈L2(Ω)成立。

在本文中，向量值函数和空间用粗体表示；例如，Wk,p(Ω)?Wk,p(Ω)d，Hk。

IMHD方程的有限元近似
接下来，考虑问题(2.1)的有限元离散化。设{Khj}j=0J为一组符合要求的、形状规则的、准均匀的Ω的划分（二维为三角形，三维为四面体），使得Ω=?K∈KhjK，hK:=diam(K)，hj:=maxK∈KhjhK，且h0>h1>?>hJ。假设Khj是通过细化Khj?1得到的。对于任意整数k≥0，记Pk(K)为K上总次数不超过k的多项式，并令Pk(K)?Pk(K)d。基于Khj，我们定义符合要求的有限元子空间Vhj。

多级算法
在本节中，我们基于斯托克斯、奥森和牛顿线性化方法设计了三种多级算法，以提高精度和效率。我们为所提方法建立了稳定性并给出了误差估计，同时推荐了合适的网格层次结构。设域Ω被有限元网格族Kh0,Kh1,…,KhJ划分，相应的空间(Vhj,Qhj,Dhj,Shj)（j=0,1,…,J）满足嵌套关系(Vh0,Qh0,Dh0,Sh0)?(Vh1,Qh1,Dh1,Sh1)???(VhJ,QhJ,DhJ,ShJ)。

数值实验
在本节中，我们通过二维和三维案例给出若干数值结果，以评估所提出的多级算法在解决稳态IMHD方程时的性能和精度。首先研究具有光滑解的二维和三维问题，以验证所提方法的收敛速率。下一个例子涉及在非凸L形区域上的奇异解，用于验证收敛速率并评估多级方法捕捉奇异解的能力。

结论
本文提出了三种基于斯托克斯、奥森和牛顿线性化方法的多级算法，用于求解稳态IMHD方程。我们对所有所提方法都建立了严格的稳定性分析及最优误差估计。这些多级策略仅在最粗的网格上求解完整的非线性IMHD系统，而通过在逐步细化的网格上求解线性化的IMHD问题，高效地获得后续近似解。大量的数值实验也验证了这些方法的有效性。

CRediT作者贡献声明
夏彦德：撰写——初稿、软件实现、形式分析；杨云波：撰写——审稿与编辑、指导、方法设计、资金获取、形式分析。