Getting started with MFEM-MOOSE

Introduction

MFEM-MOOSE is a part of the MOOSE framework which utilizes MFEM as its main Finite Element Method backend. It expands upon MOOSE's base capabilities to support functionalities including but not limited to:

  • High-order H(Div) and H(Curl) elements

  • Various assembly levels (including partial assembly targeting matrix-free algorithms)

  • GPU offloading (CUDA or HIP)

  • Low-Order-Refined solvers

Installing MFEM-MOOSE

Installation instructions for MFEM-MOOSE can be found in this page.

Solving a problem with MFEM-MOOSE

Much of the syntax of the usual MOOSE input file is preserved when creating input files for MFEM-MOOSE. For each input file block, the user can browse the syntax page for classes prefixed with MFEM or, alternatively, browse the MFEM section of the source page. Example input files may be found in the kernel tests directory, test/tests/mfem/kernels/. Here, we lay out the step-by-step process of writing a MFEM-MOOSE input file to solve a simple diffusion problem. The full input file may be found here. We roughly split the input file into five parts: Problem, Geometry, Equation System, Solver and Executioner, and Output.

Problem

First of all, we must specify that the type of problem we wish to solve is an MFEMProblem, so that we may take advantage of MFEM capabilities.

[Problem<<<{"href": "../Problem/index.html"}>>>]
  type = MFEMProblem
[]
(test/tests/mfem/kernels/diffusion.i)

Geometry - Mesh and Finite Element Spaces

Given that we wish to utilize MFEM as the backend, the mesh we import into the problem must be of MFEMMesh type. Therefore, this must be specified in the parameter "type" within the Mesh block.

[Mesh<<<{"href": "../Mesh/index.html"}>>>]
  type = MFEMMesh
  file = ../mesh/mug.e
  dim = 3
[]
(test/tests/mfem/kernels/diffusion.i)

Then, we must set up the finite element spaces our problem will make use of. They can be of MFEMScalarFESpace type, which allows for continuous H1 or discontinuous L2 elements, or MFEMVectorFESpace type, thereby allowing tangentially-continuous ND or normally-continuous RT elements.

[FESpaces<<<{"href": "../FESpaces/index.html"}>>>]
  [H1FESpace]
    type = MFEMScalarFESpace<<<{"description": "Convenience class to construct scalar finite element spaces.", "href": "../../source/mfem/fespaces/MFEMScalarFESpace.html"}>>>
    fec_type<<<{"description": "Specifies the family of FE shape functions."}>>> = H1
    fec_order<<<{"description": "Order of the FE shape function to use."}>>> = FIRST
  []
  [HCurlFESpace]
    type = MFEMVectorFESpace<<<{"description": "Convenience class to construct vector finite element spaces, abstracting away some of the mathematical complexity of specifying the dimensions.", "href": "../../source/mfem/fespaces/MFEMVectorFESpace.html"}>>>
    fec_type<<<{"description": "Specifies the family of FE shape functions."}>>> = ND
    fec_order<<<{"description": "Order of the FE shape function to use."}>>> = FIRST
  []
[]
(test/tests/mfem/kernels/diffusion.i)

In this diffusion example, besides the usual H1 space required to solve the system, we also include an ND space so that we may output the gradient of the result as an AuxVariable later on. However, if you only wish to visualize the scalar result, it is not necessary to include the ND space.

Equation System - Variables, Kernels, and Boundary Conditions

Having created the necessary finite element spaces, we may now set up the variables to be solved for. They should be of type MFEMVariable. Each variable should also be associated with a relevant finite element space.

[Variables<<<{"href": "../Variables/index.html"}>>>]
  [concentration]
    type = MFEMVariable<<<{"description": "Class for adding MFEM variables to the problem (`mfem::ParGridFunction`s).", "href": "../../source/mfem/variables/MFEMVariable.html"}>>>
    fespace<<<{"description": "The finite element space this variable is defined on."}>>> = H1FESpace
  []
[]
(test/tests/mfem/kernels/diffusion.i)

Here, if we also wish to output the gradient of the result, we can add it as an AuxVariable and set up its corresponding AuxKernel.

[AuxVariables<<<{"href": "../AuxVariables/index.html"}>>>]
  [concentration_gradient]
    type = MFEMVariable<<<{"description": "Class for adding MFEM variables to the problem (`mfem::ParGridFunction`s).", "href": "../../source/mfem/variables/MFEMVariable.html"}>>>
    fespace<<<{"description": "The finite element space this variable is defined on."}>>> = HCurlFESpace
  []
[]

