DGKernels System

DGKernels System

Overview

DGKernels are the kernels defined on internal sides. DGKernels are typically for elemental variables, i.e. variables that allows solutions to be discontinous aross element sides. DGKernels along with normal kernels allow the definition of weak forms rised from discontinous finite element methods (DFEM). DGKernels can be block restricted for calculatons with DFEM on subdomains. Internal sides are visited once during residual or Jacobian evaluations by MOOSE. DGKernels handle two pieces of residual, marked as Element and Neighbor, on an internal side and corresponding four pieces of Jacobian, marked as ElementElement, ElementNeighbor, NeighborElement and NeighborNeighbor. The normals on internal sides are pointing towards neighboring element from the current element. Typically DGKernels are irrelevant with the normal direction. When there are mesh refinement, MOOSE visits all the active internal sides, meaning that if there is a hanging node for an internal side, MOOSE visit the child internal sides. DGKernels can make use of the material properties defined on both Element and Neighbor sides. The DGKernel with interior penalty (IP) method for diffusion equations can be found at DGDiffusion. The DGKernel with upwinding scheme for hyperbolic equations can be found at DGConvection.

Extension for Hybrid Finite Element Methods

DGKernels are extended to support hybrid finite element methods (HFEM) (Raviart and Thomas, 1977).

Considering Poisson's equation of the form

2u=fΩu=gΩDun=hΩN,αu+un=cΩR,\begin{aligned} -\nabla^2 u &= f && \quad \in \Omega \\ u &= g && \quad \in \partial \Omega_D \\ \frac{\partial u}{\partial n} &= h && \quad \in \partial \Omega_N, \\ \alpha u + \frac{\partial u}{\partial n} &= c && \quad \in \partial \Omega_R, \end{aligned}

where ΩRn\Omega \subset \mathbb{R}^n is the domain, and Ω=ΩDΩNΩR\partial \Omega = \partial \Omega_D \cup \partial \Omega_N \cup \partial \Omega_R is its boundary. α\alpha is a given function on ΩR\Omega_R, which has typical constant value of 1/21/2.

The weak form with HFEM for this PDE is to find a triple (u,λ,λD)(u, \lambda, \lambda_D), in discontinuous function spaces on Ω\Omega, all internal sides Γ\Gamma and ΩD\partial \Omega_D respectively, such that, (u,λ,λD)\forall (u^\ast, \lambda^\ast, \lambda^\ast_D) in the same function spaces,

(u,u)Ω+(λ,[u])Γ+([u],λ)Γ+(u,λD)ΩD+(λD,ug)ΩD(u,h)ΩN+(u,αuc)ΩR(u,f)Ω=0.\begin{aligned} \left( \nabla u^\ast, \nabla u \right)_\Omega + \left( \lambda^\ast, [ u ] \right)_\Gamma + \left( [ u^\ast ], \lambda \right)_\Gamma +\left( u^\ast, \lambda_D \right)_{\partial \Omega_D} + \left( \lambda_D^\ast, u - g \right)_{\partial \Omega_D} -\left( u^\ast, h \right)_{\partial \Omega_N} +\left( u^\ast, \alpha u - c \right)_{\partial \Omega_R} -\left( u^\ast, f \right)_\Omega = 0. \end{aligned}

[u][u] represents the jump of uu on an internal side. It is noted that the orientation of normals on internal sides does not affect the solution of uu but flips the sign of λ\lambda. λ\lambda and λD\lambda_D are also known as Lagrange multipliers for weakly imposing the continuity of uu across internal sides on Γ\Gamma and imposing the Dirichlet boundary condition at ΩD\Omega_D. They resemble the current (un-\frac{\partial u}{\partial n}) and converge to the current when discretization error gets smaller and smaller with mesh refinement. HFEM has explicit local conservation, which can be seen if we substitute a test function of uu^\ast with constant value of one element of interest and zero elsewhere. The local conservation is evaluated with Lagrange multiplier λ\lambda and the source function ff for an element inside of the domain Ω\Omega.

