Kernels System

A "Kernel" is a piece of physics. It can represent one or more operators or terms in the weak form of a partial differential equation. With all terms on the left-hand-side, their sum is referred to as the "residual". The residual is evaluated at several integration quadrature points over the problem domain. To implement your own physics in MOOSE, you create your own kernel by subclassing the MOOSE Kernel class.

In a Kernel subclass the computeQpResidual() function must be overridden. This is where you implement your PDE weak form terms. The following member functions can optionally be overridden:

  • computeQpJacobian()

  • computeQpOffDiagJacobian()

These two functions provide extra information that can help the numerical solver(s) converge faster and better. Inside your Kernel class, you have access to several member variables for computing the residual and Jacobian values in the above mentioned functions:

  • _i, _j: indices for the current test and trial shape functions respectively.

  • _qp: current quadrature point index.

  • _u, _grad_u: value and gradient of the variable this Kernel operates on; indexed by _qp (i.e. _u[_qp]).

  • _test, _grad_test: value () and gradient () of the test functions at the q-points; indexed by _i and then _qp (i.e., _test[_i][_qp]).

  • _phi, _grad_phi: value () and gradient () of the trial functions at the q-points; indexed by _j and then _qp (i.e., _phi[_j][_qp]).

  • _q_point: XYZ coordinates of the current quadrature point.

  • _current_elem: pointer to the current element being operated on.

Custom Kernel Creation

To create a custom kernel, you can follow the pattern of the Diffusion kernel implemented and included in the MOOSE framework. Additionally, Example 2 provides a step-by-step overview of creating your own custom kernel.

The strong-form of the diffusion equations is defined on a 3-D domain as: find such that

(1) where is defined as the boundary on which the value of is fixed to a known constant , is defined as the boundary on which the flux across the boundary is fixed to a known constant , and $\hat{n} is the boundary outward normal.

The weak form is generated by multiplying by a test function () and integrating over the domain (using inner-product notation):

(-\nabla\cdot\nabla u_h, \psi_i) - (f, \psi_i) = 0\quad \forall\,\psi_i and then integrating by parts which gives the weak form:

(2) where is known as the trial function that defines the finite element discretization, , with being the basis functions.

The Jacobian, which is the derivative of Eq. 2 with respect to , is defined as:

(3)

The diffusion kernel header and implementation files are:

// This file is part of the MOOSE framework
// https://www.mooseframework.org
//
// All rights reserved, see COPYRIGHT for full restrictions
// https://github.com/idaholab/moose/blob/master/COPYRIGHT
//
// Licensed under LGPL 2.1, please see LICENSE for details
// https://www.gnu.org/licenses/lgpl-2.1.aspx

#ifndef DIFFUSION_H
#define DIFFUSION_H

#include "Kernel.h"

class Diffusion;

template <>
InputParameters validParams<Diffusion>();

/**
 * This kernel implements the Laplacian operator:
 * $\nabla u \cdot \nabla \phi_i$
 */
class Diffusion : public Kernel
{
public:
  Diffusion(const InputParameters & parameters);

protected:
  virtual Real computeQpResidual() override;

  virtual Real computeQpJacobian() override;
};

#endif /* DIFFUSION_H */
(moose/framework/include/kernels/Diffusion.h)
// This file is part of the MOOSE framework
// https://www.mooseframework.org
//
// All rights reserved, see COPYRIGHT for full restrictions
// https://github.com/idaholab/moose/blob/master/COPYRIGHT
//
// Licensed under LGPL 2.1, please see LICENSE for details
// https://www.gnu.org/licenses/lgpl-2.1.aspx

#include "Diffusion.h"

registerMooseObject("MooseApp", Diffusion);

template <>
InputParameters
validParams<Diffusion>()
{
  InputParameters params = validParams<Kernel>();
  params.addClassDescription("The Laplacian operator ($-\\nabla \\cdot \\nabla u$), with the weak "
                             "form of $(\\nabla \\phi_i, \\nabla u_h)$.");
  return params;
}

