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