FVKernels System

For an overview of MOOSE FV please see Finite Volume Design Decisions in MOOSE.

For the finite volume method (FVM), FVKernels are the base class for FVFluxKernel, FVElementalKernel. These specialized objects satisfy the following tasks:

* FVFluxKernel represents numerical fluxes evaluate on the element faces. These terms originate from applying Gauss' divergence theorem.

* FVElementalKernel represents terms that do not contain a spatial derivative so that Gauss' theorem cannot be applied. These terms include time derivatives, externally imposed source terms, and reaction terms.

Note: Currently, the FVElementalKernel category only contains kernels (subclasses) representing time derivatives. Kernels representing externally imposed sources or reaction terms will be added in the near future.

commentnote

In the documentation that follows, we will use '-' and '' to represent different sides of a face. This is purely notation. In the MOOSE code base, the '-' side is represented with an _elem suffix and the '' side is represented with a _neighbor suffix. We could just as well have chosen _left and _right, or _1 and _2, or _minus and _plus, but for consistency with previous MOOSE framework code such as discontinuous Galerkin kernels and node-face constraints, we have elected to go with the _elem and _neighbor suffixes.

FVKernels block

FVM kernels are added to simulation input files in the FVKernels block. The FVKernels block in the example below sets up a transient diffusion problem defined by the equation:

The time derivative term corresponds to the kernel named time, while the diffusion term is represented by the kernel named diff.

Listing 1: Example of the FVKernels block in a MOOSE input file.

[FVKernels]
  [./time]
    type = FVFunctorTimeKernel
    variable = v
  [../]
  [diff]
    type = FVDiffusion
    variable = v
    coeff = coeff
  []
[]
(test/tests/fvkernels/fv_simple_diffusion/transient.i)

The FVTimeKernel in the example derives from FVElementalKernel so it's a volumetric contribution to the residual, while the FVDiffusion kernel is an FVFluxKernel and it's a face contribution to the residual. The remaining MOOSE syntax is what you would expect to see in finite element kernel objects:

* variable refers to the variable that this kernel is acting on (i.e. into which equation does the residual of this term go). This must be a finite-volume variable (defined with fv = true) for all FVM kernels.

* coeff in kernel diff is a material property corresponding to the heat conduction or diffusion coefficient.

The next example shows an FVKernels block that solves the one-dimensional Burgers' equation. The Burgers' equation for speed v is given by:

Listing 2: Example of the FVKernels block in a MOOSE input file for solving one-dimensional Burgers' equation.

[FVKernels]
  [./burgers]
    type = FVBurgers1D
    variable = v
  [../]
  [./time]
    type = FVTimeKernel
    variable = v
  [../]
[]
(test/tests/fvkernels/fv_burgers/fv_burgers.i)

Note that the FVBurgers1D kernel only works for one-dimensional problems. In this example, the exact same time derivative kernels as for the diffusion example is used, but the spatial derivative term is different.

Boundary conditions are not discussed in these examples. Look at syntax files for details about boundary conditions.

