中国学术期刊网

首页 > 社会科学II辑 > >求解一维磁流体方程的无震荡中心差分格式
社会科学II辑

求解一维磁流体方程的无震荡中心差分格式

时间:2018-09-11 21:03作者:admin打印字号:

Non-oscillatory central scheme for one-dimensional
ideal MHD equations
YAN Ke-qing, WEI Wei-ping
 (School of Science, Chang'an University, Xi'an 710064, China)
Abstract: In this paper,We utilize a family of high-resolution,non-oscillatory central schemes for the approximate solution of the one-dimensional MHD equations. Solution of one-dimensional shock-tube problems is carried out using second- and third- order central schemes. A qualitative comparison reveals an excellent agreement with previous results based on upwind schemes.Central schemes,however, require little knowledge about the eigenstructure of the problem-in fact,we even aviod an explicit evalution of the corresponding Jacobians, while at the same time they eliminate the need for dimensional splitting.
Keywords: Ideal magnetohydrodynamics equations; High-resolution central schemes;Non-oscillatory reconstructions
引言
在本文中对理想的磁流体力学方程(1)-(4)的近似解我们提出二阶和三阶无震荡中心差分格式
其中,是标量,分别代表质量密度和总内能,是速度场,它的-二范数表示为,和分别代表磁场及其-二范数。最后,结合压力,内能表示为,在这儿是固定的比热比。通过螺线管型的约束增强系统,也就是说,如果初始时满足条件,那么通过(3)式它将保持不变。
由于磁流体方程的内在复杂性,对于计算磁流体方程(1)-(4)的近似解,中心格式类方法将有效替代迎风格式类。本文我们使用的中心差分格式是基于交错网格上单元平均的演化。中心差分格式的交错型版本[1]-[3]将在下面第2小节中引出。中心差分格式去除了对雅克比矩阵的特征结构详细知识的需要。时间演化上使用简单的求积公式代替了迎风格式的近似黎曼求解器模块。这个方法不仅让我们节省了高成本的雅克比矩阵的特征分解,事实上让我们完全地避免了一维空间中的雅克比矩阵的高代价评估。此外,中心差分格式避免了在多维MHD系统中特别相关的维数分裂。本文用中心差分格式解决一维磁流体方程更简便,更精确。这一点,在后文的数值模拟中我们将看到证明。
在本文中我集中关注2个原型磁流体问题。我们使用二阶和三阶无震荡中心差分格式去计算一维Brio-Wu激波管模型[1]的近似解。两个不同的MHD激波管问题在3.1和3.2小节中研究。
我们的结果发现与解决同类问题的迎风格式[4][5]的模拟结果极好的一致。结果表明用中心格式去探测和解决这类模型的间断解,保持了有效性和简单性。
2   一维中心格式
我们实施交错网格上的预估-校正中心格式近似(1)-(4)的解。此格式分两步:
  1. 基于单元平均的点值,进行无震荡分片多项式重构;
  2. 基于交错单元平均,实现这些重构多项式的演化。
下面为了第3小节的计算,描述一维中心差分格式。
2.1 无震荡重构
方程(1)-(4)系统可写成一般守恒律形式
                                   (5)
(一维磁流体方程中和的显示形式参见第3小节)。在此一般形式的开始,我们借鉴Godunov近似守恒律间断解的开创性思维[2](见图1)。

图1  Godunov型中心差分格式
我们考虑滑动平均:;这时守恒律(5)可化为
                     (6)                                      
引入小的时间步长,且在区上积分我们得到
        (7)
到目前为止,(7)是精确的。现在我们通过分片多项式近似,,实现在离散时间层上(7)的解,即
,                          (8)
在这儿,在离散单元上,其中半整数网格点是单元界面断点,是代数多项式。取样,我们得到新交错单元平均,中心在,
           (9)
表达式(9)右端分为两步计算:
第一步,先用已知的单元平均,重构无震荡分片多项式近似。我们使用二阶分片线性函数
                                (10)
和三阶分片二次多项式函数
                        (11)
