Search

Maintaining conservation laws in the fully discrete setting is critical for accurate long-time behavior of numerical simulations and requires accounting for discrete conservation properties in both space and time. This paper derives arbitrary order finite element exterior calculus spatial discretizations for the two-dimensional (2D) Navier-Stokes and drift-reduced magnetohydrodynamic equations that conserve both energy and enstrophy to machine precision when coupled with generally symplectic time-integration methods. Both continuous and discontinuous-Galerkin (DG) weak formulations can ensure conservation, but only generally symplectic time integration methods, such as the implicit midpoint method, permit exact conservation in time. Moreover, the symplectic implicit midpoint method yields an order of magnitude speedup over explicit schemes. The methods are implemented using the MFEM library and the solutions are verified for an extensive suite of 2D neutral fluid turbulence test problems. Numerical solutions are verified via comparison to a semi-analytic linear eigensolver as well as to the finite difference Global Drift Ballooning (GDB) code. However, it is found that turbulent simulations that conserve both energy and enstrophy tend to have too much power at high wavenumber and that this part of the spectrum should be controlled by reintroducing artificial dissipation. The DG formulation allows upwinding of the advection operator which dissipates enstrophy while still maintaining conservation of energy. Coupling upwinded DG with implicit symplectic integration appears to offer the best compromise of allowing mid-range wavenumbers to reach the appropriate amplitude while still controlling the high-wavenumber part of the spectrum.

Efficient exploitation of exascale architectures requires rethinking of the numerical algorithms used in many large-scale applications. These architectures favor algorithms that expose ultra fine-grain parallelism and maximize the ratio of floating point operations to energy intensive data movement. One of the few viable approaches to achieve high efficiency in the area of PDE discretizations on unstructured grids is to use matrix-free/partially-assembled high-order finite element methods, since these methods can increase the accuracy and/or lower the computational time due to reduced data motion. In this paper we provide an overview of the research and development activities in the Center for Efficient Exascale Discretizations (CEED), a co-design center in the Exascale Computing Project that is focused on the development of next-generation discretization software and algorithms to enable a wide range of finite element applications to run efficiently on future hardware. CEED is a research partnership involving more than 30 computational scientists from two US national labs and five universities, including members of the Nek5000, MFEM, MAGMA and PETSc projects. We discuss the CEED co-design activities based on targeted benchmarks, miniapps and discretization libraries and our work on performance optimizations for large-scale GPU architectures. We also provide a broad overview of research and development activities in areas such as unstructured adaptive mesh refinement algorithms, matrix-free linear solvers, high-order data visualization, and list examples of collaborations with several ECP and external applications., 22 pages, 18 figures

We study the convergence of the discontinuous Galerkin (DG) method applied to the advection-reaction equation on meshes with reentrant faces. On such meshes, the upwind numerical flux is not smooth, and so the numerical integration of the resulting face terms can only be expected to be first-order accurate. Despite this inexact integration, we prove that the DG method converges with order $\mathcal{O}(h^{p+1/2})$, which is the same rate as in the case of exact integration. Consequently, specialized quadrature rules that accurately integrate the non-smooth numerical fluxes are not required for high-order accuracy. These results are numerically corroborated on examples of linear advection and discrete ordinates transport equations., 9 pages, 2 figures

Fully implicit Runge-Kutta (IRK) methods have many desirable properties as time integration schemes in terms of accuracy and stability, but high-order IRK methods are not commonly used in practice with numerical PDEs due to the difficulty of solving the stage equations. This paper introduces a theoretical and algorithmic preconditioning framework for solving the systems of equations that arise from IRK methods applied to linear numerical PDEs (without algebraic constraints). This framework also naturally applies to discontinuous Galerkin discretizations in time. Under quite general assumptions on the spatial discretization that yield stable time integration, the preconditioned operator is proven to have condition number bounded by a small, order-one constant, independent of the spatial mesh and time-step size, and with only weak dependence on number of stages/polynomial order; for example, the preconditioned operator for 10th-order Gauss IRK has condition number less than two, independent of the spatial discretization and time step. The new method can be used with arbitrary existing preconditioners for backward Euler-type time stepping schemes, and is amenable to the use of three-term recursion Krylov methods when the underlying spatial discretization is symmetric. The new method is demonstrated to be effective on various high-order finite-difference and finite-element discretizations of linear parabolic and hyperbolic problems, demonstrating fast, scalable solution of up to 10th order accuracy. The new method consistently outperforms existing block preconditioning approaches, and in several cases, the new method can achieve 4th-order accuracy using Gauss integration with roughly half the number of preconditioner applications and wallclock time as required using standard diagonally implicit RK methods., 30 pages, accepted to SISC

