https://mooseframework.inl.gov
MooseVariableField.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "MooseTypes.h"
13 #include "MooseVariableFieldBase.h"
14 #include "SubProblem.h"
15 #include "MooseMesh.h"
16 #include "MooseVariableData.h"
17 #include "MooseFunctor.h"
18 #include "MeshChangedInterface.h"
19 
20 #include "libmesh/numeric_vector.h"
21 #include "libmesh/dof_map.h"
22 #include "libmesh/elem.h"
23 #include "libmesh/quadrature.h"
24 #include "libmesh/dense_vector.h"
25 #include "libmesh/dense_vector.h"
26 
39 template <typename OutputType>
41  public Moose::FunctorBase<typename Moose::ADType<OutputType>::type>,
43 {
44 public:
45  // type for gradient, second and divergence of template class OutputType
49 
50  // shortcut for types storing values on quadrature points
56 
57  // shape function type for the template class OutputType
59 
60  // type for gradient, second and divergence of shape functions of template class OutputType
64 
65  // shortcut for types storing shape function values on quadrature points
71 
72  // shortcut for types storing test function values on quadrature points
73  // Note: here we assume the types are the same as of shape functions.
79 
80  // DoF value type for the template class OutputType
83 
85 
86  virtual Moose::VarFieldType fieldType() const override;
87  virtual bool isArray() const override;
88  virtual bool isVector() const override;
89 
91 
92  virtual void setNodalValue(const OutputType & value, unsigned int idx = 0) = 0;
93 
95 
98  virtual void setDofValue(const OutputData & value, unsigned int index) = 0;
100 
104  virtual const ADTemplateVariableValue<OutputType> & adSln() const = 0;
105 
109  virtual const ADTemplateVariableValue<OutputType> & adSlnNeighbor() const = 0;
110 
114  virtual const ADTemplateVariableGradient<OutputType> & adGradSln() const = 0;
115 
119  virtual const ADTemplateVariableGradient<OutputType> & adGradSlnDot() const = 0;
120 
124  virtual const ADTemplateVariableCurl<OutputType> & adCurlSln() const = 0;
125 
129  virtual const ADTemplateVariableGradient<OutputType> & adGradSlnNeighbor() const = 0;
130 
135 
139  virtual const ADTemplateVariableCurl<OutputType> & adCurlSlnNeighbor() const = 0;
140 
144  virtual const ADTemplateVariableSecond<OutputType> & adSecondSln() const = 0;
145 
149  virtual const ADTemplateVariableSecond<OutputType> & adSecondSlnNeighbor() const = 0;
150 
154  virtual const ADTemplateVariableValue<OutputType> & adUDot() const = 0;
155 
159  virtual const ADTemplateVariableValue<OutputType> & adUDotDot() const = 0;
160 
164  virtual const ADTemplateVariableValue<OutputType> & adUDotNeighbor() const = 0;
165 
169  virtual const ADTemplateVariableValue<OutputType> & adUDotDotNeighbor() const = 0;
170 
174  virtual const MooseArray<ADReal> & adDofValues() const = 0;
175 
179  virtual const MooseArray<ADReal> & adDofValuesNeighbor() const = 0;
180 
184  virtual const MooseArray<ADReal> & adDofValuesDot() const = 0;
185 
187 
190  virtual const MooseArray<OutputType> & nodalValueArray() const = 0;
191  virtual const MooseArray<OutputType> & nodalValueOldArray() const = 0;
192  virtual const MooseArray<OutputType> & nodalValueOlderArray() const = 0;
194 
198  virtual const FieldVariableValue & sln() const = 0;
199 
203  virtual const FieldVariableValue & slnOld() const = 0;
204 
208  virtual const FieldVariableValue & slnNeighbor() const = 0;
209 
213  virtual const FieldVariableValue & slnOldNeighbor() const = 0;
214 
218  virtual const FieldVariableValue & slnOlder() const = 0;
219 
221  virtual const FieldVariableGradient & gradSln() const = 0;
222  virtual const FieldVariableGradient & gradSlnOld() const = 0;
223 
225  virtual const FieldVariableGradient & gradSlnNeighbor() const = 0;
226  virtual const FieldVariableGradient & gradSlnOldNeighbor() const = 0;
227 
231  virtual bool computingSecond() const = 0;
232 
236  virtual bool computingCurl() const = 0;
237 
241  virtual bool computingDiv() const = 0;
242 
246  virtual const FieldVariablePhiValue & phi() const = 0;
247 
251  virtual const FieldVariablePhiGradient & gradPhi() const = 0;
252 
256  virtual const FieldVariablePhiSecond & secondPhi() const = 0;
257 
261  virtual const FieldVariablePhiValue & curlPhi() const = 0;
262 
266  virtual const FieldVariablePhiDivergence & divPhi() const = 0;
267 
271  virtual const FieldVariablePhiValue & phiFace() const = 0;
272 
276  virtual const FieldVariablePhiGradient & gradPhiFace() const = 0;
277 
282  virtual const FieldVariablePhiSecond & secondPhiFace() const = 0;
283 
287  virtual const FieldVariablePhiValue & phiFaceNeighbor() const = 0;
288 
292  virtual const FieldVariablePhiGradient & gradPhiFaceNeighbor() const = 0;
293 
298  virtual const FieldVariablePhiSecond & secondPhiFaceNeighbor() const = 0;
299 
303  virtual const FieldVariablePhiValue & phiNeighbor() const = 0;
304 
308  virtual const FieldVariablePhiGradient & gradPhiNeighbor() const = 0;
309 
314  virtual const FieldVariablePhiSecond & secondPhiNeighbor() const = 0;
315 
319  virtual const FieldVariablePhiValue & phiLower() const = 0;
320 
324  virtual void setDofValues(const DenseVector<OutputData> & values) = 0;
325 
330  virtual void setLowerDofValues(const DenseVector<OutputData> & values) = 0;
331 
337  bool usesPhiNeighbor() const { return true; }
338 
344  bool usesGradPhiNeighbor() const { return true; }
345 
349  virtual bool usesSecondPhiNeighbor() const = 0;
350 
352 
355  virtual const DoFValue & dofValues() const = 0;
356  virtual const DoFValue & dofValuesOld() const = 0;
357  virtual const DoFValue & dofValuesOlder() const = 0;
358  virtual const DoFValue & dofValuesPreviousNL() const = 0;
359  virtual const DoFValue & dofValuesNeighbor() const = 0;
360  virtual const DoFValue & dofValuesOldNeighbor() const = 0;
361  virtual const DoFValue & dofValuesOlderNeighbor() const = 0;
362  virtual const DoFValue & dofValuesPreviousNLNeighbor() const = 0;
363  virtual const DoFValue & dofValuesDot() const = 0;
364  virtual const DoFValue & dofValuesDotNeighbor() const = 0;
365  virtual const DoFValue & dofValuesDotOld() const = 0;
366  virtual const DoFValue & dofValuesDotOldNeighbor() const = 0;
367  virtual const DoFValue & dofValuesDotDot() const = 0;
368  virtual const DoFValue & dofValuesDotDotNeighbor() const = 0;
369  virtual const DoFValue & dofValuesDotDotOld() const = 0;
370  virtual const DoFValue & dofValuesDotDotOldNeighbor() const = 0;
371  virtual const MooseArray<libMesh::Number> & dofValuesDuDotDu() const = 0;
372  virtual const MooseArray<libMesh::Number> & dofValuesDuDotDuNeighbor() const = 0;
373  virtual const MooseArray<libMesh::Number> & dofValuesDuDotDotDu() const = 0;
374  virtual const MooseArray<libMesh::Number> & dofValuesDuDotDotDuNeighbor() const = 0;
375 
376  template <bool is_ad>
378 
382  virtual const FieldVariableValue & vectorTagValue(TagID tag) const = 0;
383  virtual const DoFValue & nodalVectorTagValue(TagID tag) const = 0;
384  virtual const DoFValue & vectorTagDofValue(TagID tag) const = 0;
385  virtual const DoFValue & nodalMatrixTagValue(TagID tag) const = 0;
386  virtual const FieldVariableValue & matrixTagValue(TagID tag) const = 0;
387 
388  virtual void residualSetup() override;
389  virtual void jacobianSetup() override;
390  virtual void timestepSetup() override;
391 
393  /*
394  * Returns whether a variable is defined on a block as a functor.
395  * This makes the link between functor block restriction and the
396  * BlockRestrictable interface.
397  * @param id subdomain id we want to know whether the variable is defined on
398  * @return whether the variable is defined on this domain
399  */
400  bool hasBlocks(const SubdomainID id) const override { return BlockRestrictable::hasBlocks(id); }
401 
402 protected:
407 
410 
412  mutable ADReal _ad_real_dummy = 0;
413 };
414 
415 template <>
416 template <>
418 template <>
419 template <>
421 template <>
422 template <>
424 
425 template <typename OutputType>
426 template <bool is_ad>
429 {
430  return adDofValues();
431 }
432 
433 #define usingMooseVariableFieldMembers \
434  usingMooseVariableFieldBaseMembers; \
435  using MooseVariableField<OutputType>::_time_integrator; \
436  using MooseVariableField<OutputType>::_ad_real_dummy; \
437  using MooseVariableField<OutputType>::getSolution
438 
439 // Prevent implicit instantiation in other translation units where these classes are used
440 extern template class MooseVariableField<Real>;
441 extern template class MooseVariableField<RealVectorValue>;
442 extern template class MooseVariableField<RealEigenVector>;
MooseVariableField(const InputParameters &parameters)
typename OutputTools< typename Moose::ADType< T >::type >::VariableSecond ADTemplateVariableSecond
Definition: MooseTypes.h:609
virtual const DoFValue & dofValuesPreviousNL() const =0
VarFieldType
Definition: MooseTypes.h:721
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnNeighborDot() const =0
AD grad of time derivative neighbor solution getter.
virtual const DoFValue & dofValuesDotDotNeighbor() const =0
MooseArray< OutputDivergence > FieldVariableDivergence
virtual const MooseArray< ADReal > & adDofValuesNeighbor() const =0
Return the AD neighbor dof values.
typename OutputTools< typename Moose::ADType< T >::type >::VariableCurl ADTemplateVariableCurl
Definition: MooseTypes.h:611
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDuNeighbor() const =0
virtual const MooseArray< OutputType > & nodalValueOlderArray() const =0
virtual const DoFValue & dofValuesDotDotOldNeighbor() const =0
virtual const FieldVariablePhiValue & phi() const =0
Return the variable&#39;s elemental shape functions.
virtual const DoFValue & dofValuesDotNeighbor() const =0
virtual const ADTemplateVariableValue< OutputType > & adUDot() const =0
AD time derivative getter.
virtual const DoFValue & nodalMatrixTagValue(TagID tag) const =0
Base class template for functor objects.
Definition: MooseFunctor.h:137
virtual const ADTemplateVariableSecond< OutputType > & adSecondSln() const =0
AD second solution getter.
virtual const DoFValue & dofValuesDot() const =0
unsigned int TagID
Definition: MooseTypes.h:210
virtual const ADTemplateVariableValue< OutputType > & adSln() const =0
AD solution getter.
virtual bool computingDiv() const =0
Whether or not this variable is computing any divergence quantities.
virtual const ADTemplateVariableSecond< OutputType > & adSecondSlnNeighbor() const =0
AD second neighbor solution getter.
Class for stuff related to variables.
MooseArray< OutputSecond > FieldVariableSecond
ADReal _ad_real_dummy
A dummy ADReal variable.
virtual const DoFValue & dofValuesOlderNeighbor() const =0
OutputType type
Definition: MooseTypes.h:268
OutputType type
Definition: MooseTypes.h:257
virtual void setNodalValue(const OutputType &value, unsigned int idx=0)=0
virtual const ADTemplateVariableValue< OutputType > & adUDotDot() const =0
AD second time derivative getter.
virtual const DoFValue & vectorTagDofValue(TagID tag) const =0
virtual const DoFValue & dofValuesDotOldNeighbor() const =0
virtual const FieldVariablePhiGradient & gradPhiNeighbor() const =0
Return the gradients of the variable&#39;s shape functions on a neighboring element.
MooseArray< std::vector< OutputShapeSecond > > FieldVariablePhiSecond
typename OutputTools< typename Moose::ADType< T >::type >::VariableValue ADTemplateVariableValue
Definition: MooseTypes.h:603
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual const FieldVariablePhiDivergence & divPhi() const =0
Divergence of the shape functions.
virtual const FieldVariablePhiSecond & secondPhi() const =0
Return the rank-2 tensor of second derivatives of the variable&#39;s elemental shape functions.
virtual const ADTemplateVariableValue< OutputType > & adUDotDotNeighbor() const =0
AD neighbor second time derivative getter.
virtual const MooseArray< OutputType > & nodalValueOldArray() const =0
virtual const MooseArray< ADReal > & adDofValues() const =0
Return the AD dof values.
MooseArray< OutputType > FieldVariableValue
MooseArray< std::vector< OutputShape > > FieldVariableTestCurl
virtual const ADTemplateVariableCurl< OutputType > & adCurlSlnNeighbor() const =0
AD curl neighbor solution getter.
This class provides an interface for common operations on field variables of both FE and FV types wit...
virtual const ADTemplateVariableGradient< OutputType > & adGradSln() const =0
AD grad solution getter.
MooseArray< std::vector< OutputShapeDivergence > > FieldVariablePhiDivergence
virtual void setDofValue(const OutputData &value, unsigned int index)=0
Degree of freedom value setters.
MooseArray< std::vector< OutputShape > > FieldVariableTestValue
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnDot() const =0
AD grad of time derivative solution getter.
const MooseArray< GenericReal< is_ad > > & genericDofValues() const
virtual const MooseArray< ADReal > & adDofValuesDot() const =0
Return the AD time derivatives at dofs.
virtual const FieldVariablePhiValue & phiNeighbor() const =0
Return the variable&#39;s shape functions on a neighboring element.
virtual void jacobianSetup() override
virtual const ADTemplateVariableValue< OutputType > & adSlnNeighbor() const =0
AD neighbor solution getter.
MooseArray< std::vector< OutputShape > > FieldVariablePhiValue
MooseArray< std::vector< OutputShapeGradient > > FieldVariablePhiGradient
Moose::DOFType< OutputType >::type OutputData
MooseArray< std::vector< OutputShapeGradient > > FieldVariableTestGradient
const libMesh::NumericVector< libMesh::Number > & getSolution(const Moose::StateArg &state) const
Get the solution corresponding to the provided state.
static InputParameters validParams()
virtual const FieldVariablePhiValue & phiFaceNeighbor() const =0
Return the variable&#39;s shape functions on a neighboring element face.
MooseArray< OutputData > DoFValue
libMesh::TensorTools::DecrementRank< OutputShape >::type OutputShapeDivergence
virtual const FieldVariablePhiSecond & secondPhiFaceNeighbor() const =0
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on a neighboring ele...
MooseArray< std::vector< OutputShapeSecond > > FieldVariableTestSecond
virtual const FieldVariablePhiSecond & secondPhiFace() const =0
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on an element face...
virtual Moose::VarFieldType fieldType() const override
Field type of this variable.
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnNeighbor() const =0
AD grad neighbor solution getter.
Interface for notifications that the mesh has changed.
virtual const FieldVariableValue & vectorTagValue(TagID tag) const =0
tag values getters
virtual const DoFValue & nodalVectorTagValue(TagID tag) const =0
libMesh::TensorTools::IncrementRank< OutputGradient >::type OutputSecond
virtual void setLowerDofValues(const DenseVector< OutputData > &values)=0
Set local DOF values for a lower dimensional element and evaluate the values on quadrature points...
MooseArray< std::vector< OutputShape > > FieldVariablePhiCurl
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDu() const =0
typename OutputTools< typename Moose::ADType< T >::type >::VariableGradient ADTemplateVariableGradient
Definition: MooseTypes.h:606
virtual bool computingCurl() const =0
Whether or not this variable is computing any curl quantities.
virtual const DoFValue & dofValues() const =0
dof values getters
bool usesPhiNeighbor() const
Whether or not this variable is actually using the shape function value.
virtual bool computingSecond() const =0
Whether or not this variable is computing any second derivatives.
virtual const FieldVariableValue & slnOldNeighbor() const =0
virtual const FieldVariableGradient & gradSln() const =0
element gradients
virtual const FieldVariableGradient & gradSlnOld() const =0
libMesh::TensorTools::IncrementRank< OutputShapeGradient >::type OutputShapeSecond
virtual const FieldVariableGradient & gradSlnNeighbor() const =0
neighbor solution gradients
virtual const DoFValue & dofValuesDotDotOld() const =0
virtual const DoFValue & dofValuesOldNeighbor() const =0
virtual bool usesSecondPhiNeighbor() const =0
Whether or not this variable is actually using the shape function second derivatives.
virtual const FieldVariablePhiSecond & secondPhiNeighbor() const =0
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on a neighboring ele...
MooseArray< std::vector< OutputShapeDivergence > > FieldVariableTestDivergence
virtual const FieldVariablePhiGradient & gradPhiFaceNeighbor() const =0
Return the gradients of the variable&#39;s shape functions on a neighboring element face.
virtual const DoFValue & dofValuesOlder() const =0
virtual const FieldVariableValue & slnNeighbor() const =0
virtual const FieldVariablePhiValue & phiLower() const =0
Return the variable&#39;s shape functions on a lower-dimensional element.
virtual const FieldVariablePhiValue & phiFace() const =0
Return the variable&#39;s shape functions on an element face.
virtual const FieldVariablePhiGradient & gradPhi() const =0
Return the gradients of the variable&#39;s elemental shape functions.
virtual const FieldVariablePhiValue & curlPhi() const =0
Curl of the shape functions.
Base class for time integrators.
virtual bool isVector() const override
virtual const FieldVariableValue & sln() const =0
virtual void setDofValues(const DenseVector< OutputData > &values)=0
Set local DOF values and evaluate the values on quadrature points.
bool hasBlocks(const SubdomainID id) const override
Returns whether the functor is defined on this block.
virtual const FieldVariableValue & slnOlder() const =0
bool usesGradPhiNeighbor() const
Whether or not this variable is actually using the shape function gradient.
virtual const FieldVariableValue & matrixTagValue(TagID tag) const =0
virtual const FieldVariableValue & slnOld() const =0
const TimeIntegrator *const _time_integrator
the time integrator used for computing time derivatives
virtual const DoFValue & dofValuesDotOld() const =0
virtual const DoFValue & dofValuesNeighbor() const =0
virtual const MooseArray< OutputType > & nodalValueArray() const =0
Methods for retrieving values of variables at the nodes in a MooseArray for AuxKernelBase.
const InputParameters & parameters() const
Get the parameters of the object.
State argument for evaluating functors.
MooseArray< OutputType > FieldVariableCurl
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDotDu() const =0
libMesh::TensorTools::DecrementRank< OutputType >::type OutputDivergence
virtual void timestepSetup() override
libMesh::TensorTools::IncrementRank< OutputType >::type OutputGradient
bool hasBlocks(const SubdomainName &name) const
Test if the supplied block name is valid for this object.
virtual const ADTemplateVariableCurl< OutputType > & adCurlSln() const =0
AD curl solution getter.
virtual const DoFValue & dofValuesPreviousNLNeighbor() const =0
virtual const FieldVariableGradient & gradSlnOldNeighbor() const =0
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDotDuNeighbor() const =0
libMesh::TensorTools::IncrementRank< OutputShape >::type OutputShapeGradient
virtual const DoFValue & dofValuesDotDot() const =0
virtual const DoFValue & dofValuesOld() const =0
Moose::ShapeType< OutputType >::type OutputShape
virtual const FieldVariablePhiGradient & gradPhiFace() const =0
Return the gradients of the variable&#39;s shape functions on an element face.
virtual void residualSetup() override
virtual const ADTemplateVariableValue< OutputType > & adUDotNeighbor() const =0
AD neighbor time derivative getter.
MooseArray< OutputGradient > FieldVariableGradient
virtual bool isArray() const override