A high-order unstaggered constrained transport method for the 3D ideal magnetohydrodynamic equations based on the method of lines Bertram Taetz (joint work with C. Helzel & J.A. Rossmanith) Department of Numerical Mathematics Ruhr-University Bochum Research unit 1048: “Instabilities, turbulence and transport in cosmic magnetic fields” June 28, 2012 1 / 24
Contents Short introduction to the ideal Magnetohydrodynamic (MHD) 1 equations and the numerical challenge related to ∇ · B = 0 A constrained transport method in 3D 2 Numerical Tests 3 2 / 24
Motivation Ball mapping Observation of taken Coronal mass ejection (CME) from [Calhoun et al., 2008] taken from [SOHO,2002] 3 / 24
The ideal MHD equations in conservation form [Brackbill & Barnes, 1980] ρ ρ u � 2 | B | 2 � p + 1 ρ u ρ uu + Id − BB + ∇ · = 0 � 2 | B | 2 � E + p + 1 E u − B ( u · B ) uB − Bu B t ∇ · B = 0 the thermal pressure is related via the ideal gas law: p = ( γ − 1)( E − 1 2 || B || 2 − 1 2 ρ || u || 2 ) 4 / 24
Problems due to insufficient control of ∇ · B Figure: taken from [Rossmanith, 2006]. 5 / 24
Approaches to control ∇ · B Numerical methods for ideal MHD must in general satisfy (or at least control) some discrete version of the divergence free condition on the magnetic field. Some known methods to control ∇ · B on the discrete level: projection methods, e.g. [T ´ oth, 2000] 8-wave-formulation, [Powell, 1994] divergence cleaning, [Dedner et al.,2002] flux-distribution methods, [Torrilhon,2003],[Mishra & Tadmor, 2012] contrained transport (CT) methods, e.g. [Evans and Hawley, 1988, Rossmanith, 2006] 6 / 24
The idea of CT in 3D Consider the induction equation B t + ∇ × ( B × u ) = 0 and assume that u is a given vector valued function. Set B = ∇ × A to obtain ∇ × ( A t + ( ∇ × A ) × u ) = 0 ⇒ A t + ( ∇ × A ) × u = −∇ ψ ψ is an arbitrary scalar function taken to be ψ = 0 (Weyl gauge) in the following. Different choices of ψ represent different gauge conditions . See e.g. [C.Helzel, J.A.Rossmanith & B.Taetz, 2011] for discussions on different choices. 7 / 24
The evolution of the magnetic potential A t + N 1 ( u ) A x + N 2 ( u ) A y + N 3 ( u ) A z = 0 , with − u 2 − u 3 u 2 u 3 0 0 0 0 0 , N 2 = , N 3 = . u 1 − u 1 − u 3 u 3 N 1 = 0 0 0 0 0 u 1 u 2 − u 1 − u 2 0 0 0 0 0 This system is weakly hyperbolic , which means that we do not have a full set of linearly independent eigenvectors in all directions. 8 / 24
The evolution of the magnetic potential A t + N 1 ( u ) A x + N 2 ( u ) A y + N 3 ( u ) A z = 0 , with − u 2 − u 3 u 2 u 3 0 0 0 0 0 , N 2 = , N 3 = . u 1 − u 1 − u 3 u 3 N 1 = 0 0 0 0 0 u 1 u 2 − u 1 − u 2 0 0 0 0 0 This system is weakly hyperbolic , which means that we do not have a full set of linearly independent eigenvectors in all directions. 8 / 24
Weak hyperbolicity of the system matrix The system matrix in an arbitrary direction n ∈ S 2 is n 2 u 2 + n 3 u 3 − n 1 u 2 − n 1 u 3 n 1 u 1 + n 3 u 3 n 1 N 1 + n 2 N 2 + n 3 N 3 = − n 2 u 1 − n 2 u 3 n 1 u 1 + n 2 u 2 − n 3 u 1 − n 3 u 2 the eigenvalues are � � λ = 0 , n · u , n · u ; the eigenvectors read n 2 u 3 − n 3 u 2 u 1 ( u · n ) − n 1 � u � 2 n 1 n 3 u 1 − n 1 u 3 u 2 ( u · n ) − n 2 � u � 2 n 2 ; R = n 1 u 2 − n 2 u 1 u 3 ( u · n ) − n 3 � u � 2 n 3 the determinant of R can be written as det ( R ) = −� u � 3 cos( α ) sin 2 ( α ) , α is the angle between n and u . → Thus we do not have a full set of linearly independent eigenvectors in all directions n . 9 / 24
Weak hyperbolicity of the system matrix The system matrix in an arbitrary direction n ∈ S 2 is n 2 u 2 + n 3 u 3 − n 1 u 2 − n 1 u 3 n 1 u 1 + n 3 u 3 n 1 N 1 + n 2 N 2 + n 3 N 3 = − n 2 u 1 − n 2 u 3 n 1 u 1 + n 2 u 2 − n 3 u 1 − n 3 u 2 the eigenvalues are � � λ = 0 , n · u , n · u ; the eigenvectors read n 2 u 3 − n 3 u 2 u 1 ( u · n ) − n 1 � u � 2 n 1 n 3 u 1 − n 1 u 3 u 2 ( u · n ) − n 2 � u � 2 n 2 ; R = n 1 u 2 − n 2 u 1 u 3 ( u · n ) − n 3 � u � 2 n 3 the determinant of R can be written as det ( R ) = −� u � 3 cos( α ) sin 2 ( α ) , α is the angle between n and u . → Thus we do not have a full set of linearly independent eigenvectors in all directions n . 9 / 24
Weak hyperbolicity of the system matrix The system matrix in an arbitrary direction n ∈ S 2 is n 2 u 2 + n 3 u 3 − n 1 u 2 − n 1 u 3 n 1 u 1 + n 3 u 3 n 1 N 1 + n 2 N 2 + n 3 N 3 = − n 2 u 1 − n 2 u 3 n 1 u 1 + n 2 u 2 − n 3 u 1 − n 3 u 2 the eigenvalues are � � λ = 0 , n · u , n · u ; the eigenvectors read n 2 u 3 − n 3 u 2 u 1 ( u · n ) − n 1 � u � 2 n 1 n 3 u 1 − n 1 u 3 u 2 ( u · n ) − n 2 � u � 2 n 2 ; R = n 1 u 2 − n 2 u 1 u 3 ( u · n ) − n 3 � u � 2 n 3 the determinant of R can be written as det ( R ) = −� u � 3 cos( α ) sin 2 ( α ) , α is the angle between n and u . → Thus we do not have a full set of linearly independent eigenvectors in all directions n . 9 / 24
Weak hyperbolicity of the system matrix The system matrix in an arbitrary direction n ∈ S 2 is n 2 u 2 + n 3 u 3 − n 1 u 2 − n 1 u 3 n 1 u 1 + n 3 u 3 n 1 N 1 + n 2 N 2 + n 3 N 3 = − n 2 u 1 − n 2 u 3 n 1 u 1 + n 2 u 2 − n 3 u 1 − n 3 u 2 the eigenvalues are � � λ = 0 , n · u , n · u ; the eigenvectors read n 2 u 3 − n 3 u 2 u 1 ( u · n ) − n 1 � u � 2 n 1 n 3 u 1 − n 1 u 3 u 2 ( u · n ) − n 2 � u � 2 n 2 ; R = n 1 u 2 − n 2 u 1 u 3 ( u · n ) − n 3 � u � 2 n 3 the determinant of R can be written as det ( R ) = −� u � 3 cos( α ) sin 2 ( α ) , α is the angle between n and u . → Thus we do not have a full set of linearly independent eigenvectors in all directions n . 9 / 24
An operator splitting approach (dimensional splitting) A 1 t − u 2 A 2 x − u 3 A 3 Sub-problem 1: x = 0 , A 2 t + u 1 A 2 x = 0 , A 3 t + u 1 A 3 x = 0 , A 1 t + u 2 A 1 Sub-problem 2: y = 0 , A 2 t − u 1 A 1 y − u 3 A 3 y = 0 , A 3 t + u 2 A 3 y = 0 , A 1 t + u 3 A 1 Sub-problem 3: z = 0 , A 2 t + u 3 A 2 z = 0 , A 3 t − u 1 A 1 z − u 2 A 2 z = 0 . developed for Cartesian grids to 2 nd order of accuracy using Strang splitting. 10 / 24
A high-order unsplit spatial discretization for weakly hyperbolic systems Consider the integral form of a weakly hyperbolic system: q t + N ( x ) q x = 0 . We write the semi-discrete form for the cell-averages Q i ( t ) as: ∂ t Q i ( t ) = − 1 ∆ x ( A + ∆ Q i − 1 / 2 + A − ∆ Q i +1 / 2 + A ∆ Q i ) with � x i + 1 ˜ 2 N i ( x ) ˜ A ∆ Q i := q i,x dx, x i − 1 2 � x i − 1 2 + ǫ � � A − ∆ Q i − 1 2 + A + ∆ Q i − 1 ˜ Q ε = lim N ( x ) 2 ( t, x ) x dx. i − 1 ε → 0 2 x i − 1 2 − ε 11 / 24
A high-order unsplit spatial discretization for weakly hyperbolic systems Consider the integral form of a weakly hyperbolic system: q t + N ( x ) q x = 0 . We write the semi-discrete form for the cell-averages Q i ( t ) as: ∂ t Q i ( t ) = − 1 ∆ x ( A + ∆ Q i − 1 / 2 + A − ∆ Q i +1 / 2 + A ∆ Q i ) with � x i + 1 ˜ 2 N i ( x ) ˜ A ∆ Q i := q i,x dx, x i − 1 2 � x i − 1 2 + ǫ � � A − ∆ Q i − 1 2 + A + ∆ Q i − 1 ˜ Q ε = lim N ( x ) 2 ( t, x ) x dx. i − 1 ε → 0 2 x i − 1 2 − ε 11 / 24
A high-order unsplit spatial discretization for weakly hyperbolic systems Consider the integral form of a weakly hyperbolic system: q t + N ( x ) q x = 0 . We write the semi-discrete form for the cell-averages Q i ( t ) as: ∂ t Q i ( t ) = − 1 ∆ x ( A + ∆ Q i − 1 / 2 + A − ∆ Q i +1 / 2 + A ∆ Q i ) with � x i + 1 ˜ 2 N i ( x ) ˜ A ∆ Q i := q i,x dx, x i − 1 2 � x i − 1 2 + ǫ � � A − ∆ Q i − 1 2 + A + ∆ Q i − 1 ˜ Q ε = lim N ( x ) 2 ( t, x ) x dx. i − 1 ε → 0 2 x i − 1 2 − ε 11 / 24
Definition of fluctuations Use regularization Q ε 2 ( t, x ) with a straight-line path i − 1 � � 2 = Q − Q + 2 − Q − Ψ i − 1 2 + l , 0 ≤ l ≤ 1 i − 1 i − 1 i − 1 2 to derive � A − ∆ Q i − 1 2 + A + ∆ Q i − 1 ( Q + 2 − Q − = N 2 ) . � Ψ i − 1 i − 1 i − 1 2 2 � �� � i − 1 / 2 + N + 1 2 ( N − i − 1 / 2 ) Use generalized Rusanov-flux similar to [Castro et al., 2010]: � �� � 1 � A − ∆ Q i − 1 / 2 Q + 2 − Q − := N Ψ i − 1 / 2 − α i − 1 / 2 , Id � i − 1 i − 1 2 2 � �� � ≥| λ p i − 1 / 2 | ∀ p � �� � 1 � A + ∆ Q i − 1 / 2 Q + 2 − Q − := N Ψ i − 1 / 2 + α i − 1 / 2 Id . � i − 1 i − 1 2 2 12 / 24
Recommend
More recommend
