Line data Source code
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"
16 : #include "MooseVariableDataFV.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>
26 : class MooseVariableFV;
27 :
28 : typedef MooseVariableFV<Real> MooseVariableFVReal;
29 : class FVDirichletBCBase;
30 : class FVFluxBC;
31 :
32 : namespace libMesh
33 : {
34 : template <typename>
35 : class NumericVector;
36 : }
37 :
38 : /// This class provides variable solution values for other classes/objects to
39 : /// bind to when looping over faces or elements. It provides both
40 : /// elem and neighbor values when accessed in flux/face object calculations.
41 : /// The templating is to enable support for both vector and scalar variables.
42 : /// Gradient reconstruction (when enabled) occurs transparently within this
43 : /// class and kernels coupling to these values will naturally "see" according
44 : /// to the selected reconstruction methods.
45 : ///
46 : /// OutputType OutputShape DofValue
47 : /// ----------------------------------------------------
48 : /// Real Real Real
49 : /// RealVectorValue RealVectorValue Real
50 : /// RealEigenVector Real RealEigenVector
51 : template <typename OutputType>
52 : class MooseVariableFV : public MooseVariableField<OutputType>
53 : {
54 : public:
55 : using OutputGradient = typename MooseVariableField<OutputType>::OutputGradient;
56 : using OutputSecond = typename MooseVariableField<OutputType>::OutputSecond;
57 : using OutputDivergence = typename MooseVariableField<OutputType>::OutputDivergence;
58 :
59 : using FieldVariableValue = typename MooseVariableField<OutputType>::FieldVariableValue;
60 : using FieldVariableGradient = typename MooseVariableField<OutputType>::FieldVariableGradient;
61 : using FieldVariableSecond = typename MooseVariableField<OutputType>::FieldVariableSecond;
62 : using FieldVariableCurl = typename MooseVariableField<OutputType>::FieldVariableCurl;
63 : using FieldVariableDivergence = typename MooseVariableField<OutputType>::FieldVariableDivergence;
64 :
65 : using OutputShape = typename MooseVariableField<OutputType>::OutputShape;
66 : using OutputShapeGradient = typename MooseVariableField<OutputType>::OutputShapeGradient;
67 : using OutputShapeSecond = typename MooseVariableField<OutputType>::OutputShapeSecond;
68 : using OutputShapeDivergence = typename MooseVariableField<OutputType>::OutputShapeDivergence;
69 :
70 : using typename MooseVariableField<OutputType>::DofValue;
71 : using typename MooseVariableField<OutputType>::DofValues;
72 : using typename MooseVariableField<OutputType>::ADDofValue;
73 : using typename MooseVariableField<OutputType>::ADDofValues;
74 :
75 : using FieldVariablePhiValue = typename MooseVariableField<OutputType>::FieldVariablePhiValue;
76 : using FieldVariablePhiDivergence =
77 : typename MooseVariableField<OutputType>::FieldVariablePhiDivergence;
78 : using FieldVariablePhiGradient =
79 : typename MooseVariableField<OutputType>::FieldVariablePhiGradient;
80 : using FieldVariablePhiSecond = typename MooseVariableField<OutputType>::FieldVariablePhiSecond;
81 : using ElemQpArg = Moose::ElemQpArg;
82 : using ElemSideQpArg = Moose::ElemSideQpArg;
83 : using ElemArg = Moose::ElemArg;
84 : using FaceArg = Moose::FaceArg;
85 : using StateArg = Moose::StateArg;
86 : using NodeArg = Moose::NodeArg;
87 : using ElemPointArg = Moose::ElemPointArg;
88 : using typename MooseVariableField<OutputType>::ValueType;
89 : using typename MooseVariableField<OutputType>::DotType;
90 : using typename MooseVariableField<OutputType>::GradientType;
91 :
92 : static InputParameters validParams();
93 :
94 : MooseVariableFV(const InputParameters & parameters);
95 :
96 5027638 : virtual bool isFV() const override { return true; }
97 :
98 : // TODO: many of these functions are not relevant to FV variables but are
99 : // still called at various points from existing moose code paths. Ideally we
100 : // would figure out how to remove calls to these functions and then allow
101 : // throwing mooseError's from them instead of silently doing nothing (e.g.
102 : // reinitNodes, reinitAux, prepareLowerD, etc.).
103 :
104 15793038 : virtual void prepare() override final {}
105 147180 : virtual void prepareNeighbor() override final {}
106 : virtual void prepareAux() override final;
107 197716 : virtual void reinitNode() override final {}
108 0 : virtual void reinitNodes(const std::vector<dof_id_type> & /*nodes*/) override final {}
109 0 : virtual void reinitNodesNeighbor(const std::vector<dof_id_type> & /*nodes*/) override final {}
110 4362702 : virtual void reinitAux() override final {}
111 56664 : virtual void reinitAuxNeighbor() override final {}
112 3756 : virtual void prepareLowerD() override final {}
113 :
114 0 : virtual const dof_id_type & nodalDofIndex() const override final
115 : {
116 0 : mooseError("nodalDofIndex not supported by MooseVariableFVBase");
117 : }
118 0 : virtual const dof_id_type & nodalDofIndexNeighbor() const override final
119 : {
120 0 : mooseError("nodalDofIndexNeighbor not supported by MooseVariableFVBase");
121 : }
122 36636 : virtual std::size_t phiSize() const override final { return _phi.size(); }
123 0 : virtual std::size_t phiFaceSize() const override final { return _phi_face.size(); }
124 0 : virtual std::size_t phiNeighborSize() const override final { return _phi_neighbor.size(); }
125 0 : virtual std::size_t phiFaceNeighborSize() const override final
126 : {
127 0 : return _phi_face_neighbor.size();
128 : }
129 0 : virtual std::size_t phiLowerSize() const override final
130 : {
131 0 : mooseError("phiLowerSize not supported by MooseVariableFVBase");
132 : }
133 :
134 : virtual void computeElemValuesFace() override;
135 : virtual void computeNeighborValuesFace() override;
136 : virtual void computeNeighborValues() override;
137 3756 : virtual void computeLowerDValues() override final
138 : {
139 : // mooseError("computeLowerDValues not supported by MooseVariableFVBase");
140 3756 : }
141 0 : virtual void computeNodalNeighborValues() override final
142 : {
143 : // mooseError("computeNodalNeighborValues not supported by MooseVariableFVBase");
144 0 : }
145 0 : virtual void computeNodalValues() override final
146 : {
147 : // mooseError("computeNodalValues not supported by MooseVariableFVBase");
148 0 : }
149 : virtual const std::vector<dof_id_type> & dofIndicesLower() const override final;
150 :
151 177059 : unsigned int numberOfDofs() const override final { return _element_data->numberOfDofs(); }
152 :
153 0 : virtual unsigned int numberOfDofsNeighbor() override final
154 : {
155 0 : mooseError("numberOfDofsNeighbor not supported by MooseVariableFVBase");
156 : }
157 :
158 13693 : virtual bool isNodal() const override final { return false; }
159 :
160 0 : bool hasDoFsOnNodes() const override final { return false; }
161 :
162 344 : libMesh::FEContinuity getContinuity() const override final
163 : {
164 344 : return _element_data->getContinuity();
165 : };
166 :
167 197740 : virtual bool isNodalDefined() const override final { return false; }
168 :
169 : virtual void setNodalValue(const OutputType & value, unsigned int idx = 0) override;
170 :
171 : virtual void setDofValue(const DofValue & value, unsigned int index) override;
172 :
173 : void clearDofIndices() override;
174 :
175 : virtual void prepareIC() override;
176 :
177 664 : virtual const Elem * const & currentElem() const override { return _element_data->currentElem(); }
178 :
179 0 : virtual void getDofIndices(const Elem * elem,
180 : std::vector<dof_id_type> & dof_indices) const override
181 : {
182 0 : return _element_data->getDofIndices(elem, dof_indices);
183 : }
184 :
185 341442872 : virtual const std::vector<dof_id_type> & dofIndices() const final
186 : {
187 341442872 : return _element_data->dofIndices();
188 : }
189 298610055 : virtual const std::vector<dof_id_type> & dofIndicesNeighbor() const final
190 : {
191 298610055 : return _neighbor_data->dofIndices();
192 : }
193 :
194 0 : Moose::FV::InterpMethod faceInterpolationMethod() const { return _face_interp_method; }
195 :
196 : void clearAllDofIndices() final;
197 :
198 0 : const DofValues & nodalVectorTagValue(TagID) const override
199 : {
200 0 : mooseError("nodalVectorTagValue not implemented for finite volume variables.");
201 : }
202 0 : const DofValues & nodalMatrixTagValue(TagID) const override
203 : {
204 0 : mooseError("nodalMatrixTagValue not implemented for finite volume variables.");
205 : }
206 :
207 52 : const FieldVariableValue & vectorTagValue(TagID tag) const override
208 : {
209 52 : return _element_data->vectorTagValue(tag);
210 : }
211 52 : const DofValues & vectorTagDofValue(TagID tag) const override
212 : {
213 52 : return _element_data->vectorTagDofValue(tag);
214 : }
215 13 : const FieldVariableValue & matrixTagValue(TagID tag) const override
216 : {
217 13 : return _element_data->matrixTagValue(tag);
218 : }
219 0 : const FieldVariableValue & vectorTagValueNeighbor(TagID tag)
220 : {
221 0 : return _neighbor_data->vectorTagValue(tag);
222 : }
223 0 : const FieldVariableValue & matrixTagValueNeighbor(TagID tag)
224 : {
225 0 : return _neighbor_data->matrixTagValue(tag);
226 : }
227 :
228 0 : const FieldVariableValue & uDot() const { return _element_data->uDot(); }
229 6357 : const FieldVariableValue & sln() const override { return _element_data->sln(Moose::Current); }
230 0 : const FieldVariableValue & slnOld() const override { return _element_data->sln(Moose::Old); }
231 0 : const FieldVariableValue & slnOlder() const override { return _element_data->sln(Moose::Older); }
232 2035 : const FieldVariableGradient & gradSln() const override
233 : {
234 2035 : return _element_data->gradSln(Moose::Current);
235 : }
236 0 : const FieldVariableGradient & gradSlnOld() const override
237 : {
238 0 : return _element_data->gradSln(Moose::Old);
239 : }
240 0 : const FieldVariableValue & uDotNeighbor() const { return _neighbor_data->uDot(); }
241 78 : const FieldVariableValue & slnNeighbor() const override
242 : {
243 78 : return _neighbor_data->sln(Moose::Current);
244 : }
245 0 : const FieldVariableValue & slnOldNeighbor() const override
246 : {
247 0 : return _neighbor_data->sln(Moose::Old);
248 : }
249 26 : const FieldVariableGradient & gradSlnNeighbor() const override
250 : {
251 26 : return _neighbor_data->gradSln(Moose::Current);
252 : }
253 0 : const FieldVariableGradient & gradSlnOldNeighbor() const override
254 : {
255 0 : return _neighbor_data->gradSln(Moose::Old);
256 : }
257 :
258 0 : const VariableValue & duDotDu() const { return _element_data->duDotDu(); }
259 0 : const VariableValue & duDotDotDu() const { return _element_data->duDotDotDu(); }
260 0 : const VariableValue & duDotDuNeighbor() const { return _neighbor_data->duDotDu(); }
261 0 : const VariableValue & duDotDotDuNeighbor() const { return _neighbor_data->duDotDotDu(); }
262 :
263 : /// AD
264 8473 : const ADTemplateVariableValue<OutputType> & adSln() const override
265 : {
266 8473 : return _element_data->adSln();
267 : }
268 316 : const ADTemplateVariableGradient<OutputType> & adGradSln() const override
269 : {
270 316 : return _element_data->adGradSln();
271 : }
272 :
273 : /**
274 : * Retrieve (or potentially compute) the gradient on the provided element. Overriders of this
275 : * method *cannot* call \p getBoundaryFaceValue because that method itself may lead to a call to
276 : * \p adGradSln(const Elem * const) resulting in infinite recursion
277 : * @param elem The element for which to retrieve the gradient
278 : * @param state State argument which describes at what time / solution iteration state we want to
279 : * evaluate the variable
280 : * @param correct_skewness Whether to perform skew corrections
281 : * @return The gradient at the element centroid
282 : */
283 : virtual const VectorValue<ADReal> & adGradSln(const Elem * const elem,
284 : const StateArg & state,
285 : const bool correct_skewness = false) const;
286 :
287 : /**
288 : * Retrieve (or potentially compute) a cross-diffusion-corrected gradient on the provided face.
289 : * "Correcting" the face gradient involves weighting the gradient stencil more heavily on the
290 : * solution values on the face-neighbor cells than a linear interpolation between cell center
291 : * gradients does
292 : * @param face The face for which to retrieve the gradient.
293 : * @param state State argument which describes at what time / solution iteration state we want to
294 : * evaluate the variable
295 : * @param correct_skewness Whether to perform skew corrections
296 : */
297 : virtual VectorValue<ADReal>
298 : adGradSln(const FaceInfo & fi, const StateArg & state, const bool correct_skewness = false) const;
299 :
300 : /**
301 : * Retrieve (or potentially compute) the uncorrected gradient on the provided face. This
302 : * uncorrected gradient is a simple linear interpolation between cell gradients computed at the
303 : * centers of the two neighboring cells. "Correcting" the face gradient involves weighting the
304 : * gradient stencil more heavily on the solution values on the face-neighbor cells than the linear
305 : * interpolation process does. This is commonly known as a cross-diffusion correction. Correction
306 : * is done in \p adGradSln(const FaceInfo & fi)
307 : * @param face The face for which to retrieve the gradient
308 : * @param state State argument which describes at what time / solution iteration state we want to
309 : * evaluate the variable
310 : * @param correct_skewness Whether to perform skew corrections
311 : */
312 : virtual VectorValue<ADReal> uncorrectedAdGradSln(const FaceInfo & fi,
313 : const StateArg & state,
314 : const bool correct_skewness = false) const;
315 :
316 : /**
317 : * Retrieve the solution value at a boundary face. If we're using a one term Taylor series
318 : * expansion we'll simply return the ajacent element centroid value. If we're doing two terms,
319 : * then we will compute the gradient if necessary to help us interpolate from the element centroid
320 : * value to the face
321 : */
322 : ADReal getBoundaryFaceValue(const FaceInfo & fi,
323 : const StateArg & state,
324 : bool correct_skewness = false) const;
325 :
326 0 : const ADTemplateVariableSecond<OutputType> & adSecondSln() const override
327 : {
328 0 : return _element_data->adSecondSln();
329 : }
330 374 : const ADTemplateVariableValue<OutputType> & adUDot() const override
331 : {
332 374 : return _element_data->adUDot();
333 : }
334 13 : const ADTemplateVariableValue<OutputType> & adUDotDot() const override
335 : {
336 13 : return _element_data->adUDotDot();
337 : }
338 0 : const ADTemplateVariableGradient<OutputType> & adGradSlnDot() const override
339 : {
340 0 : return _element_data->adGradSlnDot();
341 : }
342 0 : const ADTemplateVariableCurl<OutputType> & adCurlSln() const override
343 : {
344 0 : mooseError("We don't currently implement curl for FV");
345 : }
346 :
347 : /// neighbor AD
348 5549 : const ADTemplateVariableValue<OutputType> & adSlnNeighbor() const override
349 : {
350 5549 : return _neighbor_data->adSln();
351 : }
352 290 : const ADTemplateVariableGradient<OutputType> & adGradSlnNeighbor() const override
353 : {
354 290 : return _neighbor_data->adGradSln();
355 : }
356 0 : const ADTemplateVariableSecond<OutputType> & adSecondSlnNeighbor() const override
357 : {
358 0 : return _neighbor_data->adSecondSln();
359 : }
360 0 : const ADTemplateVariableValue<OutputType> & adUDotNeighbor() const override
361 : {
362 0 : return _neighbor_data->adUDot();
363 : }
364 0 : const ADTemplateVariableValue<OutputType> & adUDotDotNeighbor() const override
365 : {
366 0 : return _neighbor_data->adUDotDot();
367 : }
368 0 : const ADTemplateVariableGradient<OutputType> & adGradSlnNeighborDot() const override
369 : {
370 0 : return _neighbor_data->adGradSlnDot();
371 : }
372 0 : const ADTemplateVariableCurl<OutputType> & adCurlSlnNeighbor() const override
373 : {
374 0 : mooseError("We don't currently implement curl for FV");
375 : }
376 :
377 : /// Initializes/computes variable values from the solution vectors for the
378 : /// current element being operated on in assembly. This
379 : /// must be called every time the current element (as defined by the assembly
380 : /// data structure changes in order for objects that bind to sln, adSln,
381 : /// uDot, etc. to be updated/correct.
382 : virtual void computeElemValues() override;
383 : /// Initializes/computes variable values from the solution vectors for the
384 : /// face represented by fi. This includes initializing data for elements on
385 : /// both sides of the face (elem and neighbor) and handles the case where
386 : /// there is no element on one side of the face gracefully - i.e. you call
387 : /// this uniformly for every face regardless of whether you are on a boundary
388 : /// or not - or whether the variable is only defined on one side of the face
389 : /// or not.
390 : virtual void computeFaceValues(const FaceInfo & fi) override;
391 :
392 : /**
393 : * Set local DOF values and evaluate the values on quadrature points
394 : */
395 : virtual void setDofValues(const DenseVector<DofValue> & values) override;
396 : virtual void setLowerDofValues(const DenseVector<DofValue> & values) override;
397 :
398 : /// Get the current value of this variable on an element
399 : /// @param[in] elem Element at which to get value
400 : /// @param[in] idx Local index of this variable's element DoFs
401 : /// @return Variable value
402 : DofValue getElementalValue(const Elem * elem, unsigned int idx = 0) const;
403 : /// Get the old value of this variable on an element
404 : /// @param[in] elem Element at which to get value
405 : /// @param[in] idx Local index of this variable's element DoFs
406 : /// @return Variable value
407 : DofValue getElementalValueOld(const Elem * elem, unsigned int idx = 0) const;
408 : /// Get the older value of this variable on an element
409 : /// @param[in] elem Element at which to get value
410 : /// @param[in] idx Local index of this variable's element DoFs
411 : /// @return Variable value
412 : DofValue getElementalValueOlder(const Elem * elem, unsigned int idx = 0) const;
413 :
414 : virtual void insert(libMesh::NumericVector<libMesh::Number> & vector) override;
415 : virtual void insertLower(libMesh::NumericVector<libMesh::Number> & vector) override;
416 : virtual void add(libMesh::NumericVector<libMesh::Number> & vector) override;
417 :
418 : const DofValues & dofValues() const override;
419 : const DofValues & dofValuesOld() const override;
420 : const DofValues & dofValuesOlder() const override;
421 : const DofValues & dofValuesPreviousNL() const override;
422 : const DofValues & dofValuesNeighbor() const override;
423 : const DofValues & dofValuesOldNeighbor() const override;
424 : const DofValues & dofValuesOlderNeighbor() const override;
425 : const DofValues & dofValuesPreviousNLNeighbor() const override;
426 : const DofValues & dofValuesDot() const override;
427 : const DofValues & dofValuesDotNeighbor() const override;
428 : const DofValues & dofValuesDotOld() const override;
429 : const DofValues & dofValuesDotOldNeighbor() const override;
430 : const DofValues & dofValuesDotDot() const override;
431 : const DofValues & dofValuesDotDotNeighbor() const override;
432 : const DofValues & dofValuesDotDotOld() const override;
433 : const DofValues & dofValuesDotDotOldNeighbor() const override;
434 : const MooseArray<libMesh::Number> & dofValuesDuDotDu() const override;
435 : const MooseArray<libMesh::Number> & dofValuesDuDotDuNeighbor() const override;
436 : const MooseArray<libMesh::Number> & dofValuesDuDotDotDu() const override;
437 : const MooseArray<libMesh::Number> & dofValuesDuDotDotDuNeighbor() const override;
438 :
439 : const ADDofValues & adDofValues() const override;
440 : const ADDofValues & adDofValuesNeighbor() const override;
441 : const ADDofValues & adDofValuesDot() const override;
442 :
443 : /// Note: const monomial is always the case - higher order solns are
444 : /// reconstructed - so this is simpler func than FE equivalent.
445 : OutputType getValue(const Elem * elem) const;
446 :
447 : /// Compute the variable gradient value at a point on an element
448 : /// @param elem The element we are computing on
449 : /// @param phi Evaluated shape functions at a point
450 : /// @return The variable gradient value
451 : typename OutputTools<OutputType>::OutputGradient getGradient(const Elem * elem) const;
452 :
453 : /// Returns true if a Dirichlet BC exists on the current face. This only
454 : /// works if the variable has been initialized on a face with
455 : /// computeFaceValues. Its return value is nonsense if initialized on a
456 : /// single element via computeElemValues.
457 0 : bool hasDirichletBC() const
458 : {
459 0 : return _element_data->hasDirichletBC() || _neighbor_data->hasDirichletBC();
460 : }
461 :
462 : std::pair<bool, const FVDirichletBCBase *> getDirichletBC(const FaceInfo & fi) const;
463 :
464 : std::pair<bool, std::vector<const FVFluxBC *>> getFluxBCs(const FaceInfo & fi) const;
465 :
466 : virtual void residualSetup() override;
467 : virtual void initialSetup() override;
468 : virtual void jacobianSetup() override;
469 : virtual void timestepSetup() override;
470 : virtual void meshChanged() override;
471 :
472 : /**
473 : * Get the solution value for the provided element and seed the derivative for the corresponding
474 : * dof index
475 : * @param elem The element to retrieive the solution value for
476 : * @param state State argument which describes at what time / solution iteration state we want to
477 : * evaluate the variable
478 : */
479 : ADReal getElemValue(const Elem * elem, const StateArg & state) const;
480 :
481 : void setActiveTags(const std::set<TagID> & vtags) override;
482 :
483 : /**
484 : * Request that quadrature point data be (pre)computed. Quadrature point data is (pre)computed by
485 : * default for this base class but derived variable classes may choose not to unless this API is
486 : * called
487 : */
488 10201 : virtual void requireQpComputations() const {}
489 :
490 : /**
491 : * Determine whether a specified face side is a Dirichlet boundary face. In the base
492 : * implementation we only inspect the face information object for whether there are Dirichlet
493 : * conditions. However, derived classes may allow discontinuities between + and - side face
494 : * values, e.g. one side may have a Dirichlet condition and the other side may perform
495 : * extrapolation to determine its value
496 : * @param fi The face information object
497 : * @param elem An element that can be used to indicate sidedness of the face
498 : * @param state The state at which to determine whether the face is a Dirichlet face or not
499 : * @return Whether the potentially sided (as indicated by \p elem) \p fi is a Dirichlet boundary
500 : * face for this variable
501 : */
502 : virtual bool isDirichletBoundaryFace(const FaceInfo & fi,
503 : const Elem * elem,
504 : const Moose::StateArg & state) const;
505 :
506 0 : bool supportsFaceArg() const override final { return true; }
507 0 : bool supportsElemSideQpArg() const override final { return true; }
508 :
509 : virtual void sizeMatrixTagData() override;
510 :
511 : protected:
512 : /**
513 : * Retrieves a Dirichlet boundary value for the provided face. Callers of this method should be
514 : * sure that \p isDirichletBoundaryFace returns true. In the base implementation we only inspect
515 : * the face information object for the Dirichlet value. However, derived
516 : * classes may allow discontinuities between + and - side face values, e.g. one side may have a
517 : * Dirichlet condition and the other side may perform extrapolation to determine its value. This
518 : * is the reason for the existence of the \p elem parameter, to indicate sidedness
519 : * @param fi The face information object
520 : * @param elem An element that can be used to indicate sidedness of the face
521 : * @param state State argument which describes at what time / solution iteration state we want to
522 : * evaluate the variable
523 : * @return The Dirichlet value on the boundary face associated with \p fi (and potentially \p
524 : * elem)
525 : */
526 : virtual ADReal getDirichletBoundaryFaceValue(const FaceInfo & fi,
527 : const Elem * elem,
528 : const Moose::StateArg & state) const;
529 :
530 : /**
531 : * Returns whether this is an extrapolated boundary face. An extrapolated boundary face is
532 : * boundary face for which is not a corresponding Dirichlet condition, e.g. we need to compute
533 : * some approximation for the boundary face value using the adjacent cell centroid information. In
534 : * the base implementation we only inspect the face information object for whether we should
535 : * perform extrapolation. However, derived classes may allow discontinuities between + and - side
536 : * face values, e.g. one side may have a Dirichlet condition and the other side may perform
537 : * extrapolation to determine its value
538 : * @param fi The face information object
539 : * @param elem An element that can be used to indicate sidedness of the face
540 : * @param state The state at which to determine whether the face is extrapolated or not
541 : * @return Whether the potentially sided (as indicated by \p elem) \p fi is an extrapolated
542 : * boundary face for this variable
543 : */
544 : bool isExtrapolatedBoundaryFace(const FaceInfo & fi,
545 : const Elem * elem,
546 : const Moose::StateArg & state) const override;
547 :
548 : private:
549 : using MooseVariableField<OutputType>::evaluate;
550 : using MooseVariableField<OutputType>::evaluateGradient;
551 : using MooseVariableField<OutputType>::evaluateDot;
552 :
553 : ValueType evaluate(const ElemArg & elem, const StateArg &) const override final;
554 : ValueType evaluate(const FaceArg & face, const StateArg &) const override final;
555 : ValueType evaluate(const NodeArg & node, const StateArg &) const override final;
556 : ValueType evaluate(const ElemPointArg & elem_point, const StateArg & state) const override final;
557 : ValueType evaluate(const ElemQpArg & elem_qp, const StateArg & state) const override final;
558 : ValueType evaluate(const ElemSideQpArg & elem_side_qp,
559 : const StateArg & state) const override final;
560 : GradientType evaluateGradient(const ElemQpArg & qp_arg, const StateArg &) const override final;
561 : GradientType evaluateGradient(const ElemArg & elem_arg, const StateArg &) const override final;
562 : GradientType evaluateGradient(const FaceArg & face, const StateArg &) const override final;
563 : DotType evaluateDot(const ElemArg & elem, const StateArg &) const override final;
564 : DotType evaluateDot(const FaceArg & face, const StateArg &) const override final;
565 : DotType evaluateDot(const ElemQpArg & elem_qp, const StateArg &) const override final;
566 :
567 : /**
568 : * Setup the boundary to Dirichlet BC map
569 : */
570 : void determineBoundaryToDirichletBCMap();
571 :
572 : /// Whether the boundary to Dirichlet cache map has been setup yet
573 : bool _dirichlet_map_setup = false;
574 :
575 : /**
576 : * Setup the boundary to Flux BC map
577 : */
578 : void determineBoundaryToFluxBCMap();
579 :
580 : /// Whether the boundary to fluxBC cache map has been setup yet
581 : bool _flux_map_setup = false;
582 :
583 : public:
584 0 : const MooseArray<OutputType> & nodalValueArray() const override
585 : {
586 0 : mooseError("Finite volume variables do not have defined values at nodes.");
587 : }
588 0 : const MooseArray<OutputType> & nodalValueOldArray() const override
589 : {
590 0 : mooseError("Finite volume variables do not have defined values at nodes.");
591 : }
592 0 : const MooseArray<OutputType> & nodalValueOlderArray() const override
593 : {
594 0 : mooseError("Finite volume variables do not have defined values at nodes.");
595 : }
596 :
597 51381 : bool computingSecond() const override final { return false; }
598 0 : bool computingCurl() const override final { return false; }
599 0 : bool computingDiv() const override final { return false; }
600 0 : bool usesSecondPhiNeighbor() const override final { return false; }
601 :
602 52236 : const FieldVariablePhiValue & phi() const override final { return _phi; }
603 51407 : const FieldVariablePhiGradient & gradPhi() const override final { return _grad_phi; }
604 0 : const FieldVariablePhiSecond & secondPhi() const override final
605 : {
606 0 : mooseError("We don't currently implement second derivatives for FV");
607 : }
608 0 : const FieldVariablePhiValue & curlPhi() const override final
609 : {
610 0 : mooseError("We don't currently implement curl for FV");
611 : }
612 0 : const FieldVariablePhiDivergence & divPhi() const override final
613 : {
614 0 : mooseError("We don't currently implement divergence for FV");
615 : }
616 :
617 39 : const FieldVariablePhiValue & phiFace() const override final { return _phi_face; }
618 36 : const FieldVariablePhiGradient & gradPhiFace() const override final { return _grad_phi_face; }
619 0 : const FieldVariablePhiSecond & secondPhiFace() const override final
620 : {
621 0 : mooseError("We don't currently implement second derivatives for FV");
622 : }
623 :
624 36 : const FieldVariablePhiValue & phiFaceNeighbor() const override final
625 : {
626 36 : return _phi_face_neighbor;
627 : }
628 36 : const FieldVariablePhiGradient & gradPhiFaceNeighbor() const override final
629 : {
630 36 : return _grad_phi_face_neighbor;
631 : }
632 0 : const FieldVariablePhiSecond & secondPhiFaceNeighbor() const override final
633 : {
634 0 : mooseError("We don't currently implement second derivatives for FV");
635 : }
636 :
637 0 : const FieldVariablePhiValue & phiNeighbor() const override final { return _phi_neighbor; }
638 0 : const FieldVariablePhiGradient & gradPhiNeighbor() const override final
639 : {
640 0 : return _grad_phi_neighbor;
641 : }
642 0 : const FieldVariablePhiSecond & secondPhiNeighbor() const override final
643 : {
644 0 : mooseError("We don't currently implement second derivatives for FV");
645 : }
646 :
647 : virtual const FieldVariablePhiValue & phiLower() const override;
648 :
649 : unsigned int oldestSolutionStateRequested() const override final;
650 :
651 : /**
652 : * Retrieves an extrapolated boundary value for the provided face. Callers of this method should
653 : * be sure that \p isExtrapolatedBoundaryFace returns true. In the base implementation we only
654 : * inspect the face information object for the extrapolated value. However, derived classes may
655 : * allow discontinuities between + and - side face values, e.g. one side may have a Dirichlet
656 : * condition and the other side may perform extrapolation to determine its value. This is the
657 : * reason for the existence of the \p elem parameter, to indicate sidedness
658 : * @param fi The face information object
659 : * @param two_term_expansion Whether to use the cell gradient in addition to the cell center value
660 : * to compute the extrapolated boundary face value. If this is false, then the cell center value
661 : * will be used
662 : * @param correct_skewness Whether to perform skew corrections. This is relevant when performing
663 : * two term expansions as the gradient evaluation may involve evaluating face values on internal
664 : * skewed faces
665 : * @param elem_side_to_extrapolate_from An element that can be used to indicate sidedness of the
666 : * face
667 : * @param state State argument which describes at what time / solution iteration state we want to
668 : * evaluate the variable
669 : * @return The extrapolated value on the boundary face associated with \p fi (and potentially \p
670 : * elem_side_to_extrapolate_from)
671 : */
672 : virtual ADReal getExtrapolatedBoundaryFaceValue(const FaceInfo & fi,
673 : bool two_term_expansion,
674 : bool correct_skewness,
675 : const Elem * elem_side_to_extrapolate_from,
676 : const StateArg & state) const;
677 :
678 : /// Function to get wether two term boundary expansion is used for the variable
679 0 : const bool & getTwoTermBoundaryExpansion() const { return _two_term_boundary_expansion; }
680 :
681 : protected:
682 : /**
683 : * clear finite volume caches
684 : */
685 : void clearCaches();
686 :
687 : usingMooseVariableFieldMembers;
688 :
689 : /// Holder for all the data associated with the "main" element
690 : std::unique_ptr<MooseVariableDataFV<OutputType>> _element_data;
691 :
692 : /// Holder for all the data associated with the neighbor element
693 : std::unique_ptr<MooseVariableDataFV<OutputType>> _neighbor_data;
694 :
695 : private:
696 : /// The current (ghosted) solution. Note that this needs to be stored as a reference to a pointer
697 : /// because the solution might not exist at the time that this variable is constructed, so we
698 : /// cannot safely dereference at that time
699 : const libMesh::NumericVector<libMesh::Number> * const & _solution;
700 :
701 : /// Shape functions
702 : const FieldVariablePhiValue & _phi;
703 : const FieldVariablePhiGradient & _grad_phi;
704 : const FieldVariablePhiValue & _phi_face;
705 : const FieldVariablePhiGradient & _grad_phi_face;
706 : const FieldVariablePhiValue & _phi_face_neighbor;
707 : const FieldVariablePhiGradient & _grad_phi_face_neighbor;
708 : const FieldVariablePhiValue & _phi_neighbor;
709 : const FieldVariablePhiGradient & _grad_phi_neighbor;
710 :
711 : /// A member used to help determine when we can return cached data as opposed to computing new
712 : /// data
713 : mutable const Elem * _prev_elem;
714 :
715 : /// Map from boundary ID to Dirichlet boundary conditions. Added to speed up Dirichlet BC lookups
716 : /// in \p getDirichletBC
717 : std::unordered_map<BoundaryID, const FVDirichletBCBase *> _boundary_id_to_dirichlet_bc;
718 :
719 : /// Map from boundary ID to flux boundary conditions. Added to enable internal separator
720 : /// boundaries.
721 : std::unordered_map<BoundaryID, std::vector<const FVFluxBC *>> _boundary_id_to_flux_bc;
722 :
723 : /**
724 : * Emit an error message for unsupported lower-d ops
725 : */
726 : [[noreturn]] void lowerDError() const;
727 :
728 : protected:
729 : /// A cache for storing gradients on elements
730 : mutable std::unordered_map<const Elem *, VectorValue<ADReal>> _elem_to_grad;
731 :
732 : /// Whether to use a two term expansion for computing boundary face values
733 : bool _two_term_boundary_expansion;
734 :
735 : /// A member to hold the cell gradient when not caching, used to return a reference (due to
736 : /// expensive ADReal copy)
737 : mutable VectorValue<ADReal> _temp_cell_gradient;
738 :
739 : /// Whether to cache cell gradients
740 : const bool _cache_cell_gradients;
741 :
742 : /// Decides if an average or skewed corrected average is used for the
743 : /// face interpolation. Other options are not taken into account here,
744 : /// but at higher, kernel-based levels.
745 : Moose::FV::InterpMethod _face_interp_method;
746 :
747 : friend void Moose::initDofIndices<>(MooseVariableFV<OutputType> &, const Elem &);
748 : };
749 :
750 : template <typename OutputType>
751 : inline const typename MooseVariableFV<OutputType>::ADDofValues &
752 0 : MooseVariableFV<OutputType>::adDofValues() const
753 : {
754 0 : return _element_data->adDofValues();
755 : }
756 :
757 : template <typename OutputType>
758 : inline const typename MooseVariableFV<OutputType>::ADDofValues &
759 0 : MooseVariableFV<OutputType>::adDofValuesNeighbor() const
760 : {
761 0 : return _neighbor_data->adDofValues();
762 : }
763 :
764 : template <typename OutputType>
765 : inline const typename MooseVariableFV<OutputType>::ADDofValues &
766 0 : MooseVariableFV<OutputType>::adDofValuesDot() const
767 : {
768 0 : return _element_data->adDofValuesDot();
769 : }
770 :
771 : template <typename OutputType>
772 : typename MooseVariableFV<OutputType>::ValueType
773 179503449 : MooseVariableFV<OutputType>::evaluate(const ElemArg & elem_arg, const StateArg & state) const
774 : {
775 179503449 : return getElemValue(elem_arg.elem, state);
776 : }
777 :
778 : template <typename OutputType>
779 : typename MooseVariableFV<OutputType>::ValueType
780 29048 : MooseVariableFV<OutputType>::evaluate(const ElemPointArg & elem_point, const StateArg & state) const
781 : {
782 0 : return (*this)(elem_point.makeElem(), state) +
783 29048 : (elem_point.point - elem_point.elem->vertex_average()) *
784 58096 : this->gradient(elem_point.makeElem(), state);
785 : }
786 :
787 : template <typename OutputType>
788 : typename MooseVariableFV<OutputType>::ValueType
789 21832 : MooseVariableFV<OutputType>::evaluate(const ElemQpArg & elem_qp, const StateArg & state) const
790 : {
791 21832 : return (*this)(ElemPointArg{elem_qp.elem,
792 : elem_qp.point,
793 21832 : _face_interp_method == Moose::FV::InterpMethod::SkewCorrectedAverage},
794 21832 : state);
795 : }
796 :
797 : template <typename OutputType>
798 : typename MooseVariableFV<OutputType>::ValueType
799 3436 : MooseVariableFV<OutputType>::evaluate(const ElemSideQpArg & elem_side_qp,
800 : const StateArg & state) const
801 : {
802 3436 : return (*this)(ElemPointArg{elem_side_qp.elem,
803 : elem_side_qp.point,
804 3436 : _face_interp_method == Moose::FV::InterpMethod::SkewCorrectedAverage},
805 3436 : state);
806 : }
807 :
808 : template <typename OutputType>
809 : typename MooseVariableFV<OutputType>::GradientType
810 6460 : MooseVariableFV<OutputType>::evaluateGradient(const ElemQpArg & qp_arg,
811 : const StateArg & state) const
812 : {
813 6460 : return adGradSln(
814 6460 : qp_arg.elem, state, _face_interp_method == Moose::FV::InterpMethod::SkewCorrectedAverage);
815 : }
816 :
817 : template <typename OutputType>
818 : typename MooseVariableFV<OutputType>::GradientType
819 4920481 : MooseVariableFV<OutputType>::evaluateGradient(const ElemArg & elem_arg,
820 : const StateArg & state) const
821 : {
822 4920481 : return adGradSln(elem_arg.elem, state, elem_arg.correct_skewness);
823 : }
824 :
825 : template <typename OutputType>
826 : typename MooseVariableFV<OutputType>::GradientType
827 4708741 : MooseVariableFV<OutputType>::evaluateGradient(const FaceArg & face, const StateArg & state) const
828 : {
829 : mooseAssert(face.fi, "We must have a non-null face information");
830 4708741 : return adGradSln(*face.fi, state, face.correct_skewness);
831 : }
832 :
833 : template <typename OutputType>
834 : void
835 280833 : MooseVariableFV<OutputType>::setActiveTags(const std::set<TagID> & vtags)
836 : {
837 280833 : _element_data->setActiveTags(vtags);
838 280833 : _neighbor_data->setActiveTags(vtags);
839 280833 : }
840 :
841 : template <typename OutputType>
842 : const std::vector<dof_id_type> &
843 102168 : MooseVariableFV<OutputType>::dofIndicesLower() const
844 : {
845 102168 : static const std::vector<dof_id_type> empty;
846 102168 : return empty;
847 : }
848 :
849 : template <typename OutputType>
850 : void
851 5245 : MooseVariableFV<OutputType>::initialSetup()
852 : {
853 5245 : determineBoundaryToDirichletBCMap();
854 5245 : determineBoundaryToFluxBCMap();
855 5245 : MooseVariableField<OutputType>::initialSetup();
856 5245 : }
857 :
858 : template <typename OutputType>
859 : void
860 624 : MooseVariableFV<OutputType>::meshChanged()
861 : {
862 624 : _prev_elem = nullptr;
863 624 : _dirichlet_map_setup = false;
864 624 : _flux_map_setup = false;
865 624 : MooseVariableField<OutputType>::meshChanged();
866 624 : }
867 :
868 : template <typename OutputType>
869 : void
870 17534 : MooseVariableFV<OutputType>::timestepSetup()
871 : {
872 17534 : _dirichlet_map_setup = false;
873 17534 : _flux_map_setup = false;
874 17534 : MooseVariableField<OutputType>::timestepSetup();
875 17534 : }
876 :
877 : template <typename OutputType>
878 : const typename MooseVariableFV<OutputType>::FieldVariablePhiValue &
879 0 : MooseVariableFV<OutputType>::phiLower() const
880 : {
881 0 : lowerDError();
882 : }
883 :
884 : template <typename OutputType>
885 : void
886 0 : MooseVariableFV<OutputType>::lowerDError() const
887 : {
888 0 : mooseError("Lower dimensional element support not implemented for finite volume variables");
889 : }
890 :
891 : // Declare all the specializations, as the template specialization declaration below must know
892 : template <>
893 : ADReal MooseVariableFV<Real>::evaluateDot(const ElemArg & elem, const StateArg & state) const;
894 : template <>
895 : ADReal MooseVariableFV<Real>::evaluateDot(const FaceArg & elem_arg, const StateArg & state) const;
896 : template <>
897 : ADReal MooseVariableFV<Real>::evaluateDot(const ElemQpArg & elem_arg, const StateArg & state) const;
898 :
899 : // Prevent implicit instantiation in other translation units where these classes are used
900 : extern template class MooseVariableFV<Real>;
|