MFEM is an open-source, lightweight, flexible and scalable C++ library for modular finite element methods that features arbitrary high-order finite element meshes and spaces, support for a wide variety of discretization approaches and emphasis on usability, portability, and high-performance computing efficiency. MFEM's goal is to provide application scientists with access to cutting-edge algorithms for high-order finite element meshing, discretizations and linear solvers, while enabling researchers to quickly and easily develop and test new algorithms in very general, fully unstructured, high-order, parallel and GPU-accelerated settings. In this paper we describe the underlying algorithms and finite element abstractions provided by MFEM, discuss the software implementation, and illustrate various applications of the library., Comment: 36 pages, 21 figures

In this paper, we develop subspace correction preconditioners for discontinuous Galerkin (DG) discretizations of elliptic problems with $hp$-refinement. These preconditioners are based on the decomposition of the DG finite element space into a conforming subspace, and a set of small nonconforming edge spaces. The conforming subspace is preconditioned using a matrix-free low-order refined technique, which in this work we extend to the $hp$-refinement context using a variational restriction approach. The condition number of the resulting linear system is independent of the granularity of the mesh $h$, and the degree of polynomial approximation $p$. The method is amenable to use with meshes of any degree of irregularity and arbitrary distribution of polynomial degrees. Numerical examples are shown on several test cases involving adaptively and randomly refined meshes, using both the symmetric interior penalty method and the second method of Bassi and Rebay (BR2)., Comment: 24 pages, 9 figures

In this paper, we develop high-order nodal discontinuous Galerkin (DG) methods for hyperbolic conservation laws that satisfy invariant domain preserving properties using a subcell flux corrections and convex limiting. These methods are based on a subcell flux corrected transport (FCT) methodology that involves blending a high-order target scheme with a robust, low-order invariant domain preserving method that is obtained using a graph viscosity technique. The new low-order discretizations are based on sparse stencils which do not increase with the polynomial degree of the high-order DG method. As a result, the accuracy of the low-order method does not degrade when used with high-order target methods. The method is applied to both scalar conservation laws, for which the discrete maximum principle is naturally enforced, and to systems of conservation laws such as the Euler equations, for which positivity of density and a minimum principle for specific entropy are enforced. Numerical results are presented on a number of benchmark test cases., Comment: 27 pages, 13 figures

We present a matrix-free flow solver for high-order finite element discretizations of the incompressible Navier-Stokes and Stokes equations with GPU acceleration. For high polynomial degrees, assembling the matrix for the linear systems resulting from the finite element discretization can be prohibitively expensive, both in terms of computational complexity and memory. For this reason, it is necessary to develop matrix-free operators and preconditioners, which can be used to efficiently solve these linear systems without access to the matrix entries themselves. The matrix-free operator evaluations utilize GPU-accelerated sum-factorization techniques to minimize memory movement and maximize throughput. The preconditioners developed in this work are based on a low-order refined methodology with parallel subspace corrections. The saddle-point Stokes system is solved using block-preconditioning techniques, which are robust in mesh size, polynomial degree, time step, and viscosity. For the incompressible Navier-Stokes equations, we make use of projection (fractional step) methods, which require Helmholtz and Poisson solves at each time step. The performance of our flow solvers is assessed on several benchmark problems in two and three spatial dimensions., 20 pages, 11 figures