Darcy Problem

Summary

Solves a 2D mixed Darcy problem. This is a saddle point problem, discretized using Raviart-Thomas finite elements (velocity $var element = document.getElementById("moose-equation-6cee8556-e2fc-4a1c-a517-b5f24ef89482");katex.render("\\vec u", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) and piecewise discontinuous polynomials (pressure $var element = document.getElementById("moose-equation-d6ff9d13-4983-4199-8244-fc4fe7de62c1");katex.render("p", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ). This example demonstrates the use of transposition in the input file for mixed problems with different trial and test variables. This example is based on MFEM Example 5.

Description

This problem solves a mixed Darcy problem with strong form:

$var element = document.getElementById("moose-equation-1e38c2c4-5d15-48f4-9195-28ac7194d785");katex.render("\\begin{split} k \\vec u + \\vec \\nabla p &= \\vec f \\\\ -\\vec \\nabla \\cdot \\vec u &= g \\\\ p = h \\,\\,\\, &\\mathrm{on}\\,\\, \\partial\\Omega \\end{split}", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

where $var element = document.getElementById("moose-equation-3134297b-8ba1-40b1-92b8-057c6c9cb73b");katex.render("k = 1", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , $var element = document.getElementById("moose-equation-ae1ea327-169e-4f38-b899-736cd7b9de72");katex.render("\\vec f = \\vec 0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , $var element = document.getElementById("moose-equation-b758230e-442d-4b59-8a0d-6f625bf8592a");katex.render("g = 0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ and $var element = document.getElementById("moose-equation-1a690237-9696-4daa-adb5-ec9a1859addb");katex.render("h = e^x \\sin(y)", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .

In this example, we solve this using the weak form

$var element = document.getElementById("moose-equation-43ee70f3-c577-4ed9-9b3b-9b68240d2a37");katex.render("\\begin{split} (k \\vec u, \\vec v)_\\Omega - (p, \\vec \\nabla \\cdot \\vec v)_\\Omega + (h, \\vec v \\cdot \\hat n)_{\\partial\\Omega} &= 0 \\,\\,\\, \\forall \\vec v \\in V \\\\ - (\\vec \\nabla \\cdot \\vec u, w)_\\Omega &= 0 \\,\\,\\, \\forall w \\in W \\end{split}", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

where

$var element = document.getElementById("moose-equation-19e1546f-0928-4a81-aa7c-b02479c98276");katex.render("\\begin{split} \\vec u, \\vec v &\\in H(\\mathrm{div})(\\Omega) \\\\ p, w &\\in L^2(\\Omega) \\end{split}", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

Example File

[Mesh<<<{"href": "../Mesh/index.html"}>>>]
  type = MFEMMesh
  file = ../mesh/star.mesh
  uniform_refine<<<{"description": "Specify the level of uniform refinement applied to the initial mesh"}>>> = 4
[]

[Problem<<<{"href": "../Problem/index.html"}>>>]
  type = MFEMProblem
[]

[Functions<<<{"href": "../Functions/index.html"}>>>]
  [exact_velocity]
    type = ParsedVectorFunction<<<{"description": "Returns a vector function based on string descriptions for each component.", "href": "../../source/functions/MooseParsedVectorFunction.html"}>>>
    expression_x<<<{"description": "x-component of function."}>>> = '-exp(x) * sin(y)'
    expression_y<<<{"description": "y-component of function."}>>> = '-exp(x) * cos(y)'
  []
  [exact_pressure]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = 'exp(x) * sin(y)'
  []
  [exact_pressure_rhs]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = '-exp(x) * sin(y)'
  []
[]

[FESpaces<<<{"href": "../FESpaces/index.html"}>>>]
  [HDivFESpace]
    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."}>>> = RT
    fec_order<<<{"description": "Order of the FE shape function to use."}>>> = SECOND
  []
  [L2FESpace]
    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."}>>> = L2
    fec_order<<<{"description": "Order of the FE shape function to use."}>>> = SECOND
  []
[]

[Variables<<<{"href": "../Variables/index.html"}>>>]
  [velocity]
    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."}>>> = HDivFESpace
  []
  [pressure]
    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."}>>> = L2FESpace
  []
[]

[BCs<<<{"href": "../BCs/index.html"}>>>]
  [flux_boundaries]
    type = MFEMVectorFEBoundaryFluxIntegratedBC<<<{"description": "Adds the boundary integrator to an MFEM problem for the linear form $(f, \\vec v \\cdot \\hat n)_{\\partial\\Omega}$", "href": "../../source/mfem/bcs/MFEMVectorFEBoundaryFluxIntegratedBC.html"}>>>
    variable<<<{"description": "Variable on which to apply the boundary condition"}>>> = velocity
    coefficient<<<{"description": "The coefficient which will be used in the integrated BC. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = exact_pressure_rhs
  []
[]

[Kernels<<<{"href": "../Kernels/index.html"}>>>]
  [VelocityMass]
    type = MFEMVectorFEMassKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the bilinear form $(k \\vec u, \\vec v)_\\Omega$ arising from the weak form of the mass operator $k \\vec u$.", "href": "../../source/mfem/kernels/MFEMVectorFEMassKernel.html"}>>>
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = velocity
  []
  [PressureGrad]
    type = MFEMVectorFEDivergenceKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the mixed bilinear form $(k \\vec \\nabla \\cdot \\vec u, v)_\\Omega$ arising from the weak form of the divergence operator $k \\vec \\nabla \\cdot \\vec u$.", "href": "../../source/mfem/kernels/MFEMVectorFEDivergenceKernel.html"}>>>
    trial_variable<<<{"description": "The trial variable this kernel is acting on and which will be solved for. If empty (default), it will be the same as the test variable."}>>> = pressure
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = velocity
    coefficient<<<{"description": "Name of property k to use. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = -1
    transpose<<<{"description": "If true, adds the transpose of the integrator to the system instead."}>>> = true
  []
  [VelocityDiv]
    type = MFEMVectorFEDivergenceKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the mixed bilinear form $(k \\vec \\nabla \\cdot \\vec u, v)_\\Omega$ arising from the weak form of the divergence operator $k \\vec \\nabla \\cdot \\vec u$.", "href": "../../source/mfem/kernels/MFEMVectorFEDivergenceKernel.html"}>>>
    trial_variable<<<{"description": "The trial variable this kernel is acting on and which will be solved for. If empty (default), it will be the same as the test variable."}>>> = velocity
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = pressure
    coefficient<<<{"description": "Name of property k to use. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = -1
  []
[]

[Solver<<<{"href": "../Solver/index.html"}>>>]
  type = MFEMSuperLU
[]

[Executioner<<<{"href": "../Executioner/index.html"}>>>]
  type = MFEMSteady
  device = cpu
[]

[Postprocessors<<<{"href": "../Postprocessors/index.html"}>>>]
  [velocity_error]
    type = MFEMVectorL2Error<<<{"description": "Computes L2 error $\\left\\Vert \\vec u_{ex} - \\vec u_{h}\\right\\Vert_{\\rm L2}$ for vector gridfunctions.", "href": "../../source/mfem/postprocessors/MFEMVectorL2Error.html"}>>>
    variable<<<{"description": "Name of the vector variable of which to find the norm of the error."}>>> = velocity
    function<<<{"description": "The analytic solution to compare against. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a numeric vector value (enclosed in curly braces)."}>>> = exact_velocity
  []
  [pressure_error]
    type = MFEML2Error<<<{"description": "Computes L2 error $\\left\\Vert u_{ex} - u_{h}\\right\\Vert_{\\rm L2}$ for gridfunctions using H1 or L2 elements.", "href": "../../source/mfem/postprocessors/MFEML2Error.html"}>>>
    variable<<<{"description": "Name of the variable of which to find the norm of the error."}>>> = pressure
    function<<<{"description": "The analytic solution to compare against. A functor is any of the following: a variable, an MFEM material property, a function, a postprocessor or a number."}>>> = exact_pressure
  []
[]

[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/Darcy
    vtk_format<<<{"description": "Select VTK data format to use, choosing between BINARY, BINARY32, and ASCII."}>>> = ASCII
  []
  [DarcyErrorCSV]
    type = CSV<<<{"description": "Output for postprocessors, vector postprocessors, and scalar variables using comma seperated values (CSV).", "href": "../../source/outputs/CSV.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/Darcy
  []
[]

(test/tests/mfem/kernels/darcy.i)

[Mesh]
  type = MFEMMesh
  file = ../mesh/star.mesh
  uniform_refine = 4
[]

[Problem]
  type = MFEMProblem
[]

[Functions]
  [exact_velocity]
    type = ParsedVectorFunction
    expression_x = '-exp(x) * sin(y)'
    expression_y = '-exp(x) * cos(y)'
  []
  [exact_pressure]
    type = ParsedFunction
    expression = 'exp(x) * sin(y)'
  []
  [exact_pressure_rhs]
    type = ParsedFunction
    expression = '-exp(x) * sin(y)'
  []
[]

[FESpaces]
  [HDivFESpace]
    type = MFEMVectorFESpace
    fec_type = RT
    fec_order = SECOND
  []
  [L2FESpace]
    type = MFEMScalarFESpace
    fec_type = L2
    fec_order = SECOND
  []
[]

[Variables]
  [velocity]
    type = MFEMVariable
    fespace = HDivFESpace
  []
  [pressure]
    type = MFEMVariable
    fespace = L2FESpace
  []
[]

[BCs]
  [flux_boundaries]
    type = MFEMVectorFEBoundaryFluxIntegratedBC
    variable = velocity
    coefficient = exact_pressure_rhs
  []
[]

[Kernels]
  [VelocityMass]
    type = MFEMVectorFEMassKernel
    variable = velocity
  []
  [PressureGrad]
    type = MFEMVectorFEDivergenceKernel
    trial_variable = pressure
    variable = velocity
    coefficient = -1
    transpose = true
  []
  [VelocityDiv]
    type = MFEMVectorFEDivergenceKernel
    trial_variable = velocity
    variable = pressure
    coefficient = -1
  []
[]

[Solver]
  type = MFEMSuperLU
[]

[Executioner]
  type = MFEMSteady
  device = cpu
[]

[Postprocessors]
  [velocity_error]
    type = MFEMVectorL2Error
    variable = velocity
    function = exact_velocity
  []
  [pressure_error]
    type = MFEML2Error
    variable = pressure
    function = exact_pressure
  []
[]

[Outputs]
  [ParaViewDataCollection]
    type = MFEMParaViewDataCollection
    file_base = OutputData/Darcy
    vtk_format = ASCII
  []
  [DarcyErrorCSV]
    type = CSV
    file_base = OutputData/Darcy
  []
[]

Summary
Description
Example File