Available Objects

  • Moose App
  • FVAdvectionResidual contribution from advection operator for finite volume method.
  • FVAnisotropicDiffusionComputes residual for anisotropic diffusion operator for finite volume method.
  • FVBodyForceDemonstrates the multiple ways that scalar values can be introduced into finite volume kernels, e.g. (controllable) constants, functions, and postprocessors.
  • FVBoundedValueConstraintThis class is used to enforce a min or max value for a finite volume variable
  • FVCoupledForceImplements a source term proportional to the value of a coupled variable.
  • FVDiffusionComputes residual for diffusion operator for finite volume method.
  • FVDivergenceComputes the residual coming from the divergence of a vector fieldthat can be represented as a functor.
  • FVFunctorTimeKernelResidual contribution from time derivative of an AD functor (default is the variable this kernel is acting upon if the 'functor' parameter is not supplied) for the finite volume method.
  • FVIntegralValueConstraintThis class is used to enforce integral of phi = volume * phi_0 with a Lagrange multiplier approach.
  • FVMassMatrixComputes a 'mass matrix', which will just be a diagonal matrix for the finite volume method, meant for use in preconditioning schemes which require one
  • FVMatAdvectionComputes the residual of advective term using finite volume method.
  • FVOrthogonalDiffusionImposes an orthogonal diffusion term.
  • FVPointValueConstraintThis class is used to enforce integral of phi = volume * phi_0 with a Lagrange multiplier approach.
  • FVReactionSimple consuming reaction term
  • FVScalarLagrangeMultiplierThis class is used to enforce integral of phi = volume * phi_0 with a Lagrange multiplier approach.
  • FVTimeKernelResidual contribution from time derivative of a variable for the finite volume method.
  • Navier Stokes App
  • CNSFVFluidEnergyHLLCImplements the fluid energy flux portion of the free-flow HLLC discretization.
  • CNSFVMassHLLCImplements the mass flux portion of the free-flow HLLC discretization.
  • CNSFVMomentumHLLCImplements the momentum flux portion of the free-flow HLLC discretization.
  • FVMatPropTimeKernelReturns a material property which should correspond to a time derivative.
  • FVPorosityTimeDerivativeA time derivative multiplied by a porosity material property
  • INSFVBodyForceBody force that contributes to the Rhie-Chow interpolation
  • INSFVEnergyAdvectionAdvects energy, e.g. rho*cp*T. A user may still override what quantity is advected, but the default is rho*cp*T
  • INSFVEnergyTimeDerivativeAdds the time derivative term to the incompressible Navier-Stokes energy equation.
  • INSFVMassAdvectionObject for advecting mass, e.g. rho
  • INSFVMeshAdvectionImplements a source/sink term for this object's variable/advected-quantity proportional to the divergence of the mesh velocity
  • INSFVMixingLengthReynoldsStressComputes the force due to the Reynolds stress term in the incompressible Reynolds-averaged Navier-Stokes equations.
  • INSFVMixingLengthScalarDiffusionComputes the turbulent diffusive flux that appears in Reynolds-averaged fluid conservation equations.
  • INSFVMomentumAdvectionObject for advecting momentum, e.g. rho*u
  • INSFVMomentumBoussinesqComputes a body force for natural convection buoyancy.
  • INSFVMomentumDiffusionImplements the Laplace form of the viscous stress in the Navier-Stokes equation.
  • INSFVMomentumFrictionImplements a basic linear or quadratic friction model as a volumetric force, for example for the X-momentum equation: and for the linear and quadratic models respectively. A linear dependence is expected for laminar flow, while a quadratic dependence is more common for turbulent flow.
  • INSFVMomentumGravityComputes a body force due to gravity in Rhie-Chow based simulations.
  • INSFVMomentumMeshAdvectionImplements a momentum source/sink term proportional to the divergence of the mesh velocity
  • INSFVMomentumPressureIntroduces the coupled pressure term into the Navier-Stokes momentum equation.
  • INSFVMomentumPressureFluxMomentum pressure term eps grad_P, as a flux kernel using the divergence theoreom, in the incompressible Navier-Stokes momentum equation.
  • INSFVMomentumTimeDerivativeAdds the time derivative term to the incompressible Navier-Stokes momentum equation.
  • INSFVPumpEffective body force for a pump that contributes to the Rhie-Chow interpolation
  • INSFVScalarFieldAdvectionAdvects an arbitrary quantity, the associated nonlinear 'variable'.
  • INSFVTKEDSourceSinkElemental kernel to compute the production and destruction terms of turbulent kinetic energy dissipation (TKED).
  • INSFVTKESourceSinkElemental kernel to compute the production and destruction terms of turbulent kinetic energy (TKE).
  • INSFVTurbulentAdvectionAdvects an arbitrary turbulent quantity, the associated nonlinear 'variable'.
  • INSFVTurbulentDiffusionComputes residual for the turbulent scaled diffusion operator for finite volume method.
  • NSFVEnergyAmbientConvectionImplements a solid-fluid ambient convection volumetric term proportional to the difference between the fluid and ambient temperatures : .
  • NSFVMixturePhaseInterfaceImplements a phase-to-phase volumetric exchange.
  • NSFVPhaseChangeSourceComputes the energy source due to solidification/melting.
  • PCNSFVDensityTimeDerivativeA time derivative kernel for which the form is eps * ddt(rho*var).
  • PCNSFVFluidEnergyHLLCImplements the fluid energy flux portion of the porous HLLC discretization.
  • PCNSFVKTComputes the residual of advective term using finite volume method.
  • PCNSFVKTDCComputes the residual of advective term using finite volume method using a deferred correction approach.
  • PCNSFVMassHLLCImplements the mass flux portion of the porous HLLC discretization.
  • PCNSFVMomentumFrictionComputes a friction force term on fluid in porous media in the Navier Stokes i-th momentum equation.
  • PCNSFVMomentumHLLCImplements the momentum flux portion of the porous HLLC discretization.
  • PINSFVEnergyAdvectionAdvects energy, e.g. rho*cp*T. A user may still override what quantity is advected, but the default is rho*cp*T
  • PINSFVEnergyAmbientConvectionImplements the solid-fluid ambient convection term in the porous media Navier Stokes energy equation.
  • PINSFVEnergyAnisotropicDiffusionAnisotropic diffusion term in the porous media incompressible Navier-Stokes equations : -div(kappa grad(T))
  • PINSFVEnergyDiffusionDiffusion term in the porous media incompressible Navier-Stokes fluid energy equations :
  • PINSFVEnergyTimeDerivativeAdds the time derivative term to the Navier-Stokes energy equation: for fluids: d(eps * rho * cp * T)/dt, for solids: (1 - eps) * d(rho * cp * T)/dtMaterial property derivatives are ignored if not provided.
  • PINSFVMassAdvectionObject for advecting mass in porous media mass equation
  • PINSFVMomentumAdvectionObject for advecting superficial momentum, e.g. rho*u_d, in the porous media momentum equation
  • PINSFVMomentumBoussinesqComputes a body force for natural convection buoyancy in porous media: eps alpha (T-T_0)
  • PINSFVMomentumDiffusionViscous diffusion term, div(mu eps grad(u_d / eps)), in the porous media incompressible Navier-Stokes momentum equation.
  • PINSFVMomentumFrictionComputes a friction force term on fluid in porous media in the Navier Stokes i-th momentum equation in Rhie-Chow (incompressible) contexts.
  • PINSFVMomentumFrictionCorrectionComputes a correction term to avoid oscillations from average pressure interpolation in regions of high changes in friction coefficients.
  • PINSFVMomentumGravityComputes a body force, due to gravity on fluid in porous media in Rhie-Chow (incompressible) contexts.
  • PINSFVMomentumPressureIntroduces the coupled pressure term into the Navier-Stokes porous media momentum equation.
  • PINSFVMomentumPressureFluxMomentum pressure term eps grad_P, as a flux kernel using the divergence theoreom, in the porous media incompressible Navier-Stokes momentum equation. This kernel is also executed on boundaries.
  • PINSFVMomentumPressurePorosityGradientIntroduces the coupled pressure times porosity gradient term into the Navier-Stokes porous media momentum equation.
  • PINSFVMomentumTimeDerivativeAdds the time derivative term: d(rho u_d) / dt to the porous media incompressible Navier-Stokes momentum equation.
  • PNSFVMomentumPressureFluxRZAdds the porous term into the radial component of the Navier-Stokes momentum equation for the problems in the RZ coordinate system when integrating by parts.
  • PNSFVMomentumPressureRZAdds the porous term into the radial component of the Navier-Stokes momentum equation for the problems in the RZ coordinate system when integrating by parts.
  • PNSFVPGradEpsilonIntroduces a -p * grad_eps term.
  • PWCNSFVMassAdvectionObject for advecting mass in porous media mass equation
  • PWCNSFVMassTimeDerivativeAdds the time derivative term to the porous weakly-compressible Navier-Stokes continuity equation.
  • WCNSFV2PMomentumAdvectionSlipComputes the slip velocity advection kernel for two-phase mixture model.
  • WCNSFV2PMomentumDriftFluxImplements the drift momentum flux source.
  • WCNSFVEnergyTimeDerivativeAdds the time derivative term to the incompressible Navier-Stokes momentum equation.
  • WCNSFVMassAdvectionObject for advecting mass, e.g. rho
  • WCNSFVMassTimeDerivativeAdds the time derivative term to the weakly-compressible Navier-Stokes continuity equation.
  • WCNSFVMixingLengthEnergyDiffusionComputes the turbulent diffusive flux that appears in Reynolds-averaged fluid energy conservation equations.
  • WCNSFVMomentumTimeDerivativeAdds the time derivative term to the incompressible Navier-Stokes momentum equation.
  • Heat Transfer App
  • FVHeatConductionTimeDerivativeAD Time derivative term of the heat equation for quasi-constant specific heat and the density .
  • Porous Flow App
  • FVPorousFlowAdvectiveFluxAdvective Darcy flux
  • FVPorousFlowDispersiveFluxAdvective Darcy flux
  • FVPorousFlowEnergyTimeDerivativeDerivative of heat energy with respect to time
  • FVPorousFlowHeatAdvectionHeat flux advected by the fluid
  • FVPorousFlowHeatConductionConductive heat flux
  • FVPorousFlowMassTimeDerivativeDerivative of fluid-component mass with respect to time

