Transient Heat Transfer (Mixed Form)

Summary

Solves a transient heat conduction problem using a mixed weak form, representing temperature on piecewise constant $var element = document.getElementById("moose-equation-e056002b-146d-40c4-aa38-9dcb4468de3a");katex.render("L^2", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ conforming finite elements, and heat fluxes on $var element = document.getElementById("moose-equation-58114d64-7f1a-45c7-a0f2-51bfb852e9e2");katex.render("H(\\mathrm{div})", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ conforming Raviart-Thomas elements.

Description

This problem solves the transient heat equation with strong form:

$var element = document.getElementById("moose-equation-b5077298-7c2f-44b1-8e81-d5b248b642b2");katex.render("\\begin{split} \\rho c_p \\frac{\\partial T}{\\partial t} - \\vec \\nabla \\cdot \\left(k\\vec\\nabla T\\right) = 0 \\,\\,\\,&\\mathrm{on}\\,\\, \\Omega \\\\ k \\vec \\nabla T \\cdot \\hat n = 0 \\,\\,\\, &\\mathrm{on}\\,\\, \\Gamma_\\mathrm{h} \\\\ T = T_0 \\,\\,\\, &\\mathrm{on}\\,\\, \\Gamma_\\mathrm{T} \\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-3cadd504-a8ce-44fb-af8a-4e30fd9430da");katex.render("\\rho c_p = 1.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-9c0954ee-d84c-489d-9696-c639e25b5f8d");katex.render("k^{-1} = 1.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}"}});$ , subject to the initial condition $var element = document.getElementById("moose-equation-0503364e-9477-47d7-b920-1d321b18bd64");katex.render("T(t=0)=0.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}"}});$ .

In this example, we solve this by introducing the heat flux $var element = document.getElementById("moose-equation-016781c1-22cd-4f0c-8a5f-69673564da39");katex.render("\\vec h = -k \\vec \\nabla T", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ to generate the mixed weak form of the transient heat equation

$var element = document.getElementById("moose-equation-1f5a7147-69fa-46fd-81b8-5a101f4796ca");katex.render("\\begin{split} \\left(\\rho c_p \\frac{\\partial T}{\\partial t}, T'\\right)_\\Omega + \\left(\\vec \\nabla \\cdot \\vec h, T'\\right)_\\Omega = 0 \\,\\,\\, \\forall T' \\in V_T \\\\ \\left(k^{-1} \\vec h, \\vec h'\\right)_\\Omega - (T, \\vec \\nabla \\cdot \\vec h')_\\Omega = (T , \\vec h' \\cdot \\hat n)_{\\partial \\Omega} \\,\\,\\, \\forall \\vec h' \\in V_h \\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-dfeba0a0-2d97-45c5-8ef0-29041739a9bb");katex.render("\\begin{split} T, T' &\\in L_2(\\Omega)\\\\ \\vec h, \\vec h' &\\in H(\\mathrm{div})(\\Omega) : \\vec h \\cdot \\hat n = 0 \\,\\,\\, \\mathrm{on}\\,\\, \\Gamma_\\mathrm{h} \\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}"}});$

