Search

The MFEM (Modular Finite Element Methods) library is a high-performance C++ library for finite element discretizations. MFEM supports numerous types of finite element methods and is the discretization engine powering many computational physics and engineering applications across a number of domains. This paper describes some of the recent research and development in MFEM, focusing on performance portability across leadership-class supercomputing facilities, including exascale supercomputers, as well as new capabilities and functionality, enabling a wider range of applications. Much of this work was undertaken as part of the Department of Energy's Exascale Computing Project (ECP) in collaboration with the Center for Efficient Exascale Discretizations (CEED)., Comment: 19 pages, 19 figures

Computational analysis with the finite element method requires geometrically accurate meshes. It is well known that high-order meshes can accurately capture curved surfaces with fewer degrees of freedom in comparison to low-order meshes. Existing techniques for high-order mesh generation typically output meshes with same polynomial order for all elements. However, high order elements away from curvilinear boundaries or interfaces increase the computational cost of the simulation without increasing geometric accuracy. In prior work, we have presented one such approach for generating body-fitted uniform-order meshes that takes a given mesh and morphs it to align with the surface of interest prescribed as the zero isocontour of a level-set function. We extend this method to generate mixed-order meshes such that curved surfaces of the domain are discretized with high-order elements, while low-order elements are used elsewhere. Numerical experiments demonstrate the robustness of the approach and show that it can be used to generate mixed-order meshes that are much more efficient than high uniform-order meshes. The proposed approach is purely algebraic, and extends to different types of elements (quadrilaterals/triangles/tetrahedron/hexahedra) in two- and three-dimensions., Comment: 14 pages, 11 figures

We present a new method for two-material Lagrangian hydrodynamics, which combines the Shifted Interface Method (SIM) with a high-order Finite Element Method. Our approach relies on an exact (or sharp) material interface representation, that is, it uses the precise location of the material interface. The interface is represented by the zero level-set of a continuous high-order finite element function that moves with the material velocity. This strategy allows to evolve curved material interfaces inside curved elements. By reformulating the original interface problem over a surrogate (approximate) interface, located in proximity of the true interface, the SIM avoids cut cells and the associated problematic issues regarding implementation, numerical stability, and matrix conditioning. Accuracy is maintained by modifying the original interface conditions using Taylor expansions. We demonstrate the performance of the proposed algorithms on established numerical benchmarks in one, two and three dimensions.

We propose a new Nitsche-type approach for weak enforcement of normal velocity boundary conditions for a Lagrangian discretization of the compressible shock-hydrodynamics equations using high-order finite elements on curved boundaries. Specifically, the variational formulation is appropriately modified to enforce free-slip wall boundary conditions, without perturbing the structure of the function spaces used to represent the solution, with a considerable simplification with respect to traditional approaches. Total energy is conserved and the resulting mass matrices are constant in time. The robustness and accuracy of the proposed method are validated with an extensive set of tests involving nontrivial curved boundaries.

We propose a method that morphs high-orger meshes such that their boundaries and interfaces coincide/align with implicitly defined geometries. Our focus is particularly on the case when the target surface is prescribed as the zero isocontour of a smooth discrete function. Common examples of this scenario include using level set functions to represent material interfaces in multimaterial configurations, and evolving geometries in shape and topology optimization. The proposed method formulates the mesh optimization problem as a variational minimization of the sum of a chosen mesh-quality metric using the Target-Matrix Optimization Paradigm (TMOP) and a penalty term that weakly forces the selected faces of the mesh to align with the target surface. The distinct features of the method are use of a source mesh to represent the level set function with sufficient accuracy, and adaptive strategies for setting the penalization weight and selecting the faces of the mesh to be fit to the target isocontour of the level set field. We demonstrate that the proposed method is robust for generating boundary- and interface-fitted meshes for curvilinear domains using different element types in 2D and 3D., Comment: 30 pages, 16 figures

In this paper we present a new GPU-oriented mesh optimization method based on high-order finite elements. Our approach relies on node movement with fixed topology, through the Target-Matrix Optimization Paradigm (TMOP) and uses a global nonlinear solve over the whole computational mesh, i.e., all mesh nodes are moved together. A key property of the method is that the mesh optimization process is recast in terms of finite element operations, which allows us to utilize recent advances in the field of GPU-accelerated high-order finite element algorithms. For example, we reduce data motion by using tensor factorization and matrix-free methods, which have superior performance characteristics compared to traditional full finite element matrix assembly and offer advantages for GPU-based HPC hardware. We describe the major mathematical components of the method along with their efficient GPU-oriented implementation. In addition, we propose an easily reproducible mesh optimization test that can serve as a performance benchmark for the mesh optimization community., Comment: 28 pages, 10 figures, 3 tables