Diffusion::Diffusion(const InputParameters & parameters) : Kernel(parameters) {}

Real
Diffusion::computeQpResidual()
{
  return _grad_u[_qp] * _grad_test[_i][_qp];
}

Real
Diffusion::computeQpJacobian()
{
  return _grad_phi[_j][_qp] * _grad_test[_i][_qp];
}
(moose/framework/src/kernels/Diffusion.C)

Before a custom physics kernel is available for use, it must be registered in your application. This is done in e.g. src/base/YourApp.C for the YourApp application.


// src/base/YourApp.C contents

#include "YourKernel.h"

...

int YourApp::registerObjects(Factory & factory)
{
  ...
  registerKernel(YourKernel);
  ...
}

...

Time Derivative Kernels

You can create a time-derivative term/kernel by subclassing TimeKernel instead of Kernel. For example, the residual contribution for a time derivative term is:

(4)

where is the finite element solution, and

(5)

because you can interchange the order of differentiation and summation.

In the equation above, is the time derivative of the th finite element coefficient of . While the exact form of this derivative depends on the time stepping scheme, without much loss of generality, we can assume the following form for the time derivative:

(6)

for some constants , which depend on and the timestepping method.

The derivative of equation Eq. 5 with respect to is then:

(7)

So that the Jacobian term for equation Eq. 5 is

(8)

where is what we call du_dot_du in MOOSE.

Therefore the computeQpResidual() function for our time-derivative term kernel looks like:

text cpp return _test[_i][_qp] * _u_dot[_qp]; ``

And the corresponding computeQpJacobian() is:


return _test[_i][_qp] * _phi[_j][_qp] * _du_dot_du[_qp];

Further Kernel Documentation