and the polynomial order of the FEs used to represent $var element = document.getElementById("moose-equation-ef3c1dd4-69d1-4f48-bdeb-b8bee7e3017d");katex.render("\\vec h", 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-3be14c10-3097-4c3a-b0bc-852775b4a4d8");katex.render("\\vec h'", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is one order higher than $var element = document.getElementById("moose-equation-f13abcca-bb03-4c7a-8102-361d382c0bbc");katex.render("T", 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-cee239b9-baec-4471-a5b4-e53c52929dc9");katex.render("T'", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .

Notably, in contrast to the primal weak form for the transient heat equation solved here, the heat flux on $var element = document.getElementById("moose-equation-1f9758c3-b798-45e7-a2c2-6b31e5aa9faa");katex.render("\\Gamma_\\mathrm{h}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is strongly imposed via a Dirichlet condition, and temperatures on boundaries $var element = document.getElementById("moose-equation-51c4498e-3956-43ea-a59f-8fe3230f9005");katex.render("\\Gamma_\\mathrm{T}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ are imposed weakly via boundary integrals.

Example File

[Mesh<<<{"href": "../Mesh/index.html"}>>>]
  type = MFEMMesh
  file = ../mesh/square.e
[]

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

[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."}>>> = FIRST
  []
  [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."}>>> = CONSTANT
  []
[]

[Variables<<<{"href": "../Variables/index.html"}>>>]
  [time_integrated_heat_flux]
    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
    time_derivative<<<{"description": "Optional name to assign to the time derivative of the variable in transient problems."}>>> = heat_flux
  []
  [temperature]
    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
  []
[]

[Kernels<<<{"href": "../Kernels/index.html"}>>>]
  [dT_dt,T']
    type = MFEMTimeDerivativeMassKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the bilinear form $(k \\dot{u}, v)_\\Omega$ arising from the weak form of the operator $k \\dot{u}$.", "href": "../../source/mfem/kernels/MFEMTimeDerivativeMassKernel.html"}>>>
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = temperature
  []
  [divh,T']
    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."}>>> = heat_flux
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = temperature
  []
  [h,h']
    type = MFEMTimeDerivativeVectorFEMassKernel<<<{"description": "Adds the domain integrator to an MFEM problem for the bilinear form $(k \\dot{\\vec u}, \\vec v)_\\Omega$ arising from the weak form of the operator $k \\dot{\\vec u}$.", "href": "../../source/mfem/kernels/MFEMTimeDerivativeVectorFEMassKernel.html"}>>>
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = time_integrated_heat_flux
  []
  [-T,div.h']
    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."}>>> = temperature
    variable<<<{"description": "Variable labelling the weak form this kernel is added to"}>>> = time_integrated_heat_flux
    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.0
    transpose<<<{"description": "If true, adds the transpose of the integrator to the system instead."}>>> = true
  []
[]

[BCs<<<{"href": "../BCs/index.html"}>>>]
  [right]
    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"}>>> = time_integrated_heat_flux
    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."}>>> = 0.0
    boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies. Defaults to all boundaries."}>>> = 2
  []
  [left]
    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"}>>> = time_integrated_heat_flux
    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."}>>> = -1.0
    boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies. Defaults to all boundaries."}>>> = 4
  []
  [topbottom]
    type = MFEMVectorDirichletBC<<<{"description": "Applies a Dirichlet condition to all components of a vector variable.", "href": "../../source/mfem/bcs/MFEMVectorDirichletBC.html"}>>>
    variable<<<{"description": "Variable on which to apply the boundary condition"}>>> = heat_flux
    vector_coefficient<<<{"description": "Vector coefficient specifying the values the variable takes on the boundary. 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)."}>>> = '0.0 0.0'
    boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies. Defaults to all boundaries."}>>> = '1 3'
  []
[]

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

[Executioner<<<{"href": "../Executioner/index.html"}>>>]
  type = MFEMTransient
  device = cpu
  assembly_level = legacy
  dt = 0.03
  start_time = 0.0
  end_time = 0.09
[]

[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/MixedHeatTransfer
    vtk_format<<<{"description": "Select VTK data format to use, choosing between BINARY, BINARY32, and ASCII."}>>> = ASCII
  []
[]

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

[Mesh]
  type = MFEMMesh
  file = ../mesh/square.e
[]

[Problem]
  type = MFEMProblem
[]

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

[Variables]
  [time_integrated_heat_flux]
    type = MFEMVariable
    fespace = HDivFESpace
    time_derivative = heat_flux
  []
  [temperature]
    type = MFEMVariable
    fespace = L2FESpace
  []
[]

[Kernels]
  [dT_dt,T']
    type = MFEMTimeDerivativeMassKernel
    variable = temperature
  []
  [divh,T']
    type = MFEMVectorFEDivergenceKernel
    trial_variable = heat_flux
    variable = temperature
  []
  [h,h']
    type = MFEMTimeDerivativeVectorFEMassKernel
    variable = time_integrated_heat_flux
  []
  [-T,div.h']
    type = MFEMVectorFEDivergenceKernel
    trial_variable = temperature
    variable = time_integrated_heat_flux
    coefficient = -1.0
    transpose = true
  []
[]

[BCs]
  [right]
    type = MFEMVectorFEBoundaryFluxIntegratedBC
    variable = time_integrated_heat_flux
    coefficient = 0.0
    boundary = 2
  []
  [left]
    type = MFEMVectorFEBoundaryFluxIntegratedBC
    variable = time_integrated_heat_flux
    coefficient = -1.0
    boundary = 4
  []
  [topbottom]
    type = MFEMVectorDirichletBC
    variable = heat_flux
    vector_coefficient = '0.0 0.0'
    boundary = '1 3'
  []
[]

[Solver]
  type = MFEMSuperLU
[]

[Executioner]
  type = MFEMTransient
  device = cpu
  assembly_level = legacy
  dt = 0.03
  start_time = 0.0
  end_time = 0.09
[]

[Outputs]
  [ParaViewDataCollection]
    type = MFEMParaViewDataCollection
    file_base = OutputData/MixedHeatTransfer
    vtk_format = ASCII
  []
[]

Summary
Description
Example File