This weak form requires a compatibility condition to have a unique solution (Raviart and Thomas, 1977). Typically we satisfy this condition by letting the order of the shape function for uu two order higher (including two) than the order for λ\lambda.

An alternative way of imposing the Robin boundary condition (αu+un=c\alpha u + \frac{\partial u}{\partial n} = c) is to replace (u,αuc)ΩR\left( u^\ast, \alpha u - c \right)_{\partial \Omega_R} with

(u,λR)ΩR+(λR,uuR)ΩR+(uR,αuRcλR)ΩR,\begin{aligned} \left( u^\ast, \lambda_R \right)_{\partial \Omega_R} + \left( \lambda_R^\ast, u - u_R \right)_{\partial \Omega_R} + \left( u_R^\ast, \alpha u_R - c - \lambda_R \right)_{\partial \Omega_R}, \end{aligned}

with Lagrange multiplier λR\lambda_R and the projected solution uRu_R on ΩR\partial \Omega_R and their corresponding test functions λR\lambda_R^\ast and uRu_R^\ast.

A variable for the Lagrange multiplier defined on all interior sides Γ\Gamma can be coupled in DGKernels with a lower-dimensional mesh derived from the main mesh for Γ\Gamma. With this extension, DGKernels can handle three pieces of residual, marked as Element, Neighbor and Lower, on an internal side and corresponding nine pieces of Jacobian, marked as ElementElement, ElementNeighbor, NeighborElement, NeighborNeighbor, PrimaryLower, SecondaryLower, LowerPrimary, LowerSecondary and LowerLower. Similarly, a variable for the Lagrange multiplier on boundary ΩD\Omega_D can be coupled in integrated boundary conditions.

The DGKernal for (λ,[u])Γ+([u],λ)Γ\left( \lambda^\ast, [ u ] \right)_\Gamma + \left( [ u^\ast ], \lambda \right)_\Gamma can be found at HFEMDiffusion. The boundary condition for (u,λD)ΩD+(λD,ug)ΩD\left( u^\ast, \lambda_D \right)_{\partial \Omega_D} + \left( \lambda_D^\ast, u - g \right)_{\partial \Omega_D} can be found at HFEMDirichletBC with the generalization of gg being either a fixed value or a variable defined on boundary. With this generalization, HFEMDirichletBC can be used along with kernels defined on the Robin boundary for (uR,αuRcλR)ΩR\left( u_R^\ast, \alpha u_R - c - \lambda_R \right)_{\partial \Omega_R} for the alternative way of imposing the Robin boundary condition above.

warningwarning

MOOSE does not support mesh adaptation with HFEM currently.

Example Input File Syntax

DGKernels are added through DGKernels input syntax.

Available Objects

  • Moose App
  • ADDGAdvectionAdds internal face advection flux contributions for discontinuous Galerkin discretizations
  • ADDGDiffusionDG kernel for diffusion operator
  • ArrayDGDiffusionImplements interior penalty method for array diffusion equations.
  • ArrayHFEMDiffusionImposes the constraints on internal sides with HFEM.
  • DGConvectionDG upwinding for the convection
  • DGDiffusionComputes residual contribution for the diffusion operator using discontinous Galerkin method.
  • HFEMDiffusionImposes the constraints on internal sides with HFEM.
  • HFEMTestJumpImposes constraints for HFEM with side-discontinuous variables.
  • HFEMTrialJumpImposes constraints for HFEM with side-discontinuous variables.
  • Rdg App
  • AEFVKernelA dgkernel for the advection equation using a cell-centered finite volume method.

Available Actions

References

  1. P. A. Raviart and J. M. Thomas. Primal hybrid finite element methods for 2nd order elliptic equations. MATHEMATICS OF COMPUTATION, 31(138):391–413, 1977.[BibTeX]