https://mooseframework.inl.gov
MooseLinearVariableFV.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 "MooseVariableField.h"
14 #include "SubProblem.h"
15 #include "MooseMesh.h"
17 
18 #include "libmesh/numeric_vector.h"
19 #include "libmesh/dof_map.h"
20 #include "libmesh/elem.h"
21 #include "libmesh/quadrature.h"
22 #include "libmesh/dense_vector.h"
23 #include "libmesh/enum_fe_family.h"
24 
25 template <typename>
27 
29 class FVDirichletBCBase;
30 class FVFluxBC;
32 
33 namespace libMesh
34 {
35 template <typename>
36 class NumericVector;
37 }
38 
42 template <typename OutputType>
43 class MooseLinearVariableFV : public MooseVariableField<OutputType>
44 {
45 public:
49 
55 
60 
63 
80 
83 
84  virtual bool isFV() const override { return true; }
85 
89  virtual bool isDirichletBoundaryFace(const FaceInfo & fi) const;
90 
95 
99  virtual bool needsGradientVectorStorage() const override { return _needs_cell_gradients; }
100 
101  virtual bool isExtrapolatedBoundaryFace(const FaceInfo & fi,
102  const Elem * elem,
103  const Moose::StateArg & state) const override;
104 
109  const VectorValue<Real> gradSln(const ElemInfo & elem_info) const;
110 
117  VectorValue<Real> gradSln(const FaceInfo & fi, const StateArg & state) const;
118 
119  virtual void initialSetup() override;
120 
128  Real getElemValue(const ElemInfo & elem_info, const StateArg & state) const;
129 
135 
136  const std::unordered_map<BoundaryID, LinearFVBoundaryCondition *> & getBoundaryConditionMap()
137  {
138  return _boundary_id_to_bc;
139  }
140 
141  virtual void prepareIC() override {}
142 
143  virtual bool isNodal() const override final { return false; }
144 
145  virtual bool hasDoFsOnNodes() const override final { return false; }
146 
147  virtual bool isNodalDefined() const override final { return false; }
148 
149  virtual bool supportsFaceArg() const override final { return true; }
150  virtual bool supportsElemSideQpArg() const override final { return false; }
151 
152  virtual const Elem * const & currentElem() const override;
153 
154  virtual bool computingSecond() const override final { return false; }
155  virtual bool computingCurl() const override final { return false; }
156  virtual bool computingDiv() const override final { return false; }
157  virtual bool usesSecondPhiNeighbor() const override final { return false; }
158 
159 protected:
161  [[noreturn]] void timeIntegratorError() const;
162 
164  [[noreturn]] void lowerDError() const;
165 
167  [[noreturn]] void nodalError() const;
168 
170  [[noreturn]] void adError() const;
171 
175  void cacheBoundaryBCMap();
176 
178 
181 
184 
186  const std::vector<std::unique_ptr<libMesh::NumericVector<libMesh::Number>>> & _grad_container;
187 
191  std::unique_ptr<MooseVariableDataLinearFV<OutputType>> _element_data;
192 
196  std::unique_ptr<MooseVariableDataLinearFV<OutputType>> _neighbor_data;
197 
200  std::unordered_map<BoundaryID, LinearFVBoundaryCondition *> _boundary_id_to_bc;
201 
203  const unsigned int _sys_num;
204 
205  friend void Moose::initDofIndices<>(MooseLinearVariableFV<OutputType> &, const Elem &);
206 
207 private:
211 
212  virtual ValueType evaluate(const ElemArg & elem, const StateArg &) const override final;
213  virtual ValueType evaluate(const FaceArg & face, const StateArg &) const override final;
214  virtual ValueType evaluate(const NodeArg & node, const StateArg &) const override final;
215  virtual ValueType evaluate(const ElemPointArg & elem_point,
216  const StateArg & state) const override final;
217  virtual ValueType evaluate(const ElemQpArg & elem_qp,
218  const StateArg & state) const override final;
219  virtual ValueType evaluate(const ElemSideQpArg & elem_side_qp,
220  const StateArg & state) const override final;
221  virtual GradientType evaluateGradient(const ElemQpArg & qp_arg,
222  const StateArg &) const override final;
223  virtual GradientType evaluateGradient(const ElemArg & elem_arg,
224  const StateArg &) const override final;
225  virtual GradientType evaluateGradient(const FaceArg & face,
226  const StateArg &) const override final;
227  virtual DotType evaluateDot(const ElemArg & elem, const StateArg &) const override final;
228 
233 
244 
245 public:
246  // *********************************************************************************
247  // *********************************************************************************
248  // The following functions are separated here because they are not essential for the
249  // solver but are necessary to interface with the auxiliary and postprocessor
250  // systems.
251  // *********************************************************************************
252  // *********************************************************************************
253 
254  virtual void setDofValue(const OutputData & /*value*/, unsigned int /*index*/) override;
255 
256  virtual void getDofIndices(const Elem * elem,
257  std::vector<dof_id_type> & dof_indices) const override;
258 
259  virtual void setDofValues(const DenseVector<OutputData> & values) override;
260 
261  virtual void clearDofIndices() override;
262 
263  virtual unsigned int numberOfDofs() const override final { return 1; }
264  virtual unsigned int numberOfDofsNeighbor() override final { return 1; }
265 
266  virtual unsigned int oldestSolutionStateRequested() const override final;
267 
268  virtual void clearAllDofIndices() override final;
269 
270  [[noreturn]] virtual const std::vector<dof_id_type> & dofIndicesLower() const override final;
271  [[noreturn]] virtual const FieldVariablePhiValue & phiLower() const override;
272 
273  // Overriding these to make sure nothing happens during residual/jacobian setup.
274  // The only time this can actually happen is when residual setup is called on the auxiliary
275  // system.
276  virtual void residualSetup() override {}
277  virtual void jacobianSetup() override {}
278 
279  virtual libMesh::FEContinuity getContinuity() const override
280  {
281  return _element_data->getContinuity();
282  };
283 
284  virtual void setNodalValue(const OutputType & value, unsigned int idx = 0) override;
285 
286  [[noreturn]] virtual const DoFValue & nodalVectorTagValue(TagID) const override;
287 
288  virtual const std::vector<dof_id_type> & dofIndices() const final;
289  virtual const std::vector<dof_id_type> & dofIndicesNeighbor() const final;
290 
291  virtual void prepare() override final {}
292  virtual void prepareNeighbor() override final {}
293  virtual void prepareAux() override final {}
294  virtual void reinitNode() override final {}
295  virtual void reinitNodes(const std::vector<dof_id_type> & /*nodes*/) override final {}
296  virtual void reinitNodesNeighbor(const std::vector<dof_id_type> & /*nodes*/) override final {}
297  virtual void reinitAux() override final {}
298  virtual void reinitAuxNeighbor() override final {}
299  virtual void prepareLowerD() override final {}
300 
301  virtual void computeElemValuesFace() override;
302  virtual void computeNeighborValuesFace() override;
303  virtual void computeNeighborValues() override;
304  virtual void computeLowerDValues() override final;
305 
306  virtual void computeNodalNeighborValues() override final;
307  virtual void computeNodalValues() override final;
308 
309  virtual void computeElemValues() override;
310  virtual void computeFaceValues(const FaceInfo & /*fi*/) override {}
311 
312  virtual void setLowerDofValues(const DenseVector<OutputData> & values) override;
313 
314  virtual void insert(libMesh::NumericVector<libMesh::Number> & vector) override;
315  virtual void insertLower(libMesh::NumericVector<libMesh::Number> & vector) override;
316  virtual void add(libMesh::NumericVector<libMesh::Number> & vector) override;
317 
318  virtual void setActiveTags(const std::set<TagID> & vtags) override;
319 
320  [[noreturn]] virtual const MooseArray<OutputType> & nodalValueArray() const override;
321  [[noreturn]] virtual const MooseArray<OutputType> & nodalValueOldArray() const override;
322  [[noreturn]] virtual const MooseArray<OutputType> & nodalValueOlderArray() const override;
323 
324  virtual const FieldVariablePhiValue & phi() const override final { return _phi; }
325  virtual const FieldVariablePhiGradient & gradPhi() const override final { return _grad_phi; }
326  [[noreturn]] virtual const FieldVariablePhiSecond & secondPhi() const override final;
327  [[noreturn]] const FieldVariablePhiValue & curlPhi() const override final;
328  [[noreturn]] const FieldVariablePhiDivergence & divPhi() const override final;
329 
330  virtual const FieldVariablePhiValue & phiFace() const override final { return _phi_face; }
331  virtual const FieldVariablePhiGradient & gradPhiFace() const override final
332  {
333  return _grad_phi_face;
334  }
335  [[noreturn]] virtual const FieldVariablePhiSecond & secondPhiFace() const override final;
336 
337  virtual const FieldVariablePhiValue & phiFaceNeighbor() const override final
338  {
339  return _phi_face_neighbor;
340  }
341  virtual const FieldVariablePhiGradient & gradPhiFaceNeighbor() const override final
342  {
344  }
345  [[noreturn]] virtual const FieldVariablePhiSecond & secondPhiFaceNeighbor() const override final;
346 
347  virtual const FieldVariablePhiValue & phiNeighbor() const override final { return _phi_neighbor; }
348  virtual const FieldVariablePhiGradient & gradPhiNeighbor() const override final
349  {
350  return _grad_phi_neighbor;
351  }
352  [[noreturn]] virtual const FieldVariablePhiSecond & secondPhiNeighbor() const override final;
353 
354  virtual const FieldVariableValue & vectorTagValue(TagID tag) const override;
355  virtual const DoFValue & vectorTagDofValue(TagID tag) const override;
356  [[noreturn]] virtual const DoFValue & nodalMatrixTagValue(TagID tag) const override;
357  virtual const FieldVariableValue & matrixTagValue(TagID tag) const override;
358 
359  virtual const FieldVariableValue & sln() const override;
360  virtual const FieldVariableValue & slnOld() const override;
361  virtual const FieldVariableValue & slnOlder() const override;
362  virtual const FieldVariableGradient & gradSln() const override;
363  virtual const FieldVariableGradient & gradSlnOld() const override;
364  virtual const FieldVariableValue & slnNeighbor() const override;
365  virtual const FieldVariableValue & slnOldNeighbor() const override;
366  virtual const FieldVariableGradient & gradSlnNeighbor() const override;
367  virtual const FieldVariableGradient & gradSlnOldNeighbor() const override;
368 
369  [[noreturn]] virtual const ADTemplateVariableSecond<OutputType> & adSecondSln() const override;
370  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> & adUDot() const override;
371  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> & adUDotDot() const override;
372  [[noreturn]] virtual const ADTemplateVariableGradient<OutputType> & adGradSlnDot() const override;
373  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> & adSlnNeighbor() const override;
374  [[noreturn]] virtual const ADTemplateVariableGradient<OutputType> &
375  adGradSlnNeighbor() const override;
376  [[noreturn]] virtual const ADTemplateVariableSecond<OutputType> &
377  adSecondSlnNeighbor() const override;
378  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> & adUDotNeighbor() const override;
379  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> &
380  adUDotDotNeighbor() const override;
381  [[noreturn]] virtual const ADTemplateVariableGradient<OutputType> &
382  adGradSlnNeighborDot() const override;
383  [[noreturn]] virtual const ADTemplateVariableValue<OutputType> & adSln() const override;
384  [[noreturn]] virtual const ADTemplateVariableGradient<OutputType> & adGradSln() const override;
385  [[noreturn]] virtual const ADTemplateVariableCurl<OutputType> & adCurlSln() const override;
386  [[noreturn]] virtual const ADTemplateVariableCurl<OutputType> &
387  adCurlSlnNeighbor() const override;
388 
389  virtual const DoFValue & dofValues() const override;
390  virtual const DoFValue & dofValuesOld() const override;
391 
392  virtual const DoFValue & dofValuesOlder() const override;
393  virtual const DoFValue & dofValuesPreviousNL() const override;
394  virtual const DoFValue & dofValuesNeighbor() const override;
395  virtual const DoFValue & dofValuesOldNeighbor() const override;
396  virtual const DoFValue & dofValuesOlderNeighbor() const override;
397  virtual const DoFValue & dofValuesPreviousNLNeighbor() const override;
398  [[noreturn]] virtual const DoFValue & dofValuesDot() const override;
399  [[noreturn]] virtual const DoFValue & dofValuesDotNeighbor() const override;
400  [[noreturn]] virtual const DoFValue & dofValuesDotOld() const override;
401  [[noreturn]] virtual const DoFValue & dofValuesDotOldNeighbor() const override;
402  [[noreturn]] virtual const DoFValue & dofValuesDotDot() const override;
403  [[noreturn]] virtual const DoFValue & dofValuesDotDotNeighbor() const override;
404  [[noreturn]] virtual const DoFValue & dofValuesDotDotOld() const override;
405  [[noreturn]] virtual const DoFValue & dofValuesDotDotOldNeighbor() const override;
406  [[noreturn]] virtual const MooseArray<libMesh::Number> & dofValuesDuDotDu() const override;
407  [[noreturn]] virtual const MooseArray<libMesh::Number> &
408  dofValuesDuDotDuNeighbor() const override;
409  [[noreturn]] virtual const MooseArray<libMesh::Number> & dofValuesDuDotDotDu() const override;
410  [[noreturn]] virtual const MooseArray<libMesh::Number> &
411  dofValuesDuDotDotDuNeighbor() const override;
412 
413  [[noreturn]] virtual const MooseArray<ADReal> & adDofValues() const override;
414  [[noreturn]] virtual const MooseArray<ADReal> & adDofValuesNeighbor() const override;
415  [[noreturn]] virtual const MooseArray<ADReal> & adDofValuesDot() const override;
416  [[noreturn]] virtual const dof_id_type & nodalDofIndex() const override final;
417  [[noreturn]] virtual const dof_id_type & nodalDofIndexNeighbor() const override final;
418 
419  virtual std::size_t phiSize() const override final { return _phi.size(); }
420  virtual std::size_t phiFaceSize() const override final { return _phi_face.size(); }
421  virtual std::size_t phiNeighborSize() const override final { return _phi_neighbor.size(); }
422  virtual std::size_t phiFaceNeighborSize() const override final
423  {
424  return _phi_face_neighbor.size();
425  }
426  [[noreturn]] virtual std::size_t phiLowerSize() const override final;
427 };
428 
429 template <typename OutputType>
430 typename MooseLinearVariableFV<OutputType>::ValueType
431 MooseLinearVariableFV<OutputType>::evaluate(const ElemArg & elem_arg, const StateArg & state) const
432 {
433  const auto & elem_info = this->_mesh.elemInfo(elem_arg.elem->id());
434  return getElemValue(elem_info, state);
435 }
436 
437 template <typename OutputType>
440  const StateArg & state) const
441 {
442  const auto & elem_info = this->_mesh.elemInfo(elem_point.elem->id());
443  return getElemValue(elem_info, state);
444 }
445 
446 template <typename OutputType>
449 {
450  const auto & elem_info = this->_mesh.elemInfo(elem_qp.elem->id());
451  return getElemValue(elem_info, state);
452 }
453 
454 template <typename OutputType>
457  const StateArg & state) const
458 {
459  return (*this)(ElemPointArg{elem_side_qp.elem, elem_side_qp.point, false}, state);
460 }
461 
462 template <typename OutputType>
465  const StateArg & /*state*/) const
466 {
467  const auto & elem_info = this->_mesh.elemInfo(qp_arg.elem->id());
468  return gradSln(elem_info);
469 }
470 
471 template <typename OutputType>
474  const StateArg & /*state*/) const
475 {
476  const auto & elem_info = this->_mesh.elemInfo(elem_arg.elem->id());
477  return gradSln(elem_info);
478 }
479 
480 template <typename OutputType>
483  const StateArg & state) const
484 {
485  mooseAssert(face.fi, "We must have a non-null face information");
486  return gradSln(*face.fi, state);
487 }
488 
489 template <typename OutputType>
490 void
492 {
493  mooseError("MooseLinearVariableFV does not support time integration at the moment! The variable "
494  "which is causing the issue: ",
495  this->name());
496 }
497 
498 template <typename OutputType>
499 void
501 {
502  mooseError("Lower dimensional element support not implemented for finite volume variables!The "
503  "variable which is causing the issue: ",
504  this->name());
505 }
506 
507 template <typename OutputType>
508 void
510 {
511  mooseError("FV variables don't support nodal variable treatment! The variable which is causing "
512  "the issue: ",
513  this->name());
514 }
515 
516 template <typename OutputType>
517 void
519 {
520  mooseError("Linear FV variable does not support automatic differentiation, the variable which is "
521  "attempting it is: ",
522  this->name());
523 }
524 
525 // Declare all the specializations, as the template specialization declarations below must know
526 template <>
527 ADReal MooseLinearVariableFV<Real>::evaluateDot(const ElemArg & elem, const StateArg & state) const;
528 
529 // Prevent implicit instantiation in other translation units where these classes are used
530 extern template class MooseLinearVariableFV<Real>;
virtual const FieldVariablePhiGradient & gradPhiFace() const override final
Return the gradients of the variable&#39;s shape functions on an element face.
virtual const DoFValue & vectorTagDofValue(TagID tag) const override
virtual void clearAllDofIndices() override final
typename OutputTools< typename Moose::ADType< T >::type >::VariableSecond ADTemplateVariableSecond
Definition: MooseTypes.h:609
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDotDuNeighbor() const override
virtual void setActiveTags(const std::set< TagID > &vtags) override
Set the active vector tags.
typename OutputTools< typename Moose::ADType< T >::type >::VariableCurl ADTemplateVariableCurl
Definition: MooseTypes.h:611
MooseLinearVariableFV< Real > MooseLinearVariableFVReal
virtual bool isExtrapolatedBoundaryFace(const FaceInfo &fi, const Elem *elem, const Moose::StateArg &state) const override
Returns whether this (sided) face is an extrapolated boundary face for this functor.
virtual const DoFValue & dofValuesOldNeighbor() const override
virtual void prepare() override final
Prepare the elemental degrees of freedom.
virtual void computeElemValuesFace() override
Compute values at facial quadrature points.
Base class for boundary conditions for linear FV systems.
virtual const MooseArray< OutputType > & nodalValueOlderArray() const override
virtual void reinitAuxNeighbor() override final
virtual const ADTemplateVariableCurl< OutputType > & adCurlSln() const override
AD curl solution getter.
virtual DotType evaluateDot(const ElemArg &elem, const StateArg &) const override final
Evaluate the functor time derivative with a given element.
void adError() const
Throw an error when somebody wants to use this variable with automatic differentiation.
std::unordered_map< BoundaryID, LinearFVBoundaryCondition * > _boundary_id_to_bc
Map for easily accessing the boundary conditions based on the boundary IDs.
virtual const FieldVariableValue & sln() const override
virtual const DoFValue & dofValuesOlder() const override
virtual bool isDirichletBoundaryFace(const FaceInfo &fi) const
If the variable has a dirichlet boundary condition at face described by fi .
unsigned int TagID
Definition: MooseTypes.h:210
const FieldVariablePhiGradient & _grad_phi_face
Class for stuff related to variables.
const FieldVariablePhiGradient & _grad_phi_face_neighbor
virtual const MooseArray< OutputType > & nodalValueArray() const override
Methods for retrieving values of variables at the nodes in a MooseArray for AuxKernelBase.
const FieldVariablePhiGradient & _grad_phi_neighbor
virtual const Elem *const & currentElem() const override
Current element this variable is evaluated at.
virtual void setDofValues(const DenseVector< OutputData > &values) override
Set local DOF values and evaluate the values on quadrature points.
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnDot() const override
AD grad of time derivative solution getter.
virtual const DoFValue & dofValuesDotOld() const override
virtual const FieldVariablePhiSecond & secondPhiFace() const override final
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on an element face...
virtual bool isNodal() const override final
Is this variable nodal.
typename MooseVariableField< OutputType >::FieldVariablePhiSecond FieldVariablePhiSecond
const FieldVariablePhiDivergence & divPhi() const override final
Divergence of the shape functions.
virtual const FieldVariablePhiSecond & secondPhiFaceNeighbor() const override final
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on a neighboring ele...
virtual const ADTemplateVariableValue< OutputType > & adSln() const override
AD solution getter.
virtual const FieldVariableGradient & gradSlnOldNeighbor() const override
virtual void residualSetup() override
virtual const DoFValue & nodalMatrixTagValue(TagID tag) const override
Base class for finite volume Dirichlet boundaray conditions.
const unsigned int _sys_num
Cache the number of the system this variable belongs to.
virtual const DoFValue & dofValues() const override
dof values getters
virtual std::size_t phiFaceSize() const override final
Return phiFace size.
const std::string & name() const override
Get the variable name.
typename OutputTools< typename Moose::ADType< T >::type >::VariableValue ADTemplateVariableValue
Definition: MooseTypes.h:603
virtual const FieldVariablePhiValue & phiFace() const override final
Return the variable&#39;s shape functions on an element face.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual const MooseArray< ADReal > & adDofValues() const override
Return the AD dof values.
virtual void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices) const override
virtual const DoFValue & dofValuesPreviousNL() const override
virtual const FieldVariablePhiValue & phi() const override final
Return the variable&#39;s elemental shape functions.
virtual const FieldVariablePhiGradient & gradPhi() const override final
Return the gradients of the variable&#39;s elemental shape functions.
RealVectorValue _cell_gradient
Temporary storage for the cell gradient to avoid unnecessary allocations.
virtual void prepareNeighbor() override final
Prepare the neighbor element degrees of freedom.
virtual void setNodalValue(const OutputType &value, unsigned int idx=0) override
LinearFVBoundaryCondition * getBoundaryCondition(const BoundaryID bd_id) const
Get the boundary condition object which corresponds to the given boundary ID.
virtual const std::vector< dof_id_type > & dofIndicesNeighbor() const final
Get neighbor DOF indices for currently selected element.
virtual const DoFValue & dofValuesDotNeighbor() const override
virtual const FieldVariablePhiSecond & secondPhi() const override final
Return the rank-2 tensor of second derivatives of the variable&#39;s elemental shape functions.
A structure that is used to evaluate Moose functors at an arbitrary physical point contained within a...
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
virtual void initialSetup() override
Gets called at the beginning of the simulation before this object is asked to do its job...
const FieldVariablePhiValue & _phi_face_neighbor
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
virtual bool supportsElemSideQpArg() const override final
Whether this functor supports evaluation with ElemSideQpArg.
virtual void clearDofIndices() override
Clear out the dof indices.
virtual const std::vector< dof_id_type > & dofIndices() const final
Get local DoF indices.
virtual void computeElemValues() override
Compute values at interior quadrature points.
typename FunctorReturnType< Moose::ADType< OutputType >::type, FunctorEvaluationKind::Gradient >::type GradientType
This rigmarole makes it so that a user can create functors that return containers (std::vector...
Definition: MooseFunctor.h:149
virtual const ADTemplateVariableCurl< OutputType > & adCurlSlnNeighbor() const override
AD curl neighbor solution getter.
void lowerDError() const
Throw and error when somebody requests lower-dimensional data from this variable. ...
Moose::DOFType< Real >::type OutputData
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDuNeighbor() const override
virtual const FieldVariablePhiSecond & secondPhiNeighbor() const override final
Return the rank-2 tensor of second derivatives of the variable&#39;s shape functions on a neighboring ele...
virtual const ADTemplateVariableValue< OutputType > & adUDotDot() const override
AD second time derivative getter.
virtual void prepareLowerD() override final
Prepare a lower dimensional element&#39;s degrees of freedom.
libMesh::TensorTools::DecrementRank< OutputShape >::type OutputShapeDivergence
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnNeighbor() const override
AD grad neighbor solution getter.
This data structure is used to store geometric and variable related metadata about each cell face in ...
Definition: FaceInfo.h:36
virtual bool supportsFaceArg() const override final
Whether this functor supports evaluation with FaceArg.
virtual bool isFV() const override
virtual const FieldVariableValue & slnOld() const override
typename MooseVariableField< OutputType >::FieldVariablePhiValue FieldVariablePhiValue
void timeIntegratorError() const
Throw an error when somebody requests time-related data from this variable.
dof_id_type id() const
virtual std::size_t phiSize() const override final
Return phi size.
A structure defining a "face" evaluation calling argument for Moose functors.
const std::vector< std::unique_ptr< libMesh::NumericVector< libMesh::Number > > > & _grad_container
Pointer to the cell gradients which are stored on the linear system.
typename MooseVariableField< OutputType >::FieldVariablePhiGradient FieldVariablePhiGradient
virtual const MooseArray< ADReal > & adDofValuesNeighbor() const override
Return the AD neighbor dof values.
virtual void insertLower(libMesh::NumericVector< libMesh::Number > &vector) override
Insert the currently cached degree of freedom values for a lower-dimensional element into the provide...
const FaceInfo * fi
a face information object which defines our location in space
virtual const FieldVariablePhiValue & phiLower() const override
Return the variable&#39;s shape functions on a lower-dimensional element.
boundary_id_type BoundaryID
Point point
The physical location of the quadrature point.
virtual void reinitNode() override final
const FieldVariablePhiGradient & _grad_phi
virtual const FieldVariableValue & matrixTagValue(TagID tag) const override
const libMesh::Elem * elem
typename OutputTools< typename Moose::ADType< T >::type >::VariableGradient ADTemplateVariableGradient
Definition: MooseTypes.h:606
const libMesh::Elem * elem
virtual void reinitNodes(const std::vector< dof_id_type > &) override final
Provides an interface for computing residual contributions from finite volume numerical fluxes comput...
Definition: FVFluxBC.h:23
std::unique_ptr< MooseVariableDataLinearFV< OutputType > > _element_data
Holder for all the data associated with the "main" element.
virtual const ADTemplateVariableGradient< OutputType > & adGradSlnNeighborDot() const override
AD grad of time derivative neighbor solution getter.
virtual unsigned int numberOfDofsNeighbor() override final
A structure that is used to evaluate Moose functors logically at an element/cell center.
typename MooseVariableField< OutputType >::FieldVariableGradient FieldVariableGradient
MooseLinearVariableFV(const InputParameters &parameters)
Argument for requesting functor evaluation at a quadrature point location in an element.
virtual const dof_id_type & nodalDofIndex() const override final
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDu() const override
const libMesh::Elem * elem
The element.
virtual const ADTemplateVariableSecond< OutputType > & adSecondSln() const override
AD second solution getter.
const FieldVariablePhiValue & curlPhi() const override final
Curl of the shape functions.
virtual const DoFValue & dofValuesDotDotNeighbor() const override
void nodalError() const
Throw an error when somebody wants to use this variable as a nodal variable.
virtual void jacobianSetup() override
virtual void prepareAux() override final
virtual void setDofValue(const OutputData &, unsigned int) override
Degree of freedom value setters.
virtual void reinitNodesNeighbor(const std::vector< dof_id_type > &) override final
virtual bool hasDoFsOnNodes() const override final
Does this variable have DoFs on nodes.
virtual const std::vector< dof_id_type > & dofIndicesLower() const override final
Get dof indices for the current lower dimensional element (this is meaningful when performing mortar ...
virtual bool computingSecond() const override final
Whether or not this variable is computing any second derivatives.
virtual const FieldVariablePhiGradient & gradPhiFaceNeighbor() const override final
Return the gradients of the variable&#39;s shape functions on a neighboring element face.
virtual const FieldVariableValue & slnNeighbor() const override
virtual bool isNodalDefined() const override final
Is this variable defined at nodes.
virtual const FieldVariableGradient & gradSlnNeighbor() const override
neighbor solution gradients
virtual const FieldVariableGradient & gradSlnOld() const override
virtual const FieldVariablePhiValue & phiNeighbor() const override final
Return the variable&#39;s shape functions on a neighboring element.
virtual const DoFValue & dofValuesOlderNeighbor() const override
virtual const MooseArray< OutputType > & nodalValueOldArray() const override
MooseMesh & _mesh
mesh the variable is active in
virtual const DoFValue & dofValuesDotDotOldNeighbor() const override
virtual void computeNodalNeighborValues() override final
Compute nodal values of this variable in the neighbor.
virtual void reinitAux() override final
virtual void add(libMesh::NumericVector< libMesh::Number > &vector) override
Add the currently cached degree of freedom values into the provided vector.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void computeNeighborValuesFace() override
Compute values at facial quadrature points for the neighbor.
virtual const DoFValue & dofValuesNeighbor() const override
virtual bool usesSecondPhiNeighbor() const override final
Whether or not this variable is actually using the shape function second derivatives.
virtual const FieldVariableValue & vectorTagValue(TagID tag) const override
tag values getters
const FieldVariablePhiValue & _phi_face
std::unique_ptr< MooseVariableDataLinearFV< OutputType > > _neighbor_data
Holder for all the data associated with the "neighbor" element.
virtual void computeNodalValues() override final
Compute nodal values of this variable.
virtual void setLowerDofValues(const DenseVector< OutputData > &values) override
Set local DOF values for a lower dimensional element and evaluate the values on quadrature points...
const std::unordered_map< BoundaryID, LinearFVBoundaryCondition * > & getBoundaryConditionMap()
static InputParameters validParams()
virtual const DoFValue & dofValuesOld() const override
virtual std::size_t phiFaceNeighborSize() const override final
Return phiFaceNeighbor size.
virtual const dof_id_type & nodalDofIndexNeighbor() const override final
typename MooseVariableField< OutputType >::FieldVariablePhiDivergence FieldVariablePhiDivergence
Real getElemValue(const ElemInfo &elem_info, const StateArg &state) const
Get the solution value for the provided element and seed the derivative for the corresponding dof ind...
virtual void computeLowerDValues() override final
compute values at quadrature points on the lower dimensional element
void cacheBoundaryBCMap()
Setup the boundary to Dirichlet BC map.
virtual std::size_t phiLowerSize() const override final
Return the number of shape functions on the lower dimensional element for this variable.
virtual const DoFValue & dofValuesDot() const override
const FieldVariablePhiValue & _phi
Shape functions, only used when we are postprocessing or using this variable in an auxiliary system...
typename MooseVariableField< OutputType >::FieldVariableValue FieldVariableValue
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual const DoFValue & dofValuesPreviousNLNeighbor() const override
virtual const DoFValue & dofValuesDotOldNeighbor() const override
virtual void insert(libMesh::NumericVector< libMesh::Number > &vector) override
Insert the currently cached degree of freedom values into the provided vector.
virtual const ADTemplateVariableValue< OutputType > & adUDotNeighbor() const override
AD neighbor time derivative getter.
const InputParameters & parameters() const
Get the parameters of the object.
State argument for evaluating functors.
virtual std::size_t phiNeighborSize() const override final
Return phiNeighbor size.
virtual void computeFaceValues(const FaceInfo &) override
Compute values at face quadrature points for the element+neighbor (both sides of the face)...
virtual const ADTemplateVariableValue< OutputType > & adSlnNeighbor() const override
AD neighbor solution getter.
virtual void prepareIC() override
Prepare the initial condition.
virtual bool needsGradientVectorStorage() const override
Check if cell gradient computations were requested for this variable.
virtual bool computingCurl() const override final
Whether or not this variable is computing any curl quantities.
libMesh::TensorTools::DecrementRank< Real >::type OutputDivergence
virtual const DoFValue & dofValuesDotDot() const override
virtual const ADTemplateVariableSecond< OutputType > & adSecondSlnNeighbor() const override
AD second neighbor solution getter.
virtual void computeNeighborValues() override
Compute values at quadrature points for the neighbor.
Class used for caching additional information for elements such as the volume and centroid...
Definition: ElemInfo.h:25
virtual libMesh::FEContinuity getContinuity() const override
Return the continuity of this variable.
const libMesh::NumericVector< libMesh::Number > *const & _solution
The current (ghosted) solution.
Real Number
virtual const FieldVariablePhiGradient & gradPhiNeighbor() const override final
Return the gradients of the variable&#39;s shape functions on a neighboring element.
void computeCellGradients()
Switch to request cell gradient computations.
virtual const MooseArray< libMesh::Number > & dofValuesDuDotDotDu() const override
virtual unsigned int oldestSolutionStateRequested() const override final
The oldest solution state that is requested for this variable (0 = current, 1 = old, 2 = older, etc).
virtual const DoFValue & dofValuesDotDotOld() const override
const libMesh::Elem * elem
The element.
virtual const MooseArray< ADReal > & adDofValuesDot() const override
Return the AD time derivatives at dofs.
typename MooseVariableField< OutputType >::DoFValue DoFValue
const FieldVariablePhiValue & _phi_neighbor
This class provides variable solution interface for linear finite volume problems.
Definition: FVUtils.h:24
virtual unsigned int numberOfDofs() const override final
Get the number of local DoFs.
virtual const FieldVariableValue & slnOlder() const override
virtual const ADTemplateVariableValue< OutputType > & adUDot() const override
AD time derivative getter.
virtual bool computingDiv() const override final
Whether or not this variable is computing any divergence quantities.
virtual const FieldVariableGradient & gradSln() const override
element gradients
virtual const ADTemplateVariableGradient< OutputType > & adGradSln() const override
AD grad solution getter.
virtual GradientType evaluateGradient(const ElemQpArg &qp_arg, const StateArg &) const override final
Argument for requesting functor evaluation at quadrature point locations on an element side...
virtual const FieldVariablePhiValue & phiFaceNeighbor() const override final
Return the variable&#39;s shape functions on a neighboring element face.
const ElemInfo & elemInfo(const dof_id_type id) const
Accessor for the elemInfo object for a given element ID.
Definition: MooseMesh.C:3866
Moose::ShapeType< Real >::type OutputShape
virtual const DoFValue & nodalVectorTagValue(TagID) const override
uint8_t dof_id_type
bool _needs_cell_gradients
Boolean to check if this variable needs gradient computations.
virtual const FieldVariableValue & slnOldNeighbor() const override
virtual ValueType evaluate(const ElemArg &elem, const StateArg &) const override final
Evaluate the functor with a given element.
virtual const ADTemplateVariableValue< OutputType > & adUDotDotNeighbor() const override
AD neighbor second time derivative getter.