With the introduction of advanced heterogeneous computing architectures based on GPU accelerators, large-scale production codes have had to rethink their numerical algorithms and incorporate new programming models and memory management strategies in order to run efficiently on the latest supercomputers. In this work we discuss our co-design strategy to address these challenges and achieve performance and portability with MARBL, a next-generation multi-physics code in development at Lawrence Livermore National Laboratory. We present a two-fold approach, wherein new hardware is used to motivate both new algorithms and new abstraction layers, resulting in a single source application code suitable for a variety of platforms. Focusing on MARBL's ALE hydrodynamics package, we demonstrate scalability on different platforms and highlight that many of our innovations have been contributed back to open-source software libraries, such as MFEM (finite element algorithms) and RAJA (kernel abstractions)., Comment: Submitted to The International Journal of High Performance Computing Applications

In this paper we describe the research and development activities in the Center for Efficient Exascale Discretization within the US Exascale Computing Project, targeting state-of-the-art high-order finite-element algorithms for high-order applications on GPU-accelerated platforms. We discuss the GPU developments in several components of the CEED software stack, including the libCEED, MAGMA, MFEM, libParanumal, and Nek projects. We report performance and capability improvements in several CEED-enabled applications on both NVIDIA and AMD GPU systems.

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., Comment: 22 pages, 18 figures

We propose a new approach for controlling the characteristics of certain mesh faces during optimization of high-order curved meshes. The practical goals are tangential relaxation along initially aligned curved boundaries and internal surfaces, and mesh fitting to initially non-aligned surfaces. The distinct feature of the method is that it utilizes discrete finite element functions (for example level set functions) to define implicit surfaces, which are used to adapt the positions of certain mesh nodes. The algorithm does not require CAD descriptions or analytic parametrizations, and can be beneficial in computations with dynamically changing geometry, for example shape optimization and moving mesh multimaterial simulations. The main advantage of this approach is that it completely avoids geometric operations (e.g., surface projections), and all calculations can be performed through finite element operations., Comment: 11 pages, 13 figures

Large-scale finite element simulations of complex physical systems governed by partial differential equations (PDE) crucially depend on adaptive mesh refinement (AMR) to allocate computational budget to regions where higher resolution is required. Existing scalable AMR methods make heuristic refinement decisions based on instantaneous error estimation and thus do not aim for long-term optimality over an entire simulation. We propose a novel formulation of AMR as a Markov decision process and apply deep reinforcement learning (RL) to train refinement policies directly from simulation. AMR poses a new problem for RL as both the state dimension and available action set changes at every step, which we solve by proposing new policy architectures with differing generality and inductive bias. The model sizes of these policy architectures are independent of the mesh size and hence can be deployed on larger simulations than those used at train time. We demonstrate in comprehensive experiments on static function estimation and time-dependent equations that RL policies can be trained on problems without using ground truth solutions, are competitive with a widely-used error estimator, and generalize to larger, more complex, and unseen test problems., Comment: AISTATS 2023. 18 pages, 15 figures

We present an $hr$-adaptivity framework for optimization of high-order meshes. This work extends the $r$-adaptivity method for mesh optimization by Dobrev et al., where we utilized the Target-Matrix Optimization Paradigm (TMOP) to minimize a functional that depends on each element's current and target geometric parameters: element aspect-ratio, size, skew, and orientation. Since fixed mesh topology limits the ability to achieve the target size and aspect-ratio at each position, in this paper we augment the $r$-adaptivity framework with nonconforming adaptive mesh refinement to further reduce the error with respect to the target geometric parameters. The proposed formulation, referred to as $hr$-adaptivity, introduces TMOP-based quality estimators to satisfy the aspect-ratio-target via anisotropic refinements and size-target via isotropic refinements in each element of the mesh. The methodology presented is purely algebraic, extends to both simplices and hexahedra/quadrilaterals of any order, and supports nonconforming isotropic and anisotropic refinements in 2D and 3D. Using a problem with a known exact solution, we demonstrate the effectiveness of $hr$-adaptivity over both $r-$ and $h$-adaptivity in obtaining similar accuracy in the solution with significantly fewer degrees of freedom. We also present several examples that show that $hr$-adaptivity can help satisfy geometric targets even when $r$-adaptivity fails to do so, due to the topology of the initial mesh., Comment: 29 pages, 15 figures

