MHD现代算法介绍 张嫚 2011.03.24

主要内容 1 MHD算法简介 2 MHD算法存在的主要问题 3 MHD在太阳风模拟中应用

1 MHD算法简介

理想MHD方程

Classical schemes 迎风差分格式 Lax单步差分格式

Lax-Wendroff scheme MacCormack scheme Lax-Wendroff 和MacCormack都是二阶精度，耗散效应很弱，存在色散效应，激波附近会产生振荡 Others like Warming-Beam, Rusanov,etc.

Riemann Solvers 从流动的物理特性出发，采用间断分解概念求解流动方程组，在相邻的网格单元之间形成一个局部Riemann问题，转化为求解黎曼问题，在此基础上发展了一系列Riemann solvers

Godunov Riemann solvers 流动量在每个网格内以常数分布求解过程物理概念清晰，能较好反应流动特性，计算效果较好，但是该算法只有一阶精度。 发展高阶精度的Godunov 是必要的，MUSCUL和PPM差分算法是这方面的典型

要获得黎曼精确解是非常复杂的，并且很大程度上也没有必要，因此数值方法上很少涉及到精确解，广泛应用的是近似的黎曼解，下面介绍几种常用的近似黎曼解

Roe (linearized) type solvers 假定间断面两侧流动量常数分布，构造一个合理 的线性近似系数矩阵 来代替非线性Jiacobian 系数矩阵 把非线性矢量方程转化为线性矢量方 程，求出近似解

HLL (non linear) type solvers These non-linear solvers are based on approximating the wave structure of the full Riemann problem by a simplified wave structure

HLL Riemann solver .

这种方法可以保持正压，保持正压在low-beta是 非常重要的，但是具有较强的耗散性，在此基础上 发展的HLLC和HLLE克服这一缺点

HLLC Riemann solver

The HLLC flux based on the approximate Riemann solver is a modification to account for the shortcoming of the HLL flux offset the influence of intermediate waves

HLLD Riemann solver High resolution comparable to Roe scheme HLLE, HLLEM .etc

2 MHD算法存在的主要问题

保正性 MHD模拟的困难之一就是保持压力为正，由于MHD总能量包含三个分量：内能，磁能和动能，在low-beta区域，磁能相比于动能很大，因此，在计算总能量，动能和磁能时产生的离散误差很大程度上很导致负压的产生，出现非物理解。为了保持正压性，我们在局部区域用熵密度方程或者压强方程来代替能量方程求解MHD方程，为此我们需要用开关函数来判断可能出现负压的地方参考（Dinshaw S. Balsara* and Daniel Spicer. Jcp)给出开关函数

It consists of 3 sub-switches which includes or excludes some flow conditions. Switch 1 (SW1 for short) is designed to indicate whether the pressure is very low. So SW1 is on when and it is switched off in other situations,

Switch 2 the second switch SW2 is intended to exclude the regions of strong magnetosonic shocks. This is accomplished by ensuring that the local minima and maxima pressure are not too far apart. It is on when

The third switch SW3 is used to exclude the regions where mildly compressive motions(also strong compressive motions) may take place. It is on when

can be adjusted in different situations. Typically, we use and the pressure equation is used only when all of the three switches are on

保持磁场散度为零 MHD不同于Euler方程的就是存在磁场效应，如何处理磁场使之散度为零也是求解MHD的一个重要方面，有很多种方法可以用来处理磁场 Projection method， 8-wave formulation， Constrained Transport (CT) method， Central difference Diffusive control，GLM-MHD， Vector potential Hyperbolic divergence cleaning Splitting-based schemes

Projection method 在每个计算一个时间步之后对磁场进行更新

8-wave formulation 在方程右边加入Powell源项 The eight-wave scheme is computationally inexpensive, easy to add an existing code, and quite robust,

Constrained Transport (CT) method 采用交错网格方法，计算面上的磁场，之后将其投影到网格中心处 这种方法是严格守恒的，可以通过不同的构造电场磁场方法达到高阶精度，但是延伸到一般的网格上比较困难，不能直接应用于AMR

Central difference CD算法是在简化CT算法基础上发展的，在时间和空间上都是二阶精度

Diffusive control Modifying the magnetic equation by adding a source term proportional to the gradient of divergence so that the divergence satisfies an advection-diffusion equation This technique, referred as diffusive control of divergence , has the same advantages and disadvantages as eight-wave approach, but appears to keep the level of divergence lower than the eight-wave approach

GLM-MHD 采用9个方程来控制divergence hyperbolic correction

Hyperbolic divergence cleaning 显式的抛物型方法中 的取值受到稳定性条件的限制，结果不是很满意，为了克服这一缺点，将抛物型和双曲型方法结合起来 Vector potential use vector potential A as the primary variable instead of B，and the divergence free property is automatically ensured

3 MHD在太阳风模拟中应用 ——球坐标下太阳风模拟

目前太阳风模拟大多是在直角坐标系下进行，计算量比较大，直接应用于AMR比较困难，鉴于计算区域为球面网格划分，所以考虑在球坐标下实现太阳风模拟，达到稳定所需时间比较短，关于磁场处理部分，用一种可以直接应用于AMR的CT方法实现，不用求解Possion方程，大大减少了计算量

模拟区域 方程形式 考虑太阳自传，Powell 源项，加速加热，磁场分裂 直角坐标下 球坐标下

数值方法 有限体积法，将流场和磁场分开处理，磁场部分用CT方法，采用一种新 central scheme 的方法，不需要求解黎曼问题，也不需要求解特征值和特征向量，该方法适应于任何坐标，半离散形式，时间上用两步 Runge-kutta 法达到二阶精度

符号说明

流场重构 用线性重构达到二阶精度 流场离散形式

磁场重构 满足 二阶精度，其余方向上的重构类似

磁场离散形式 CT方法在面上对磁场进行离散，左边均为面上的平均值，之后有面上的值可以得到体积平均的磁场

计算区域考虑整个球面时，在极区离散方程会存在奇异性，需要特别处理

流场流通量计算 每个小单元被分割为27部分，以边和脚点为中心的矩形单元，以面为中心的非矩形单元，剩余的内部区域 分割依据：在 时间内，由边界上产生的波并没有影响到内部区域这样每个小单元的边界上都是光滑的， 边界上流通量的计算直接求积方法

每个小单元上面27个区域同时进行更新，之后把结果投影到单元中心上，得到新时刻体平均上的流体变量的值

电场计算 类似于流场，用周围四点的值投影到计算区域

边沿速度的确定

源项处理 不能简单用体心或者几何中心的值来近似体平，那样不能达到二阶精度，由于三角函数的存在，源项众不同部分取单元中不同位置的值，以此可以写出不同源项的表达式，这里不再一一列举

当源项中包含有Powell源项时，需要考虑divergence B在极区的奇异性，另外就是用CT处理磁场时，需要计算divergence B在面上的表示

初始输入 用太阳表面参数进行无量钢化，利用太阳风Parker解偶极场或者螺旋磁场作为初值 偶极场 螺旋磁场

边界条件 在半径方向上采用固定边界条件

Thankyou