Implicit adaptive mesh refinement for 2D reduced resistive magnetohydrodynamics
Introduction
The magnetohydrodynamics (MHD) model is useful for studying the macroscopic behavior of fully ionized gases (plasmas). Plasmas exhibit a wide range of complex behavior, and are intrinsically multiscale both temporally and spatially. While MHD provides a tractable model for the macroscopic description of plasmas, it still presents formidable challenges for the numerical modeler. In particular, MHD (even in its simplest form) supports multiple time scales (which manifest in the form of waves) and multiple length scales (which manifest in the form of microscopic layer formation, often with macroscopic relevance).
Algorithmically, the multiscale nature of MHD needs to be addressed separately in time and space. Spatially, the dynamic formation of thin layers requires grid adaptation that can respond dynamically. While there are many options available for dynamic grid adaptation depending on the spatial representation of choice (e.g., r-refinement, h-refinement, p-refinement), our focus here is on h-refinement in the finite-volume context via structured adaptive mesh refinement (SAMR) [8], [7]. A SAMR grid is organized as a hierarchy of nested refinement levels, with each level comprised of a union of rectangular patches. As the locally refined grid evolves to follow important features in the solution, these levels are created and destroyed as needed, and the solution is transferred from the old grid to the new grid to continue the simulation. While all approaches to adaptive mesh refinement have different strengths and weaknesses, we choose SAMR for a specific set of advantages it provides us. Firstly, SAMR employs uniform-grid stencils in large parts of the domain, thus providing higher accuracy than stencils with the same structure on non-uniform grids of comparable resolution [5]. Secondly, uniform-grid patches on all levels allow us to reuse single-grid smoothers and solvers developed for uniform-grid applications, such as geometric multigrid. This is important, since it is well known that exploiting geometric information when available results in better algorithmic performance [11], [62]. Finally, upon regridding, SAMR only requires solver-setup operations in newly created patches, thus minimizing setup overhead. In this application, we employ the SAMR package SAMRAI [32] to handle grid and data management operations. SAMRAI is an obvious choice, since it provides interfaces [46] to state-of-the-art nonlinear solver packages such as PETSc [3].
Temporally, our interest is on applications where fast time scales are parasitic to a slower dynamical time scale of interest (such applications arise, for instance, in fusion [55], [53] and space plasmas [28], [49]). Accordingly, it is of interest to step over such fast time scales in order to resolve those of dynamical interest, while preserving the temporal accuracy of the approach. Explicit time-integration methods are subject to stability constraints that arise from the fastest time scales, and are inappropriate for this purpose because they force the modeler to follow the fastest time scale supported. Fully-implicit time integration methods allow stepping over fast time scales, since time steps are generally constrained only by accuracy, not stability [37]. However, they require the solution of large-scale systems of nonlinear equations at each time step, and fast, robust solution methods are necessary for implicit methods to be practical. Fortunately, Newton–Krylov methods [12] have provided such robust solvers in a variety of contexts [38], including MHD [16], [15], [14], provided effective preconditioning is used. In Refs. [16], [15], [14], the key for algorithmic performance was the use of multigrid methods in the preconditioner stage.
Patch-based refinement in the context of MHD has been explored by many previous studies in the literature, both in the context of finite-volumes (see e.g. [4], [35], [55], [54], [61], [27]) and finite (and spectral) elements (e.g., [58], [40], [51]). In the finite-volume context, these studies have focused on various aspects of both the temporal and spatial discretization of the MHD equations on AMR grids. Spatially, authors have explored both staggered [4] and cell-centered [35], [55], [54] representations, with special emphasis on the preservation of conserved quantities and the solenoidal property of the magnetic field. An interesting study comparing the accuracy of finite-volumes/differences vs. spectral elements in an MHD-AMR context can be found in Ref. [45]. Temporally, most AMR implementations have relied on explicit methods, albeit with some flavor of time step subcycling for better performance (see e.g. [4], [35]). However, a number of authors have explored more advanced time stepping algorithms, such as partially implicit [54] (where hyperbolic terms are treated explicitly, and diffusive terms implicitly), implicit/explicit [61] (where some blocks are treated explicitly, while others are treated with a linearly implicit method), and fully implicit [27] (albeit using unscalable direct solvers).
The focus of this study is to merge the SAMR dynamic adaptive-grid approach with efficient, scalable, fully-implicit time integration, in the context of MHD. For simplicity, we focus our attention on the 2D reduced resistive MHD model [59], [20], [30], which is rigorously valid in the presence of a large guide magnetic field. The reduced resistive MHD model has the advantage of simplicity while maintaining a truly multiscale character, both temporally (because it supports the fast Alfvén wave) and spatially (because it develops thin layers). Furthermore, mature fully-implicit technology is available [16], which will be reused for this study.
The advantages of fully-implicit SAMR are obvious, as it enables dynamic refinement while decoupling the time integrator from the small explicit time step stability limits (which scale with the mesh size) that would arise in the patches of finest resolution. Key to the effectiveness and scalability of the proposed approach is to generalize the multigrid treatment proposed in Refs. [16], [15] to SAMR grids. This can be achieved with fast multilevel methods that exploit the structure of the mesh, such as the Fast Adaptive Composite grid (FAC) method [42], as has been already demonstrated in the context of 2D radiation diffusion [47].
Preliminary results on combining implicit time integration with SAMR for resistive MHD were first reported in [48]. Here we expand the study to include considerations of accuracy and provide details of our treatment of discretization at coarse–fine interfaces. Section 2 describes the mathematical model and its numerical discretization. Section 3 introduces the nonlinear solver of choice, Jacobian-free Newton–Krylov methods, and the preconditioning approach to make it efficient. The specifics of the coarse–fine interface treatment for this application are provided in Section 4. Finally, numerical results focusing on performance and accuracy aspects of the solver are presented in Section 5, and we conclude in Section 6.
Section snippets
Numerical model: Current–vorticity formulation of reduced MHD
In the 2D reduced MHD (RMHD) formalism, the magnetic field component in the ignorable direction is much larger than the magnitude of the in-plane magnetic field . As a result, constant and the velocity is nearly incompressible (), and the general MHD formalism reduces to [59], [20], [30]:where is the velocity stream function (), is the z-component of the fluid vorticity (), is the flux function (which
Nonlinear solution algorithm
Our general approach to the solution of (6) is via preconditioned Jacobian-free Newton–Krylov methods (JFNK). These methods have demonstrated their effectiveness in many similar applications [38], including 2D reduced resistive MHD [16] and Hall MHD [15], and 3D resistive MHD [14]. Our approach generalizes that of [16] in two fundamental ways: firstly, we have adapted the preconditioning strategy to deal with the formulation instead of the formulation; and secondly, we have generalized
Adaptive mesh refinement
The previous discussion has considered the generalization of the physics-based preconditioner proposed in Ref. [16] to the application at hand, without regard to the specifics of the spatial discretization employed. In what follows, we describe the AMR-specific details of our treatment of the MHD equations, with particular emphasis on (1) the spatial discretization at coarse–fine interfaces, (2) the generalization of multilevel solvers for SAMR grids, and (3) regridding and its impact on time
Numerical results
This section introduces several challenging test cases with the goal of demonstrating two main aspects of our implicit AMR implementation: algorithmic performance, and its accuracy properties. In regards to performance, we demonstrate that the convergence properties of the iterative approach are essentially independent of the number of grid levels present in the simulation (for equivalent fine-level resolution), and that it has good scaling properties with respect to the total number of
Conclusions
We have described the implementation of an algorithmically scalable, fully implicit, SAMR simulation tool for 2D reduced resistive MHD. The tool employs Jacobian-free Newton–Krylov methods as the solver engine. We use the reduced MHD solver developed in [16] as a starting point, albeit with several modifications. Following [58], we have reformulated the original problem (in terms of flux, vorticity, and streamfunction) using the current as a dynamic variable instead of the flux. This avoids
Acknowledgements
The authors acknowledge useful discussions with D.A. Knoll. The work was funded by the Los Alamos Directed Research and Development program at Los Alamos National Laboratory, operated for DOE under contract No. DE-AC52-06NA25396, and by the DOE Office of ASCR program in Applied Mathematical Sciences.
References (62)
- et al.
A conservative adaptive projection method for the variable density incompressible Navier–Stokes equations
J. Comput. Phys.
(1998) Divergence-free adaptive mesh refinement for magnetohydrodynamics
J. Comput. Phys.
(2001)- et al.
Local adaptive mesh refinement for shock hydrodynamics
J. Comput. Phys.
(1989) - et al.
Adaptive mesh refinement for hyperbolic partial differential equations
J. Comput. Phys.
(1984) A non-staggered, conservative, , finite-volume scheme for 3D implicit extended magnetohydrodynamics in curvilinear geometries
Comput. Phys. Comm.
(2004)- et al.
A 2D high- Hall MHD implicit nonlinear solver
J. Comput. Phys.
(2003) - et al.
An implicit, nonlinear reduced resistive MHD solver
J. Comput. Phys.
(2002) - et al.
Anisotropic grid adaptation for Navier–Stokes equations
J. Comput. Phys.
(2003) Heuristic stability theory for finite-difference equations
J. Comput. Phys
(1968)Adaptive mesh refinement for conservative systems: multi-dimensional efficiency evaluation
Comput. Phys. Comm.
(2003)