In this paper we propose tools for high-order mesh optimization and demonstrate their benefits in the context of multi-material Arbitrary Lagrangian-Eulerian (ALE) compressible shock hydrodynamic applications. The mesh optimization process is driven by information provided by the simulation which uses the optimized mesh, such as shock positions, material regions, known error estimates, etc. These simulation features are usually represented discretely, for instance, as finite element functions on the Lagrangian mesh. The discrete nature of the input is critical for the practical applicability of the algorithms we propose and distinguishes this work from approaches that strictly require analytical information. Our methods are based on node movement through a high-order extension of the Target-Matrix Optimization Paradigm (TMOP) of Knupp. The proposed formulation is fully algebraic and relies only on local Jacobian matrices, so it is applicable to all types of mesh elements, in 2D and 3D, and any order of the mesh. We discuss the notions of constructing adaptive target matrices and obtaining their derivatives, reconstructing discrete data in intermediate meshes, node limiting that enables improvement of global mesh quality while preserving space-dependent local mesh features, and appropriate normalization of the objective function. The adaptivity methods are combined with automatic ALE triggers that can provide robustness of the mesh evolution and avoid excessive remap procedures. The benefits of the new high-order TMOP technology are illustrated on several simulations performed in the high-order ALE application BLAST., Comment: 26 pages, 10 figures. Submitted to Computers & Fluids

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

We present a new flux-fixup approach for arbitrarily high-order discontinuous Galerkin discretizations of the SN transport equation. This approach is sweep-compatible: as the transport sweep is performed, a local quadratic programming (QP) problem is solved in each spatial cell to ensure that the solution satisfies certain physical constraints, including local particle balance. The constraints can be chosen in two ways, leading to two variants of the method: QP Zero (QPZ) and QP Maximum Principle (QPMP). The coefficients of the solution are constrained to be nonnegative in QPZ, and they are constrained by an approximate discrete maximum principle in QPMP. There are two primary takeaways in this paper. First, it is shown that the QPMP method, when used with the positive Bernstein basis, eliminates negativities, preserves high-order accuracy for smooth problems, and significantly dampens unphysical oscillations in the solution. The second takeaway is that the Variable Eddington Factor (VEF) method can be used to accelerate the convergence of source iteration with fixup for problems with optically thick regions. Regardless of whether a fixup is applied, source iteration converges slowly when optically thick regions are present and acceleration is needed. When VEF is combined with fixed-up transport sweeps, the result is a scheme that produces a nonnegative solution, converges independently of the mean free path, and, in the case of the QPMP fixup, adheres to an approximate discrete maximum principle., Comment: Revised draft submitted to Journal of Computational Physics

We describe a framework for controlling and improving the quality of high-order finite element meshes based on extensions of the Target-Matrix Optimization Paradigm (TMOP) of Knupp. This approach allows high-order applications to have a very precise control over local mesh quality, while still improving the mesh globally. We address the adaption of various TMOP components to the settings of general isoparametric element mappings, including the mesh quality metric in 2D and 3D, the selection of sample points and the solution of the resulting mesh optimization problem. We also investigate additional practical concerns, such as tangential relaxation and restricting the deviation from the original mesh. The benefits of the new high-order TMOP algorithms are illustrated on a number of test problems and examples from a high-order arbitrary Eulerian-Lagrangian (ALE) application. Our implementation is freely available in an open-source library form.

We propose an unified algebraic approach for static condensation and hybridization, two popular techniques in finite element discretizations. The algebraic approach is supported by the construction of scalable solvers for problems involving H(div)-spaces discretized by conforming (Raviart-Thomas) elements of arbitrary order. We illustrate through numerical experiments the relative performance of the two (in some sense dual) techniques in comparison with a state-of-the-art parallel solver, ADS, available in the software hypre and MFEM. The superior performance of the hybridization technique is clearly demonstrated with increased benefit for higher order elements., Comment: 18 pages, 9 figures