Implicit adaptive mesh refinement for 2D reduced resistive magnetohydrodynamics

https://doi.org/10.1016/j.jcp.2008.06.029Get rights and content

Abstract

An implicit structured adaptive mesh refinement (SAMR) solver for 2D reduced magnetohydrodynamics (MHD) is described. The time-implicit discretization is able to step over fast normal modes, while the spatial adaptivity resolves thin, dynamically evolving features. A Jacobian-free Newton–Krylov method is used for the nonlinear solver engine. For preconditioning, we have extended the optimal “physics-based” approach developed in [L. Chacón, D.A. Knoll, J.M. Finn, An implicit, nonlinear reduced resistive MHD solver, J. Comput. Phys. 178 (2002) 15–36] (which employed multigrid solver technology in the preconditioner for scalability) to SAMR grids using the well-known Fast Adaptive Composite grid (FAC) method [S. McCormick, Multilevel Adaptive Methods for Partial Differential Equations, SIAM, Philadelphia, PA, 1989]. A grid convergence study demonstrates that the solver performance is independent of the number of grid levels and only depends on the finest resolution considered, and that it scales well with grid refinement. The study of error generation and propagation in our SAMR implementation demonstrates that high-order (cubic) interpolation during regridding, combined with a robustly damping second-order temporal scheme such as BDF2, is required to minimize impact of grid errors at coarse-fine interfaces on the overall error of the computation for this MHD application. We also demonstrate that our implementation features the desired property that the overall numerical error is dependent only on the finest resolution level considered, and not on the base-grid resolution or on the number of refinement levels present during the simulation. We demonstrate the effectiveness of the tool on several challenging problems.

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 Bz is much larger than the magnitude of the in-plane magnetic field Bp. As a result, Bz constant and the velocity v is nearly incompressible (·v0), and the general MHD formalism reduces to [59], [20], [30]:t+v·-ημ0ΔΨ+E0=0,ρ(t+v·-νΔ)U=1μ0B·J,ΔΦ=U,where Φ is the velocity stream function (v=zˆ×Φ), U is the z-component of the fluid vorticity (U=zˆ·×v), Ψ 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 J-U formulation instead of the Ψ-U 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)

  • D.A. Knoll et al.

    On balanced approximations for the time integration of multiple time scale systems

    J. Comput. Phys.

    (2003)
  • D.A. Knoll et al.

    Jacobian-free Newton-Krylov methods: a survey of approaches and applications

    J. Comput. Phys.

    (2004)
  • S. Lankalapalli et al.

    An adaptive finite element method for magnetohydrodynamics

    J. Comput. Phys.

    (2007)
  • R. Samtaney et al.

    3D adaptive mesh refinement simulations of pellet injection in tokamaks

    Comput. Phys. Comm.

    (2004)
  • M. Shashkov et al.

    Support-operator finite-difference algorithms for general elliptic problems

    J. Comput. Phys.

    (1995)
  • H. Strauss et al.

    An adaptive finite element method for magnetohydrodynamics

    J. Comput. Phys.

    (1998)
  • H.R. Strauss et al.

    An adaptive finite element method for magnetohydrodynamics

    J. Comput. Phys.

    (1998)
  • G. Tóth et al.

    A parallel explicit/implicit time stepping scheme on block-adaptive grids

    J. Comput. Phys.

    (2006)
  • M.J. Aftosmis, M. Berger, Multilevel error estimation and adaptive h-refinement for cartesian meshes with embedded...
  • S. Balay, W.D. Gropp, L. Curfman-McInnes, B.F. Smith, PETSc Users Manual, Tech. Report ANL-95/11 – Revision 2.1.6,...
  • M. Berger, M. Aftosmis, J. Melton, Accuracy, adaptive methods and complex geometry,...
  • M. Berger et al.

    An algorithm for point clustering and grid generation

    IEEE Trans. Syst. Man Cybernet.

    (1991)
  • M. Berzins et al.

    On spatial adaptivity and interpolation when using the method of lines

    Appl. Numer. Math.

    (1997)
  • J.G. Blom et al.

    Algorithm 759: Vlugr3: a vectorizable adaptive-grid solver for PDEs in 3d; part ii. code description

    ACM Trans. Math. Software

    (1996)
  • A. Brandt, Multigrid techniques: 1984 guide, with applications to fluid dynamics, gmd studien nr. 85, GMD, GMD-AIW,...
  • P.N. Brown et al.

    Hybrid Krylov methods for nonlinear systems of equations

    SIAM J. Sci. Statist. Comput.

    (1990)
  • L. Chacón

    An optimal, parallel, fully implicit Newton–Krylov solver for three-dimensional visco-resistive magnetohydrodynamics

    Phys. Plasmas

    (2008)
  • J. Crank et al.

    A practical method for numerical evaluation of solutions of partial differential equations of the heat conduction type

    Proc. Cambridge Phil. Soc.

    (1947)
  • C.F. Curtiss et al.

    Integration of stiff equations

    Proc. Nat. Acad. Sci.

    (1952)
  • R.S. Dembo et al.

    Inexact Newton methods

    SIAM J. Numer. Anal.

    (1982)
  • J.F. Drake et al.

    Nonlinear reduced fluid equations for toroidal plasmas

    Phys. Fluids

    (1984)
  • Cited by (0)

    View full text