变密度TIMHD模型下线性和守恒电荷有限元方法的收敛性分析

《Communications in Nonlinear Science and Numerical Simulation》：Convergence analysis on linear and charge conservative finite element methods for a variable density TIMHD model

　　金淼人|崔霞|尹莉 中国北京应用物理与计算数学研究所，邮编100088

摘要
我们研究了一种变密度热耦合不可压缩无感应磁流体力学模型的数值模拟方法。首先，我们基于隐显混合技术以及后处理速度方法，构建了一种线性且满足电荷守恒的有限元方案。由于每个时间步仅需求解较小的线性子问题，该方案能够高效实现。接着，我们证明了该方案具有电荷守恒性。通过引入适当的有限元投影算子以及一种新的归纳假设论证方法，我们克服了系统强耦合与非线性以及数值方案解耦结构所带来的困难，从而得到了流体密度、速度、电流密度、温度和电势变量的L∞(L2)收敛性估计。基于这些收敛性分析结果，我们确定了该方案的条件稳定性。此外，为获得无条件稳定性，我们提出了另一种线性且满足电荷守恒的有限元方案，该方法通过同时求解速度、压力、电流密度和电势变量来体现其特性，并展示了相应的性质。最后，通过若干数值实验和比较，验证了这些方案的精确度与电荷守恒性，同时体现了高度解耦方案在效率上的优势。

引言
TIMHD模型描述了在外部电磁场作用下且存在热耦合的导电流体的动力学行为，其中由电流产生的磁场与外部磁场相比可以忽略不计[1]，[2]，[3]，[4]，[5]，[6]。温度梯度会导致密度变化，进而驱动流体运动。在某些情况下，这种密度变化相对较小，因此在布辛涅斯克近似下，模拟中通常将密度视为恒定参数。然而，在许多与TIMHD相关的科学与工程应用中，较大的温度差异会导致明显的密度变化，此时布辛涅斯克近似不再适用[7]，[8]，[9]，此时就需要使用变密度TIMHD模型。

设Ω为Rd（d=2,3）中具有利普希茨边界?Ω的有界域，我们考虑以下变密度TIMHD模型：
ρt+?·(ρu)=0，在Ω内
T:=Ω×(0,T]，
ρ(ut+u·?u)?μΔu+?p?j×B0?ρg=f，在ΩT内
?·u=0，在ΩT内
ρ(θt+u·?θ)?ηΔθ=ω，在ΩT内
κ?1j+???u×B0=χ，在ΩT内
?·j=0，在ΩT内
初始条件和边界条件为：
u|t=0=u0，θ|t=0=θ0，在Ω内
ρ|t=0=ρ0，在Ω￣内
u=0，θ=0，j·n=0，在?ΩT:=?Ω×(0,T]上

在系统(1.1)–(1.2)中，未知量包括流体密度ρ、速度u、压力p、电流密度j、温度θ和电势?。各物理参数的定义如下：μ=Re?1（Re为雷诺数）为运动粘度，B0为磁感应强度，g为重力加速度，f为驱动力函数，η为热导率，ω为热源强度，κ为电导率，χ为电源强度。n表示?Ω上的单位外法向量。此外，我们假设域内不存在真空状态，且存在两个正常数ρmin0和ρmax0，使得ρmin0≤ρ0(x)≤ρmax0，在Ω￣内成立。

当d=2时，我们将u表示为(u1,u2,0)，j表示为(j1,j2,0)，B0表示为(0,0,b3)，这样该模型同时适用于d=2和d=3的情况。

与不可压缩流动相比，由于出现了密度输运方程，变密度不可压缩流动本质上更为复杂，这导致了变量之间的更强耦合，给数值模拟带来了更大挑战。近年来，人们越来越关注变密度不可压缩纳维-斯托克斯方程的数值方法，相关研究也取得了一定进展，可参考[10]，[11]，[12]，[13]，[14]，[15]，[16]，[17]，[18]，[19]，[20]，[21]，[22]，[23]。然而，对于变密度不可压缩磁流体动力学系统而言，构建线性且无条件稳定的数值方案仍然面临巨大挑战，现有研究还相对较少。在[24]中，作者提出了一阶向后欧拉半离散方案以及基于隐显混合技术分别计算密度的完全离散有限元方案，并对半离散方案进行了严格的收敛性分析。在[25]中，作者提出了另一种向后欧拉半离散方案、完全离散有限元方案，以及利用隐显混合技术分离密度并采用Gauge-Uzawa方法解耦速度和压力的二阶向后差分公式半离散方案，同时还证明了向后欧拉半离散方案的收敛性。在[26]中，作者通过应用Gauge-Uzawa方法和零能量贡献方法，开发出了完全解耦的向后欧拉半离散方案和完全离散有限元方案，并证明了半离散方案的收敛性。在[27]中，作者通过采用后处理速度以及隐显混合技术分离密度这一辅助变量，构建了向后欧拉半离散方案和完全离散方案，并对低阶空间中的完全离散有限元方案进行了收敛性分析。在[28]中，作者针对更复杂的模型，提出了一种向后欧拉-有限元完全离散方案，该方案利用不连续伽辽金方法近似密度方程，并通过隐显混合技术实现其单独计算，同时证明了随着空间和时间步长趋于零，有限元解会趋近于问题的弱解。在[29]中，对于类似复杂的模型，作者提出了一种更高效的向后欧拉-有限元完全离散方案，该方案通过隐显混合技术分离密度模拟，利用线性外推技术解耦速度和压力，同时还采用了流线扩散方法提升方案的数值稳定性，并给出了其收敛性估计，不过速度变量的L∞(L2)收敛性在空间上并非最优。然而，据作者所知，目前尚未有针对变密度热耦合无感应系统(1.1)–(1.1)的有效数值研究。而对于密度恒定的磁流体动力学和TIMHD模型，我们可以通过一些研究来了解其数值模拟和分析方法[1]，[2]，[3]，[4]，[5]，[6]，[30]，[31]，[32]，[39]。

