Kokkos AuxKernels System
Before reading this documentation, consider reading the following materials first for a better understanding of this documentation:
AuxKernels System to understand the MOOSE auxiliary kernel system,
Getting Started with Kokkos-MOOSE to understand the programming practices for Kokkos-MOOSE,
Kokkos Kernels System to understand the common design pattern of objects in Kokkos-MOOSE.
Kokkos-MOOSE auxiliary kernels do not support boundary-restricted elemental kernels yet.
The Kokkos-MOOSE auxiliary kernels can be created by inheriting Moose::Kokkos::AuxKernel and registering them with registerKokkosAuxKernel(). There is no vector or array version at this moment. The signature of the computeValue() hook method is now as follows:
KOKKOS_FUNCTION Real computeValue(const unsigned int qp, AssemblyDatum & datum) const;
To reiterate, the hook method no longer uses virtual dispatch. Instead, it should be statically defined in the derived class as inlined public method.
Like the original MOOSE, the single interface is used for both elemental and nodal kernels, and whether the object is elemental or nodal can be queried with the isNodal() method on both CPU and GPU. Depending on whether the kernel is nodal or elemental, some data that can be queried through datum may be invalid. For instance, datum.node() will throw an invalid index in elemental kernels, and vice versa. Therefore, if you need to treat elemental and nodal variables separately, it is important to properly separate their code paths. For nodal kernels, the index qp is simply zero but still passed as an argument to have a unified interface, as explained in nodal boundary conditions.
The following example shows the implementation of KokkosVariableTimeIntegrationAux:
Listing 1: The KokkosVariableTimeIntegrationAux header file.
#pragma once
#include "KokkosAuxKernel.h"
/**
* An AuxKernel that can be used to integrate a field variable in time
* using a variety of different integration methods. The result is
* stored in another field variable.
*/
class KokkosVariableTimeIntegrationAux : public Moose::Kokkos::AuxKernel
{
public:
static InputParameters validParams();
KokkosVariableTimeIntegrationAux(const InputParameters & parameters);
KOKKOS_FUNCTION Real computeValue(const unsigned int qp, AssemblyDatum & datum) const;
protected:
KOKKOS_FUNCTION Real getIntegralValue(const unsigned int qp, AssemblyDatum & datum) const;
Moose::Kokkos::Array<Moose::Kokkos::VariableValue> _coupled_vars;
const Real _coef;
const unsigned int _order;
Moose::Kokkos::Array<Real> _integration_coef;
/// The old variable value (zero if order == 3)
const Moose::Kokkos::VariableValue _u_old;
/// The older variable value (zero if order != 3)
const Moose::Kokkos::VariableValue _u_older;
};
KOKKOS_FUNCTION inline Real
KokkosVariableTimeIntegrationAux::computeValue(const unsigned int qp, AssemblyDatum & datum) const
{
Real integral = getIntegralValue(qp, datum);
if (_order == 3)
return _u_older(datum, qp) + _coef * integral;
return _u_old(datum, qp) + _coef * integral;
}
KOKKOS_FUNCTION inline Real
KokkosVariableTimeIntegrationAux::getIntegralValue(const unsigned int qp,
AssemblyDatum & datum) const
{
Real integral_value = 0.0;
for (unsigned int i = 0; i < _order; ++i)
integral_value += _integration_coef[i] * _coupled_vars[i](datum, qp) * _dt;
/**
* Subsequent timesteps may be unequal, so the standard Simpson rule
* cannot be used. Use a different set of coefficients here.
* J. McNAMEE, "A PROGRAM TO INTEGRATE A FUNCTION TABULATED AT
* UNEQUAL INTERVALS," Internation Journal for Numerical Methods in
* Engineering, Vol. 17, 217-279. (1981).
*/
if (_order == 3 && _dt != _dt_old)
{
Real x0 = 0.0;
Real x1 = _dt_old;
Real x2 = _dt + _dt_old;
Real y0 = _coupled_vars[2](datum, qp);
Real y1 = _coupled_vars[1](datum, qp);
Real y2 = _coupled_vars[0](datum, qp);
Real term1 = (x2 - x0) * (y0 + (x2 - x0) * (y1 - y0) / (2.0 * (x1 - x0)));
Real term2 = (2.0 * x2 * x2 - x0 * x2 - x0 * x0 + 3.0 * x0 * x1 - 3.0 * x1 * x2) / 6.0;
Real term3 = (y2 - y1) / (x2 - x1) - (y1 - y0) / (x1 - x0);
integral_value = term1 + term2 * term3;
}
return integral_value;
}
(framework/include/kokkos/auxkernels/KokkosVariableTimeIntegrationAux.h)Listing 2: The KokkosVariableTimeIntegrationAux source file.
#include "KokkosVariableTimeIntegrationAux.h"
registerKokkosAuxKernel("MooseApp", KokkosVariableTimeIntegrationAux);
InputParameters
KokkosVariableTimeIntegrationAux::validParams()
{
InputParameters params = AuxKernel::validParams();
params.addClassDescription("Integrates a field variable in time.");
params.addRequiredCoupledVar("variable_to_integrate", "The variable to be integrated");
params.addParam<Real>("coefficient", 1.0, "A simple coefficient multiplying the integral");
params.addParam<unsigned int>(
"order", 2, "The order of global truncation error: midpoint=1, trapezoidal=2, Simpson=3");
return params;
}
KokkosVariableTimeIntegrationAux::KokkosVariableTimeIntegrationAux(
const InputParameters & parameters)
: AuxKernel(parameters),
_coef(getParam<Real>("coefficient")),
_order(getParam<unsigned int>("order")),
_u_old(_order != 3 ? uOld() : kokkosZeroValue()),
_u_older(_order == 3 ? uOlder() : kokkosZeroValue())
{
_integration_coef.create(_order);
_coupled_vars.create(_order);
switch (_order)
{
case 1:
_integration_coef[0] = 1.0;
_coupled_vars[0] = kokkosCoupledValue("variable_to_integrate");
break;
case 2:
_integration_coef[0] = 0.5;
_integration_coef[1] = 0.5;
_coupled_vars[0] = kokkosCoupledValue("variable_to_integrate");
_coupled_vars[1] = kokkosCoupledValueOld("variable_to_integrate");
break;
case 3:
_integration_coef[0] = 1.0 / 3.0;
_integration_coef[1] = 4.0 / 3.0;
_integration_coef[2] = 1.0 / 3.0;
_coupled_vars[0] = kokkosCoupledValue("variable_to_integrate");
_coupled_vars[1] = kokkosCoupledValueOld("variable_to_integrate");
_coupled_vars[2] = kokkosCoupledValueOlder("variable_to_integrate");
break;
default:
mooseError("KokkosVariableTimeIntegrationAux: unknown time integration order specified");
}
_integration_coef.copyToDevice();
_coupled_vars.copyToDevice();
}
(framework/src/kokkos/auxkernels/KokkosVariableTimeIntegrationAux.K)Advanced Topics
The computational loops are defined in computeElementInternal() and computeNodeInternal() for elemental and nodal variables, respectively, and they call the hook method computeValue(). Depending on the type of operations, however, you may be able to write more efficient kernels by customizing the loops directly. For example, what KokkosCopyValueAux does is simply copy one auxiliary variable to another, and the shape function family and order are always identical between the two variables. In this case, you can directly loop over DOFs instead of quadrature points. See the following source code of KokkosCopyValueAux, where computeElementInternal() and computeNodeInternal() are redefined instead of computeValue():
Listing 3: The KokkosCopyValueAux header file.
#pragma once
#include "KokkosAuxKernel.h"
/**
* Copies one variable onto an auxiliary variable
*/
class KokkosCopyValueAux : public Moose::Kokkos::AuxKernel
{
public:
static InputParameters validParams();
KokkosCopyValueAux(const InputParameters & parameters);
template <typename Derived>
KOKKOS_FUNCTION void computeElementInternal(const Derived & auxkernel,
AssemblyDatum & datum) const;
template <typename Derived>
KOKKOS_FUNCTION void computeNodeInternal(const Derived & auxkernel, AssemblyDatum & datum) const;
protected:
/// Variable used to specify state being copied
unsigned short _state;
/// The variable value to copy from
const Moose::Kokkos::VariableValue _v;
/// A reference to the variable to copy from
const MooseVariable & _source_variable;
};
template <typename Derived>
KOKKOS_FUNCTION void
KokkosCopyValueAux::computeElementInternal(const Derived & /* auxkernel */,
AssemblyDatum & datum) const
{
auto & sys = kokkosSystem(_kokkos_var.sys());
auto var = _kokkos_var.var();
auto tag = _kokkos_var.tag();
auto elem = datum.elem().id;
for (unsigned int i = 0; i < datum.n_dofs(); ++i)
sys.getVectorDofValue(sys.getElemLocalDofIndex(elem, i, var), tag) = _v(datum, i);
}
template <typename Derived>
KOKKOS_FUNCTION void
KokkosCopyValueAux::computeNodeInternal(const Derived & /* auxkernel */,
AssemblyDatum & datum) const
{
auto & sys = kokkosSystem(_kokkos_var.sys());
auto var = _kokkos_var.var();
auto tag = _kokkos_var.tag();
auto node = datum.node();
sys.getVectorDofValue(sys.getNodeLocalDofIndex(node, 0, var), tag) = _v(datum, 0);
}
(framework/include/kokkos/auxkernels/KokkosCopyValueAux.h)Listing 4: The KokkosCopyValueAux source file.
#include "KokkosCopyValueAux.h"
registerKokkosAuxKernel("MooseApp", KokkosCopyValueAux);
InputParameters
KokkosCopyValueAux::validParams()
{
InputParameters params = AuxKernel::validParams();
params.addClassDescription("Returns the specified variable as an auxiliary variable with a "
"simple copy of the variable values.");
params.addRequiredCoupledVar("source", "Variable to take the value of.");
MooseEnum stateEnum("CURRENT=0 OLD=1 OLDER=2", "CURRENT");
params.addParam<MooseEnum>(
"state",
stateEnum,
"This parameter specifies the state being copied. CURRENT=0 OLD=1 OLDER=2. Copying an older "
"state allows access to previous solution information if necessary.");
return params;
}
KokkosCopyValueAux::KokkosCopyValueAux(const InputParameters & parameters)
: AuxKernel(parameters),
_state(getParam<MooseEnum>("state")),
_v(_state == 0 ? kokkosCoupledDofValue("source")
: _state == 1 ? kokkosCoupledDofValueOld("source")
: kokkosCoupledDofValueOlder("source")),
_source_variable(*getVar("source", 0))
{
if (_var.feType().family != _source_variable.feType().family)
paramError("variable",
"Source (" + Moose::stringify(_source_variable.feType().family) + ") and target (" +
Moose::stringify(_var.feType().family) +
") variable have different families. You may use a ProjectionAux "
"instead of the CopyValueAux");
if (_var.feType().order != _source_variable.feType().order)
paramError("variable",
"Source (" + Moose::stringify(_source_variable.feType().order) + ") and target (" +
Moose::stringify(_var.feType().order) +
") variable are of different orders. You may use a ProjectionAux "
"instead of the CopyValueAux");
}
(framework/src/kokkos/auxkernels/KokkosCopyValueAux.K)Available Objects
- Moose App
- ADDivergenceAuxComputes the divergence of a vector of functors.
- ADFunctorElementalAuxEvaluates a functor (variable, function or functor material property) on the current element, quadrature point, or node.
- ADFunctorElementalGradientAuxEvaluates the gradient of a functor (variable, function or functor material property) on the current element or quadrature point.
- ADFunctorMaterialRealAuxOutputs element volume-averaged material properties
- ADFunctorMaterialRealVectorValueAuxCapture a component of a vector material property in an auxiliary variable.
- ADFunctorVectorElementalAuxEvaluates a vector functor (material property usually) on the current element.For finite volume, this evaluates the vector functor at the centroid.
- ADMaterialRankFourTensorAuxAccess a component of a RankFourTensor for automatic material property output
- ADMaterialRankTwoTensorAuxAccess a component of a RankTwoTensor for automatic material property output
- ADMaterialRateRealAuxOutputs element material properties rate of change
- ADMaterialRealAuxOutputs element volume-averaged material properties
- ADMaterialRealTensorValueAuxObject for extracting a component of a rank two tensor material property to populate an auxiliary variable.
- ADMaterialRealVectorValueAuxCapture a component of a vector material property in an auxiliary variable.
- ADMaterialStdVectorAuxExtracts a component of a material type std::vector<Real> to an aux variable. If the std::vector is not of sufficient size then zero is returned
- ADMaterialSymmetricRankFourTensorAuxAccess a component of a RankTwoTensor for automatic material property output
- ADMaterialSymmetricRankTwoTensorAuxCapture a component of a vector material property in an auxiliary variable.
- ADProjectedStatefulMaterialRankFourTensorAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ADProjectedStatefulMaterialRankTwoTensorAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ADProjectedStatefulMaterialRealAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ADProjectedStatefulMaterialRealVectorValueAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ADVectorMaterialRealVectorValueAuxConverts a vector-quantity material property into a vector auxiliary variable
- AdvectiveFluxAuxCompute components of flux vector for advection problems .
- ArrayParsedAuxSets field array variable values to the evaluation of a parsed expression.
- ArrayVarReductionAuxTakes an array variable and performs a reduction operation on it (max, min, sum, average) and stores as a standard variable.
- ArrayVariableComponentCopy a component of an array variable.
- BuildArrayVariableAuxCombines multiple standard variables into an array variable.
- ConstantAuxCreates a constant field in the domain.
- ContainsPointAuxComputes a binary field where the field is 1 in the elements that contain the point and 0 everywhere else
- CopyValueAuxReturns the specified variable as an auxiliary variable with a simple copy of the variable values.
- DebugResidualAuxPopulate an auxiliary variable with the residual contribution of a variable.
- DiffusionFluxAuxCompute components of flux vector for diffusion problems .
- DivergenceAuxComputes the divergence of a vector of functors.
- ElemExtraIDAuxPuts element extra IDs into an aux variable.
- ElementAdaptivityLevelAuxstores the element hierarchy in a aux variable
- ElementH1ErrorFunctionAuxComputes the H1 or W^{1,p} error between an exact function and a coupled variable.
- ElementIntegerAuxCreates a field showing the element integer.
- ElementL2ErrorFunctionAuxA class for computing the element-wise L^2 (Euclidean) error between a function and a coupled variable.
- ElementLengthAuxCompute the element size using Elem::hmin() or Elem::hmax() from libMesh.
- ElementLpNormAuxCompute an elemental field variable (single value per element) equal to the Lp-norm of a coupled Variable.
- ElementQualityAuxGenerates a field containing the quality metric for each element. Useful for visualizing mesh quality.
- ElementUOAuxAux Kernel to display generic spatial (elemental) information from a UserObject that satisfies the underlying ElementUOProvider interface.
- ExtraElementIDAuxPuts element extra IDs into an aux variable.
- ForcingFunctionAuxAuxiliary Kernel that adds a forcing function to the value of an AuxVariable from the previous time step.
- FunctionArrayAuxAuxiliary Kernel that creates and updates an array field variable by sampling functions through space and time.
- FunctionAuxAuxiliary Kernel that creates and updates a field variable by sampling a function through space and time.
- FunctorAuxEvaluates a functor (variable, function or functor material property) on the current element, quadrature point, or node.
- FunctorCoordinatesFunctionAuxAuxiliary Kernel that creates and updates a field variable by sampling a function with functors (variables, functions, others) as the coordinates.
- FunctorElementalAuxEvaluates a functor (variable, function or functor material property) on the current element, quadrature point, or node.
- FunctorElementalGradientAuxEvaluates the gradient of a functor (variable, function or functor material property) on the current element or quadrature point.
- FunctorMaterialRealAuxOutputs element volume-averaged material properties
- FunctorMaterialRealVectorValueAuxCapture a component of a vector material property in an auxiliary variable.
- FunctorVectorElementalAuxEvaluates a vector functor (material property usually) on the current element.For finite volume, this evaluates the vector functor at the centroid.
- GapValueAuxReturn the nearest value of a variable on a boundary from across a gap.
- GhostingAuxColors the elements ghosted to the chosen PID.
- GhostingFromUOAuxColors the elements ghosted to the chosen PID.
- HardwareIDAuxCreates a field showing the assignment of partitions to physical nodes in the cluster.
- InterfaceValueUserObjectAuxGet stored value from the specified InterfaceQpUserObjectBase.
- KokkosConstantAuxCreates a constant field in the domain.
- KokkosCopyValueAuxReturns the specified variable as an auxiliary variable with a simple copy of the variable values.
- KokkosFunctionAuxAuxiliary Kernel that creates and updates a field variable by sampling a function through space and time.
- KokkosMaterialRealAuxOutputs element volume-averaged material properties
- KokkosTagVectorAuxExtract DOF values from a tagged vector into an AuxVariable
- KokkosVariableTimeIntegrationAuxIntegrates a field variable in time.
- MFEMCrossProductAuxProjects s(x) * (U x V) onto a vector MFEM auxvariable
- MFEMCurlAuxCalculates the curl of an H(curl) conforming ND source variable and stores the result on an H(div) conforming RT result auxvariable
- MFEMDivAuxCalculates the divergence of an H(div) conforming RT source variable and stores the result on an L2 conforming result auxvariable
- MFEMDotProductAuxProject s(x) * (U . V) onto a scalar MFEM auxvariable.
- MFEMGradAuxCalculates the gradient of an H1 conforming source variable and stores the result on an H(curl) conforming ND result auxvariable
- MFEMScalarProjectionAuxProjects a scalar coefficient onto a scalar MFEM auxvariable
- MFEMScalarTimeAverageAuxCalculates a running time average of a scalar MFEM variable projected onto an auxvariable
- MFEMSumAuxCalculates the sum of an arbitrary number of variables sharing an FE space, each optionally scaled by a real constant, and stores the result in an auxiliary variable.
- MFEMVectorProjectionAuxProjects a vector coefficient onto a vector MFEM auxvariable.
- MaterialRankFourTensorAuxAccess a component of a RankFourTensor for automatic material property output
- MaterialRankTwoTensorAuxAccess a component of a RankTwoTensor for automatic material property output
- MaterialRateRealAuxOutputs element material properties rate of change
- MaterialRealAuxOutputs element volume-averaged material properties
- MaterialRealDenseMatrixAuxPopulate an auxiliary variable with an entry from a dense matrix material property.
- MaterialRealTensorValueAuxObject for extracting a component of a rank two tensor material property to populate an auxiliary variable.
- MaterialRealVectorValueAuxCapture a component of a vector material property in an auxiliary variable.
- MaterialStdVectorAuxExtracts a component of a material type std::vector<Real> to an aux variable. If the std::vector is not of sufficient size then zero is returned
- MaterialStdVectorRealGradientAuxExtracts a component of a material's std::vector<RealGradient> to an aux variable. If the std::vector is not of sufficient size then zero is returned
- MaterialSymmetricRankFourTensorAuxAccess a component of a RankTwoTensor for automatic material property output
- MaterialSymmetricRankTwoTensorAuxCapture a component of a vector material property in an auxiliary variable.
- MeshDivisionAuxReturns the value of the mesh division index for each element / node
- NearestNodeDistanceAuxStores the distance between a block and boundary or between two boundaries.
- NearestNodeValueAuxRetrieves a field value from the closest node on the paired boundary and stores it on this boundary or block.
- NodalPatchRecoveryAuxThis Auxkernel solves a least squares problem at each node to fit a value from quantities defined on quadrature points.
- NormalizationAuxNormalizes a variable based on a Postprocessor value.
- ParsedAuxSets a field variable value to the evaluation of a parsed expression.
- ParsedVectorAuxSets a field vector variable value to the evaluation of a parsed expression.
- PenetrationAuxAuxiliary Kernel for computing several geometry related quantities between two contacting bodies.
- ProcessorIDAuxCreates a field showing the processors and partitioning.
- ProjectedMaterialPropertyNodalPatchRecoveryAuxThis Auxkernel solves a least squares problem at each node to fit a value from quantities defined on quadrature points.
- ProjectedStatefulMaterialRankFourTensorAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ProjectedStatefulMaterialRankTwoTensorAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ProjectedStatefulMaterialRealAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ProjectedStatefulMaterialRealVectorValueAuxPicks a component of an indexable material property to get projected on an elemental Auxvariable. For use by ProjectedStatefulMaterialStorageAction.
- ProjectionAuxReturns the specified variable as an auxiliary variable with a projection of the source variable. If they are the same type, this amounts to a simple copy.
- QuotientAuxDivides two coupled variables.
- SecondTimeDerivativeAuxReturns the second order time derivative of the specified variable as an auxiliary variable.
- SelfAuxReturns the specified variable as an auxiliary variable with a projection of the source variable. If they are the same type, this amounts to a simple copy.
- SolutionAuxCreates fields by using information from a SolutionUserObject.
- SpatialUserObjectAuxPopulates an auxiliary variable with a spatial value returned from a UserObject spatialValue method.
- TagMatrixAuxCouple the diagonal of a tag matrix, and return its nodal value
- TagVectorArrayVariableAuxCouple a tagged vector, and return its evaluations at degree of freedom indices corresponding to the coupled array variable.
- TagVectorArrayVariableValueAuxCouple a tagged vector, and return its array value.
- TagVectorAuxExtract DOF values from a tagged vector into an AuxVariable
- TimeDerivativeAuxReturns the time derivative of the specified variable/functor as an auxiliary variable.
- VariableGradientComponentCreates a field consisting of one component of the gradient of a coupled variable.
- VariableTimeIntegrationAuxIntegrates a field variable in time.
- VectorFunctionAuxAuxiliary Kernel that creates and updates a vector field variable by sampling a Function object, via the vectorValue method, through space and time.
- VectorMagnitudeAuxCreates a field representing the magnitude of three coupled variables using an Euclidean norm.
- VectorMaterialRealVectorValueAuxConverts a vector-quantity material property into a vector auxiliary variable
- VectorPostprocessorVisualizationAuxRead values from a VectorPostprocessor that is producing vectors that are 'number of processors' * in length. Puts the value for each processor into an elemental auxiliary field.
- VectorVariableComponentAuxCreates a field consisting of one component of a coupled vector variable.
- VectorVariableMagnitudeAuxCreates a field consisting of the magnitude of a coupled vector variable.
- VolumeAuxAuxiliary Kernel that samples volumes.
- WeightedGapAuxReturns the specified variable as an auxiliary variable with the same value.
- Solid Mechanics App
- ADKineticEnergyAuxCompute the kinetic energy of continuum-based finite elements
- ADRankFourAuxAccess a component of a RankFourTensor
- ADRankTwoAuxAccess a component of a RankTwoTensor
- ADRankTwoScalarAuxCompute a scalar property of a RankTwoTensor
- AccumulateAuxAccumulates one or more variables and this auxiliary variable into this auxiliary variable
- CylindricalRankTwoAuxTakes RankTwoTensor material and outputs component in cylindrical coordinates
- DomainIntegralQFunctionComputes the q-function for a segment along the crack front, used in the calculation of the J-integral
- DomainIntegralTopologicalQFunctionDetermines if a node is within the ring of the crack front defintion; this object is normally created by the DomainIntegralAction.
- ElasticEnergyAuxCompute the local elastic energy
- GlobalDisplacementAuxAuxKernel to visualize the displacements generated by the global strain tensor
- KineticEnergyAuxCompute the kinetic energy of continuum-based finite elements
- NewmarkAccelAuxComputes the current acceleration using the Newmark method.
- NewmarkVelAuxCalculates the current velocity using Newmark method.
- RadialDisplacementCylinderAuxCompute the radial component of the displacement vector for cylindrical models.
- RadialDisplacementSphereAuxCompute the radial component of the displacement vector for spherical models.
- RankFourAuxAccess a component of a RankFourTensor
- RankTwoAuxAccess a component of a RankTwoTensor
- RankTwoScalarAuxCompute a scalar property of a RankTwoTensor
- ReactionForceAuxExtract the value of the residual from an appropriately formed tag vector and save those values as reaction forces in an AuxVariable
- RotationAngleCompute the field of angular rotations of points around an axis defined by an origin point and a direction vector
- ShellLocalCoordinatesAuxThis AuxKernel stores a specific component of a shell element's local coordinate vector in an auxiliary variable.
- ShellResultantsAuxComputes the local forces, bending moments and shear forces acting on shell elements
- TestNewmarkTIAssigns the velocity/acceleration calculated by time integrator to the velocity/acceleration auxvariable.
- Geochemistry App
- GeochemistryQuantityAuxExtracts information from the Reactor and records it in the AuxVariable
- NodalVoidVolumeAuxExtracts information from the NodalVoidVolume UserObject and records it in the AuxVariable
- Fsi App
- WaveHeightAuxKernelCalculates the wave heights given pressures.
- Electromagnetics App
- ADCurrentDensityCalculates the current density vector field (in A/m^2) when given electrostatic potential (electrostatic = true, default) or electric field.
- AzimuthMagneticTimeDerivRZComputes the time derivative of the azimuthal component of the magnetic field assuming cylindrical electric field. The electric field can be supplied as a vector or scalar components.
- CurrentDensityCalculates the current density vector field (in A/m^2) when given electrostatic potential (electrostatic = true, default) or electric field.
- EMJouleHeatingHeatGeneratedAuxComputes the heating due to the electic field in the form of
- PotentialToFieldAuxAn AuxKernel that calculates the electrostatic electric field given the electrostatic potential.
- SourceCurrentHeatingComputes the heating due to the electic field in the form of
- Thermal Hydraulics App
- ADConvectiveHeatFlux1PhaseAuxComputes convective heat flux for 1-phase flow.
- ADVectorVelocityComponentAuxComputes the velocity from the 1D phase-fraction and area weighted momentum and density variables.
- ConvectiveHeatFlux1PhaseAuxComputes convective heat flux for 1-phase flow.
- FlowModelGasMixAuxComputes various quantities for FlowModelGasMix.
- MachNumberAuxComputes Mach number.
- PrandtlNumberAuxComputes the Prandtl number for the fluid in the simulation domain
- ReynoldsNumberAuxComputes the Reynolds number.
- ShaftConnectedCompressor1PhaseAuxComputes various quantities for a ShaftConnectedCompressor1Phase.
- ShaftConnectedPump1PhaseAuxComputes various quantities for a ShaftConnectedPump1Phase.
- ShaftConnectedTurbine1PhaseAuxComputes various quantities for a ShaftConnectedTurbine1Phase.
- SimpleTurbinePowerFieldAuxComputes turbine power for 1-phase flow for a simple on/off turbine
- SoundSpeedAuxComputes the speed of sound.
- SpecificTotalEnthalpyAuxCompute the specific total enthalpy
- SumAuxSum of nonlinear or auxiliary variables
- THMSpecificInternalEnergyAuxComputed the specific internal energy.
- THMSpecificVolumeAuxComputes the specific volume for the phase.
- VariableValueTransferAuxRetrieves a field value from the closest node on the paired boundary and stores it on this boundary or block.
- VectorVelocityComponentAuxComputes the component of a vector-valued velocity field given by its magnitude and direction.
- VolumeJunction1PhaseAuxComputes various quantities for a VolumeJunction1Phase.
- WeightedAverageAuxWeighted average of variables using other variables as weights
- Sub Channel App
- BlockedMassFlowRateAuxComputes inlet mass flow rate BCs, from specified mass flux and cross-sectional area and applies blocked inlet conditions
- FlatMassFlowRateAuxComputes a uniform mass flow rate at the inlet
- MassFlowRateAuxComputes mass flow rate from specified mass flux and subchannel cross-sectional area. Can read either PostprocessorValue or Real
- QuadPowerAuxComputes axial power rate (W/m) that goes into the subchannel cells or is assigned to the fuel pins, in a quadrilateral lattice arrangement
- RZPinQPrimeAuxAxial heat rate on pin surface for a 2D-RZ axi-symmetric fuel pin model
- SCMBlockedMassFlowRateAuxComputes inlet mass flow rate BCs, from specified mass flux and cross-sectional area and applies blocked inlet conditions
- SCMFlatMassFlowRateAuxComputes a uniform mass flow rate at the inlet
- SCMMassFlowRateAuxComputes mass flow rate from specified mass flux and subchannel cross-sectional area. Can read either PostprocessorValue or Real
- SCMQuadPowerAuxComputes axial power rate (W/m) that goes into the subchannel cells or is assigned to the fuel pins, in a quadrilateral lattice arrangement
- SCMRZPinQPrimeAuxAxial heat rate on pin surface for a 2D-RZ axi-symmetric fuel pin model
- SCMTriPowerAuxComputes axial power rate (W/m) that goes into the subchannel cells or is assigned to the fuel pins, in a triangular lattice arrangement
- TriPowerAuxComputes axial power rate (W/m) that goes into the subchannel cells or is assigned to the fuel pins, in a triangular lattice arrangement
- Chemical Reactions App
- AqueousEquilibriumRxnAuxConcentration of secondary equilibrium species
- EquilibriumConstantAuxEquilibrium constant for a given equilibrium species (in form log10(Keq))
- KineticDisPreConcAuxConcentration of secondary kinetic species
- KineticDisPreRateAuxKinetic rate of secondary kinetic species
- PHAuxpH of solution
- TotalConcentrationAuxTotal concentration of primary species (including stoichiometric contribution to secondary equilibrium species)
- Navier Stokes App
- CourantComputes |u| dt / h_min.
- EnthalpyAuxThis AuxKernel computes the specific enthalpy of the fluidfrom the total energy and the pressure.
- HasPorosityJumpFaceShows whether an element has any attached porosity jump faces
- INSCourantComputes h_min / |u|.
- INSFVMixingLengthTurbulentViscosityAuxComputes the turbulent viscosity for the mixing length model.
- INSQCriterionAuxThis class computes the Q criterion, a scalar whichaids in vortex identification in turbulent flows
- INSStressComponentAuxThis class computes the stress component based on pressure and velocity for incompressible Navier-Stokes
- InternalEnergyAuxThis AuxKernel computes the internal energy based on the equation of state / fluid properties and the local pressure and density.
- NSInternalEnergyAuxAuxiliary kernel for computing the internal energy of the fluid.
- NSLiquidFractionAuxComputes liquid fraction given the temperature.
- NSMachAuxAuxiliary kernel for computing the Mach number assuming an ideal gas.
- NSPressureAuxNodal auxiliary variable, for computing pressure at the nodes.
- NSSpecificTotalEnthalpyAuxNodal auxiliary variable, for computing enthalpy at the nodes.
- NSTemperatureAuxTemperature is an auxiliary value computed from the total energy based on the FluidProperties.
- NSVelocityAuxVelocity auxiliary value.
- PecletNumberFunctorAuxComputes the Peclet number: u*L/alpha.
- RANSYPlusAuxCalculates non-dimensional wall distance (y+) value.
- ReynoldsNumberFunctorAuxComputes rho*u*L/mu.
- SpecificInternalEnergyAuxThis AuxKernel computes the specific internal energy based from the total and the kinetic energy.
- SpecificVolumeAuxThis auxkernel computes the specific volume of the fluid.
- TurbulentConductivityAuxCalculates the turbulent heat conductivity.
- WallDistanceMixingLengthAuxComputes the turbulent mixing length by assuming that it is proportional to the distance from the nearest wall. The mixinglength is capped at a distance proportional to inputted parameter delta.
- WallFunctionWallShearStressAuxCalculates the wall shear stress based on algebraic standard velocity wall functions.
- WallFunctionYPlusAuxCalculates y+ value according to the algebraic velocity standard wall function.
- kEpsilonViscosityAuxCalculates the turbulent viscosity according to the k-epsilon model.
- Stochastic Tools App
- SurrogateModelArrayAuxKernelSets a value of a variable based on a surrogate model.
- SurrogateModelAuxKernelSets a value of a variable based on a surrogate model.
- Fluid Properties App
- FluidDensityAuxComputes density from pressure and temperature
- PressureAuxComputes pressure given specific volume and specific internal energy
- SaturationTemperatureAuxComputes saturation temperature from pressure and 2-phase fluid properties object
- SpecificEnthalpyAuxComputes specific enthalpy from pressure and temperature
- StagnationPressureAuxComputes stagnation pressure from specific volume, specific internal energy, and velocity
- StagnationTemperatureAuxComputes stagnation temperature from specific volume, specific internal energy, and velocity
- TemperatureAuxComputes temperature given specific volume and specific internal energy
- Phase Field App
- BndsCalcAuxCalculate location of grain boundaries in a polycrystalline sample
- CrossTermGradientFreeEnergyFree energy contribution from the cross terms in ACMultiInterface
- DiscreteNucleationAuxProject the DiscreteNucleationMap state onto an AuxVariable
- EBSDReaderAvgDataAuxOutputs the requested EBSD reader average data, for a given phase if specified, for the grain at the local node/element.
- EBSDReaderPointDataAuxOutputs the requested EBSD reader point data.
- EulerAngleProvider2RGBAuxOutput RGB representation of crystal orientation from user object to an AuxVariable. The entire domain must have the same crystal structure.
- EulerAngleVariables2RGBAuxOutputs one color or a scalar for the RGB-encoding of the local Euler angles for the grain orientation
- FeatureFloodCountAuxFeature detection by connectivity analysis
- GrainAdvectionAuxCalculates the advection velocity of grain due to rigid body translation and rotation
- GrainBoundaryVelocityCompute the velocity of grain boundaries.
- KKSGlobalFreeEnergyTotal free energy in KKS system, including chemical, barrier and gradient terms
- KKSMultiFreeEnergyTotal free energy in multi-phase KKS system, including chemical, barrier and gradient terms
- LinearizedInterfaceAuxCalculates the order parameter from the linearized interface function
- OutputEulerAnglesOutput Euler angles from user object to an AuxVariable.
- PFCEnergyDensityComputes the crystal free energy density
- PFCRFFEnergyDensityComputes the crystal free energy density for the RFF form of the phase field crystal model
- SolutionAuxMisorientationBoundaryCalculate location of grain boundaries by using information from a SolutionUserObject.
- TotalFreeEnergyTotal free energy (both the bulk and gradient parts), where the bulk free energy has been defined in a material
- Functional Expansion Tools App
- FunctionSeriesToAuxAuxKernel to convert a functional expansion (Functions object, type = FunctionSeries) to an AuxVariable
- Peridynamics App
- BoundaryOffsetPDClass for output offset of PD boundary nodes compared to initial FE mesh
- NodalRankTwoPDClass for computing and outputing components and scalar quantities of nodal rank two strain and stress tensors for bond-based and ordinary state-based peridynamic models
- NodalVolumePDClass for output nodal area(2D) or nodal volume(3D)
- RankTwoBasedFailureCriteriaNOSPDClass for rank two tensor based failure criteria in non-ordinary state-based model
- StretchBasedFailureCriterionPDClass for bond stretch failure criterion in bond-based model and ordinary state-based model
- XFEMApp
- CutSubdomainIDAuxFill the elemental variable with CutSubdomainID
- MeshCutLevelSetAuxCalculates signed distance from interface defined by InterfaceMeshCutUserObject.
- XFEMCutPlaneAuxComputes the normal and origin of a cutting plane for each partial element.
- XFEMMarkerAuxIdentify the crack tip elements.
- XFEMVolFracAuxComputes the volume fraction of the physical material in each partial element.
- Contact App
- CohesiveZoneMortarUserObjectAuxPopulates an auxiliary variable with mortar cohesive zone model quantities.
- ContactPressureAuxComputes the contact pressure from the contact force and nodal area
- MortarArchardsLawAuxReturns the weighted gap velocity at a node. This quantity is useful for mortar contact, particularly when dual basis functions are used in contact mechanics
- MortarFrictionalPressureVectorAuxThis class creates an auxiliary vector for outputting the mortar frictional pressure vector.
- MortarFrictionalStateAuxThis class creates discrete states for nodes into frictional contact, including contact/no-contact and stick/slip.
- MortarPressureComponentAuxThis class transforms the Cartesian Lagrange multiplier vector to local coordinates and outputs each individual component along the normal or tangential direction.
- MortarUserObjectAuxPopulates an auxiliary variable with a contact quantities from penalty mortar contact.
- PenaltyMortarUserObjectAuxPopulates an auxiliary variable with a contact quantities from penalty mortar contact.
- WeightedGapVelAuxReturns the weighted gap velocity at a node. This quantity is useful for mortar contact, particularly when dual basis functions are used in contact mechanics
- Heat Transfer App
- JouleHeatingHeatGeneratedAuxCompute heat generated from Joule heating.
- Porous Flow App
- ADPorousFlowDarcyVelocityComponentDarcy velocity (in m3.s-1.m-2, or m.s-1) -(k_ij * krel /mu (nabla_j P - w_j)), where k_ij is the permeability tensor, krel is the relative permeability, mu is the fluid viscosity, P is the fluid pressure, and w_j is the fluid weight.
- ADPorousFlowPropertyAuxAuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.
- PorousFlowDarcyVelocityComponentDarcy velocity (in m3.s-1.m-2, or m.s-1) -(k_ij * krel /mu (nabla_j P - w_j)), where k_ij is the permeability tensor, krel is the relative permeability, mu is the fluid viscosity, P is the fluid pressure, and w_j is the fluid weight.
- PorousFlowDarcyVelocityComponentLowerDimensionalDarcy velocity on a lower-dimensional element embedded in a higher-dimensional mesh. Units m3.s-1.m-2, or m.s-1. Darcy velocity = -(k_ij * krel /(mu * a) (nabla_j P - w_j)), where k_ij is the permeability tensor, krel is the relative permeability, mu is the fluid viscosity, P is the fluid pressure, a is the fracture aperture and w_j is the fluid weight. The difference between this AuxKernel and PorousFlowDarcyVelocity is that this one projects gravity along the element's tangent direction. NOTE! For a meaningful answer, your permeability tensor must NOT contain terms that rotate tangential vectors to non-tangential vectors.
- PorousFlowElementLengthAuxKernel to compute the 'length' of elements along a given direction. A plane is constructed through the element's centroid, with normal equal to the direction given. The average of the distance of the nodal positions to this plane is the 'length' returned. The Variable for this AuxKernel must be an elemental Variable
- PorousFlowElementNormalAuxKernel to compute components of the element normal. This is mostly designed for 2D elements living in 3D space, however, the 1D-element and 3D-element cases are handled as special cases. The Variable for this AuxKernel must be an elemental Variable
- PorousFlowPropertyAuxAuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.
- Misc App
- CoupledDirectionalMeshHeightInterpolationScales a variable based on position relative to the model bounds in a specified direction