[AuxKernels<<<{"href": "../AuxKernels/index.html"}>>>]
  [grad]
    type = MFEMGradAux<<<{"description": "Calculates the gradient of an H1 conforming source variable and stores the result on an H(curl) conforming ND result auxvariable", "href": "../../source/mfem/auxkernels/MFEMGradAux.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = concentration_gradient
    source<<<{"description": "Scalar H1 MFEMVariable to take the gradient of."}>>> = concentration
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = TIMESTEP_END
  []
[]
(test/tests/mfem/kernels/diffusion.i)

Then, within the Kernels block, we specify the weak forms to be added to our equation system. Typically, one would pick the MFEM integrators they wish to implement by checking the linear form integrators page and the bilinear form integrators page. Note that not all linear and bilinear forms that are available in MFEM have been implemented in MFEM-MOOSE. Should you wish to implement an integrator that is not yet available, please submit a feature request.

If the integrator you wish to use is available, you can specify it in the type parameter simply by taking its MFEM integrator name, swapping the word Integrator for Kernel, and prepending MFEM to the beginning of its name. The table below shows a few examples of this naming convention:

MFEM nameMFEM-MOOSE name
DiffusionIntegratorMFEMDiffusionKernel
CurlCurlIntegratorMFEMCurlCurlKernel
VectorDomainLFIntegratorMFEMVectorDomainLFKernel

Putting this together, our Kernels block might look as follows:

[Kernels<<<{"href": "../Kernels/index.html"}>>>]
  [diff]
    type = MFEMDiffusionKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the bilinear form $(k\\vec\\nabla u, \\vec\\nabla v)_\\Omega$ arising from the weak form of the Laplacian operator $- \\vec\\nabla \\cdot \\left( k \\vec \\nabla u \\right)$.", "href": "../../source/mfem/kernels/MFEMDiffusionKernel.html"}>>>
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = concentration
    coefficient<<<{"description": "Name of property for diffusion coefficient k. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = diffusivity
  []
[]
(test/tests/mfem/kernels/diffusion.i)

Now we set up boundary conditions. Here, we choose scalar Dirichlet boundary conditions, which correspond to the MFEMScalarDirichletBC type.

[BCs<<<{"href": "../BCs/index.html"}>>>]
  [bottom]
    type = MFEMScalarDirichletBC<<<{"description": "Applies a Dirichlet condition to a scalar variable.", "href": "../../source/mfem/bcs/MFEMScalarDirichletBC.html"}>>>
    variable<<<{"description": "Variable on which to apply the boundary condition"}>>> = concentration
    boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies. Defaults to all boundaries."}>>> = 'bottom'
    coefficient<<<{"description": "The coefficient setting the values on the essential boundary. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = 1.0
  []
  [top]
    type = MFEMScalarDirichletBC<<<{"description": "Applies a Dirichlet condition to a scalar variable.", "href": "../../source/mfem/bcs/MFEMScalarDirichletBC.html"}>>>
    variable<<<{"description": "Variable on which to apply the boundary condition"}>>> = concentration
    boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies. Defaults to all boundaries."}>>> = 'top'
  []
[]
(test/tests/mfem/kernels/diffusion.i)

Solver and Executioner

With the equation system set up, we specify how it is to be solved. Firstly, we choose a preconditioner and solver. For problems with high polynomial order, setting "low_order_refined" to true may greatly increase performance, as explained here.

While in principle any solver may be used as the main solver or preconditioner, the main limitation to keep in mind is that Hypre solvers may only be preconditioned by other Hypre solvers. Furthermore, when a Hypre solver has its low_order_refined parameter set to true, it ceases to be considered a Hypre solver for preconditioning purposes.

[Preconditioner<<<{"href": "../Preconditioner/index.html"}>>>]
  [boomeramg]
    type = MFEMHypreBoomerAMG<<<{"description": "Hypre BoomerAMG solver and preconditioner for the iterative solution of MFEM equation systems.", "href": "../../source/mfem/solvers/MFEMHypreBoomerAMG.html"}>>>
  []
[]

[Solver<<<{"href": "../Solver/index.html"}>>>]
  type = MFEMHypreGMRES
  preconditioner = boomeramg
  l_tol = 1e-16
  l_max_its = 1000
[]
(test/tests/mfem/kernels/diffusion.i)

Static and time-dependent executioners may be implemented respectively with the MFEMSteady and MFEMTransient types. If MFEM-MOOSE has been built with GPU offloading capabilities, it is possible to set "device" to cuda or hip to make use of GPU acceleration. For GPU runs, it is advisable to choose an "assembly_level" other than legacy, otherwise the matrix assembly step will not be offloaded. The options for "assembly_level" are legacy, full, element, partial, and none (the latter is only available if MFEM-MOOSE has been built with libCEED support).

[Executioner<<<{"href": "../Executioner/index.html"}>>>]
  type = MFEMSteady
  device = cpu
[]
(test/tests/mfem/kernels/diffusion.i)

Output

Finally, we can choose a data collection type for our result outputs.

[Outputs<<<{"href": "../Outputs/index.html"}>>>]
  [ParaViewDataCollection]
    type = MFEMParaViewDataCollection<<<{"description": "Output for controlling export to an mfem::ParaViewDataCollection.", "href": "../../source/mfem/outputs/MFEMParaViewDataCollection.html"}>>>
    file_base<<<{"description": "The desired solution output name without an extension. If not provided, MOOSE sets it with Outputs/file_base when available. Otherwise, MOOSE uses input file name and this object name for a master input or uses master file_base, the subapp name and this object name for a subapp input to set it."}>>> = OutputData/Diffusion
    vtk_format<<<{"description": "Select VTK data format to use, choosing between BINARY, BINARY32, and ASCII."}>>> = ASCII
  []
[]
(test/tests/mfem/kernels/diffusion.i)