Divergence-Free WENO Reconstruction-Based Finite Volume Scheme for Solving Ideal MHD Equations on Triangular Meshes

Here B is the megnetic field, and E is the electric field defined by E = -u × B for ideal MHD. u is the velocity. J = ∇ × B is the current density. The induction equation ensures that the magnetic field remains divergence-free if it is divergence-free initially. In numerical simulations, maintaining discrete divergence-free is also important. Previous studies [11,6] have shown that a divergence error on the order of numerical truncation error introduced by the numerical scheme can lead to spurious solutions and the production of negative pressures.

catalogued in detail in subsection 3.1. Numerical tests are given in Section 4 to demonstrate the accuracy and non-oscillatory properties of the proposed scheme by computing smooth solution and shock wave related problems. We draw conclusions in Section 5.

Ideal MHD governing equations in the conservation form can be expressed as

where U = (ρ, ρu x , ρu y , ρu z , ε, B x , B y , B z ) T , (

and

Here p = p gas + B • B/2 is the total pressure, p gas is the gas pressure that satisfies the following equation of state

with u = (u x , u y , u z ) T and B = (B x , B y , B z ) T . For a 2D ideal MHD problem, we have

We employ the CT approach and the Godunov type finite volume scheme to solve Eq. (2.1). To this end, the physical domain Ω is partitioned into a collection of N triangular cells K i so that Ω = N i=1 K i and we define

We also collect cell edges L j to form

where N E is the total number of edges in the partition. For every cell edge L j , we uniquely identify an edge unit normal n j and tangent ζ j . Here ζ j is obtained by rotating n j 90 degrees in the counterclockwise direction. For simplicity, we assume that there are no hanging nodes in the partition T h . Let the edges of cell K i be denoted as ∂K i, , = 1, 2, 3. For convenience in discussion, we define a mapping between the local cell edge index of cell K i and the global edge index j such that = i (j) and j = -1 i ( ) .

(2.7)

We also define the mesh parameter h to be

We place the magnetic field variables B x and B y at the cell edges to maintain the global divergence-free evolution of the magnetic field; the z-component of the electric field E z at the cell vertices; and the conservative variables ρ, ρu and ε and B z on the cells. B x and B y are always initialized to be divergence-free. The Godunov type finite volume scheme is utilized to evolve ρ, ρu, ε and B z on the cells and the normal component of the magnetic field within the xy-plane on the cell edges. To evaluate E z at cell vertices, the flux-interpolated approach introduced by Balsara and Spicer [3] is further developed here.

For convenience in discussion, we introduce notations

where

where

(2.12)

Thus solving Eq. (2.1) is equivalent to solving equations (2.9) and (2.11) together.

Taking the cell K i , i = 1, • • • , N , in partition (2.5) as a discrete control volume, the semidiscrete finite volume method for solving Eq. (2.9) is formulated by integrating (2.9) over the cell K i :

where

, and n i is the outward unit normal of the boundary of the cell K i . |K i | is a shorthand notation for the area of K i .

To solve Eq. (2.13) numerically, we evaluate the flux integral by Gaussian quadrature rule with the exact value of (F H , G H )•n i being replaced by the Lax-Friedrichs flux F * (x, y, t) given by

Here α is taken as an upper bound for the eigenvalues of the Jacobian in the n i direction; U - (or U H,-) and U + (or U H,+ ) are the numerical values of U (or U H ) inside the triangle and outside the triangle at the Gaussian point. To this end, we obtain the following semi-discrete finite volume scheme for solving Eq. (2.9)

where U H h,k,i (t) is the approximate cell average of the k th component of U H on the cell K i .

The 2D constrained transport scheme developed in the present paper is based upon cell edge-length-averaged magnetic field located at the edges of grid cells. On every cell edge L j ∈ E h , we solve Eq. (2.11) to evolve the normal component of B xy with respect to the defined cell edge unit normal n j . Denote the normal and tangential contribution of B xy in directions given by n j and ζ j to be B n and B ζ respectively. We rewrite Eq. (2.11) by B n and B ζ to obtain

Here u n and u ζ are the components of velocity u in the n j and ζ j directions respectively. Let B n,j be the edge-length-averaged B n on the edge L j defined by

where |L j | is a shorthand notation for the length of the edge L j . Integrating Eq. (2.16) along the cell edge L j , the semi-discrete finite volume scheme to evolve B n,j numerically on L j can