Several specialized kernel types exist in MOOSE each with useful functionality. Details for each are in the sections below.

  • AnisoHeatConduction
  • ConsistentHeatCapacityTimeDerivativeTime derivative term of the heat equation with the heat capacity as an argument.
  • ConsistentSpecificHeatTimeDerivativeTime derivative term with of the heat equation with the specific heat capacity and density as arguments.
  • HeatCapacityConductionTimeDerivativeTime derivative term of the heat equation with the heat capacity as an argument.
  • HeatConductionComputes residual/Jacobian contribution for term.
  • HeatConductionTimeDerivativeTime derivative term of the heat equation for quasi-constant specific heat and the density .
  • HeatSourceDemonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form .
  • HomogenizedHeatConduction
  • JouleHeatingSourceDemonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form .
  • SpecificHeatConductionTimeDerivativeTime derivative term of the heat equation with the specific heat and the density as arguments.
  • CosseratStressDivergenceTensorsStress divergence kernel for the Cartesian coordinate system
  • DynamicStressDivergenceTensorsResidual due to stress related Rayleigh damping and HHT time integration terms
  • GeneralizedPlaneStrainOffDiagGeneralized Plane Strain kernel to provide contribution of the out-of-plane strain to other kernels
  • GravityApply gravity. Value is in units of acceleration.
  • InertialForceCalculates the residual for the interial force () and the contribution of mass dependent Rayleigh damping and HHT time integration scheme ($\eta \cdot M \cdot ((1+\alpha)velq2-\alpha \cdot vel-old) )
  • InertialForceBeamCalculates the residual for the interial force/moment and the contribution of mass dependent Rayleigh damping and HHT time integration scheme.
  • InertialTorqueKernel for interial torque: density * displacement x acceleration
  • MomentBalancing
  • OutOfPlanePressureApply pressure in the out-of-plane direction in 2D plane stress or generalized plane strain models
  • PhaseFieldFractureMechanicsOffDiagStress divergence kernel for phase-field fracture: Computes off diagonal damage dependent Jacobian components. To be used with StressDivergenceTensors or DynamicStressDivergenceTensors.
  • PlasticHeatEnergyPlastic heat energy density = coeff * stress * plastic_strain_rate
  • PoroMechanicsCouplingAdds , where the subscript is the component.
  • StressDivergenceBeamQuasi-static and dynamic stress divergence kernel for Beam element
  • StressDivergenceRSphericalTensorsCalculate stress divergence for an spherically symmetric 1D problem in polar coordinates.
  • StressDivergenceRZTensorsCalculate stress divergence for an axisymmetric problem in cylinderical coordinates.
  • StressDivergenceTensorsStress divergence kernel for the Cartesian coordinate system
  • StressDivergenceTensorsTrussKernel for truss element
  • WeakPlaneStressPlane stress kernel to provide out-of-plane strain contribution
  • ACGBPolyGrain-Boundary model concentration dependent residual
  • ACGrGrElasticDrivingForceAdds elastic energy contribution to the Allen-Cahn equation
  • ACGrGrMultiMulti-phase poly-crystaline Allen-Cahn Kernel
  • ACGrGrPolyGrain-Boundary model poly-crystaline interface Allen-Cahn Kernel
  • ACInterfaceGradient energy Allen-Cahn Kernel
  • ACInterfaceKobayashi1Anisotropic gradient energy Allen-Cahn Kernel Part 1
  • ACInterfaceKobayashi2Anisotropic Gradient energy Allen-Cahn Kernel Part 2
  • ACInterfaceStressInterface stress driving force Allen-Cahn Kernel
  • ACMultiInterfaceGradient energy Allen-Cahn Kernel with cross terms
  • ACSEDGPolyStored Energy contribution to grain growth
  • ACSwitchingKernel for Allen-Cahn equation that adds derivatives of switching functions and energies
  • AllenCahnAllen-Cahn Kernel that uses a DerivativeMaterial Free Energy
  • CHBulkPFCTradCahn-Hilliard kernel for a polynomial phase field crystal free energy.
  • CHCpldPFCTradSplit with a variable that holds the Laplacian of a phase field variable.
  • CHInterfaceGradient energy Cahn-Hilliard Kernel with a scalar (isotropic) mobility
  • CHInterfaceAnisoGradient energy Cahn-Hilliard Kernel with a tensor (anisotropic) mobility
  • CHMathSimple demonstration Cahn-Hilliard Kernel using an algebraic double-well potential
  • CHPFCRFFCahn-Hilliard residual for the RFF form of the phase field crystal model
  • CHSplitChemicalPotentialChemical potential kernel in Split Cahn-Hilliard that solves chemical potential in a weak form
  • CHSplitConcentrationConcentration kernel in Split Cahn-Hilliard that solves chemical potential in a weak form
  • CHSplitFluxComputes flux as nodal variable
  • CahnHilliardCahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy and a scalar (isotropic) mobility
  • CahnHilliardAnisoCahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy and a tensor (anisotropic) mobility
  • CoefCoupledTimeDerivativeScaled time derivative Kernel that acts on a coupled variable
  • CoefReactionImplements the residual term (p*u, test)
  • ConservedLangevinNoiseSource term for noise from a ConservedNoise userobject
  • CoupledAllenCahnCoupled Allen-Cahn Kernel that uses a DerivativeMaterial Free Energy
  • CoupledMaterialDerivativeKernel that implements the first derivative of a function material property with respect to a coupled variable.
  • CoupledSusceptibilityTimeDerivativeA modified coupled time derivative Kernel that multiplies the time derivative of a coupled variable by a generalized susceptibility
  • CoupledSwitchingTimeDerivativeCoupled time derivative Kernel that multiplies the time derivative by $\frac{dh_\alpha}{d\eta_i} F_\alpha + \frac{dh_\beta}{d\eta_i} F_\beta + \dots)
  • GradientComponentSet the kernel variable to a specified component of the gradient of a coupled variable.
  • HHPFCRFFReaction type kernel for the RFF phase fit crystal model
  • KKSACBulkCKKS model kernel (part 2 of 2) for the Bulk Allen-Cahn. This includes all terms dependent on chemical potential.
  • KKSACBulkFKKS model kernel (part 1 of 2) for the Bulk Allen-Cahn. This includes all terms NOT dependent on chemical potential.
  • KKSCHBulkKKS model kernel for the Bulk Cahn-Hilliard term. This operates on the concentration 'c' as the non-linear variable
  • KKSMultiACBulkCMulti-phase KKS model kernel (part 2 of 2) for the Bulk Allen-Cahn. This includes all terms dependent on chemical potential.
  • KKSMultiACBulkFKKS model kernel (part 1 of 2) for the Bulk Allen-Cahn. This includes all terms NOT dependent on chemical potential.
  • KKSMultiPhaseConcentrationKKS multi-phase model kernel to enforce . The non-linear variable of this kernel is , the final phase concentration in the list.
  • KKSPhaseChemicalPotentialKKS model kernel to enforce the pointwise equality of phase chemical potentials dFa/dca = dFb/dcb. The non-linear variable of this kernel is ca.
  • KKSPhaseConcentrationKKS model kernel to enforce the decomposition of concentration into phase concentration (1-h(eta))ca + h(eta)cb - c = 0. The non-linear variable of this kernel is cb.
  • KKSSplitCHCResKKS model kernel for the split Bulk Cahn-Hilliard term. This operates on the chemical potential 'c' as the non-linear variable
  • LangevinNoiseSource term for non-conserved Langevin noise
  • LaplacianSplitSplit with a variable that holds the Laplacian of a phase field variable.
  • MaskedBodyForceKernel that defines a body force modified by a material mask
  • MatAnisoDiffusionDiffusion equation Kernel that takes an anisotropic Diffusivity from a material property
  • MatDiffusionDiffusion equation Kernel that takes an isotropic Diffusivity from a material property
  • MatGradSquareCoupledGradient square of a coupled variable.
  • MatReactionKernel to add -L*v, where L=reaction rate, v=variable
  • MultiGrainRigidBodyMotionAdds rigid mody motion to grains
  • SimpleACInterfaceGradient energy for Allen-Cahn Kernel with constant Mobility and Interfacial parameter
  • SimpleCHInterfaceGradient energy for Cahn-Hilliard equation with constant Mobility and Interfacial parameter
  • SimpleCoupledACInterfaceGradient energy for Allen-Cahn Kernel with constant Mobility and Interfacial parameter for a coupled order parameter variable.
  • SimpleSplitCHWResGradient energy for split Cahn-Hilliard equation with constant Mobility for a coupled order parameter variable.
  • SingleGrainRigidBodyMotionAdds rigid mody motion to a single grain
  • SoretDiffusionAdd Soret effect to Split formulation Cahn-Hilliard Kernel
  • SplitCHMathSimple demonstration split formulation Cahn-Hilliard Kernel using an algebraic double-well potential
  • SplitCHParsedSplit formulation Cahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy
  • SplitCHWResSplit formulation Cahn-Hilliard Kernel for the chemical potential variable with a scalar (isotropic) mobility
  • SplitCHWResAnisoSplit formulation Cahn-Hilliard Kernel for the chemical potential variable with a tensor (anisotropic) mobility
  • SusceptibilityTimeDerivativeA modified time derivative Kernel that multiplies the time derivative of a variable by a generalized susceptibility
  • SwitchingFunctionConstraintEtaLagrange multiplier kernel to constrain the sum of all switching functions in a multiphase system. This kernel acts on a non-conserved order parameter eta_i.
  • SwitchingFunctionConstraintLagrangeLagrange multiplier kernel to constrain the sum of all switching functions in a multiphase system. This kernel acts on the lagrange multiplier variable.
  • SwitchingFunctionPenaltyPenalty kernel to constrain the sum of all switching functions in a multiphase system.
  • CoefDiffusionKernel for diffusion with diffusivity = coef + function
  • ThermoDiffusionKernel for thermo-diffusion (Soret effect, thermophoresis, etc.)
  • ArrheniusDiffusionDiffusion with Arrhenius coefficient
  • CompositeHeatConductionCompute thermal conductivity
  • ConstituentDiffusion
  • ConstitutiveHeatConductionThe Laplacian operator (), with the weak form of .
  • ConstitutiveHeatConductionTimeDerivativeTime derivative term of the heat equation for quasi-constant specific heat and the density .
  • Decay
  • DiffusionLimitedReactionCalculates losses due to diffusion limited reaction
  • FissionRateHeatSource
  • HZrHSource
  • HydrideSourceAdd source (sink) term for precipitation (dissolution) of hydrogen as hydride
  • HydrogenDiffusionCalculates the diffusion of hydrogen in solid solution due to Fick's law and the Soret effect
  • HydrogenPrecipitationCalculates the precipitation of hydrogen in solid solution from McMinn's TSSp equilibrium and the dissolution hydrides to solid solution hydrogen from Marino's kinetics
  • HydrogenSourceAdd source (sink) term for dissolved hydrogen from hydride dissolution (precipitation)
  • HydrogenTimeDerivativeTime derivative for species where the volume fraction of the phase is time-dependent.
  • IsotropicDiffusionIsotropic diffusion that uses arbitrary diffusivity
  • MOXActinideRedistributionMOX kernel used to simulate actinide redistribution.
  • MOXActinideRedistributionEnhancementMOX kernel used to simulate actinide redistribution enhanced by porosity.
  • MOXOxygenDiffusionMOX oxygen diffusion kernel.
  • MOXPoreContinuityMOX kernel used to simulate pore migration.
  • MOXPoreDiffusionMOX porosity diffusion kernel used with kernel MOXPoreContinuity.
  • NeutronHeatSourceCompute heat generation due to fission.
  • OxideEnergyDepositionComputes the amount of energy released from the zirconium oxide reaction and applies it to the cladding.
  • OxygenDiffusion
  • ZirconiumDiffusionCalculates the amount of zirconium that is transported across the mesh
  • AnisotropicDiffusionAnisotropic diffusion kernel with weak form given by .
  • BodyForceDemonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form .
  • CoefTimeDerivativeThe time derivative operator with the weak form of .
  • ConservativeAdvectionConservative form of which in its weak form is given by: .
  • CoupledForceImplements a source term proportional to the value of a coupled variable. Weak form: .
  • CoupledTimeDerivativeTime derivative Kernel that acts on a coupled variable. Weak form: .
  • DiffusionThe Laplacian operator (), with the weak form of .
  • MassEigenKernelAn eigenkernel with weak form where is the eigenvalue.
  • MassLumpedTimeDerivativeLumped formulation of the time derivative . Its corresponding weak form is where denotes the time derivative of the solution coefficient associated with node .
  • MaterialDerivativeRankFourTestKernelClass used for testing derivatives of a rank four tensor material property.
  • MaterialDerivativeRankTwoTestKernelClass used for testing derivatives of a rank two tensor material property.
  • MaterialDerivativeTestKernelClass used for testing derivatives of a scalar material property.
  • NullKernelKernel that sets a zero residual.
  • ReactionImplements a simple consuming reaction term with weak form .
  • TimeDerivativeThe time derivative operator with the weak form of .
  • UserForcingFunctionDemonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form .
  • VectorBodyForceDemonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form .
  • VectorDiffusionThe Laplacian operator (), with the weak form of .
  • CrackTipEnrichmentStressDivergenceTensorsEnrich stress divergence kernel for small-strain simulations