FVKernel source code: FVDiffusion example

First, FVFluxKernels are discussed. FVFluxKernels are used to calculate numerical flux contributions from face (surface integral) terms to the residual. The residual contribution is implemented by overriding the computeQpResidual function.

In the FVM, one solves for the averages of the variables over each element. The values of the variables on the faces are unknown and must be computed from the cell average values. This interpolation/reconstruction determines the accuracy of the FVM. The discussion is based on the example of FVDiffusion that discretizes the diffusion term using a central difference approximation.

Listing 3: Example source code for a finite volume kernel discretizing the diffusion term using a central difference.


#include "FVDiffusion.h"

registerMooseObject("MooseApp", FVDiffusion);

InputParameters
FVDiffusion::validParams()
{
  InputParameters params = FVFluxKernel::validParams();
  params.addClassDescription("Computes residual for diffusion operator for finite volume method.");
  params.addRequiredParam<MooseFunctorName>("coeff", "diffusion coefficient");
  MooseEnum coeff_interp_method("average harmonic", "harmonic");
  params.addParam<MooseEnum>(
      "coeff_interp_method",
      coeff_interp_method,
      "Switch that can select face interpolation method for diffusion coefficients.");
  params.set<unsigned short>("ghost_layers") = 2;

  return params;
}

FVDiffusion::FVDiffusion(const InputParameters & params)
  : FVFluxKernel(params),
    _coeff(getFunctor<ADReal>("coeff")),
    _coeff_interp_method(
        Moose::FV::selectInterpolationMethod(getParam<MooseEnum>("coeff_interp_method")))
{
  if ((_var.faceInterpolationMethod() == Moose::FV::InterpMethod::SkewCorrectedAverage) &&
      (_tid == 0))
    adjustRMGhostLayers(std::max((unsigned short)(3), _pars.get<unsigned short>("ghost_layers")));
}