本文旨在为(1.1)–(1.1)模型开发简单、线性且高效的数值方案。为此，我们采用了多种技术组合。首先，由于直接从系统(1.1)–(1.1)出发分析方案较为困难，我们引入了一个新变量σ=ρ，从而得到一个更易于处理的等价系统。其次，我们采用有限元方法进行空间离散化，采用向后欧拉方法进行时间离散化。第三，我们通过后处理速度来确保稳定性，并解耦速度与这个新变量，同时运用隐显混合技术解耦电流密度、速度和温度，还对对流项进行了线性化处理。我们设计了两种方案：第一种方案的耦合程度较低，但其稳定性分析依赖于收敛性结果，因此是条件稳定的；第二种方案的耦合程度更高，且具有无条件稳定性。我们通过适当的有限元投影以及新的归纳假设论证方法，证明了这两种方案都满足电荷守恒性，分析了它们的稳定性，并得到了ρ、u、j、θ和?的L∞(L2)误差估计。数值实验和比较结果验证了理论预测，同时也展示了第一种方案的高效率。

本文结构如下：第2节介绍一些符号和基本概念；第3节提出完全离散方案I并证明其电荷守恒性；第4节分析其收敛性和稳定性；第5节开发完全离散方案II并确定其相应性质；第6节通过数值测试验证两种方案的性能；最后在第7节给出结论。

预备知识
对于实数q≥1和非负整数m，我们用Lq(Ω)和Wm,q(Ω)表示标准的勒贝格空间和索伯列夫空间，它们分别带有∥·∥Lq和∥·∥Wm,q范数。L2(Ω)中的内积和范数分别记为(·,·)和‖·‖。在希尔伯特空间情况下，Wm,2(Ω)简记为Hm(Ω)，其范数为∥·∥Hm。在全文中，粗体字母表示向量值函数空间，例如Lq(Ω)表示[d维]向量空间[Lq(Ω)]d。此外，我们还将……

完全离散方案I
在本节中，我们为系统(2.3)–(2.4)提出了一种有限元方案，并证明其具有电荷守恒性。

收敛性与稳定性分析
在本节中，我们对完全离散方案(3.1)–(3.4)的收敛性和稳定性进行分析。首先在4.1小节推导σ、ρ、u和θ的误差估计，然后在4.2小节和4.3小节分别得到j以及?的误差估计，最后在4.4小节阐述稳定性问题。

为便于分析，我们引入了几种有限元投影算子。具体而言，给定σ∈F，其L2投影记为Sσ∈Fh：(σ?Sσ,rh)=0，?rh∈Fh。

完全离散方案II
为了获得无条件稳定的算法，我们同时求解耦合的uhn+1、phn+1、jhn+1和?hn+1，由此得到以下方案：
步骤1：找到σhn+1∈Fh，使得对所有rh∈Fh都有(σhn+1?σhnτ,rh)+(P0huhn·?σhn+1,rh)=0。
步骤2：找到(uhn+1,phn+1,jhn+1,?hn+1)∈(Xh,Mh,Wh,Qh)，使得对所有(vh,qh,ch,ψh)∈(Xh,Mh,Wh,Qh)都有(σhn+1σhn+1uhn+1?σhnuhnτ,vh)+μ(?uhn+1,?vh)+(ρhn+1uhn·?uhn+1,vh)+12(uhn+1?·(ρhn+1uhn),vh)?(?·vh,phn+1)?(jhn+1×B0,vh)?(ρhn+1g,vh)=(fn+1,vh)，且(?·uhn+1,qh)=0，κ?1

数值结果
在本节中，我们通过若干数值测试，验证方案I和II的预测精度和电荷守恒性，并比较它们的精度和效率。

结论
在本文中，我们首先为变密度TIMHD系统构建了一种线性且满足电荷守恒性的有限元方案。在密度方程的离散化过程中使用了后处理速度，同时运用隐显混合技术来处理非线性项以及σ、u、θ和j之间的耦合关系，这使得该方案的实现十分简便。在理论分析方面，我们首先证明了该方案具有电荷守恒性。随后，通过运用多种有限元投影算子并提出一种新的……