其中,、和是在点一阶和二阶导数的近似点值, 是已给单元平均的重构。在约束格式精确性,这几个数值导数的近似是有效的。注意到点值的重构过程和几个给定单元平均的数值导数是高分辨率、无震荡中心格式的核心。特别地,这些重构应该满足以下三个基本特性:
(1)单元平均的守恒:
(2)精确性:是-阶精确方法,这里是完全光滑的。
(3)对于不同重构的不同方法,具有无震荡性。
2.1.1 二阶重构
为了第3小节的二阶结果,(10)式中数值导数如下:
                 (12)
注:,
其中,代表限制器[6],它使得重构过程在满足极大值原理时无震荡。此外,对参数的限制,这里基于的重构是总变差减少格式
                        (13)
因此,中心格式是TVD格式。
2.1.2  三阶重构

多项式(14)满足性质(1)和(2),且保证了与的保形性,即在单元上无新的极值产生。然而,当

时界面跳转导致开始阶段有伪震荡。为了阻止界面跳转,我们寻求三阶限制器[2]。限制器计算基于单元数量

限制器通过下式给出

这时新的限制器可被写成限制器凸组合

此式仍满足性质(1)和(2)及上面提到的保形性,避免了伪界面极值,相当于。也就是说极值递减性质和保形性一起刻画了该格式的无震荡性,相当于性质(3)。
注:
  1. (15)式中和值不需要明确的计算,由于它们只进入限制器的计算,序列是单调的,它们分别在递增情况和递减情况下,与界面值一致。
  2. 序列的单调性也保证了差分和是阶的,即是三阶限制器,且。
  3. 对磁流体方程,值在极值单元()为引起所谓的削波现象需要设置为0,即避免一阶版本开始阶段的伪震荡[2]
  4. 多项式(18)和(11)一样用点值和数值导数给出:
一旦由分片多项式实现,它在区间上精确的积分去计算(9)式右端项的交错单元平均
    (22)
2.2 交错演化
让我们回到第二步中心格式的构造。这时沿着单元中心我们跟随重构点值的演化,通过
                (23)
控制。
如果是雅克比矩阵的特征值,在点的界面间断传播速度不快于。因此,只要CFL条件满足,中间单元值保持非间断性。因此,演化中(11)式右端通量积分被积函数是光滑的,其计算按所要求的精度选用适合的求积公式。特别地,二阶N-T格式[3]使用中点公式

三阶L-T格式[2]使用辛普森公式

这些求积公式须求中间点值,自然的方法应用泰勒展开式和差分方程(23), 。二阶情况的预估格式可写为
                          (24)