ADReal
FVDiffusion::computeQpResidual()
{
  using namespace Moose::FV;
  const auto state = determineState();

  auto dudn = gradUDotNormal(state);
  ADReal coeff;

  // If we are on internal faces, we interpolate the diffusivity as usual
  if (_var.isInternalFace(*_face_info))
  {
    const ADReal coeff_elem = _coeff(elemArg(), state);
    const ADReal coeff_neighbor = _coeff(neighborArg(), state);
    // If the diffusion coefficients are zero, then we can early return 0 (and avoid warnings if we
    // have a harmonic interpolation)
    if (!coeff_elem.value() && !coeff_neighbor.value())
      return 0;

    interpolate(_coeff_interp_method, coeff, coeff_elem, coeff_neighbor, *_face_info, true);
  }
  // Else we just use the boundary values (which depend on how the diffusion
  // coefficient is constructed)
  else
  {
    const auto face = singleSidedFaceArg();
    coeff = _coeff(face, state);
  }

  return -1 * coeff * dudn;
}
(framework/src/fvkernels/FVDiffusion.C)

The kernel FVDiffusion discretizes the diffusion term . Integrating over the extend of an element and using Gauss' theorem leads to:

The term in parenthesis in the surface integral on the right hand side must be implemented in the FVKernel. However, there is one more step before we can implement the kernel. We must determine how the values of and depend on the values of and on the '+' and '-' side of the face . In this example, the following approximation is used:

This is a central difference approximation of the gradient on the face that neglects cross diffusion terms.

Now, the implementation of this numerical flux into FVDiffusion::computeQpResidual is discussed.

* the kernel provides the '-' and '+' values of the variable as _u_elem[_qp] and _u_neighbor[_qp]