其是固定网格比,是近似数值导数—空间导数的近似。在我们计算精确性的约束下,差分方法对数值导数是可利用的,我们应保持其全部二阶精度。在第一种方法,我们用链式法则。用表示被计算的雅克比矩阵,表示重构数值导数,, 我们有
                 (25
另外,我们对网格函数的数值差分能直接使用Min-Mod限制器。
                   (26)
注:预估式(24)和(26)完全避免了雅克比矩阵的计算。
类似地,三阶预估式
                    (27)
(27)式右端的表达式计算可用链式法则

用(19)式三阶精确数值导数。(27)中三阶精确预估中值给出两种肯那个的计算方法如下
           (28)
         (29)
在这里是限制器的上标,表明其使用网格函数计算。中间值通过R-K方法近似。这个方法对中心格式的时间演化还提供了另外的收获:基于高阶R-K方法自然连续推广,我们在上下文提到的有效的演化过程明显地减少了每个网格点的计算数量。
结合预估点值,我们计算(8)式右端的积分项,用二阶中点公式近似
                       (30)
和用三阶辛普森公式计算
           (31)
用 表示这些通量积分的近似值,对二阶和三阶中心格式的校正步骤如下:
                   (32)
3   一维数值算例
方程(1)-(4)可写成守恒律形式
                                                        (33)
             (34)
其中是总压力(静压力和磁压力的总和)。
在这一小节,我们展示了一维磁流体方程(5)-(32),(33),(34)的数值模拟。使用基于Min-Mod限制器的二阶N-T中心格式的(24),(25),(30),(32)和基于无震荡限制器的三阶L-T中心格式的(27),(28),(31),(32)。我们使用空间离散一致网格,在CFL条件约束下采用动态时间步长,其中是的雅克比矩阵的特征值。
3.1  Brio-Wu 激波管问题
     首先考虑由两个初始平衡状态组成的激波管一维黎曼问题,其初始条件为:
                (35)
且,。计算区间为,选取800个网格点,图2-3为t=0.2时的密度、速度场分量、磁感应强度场及压力基于不同的数值导数的数值计算结果。
我们注意到这个问题的流体力学数据和空气动力学Sod激波管问题是一样的。然而磁流体力学波的多样性给给高阶方法带来了相当大的挑战,如我们结果包括向左移动的快速稀疏波,和由中间激波引起的慢速复合波,以及向右移动的一个接触间断、一个慢激波和一个快速稀疏波.注意如果不等于0,则这个问题的解将不唯一。
图2-3展示了二阶和三阶格式的不同结果,这个结果堪比二阶迎风格式[4],我们的数值结果表明了中心格式避免了除最大速度的估计外任何的雅克比矩阵特征信息,却仍然捕获了间断解的主要特征。观察图3的三阶,无震荡格式的拖尾边缘没有二阶格式明显,这是由高阶多项式重构造成的[5]

  1.                      (b)

(c)                   (d)                     (e)
      图2  选择800个网格点,t=0.2时,用二阶中心格式(24),(25),(31),(32),
      Brio-Wu 激波管问题的数值计算结果
 
  1.                      (b)   
  
(c)                   (d)                     (e)
    图3  选择800个网格点,t=0.2时,用三阶中心格式(27),(28),(31),(32),
 Brio-Wu 激波管问题的数值计算结果
3.2  Brio-Wu高马赫数激波管问题
Brio和Wu[4]将该激波管问题用来考核在解决高马赫数问题时数值格式的鲁棒性[]。其初始平衡状态为:

且。向右移动激波的马赫数为15.5。如果等离子体压力由静压力和磁压力的总和所代替,那此问题变成一个标准的磁流体黎曼问题。计算区间为[-1,1],选取200个网格点,CFL条件数取0.4,图4为t=0.012时的密度、速度场、磁感应强度场及压力的数值计算结果.结果包括一个向左移动的快速稀疏波,紧接着是一个正切方向的间断和向右移动的马赫数为15.5的一个快速激波.穿过正切方向的间断,密度、磁场和压力都会变化,但是流体速度和总压力是连续的。相比较二阶迎风格式[5],我们的数值结果表明了中心格式的鲁棒性。

(a)                         (b)

(c)                     (d)
图4  选择200个网格点,t=0.012时,用三阶中心格式(27),(28),(31),(32),
Brio-Wu 高马赫数激波管问题的数值计算结果
4   结论
本文结果表明了用中心格式精确计算理想磁流体方程间断解的有效性。我们的数值算例结果与Jiang和Wu[5]的结果极好的一致。且我们的方法不需要特征分解和维数分裂。在时间层的处理上引入了NCE的RK处理方法,从而使该格式达到很好的计算效果。两个典型的一维磁流体方程数值算例结果令人满意,具有高精度、无震荡、高分辨率的特性。
参考文献
[1] Jiang G S, Tadmor E. Nonoscillatory Central Schemes for Multidimensional Hyperbolic Conservation Laws ,SIAM [J]. Scientific Computation , 1998, 19(6): 1892-1917.
[2] Liu X D, Tadmor E. Third Order Nonoscillatory Central Scheme for Hyperbolic Conservation Laws,Numberische Mathematik, 1998,79:397-425.
[3] Nessayahu H, Tadmor E. Non-oscillatory Central Differencing for Hyperbolic Conservation Laws[J] . Journal of Computational Physics, 1990, 87(2): 408 -463 .
[4] Brio M, Wu C C. An Upwind Differencing Scheme for the Equations of Magnetohydrodynamics [J]. Journal of Computational Physics, 1988, 75(2): 400 -422.
[5] Jiang G S, Wu C C. A High -order WENO Finite Difference Scheme for the Equations of Ideal Magnetohydrodynamics [J]. Journal of Computational Physics, 1999, 150: 561-594.
[6] Van Leer B. Towards the Ultimate Conservative Difference Scheme V. A Second Order Sequel to Godunov’s method [J]. Journal of Computational Physics, 1979, 32: 101-136.
 

上一篇:以“立德树人”为目标的高校艺术教育普及初探
下一篇:高职高专机械设计课程设计课程教学改革的探索和实践