* the values of the material properties on the '-' side of the face is obtained by _coeff_elem(getADMaterialProperty<Real>("coeff")) while the '+' side value is obtained by calling getNeighborADMaterialProperty<Real>("coeff").

* geometric information about the '-' and '+' adjacent elements is available from the face_info object.

The implementation is then straight forward. The first line of the code computes dudn which corresponds to the term:

while the second line computes k corresponding to:

The minus sign originates from the minus sign in the original expression. Flow from '-' to '+ is defined to be positive.

FVKernel source code: FVMatAdvection example

In this example the advection term:

is discretized using upwinding. The velocity is denoted by and represents a passive scalar quantity advected by the flow. Upwinding is a strategy that approximates the value of a variable on a face by taking the value from the upwind element (i.e. the element where the flow originates from).

Listing 4: Example source code for a finite volume kernel discretizing advection of a passive scalar.


#include "FVDiffusion.h"

registerMooseObject("MooseApp", FVDiffusion);

InputParameters
FVDiffusion::validParams()
{
  InputParameters params = FVFluxKernel::validParams();
  params.addClassDescription("Computes residual for diffusion operator for finite volume method.");
  params.addRequiredParam<MooseFunctorName>("coeff", "diffusion coefficient");
  MooseEnum coeff_interp_method("average harmonic", "harmonic");
  params.addParam<MooseEnum>(
      "coeff_interp_method",
      coeff_interp_method,
      "Switch that can select face interpolation method for diffusion coefficients.");
  params.set<unsigned short>("ghost_layers") = 2;

  return params;
}

FVDiffusion::FVDiffusion(const InputParameters & params)
  : FVFluxKernel(params),
    _coeff(getFunctor<ADReal>("coeff")),
    _coeff_interp_method(
        Moose::FV::selectInterpolationMethod(getParam<MooseEnum>("coeff_interp_method")))
{
  if ((_var.faceInterpolationMethod() == Moose::FV::InterpMethod::SkewCorrectedAverage) &&
      (_tid == 0))
    adjustRMGhostLayers(std::max((unsigned short)(3), _pars.get<unsigned short>("ghost_layers")));
}

ADReal
FVDiffusion::computeQpResidual()
{
  using namespace Moose::FV;
  const auto state = determineState();

  auto dudn = gradUDotNormal(state);
  ADReal coeff;

  // If we are on internal faces, we interpolate the diffusivity as usual
  if (_var.isInternalFace(*_face_info))
  {
    const ADReal coeff_elem = _coeff(elemArg(), state);
    const ADReal coeff_neighbor = _coeff(neighborArg(), state);
    // If the diffusion coefficients are zero, then we can early return 0 (and avoid warnings if we
    // have a harmonic interpolation)
    if (!coeff_elem.value() && !coeff_neighbor.value())
      return 0;

    interpolate(_coeff_interp_method, coeff, coeff_elem, coeff_neighbor, *_face_info, true);
  }
  // Else we just use the boundary values (which depend on how the diffusion
  // coefficient is constructed)
  else
  {
    const auto face = singleSidedFaceArg();
    coeff = _coeff(face, state);
  }

  return -1 * coeff * dudn;
}
(framework/src/fvkernels/FVDiffusion.C)

Integrating the advection term over the element and using Gauss' theorem leads to:

This term in parenthesis on the right hand side is approximated using upwinding:

where

and if and otherwise. By convention, the normal points from the '-' side to the '+' side.

The implementation is straight forward. In the constructor the '-' and '' velocities are obtained as RealVectorValue material properties. The average is computed and stored in variable v_avg. The direction of the flow is determined using the inner product of v_avg * _normal and the residual is then computed using either the '-' value of given by _u_elem[_qp] or the '' value given by _u_neighbor [_qp].

FVKernel source code: FVTimeKernel

This example demonstrates source code for an FVElementalKernel. FVElementalKernel are volumetric terms. In this case, the kernel is FVTimeKernel.

Listing 5: Example source code for the finite volume time kernel.

FVTimeKernel::computeQpResidual()
{
  return _u_dot[_qp];
}
(framework/src/fvkernels/FVTimeKernel.C)

This kernel implements the term:

The implementation is identical to the implementation of FEM kernels except that the FVM does not require multiplication by the test function.