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 "Problem.h"
13 : #include "DiracKernelInfo.h"
14 : #include "GeometricSearchData.h"
15 : #include "MooseTypes.h"
16 : #include "VectorTag.h"
17 : #include "MooseError.h"
18 : #include "FunctorMaterialProperty.h"
19 : #include "RawValueFunctor.h"
20 : #include "ADWrapperFunctor.h"
21 :
22 : #include "libmesh/coupling_matrix.h"
23 : #include "libmesh/parameters.h"
24 :
25 : #include <memory>
26 : #include <unordered_map>
27 : #include <map>
28 : #include <vector>
29 :
30 : namespace libMesh
31 : {
32 : template <typename>
33 : class VectorValue;
34 : typedef VectorValue<Real> RealVectorValue;
35 : class GhostingFunctor;
36 : }
37 :
38 : class MooseMesh;
39 : class Factory;
40 : class Assembly;
41 : class MooseVariableFieldBase;
42 : class MooseVariableScalar;
43 : template <typename>
44 : class MooseVariableFE;
45 : typedef MooseVariableFE<Real> MooseVariable;
46 : typedef MooseVariableFE<RealVectorValue> VectorMooseVariable;
47 : typedef MooseVariableFE<RealEigenVector> ArrayMooseVariable;
48 : class RestartableDataValue;
49 : class SystemBase;
50 : class LineSearch;
51 : class FaceInfo;
52 : class MooseObjectName;
53 : namespace Moose
54 : {
55 : class FunctorEnvelopeBase;
56 : }
57 :
58 : // libMesh forward declarations
59 : namespace libMesh
60 : {
61 : class EquationSystems;
62 : class DofMap;
63 : class CouplingMatrix;
64 : template <typename T>
65 : class SparseMatrix;
66 : template <typename T>
67 : class NumericVector;
68 : class System;
69 : } // namespace libMesh
70 :
71 : using libMesh::CouplingMatrix;
72 : using libMesh::EquationSystems;
73 :
74 : /**
75 : * Generic class for solving transient nonlinear problems
76 : *
77 : */
78 : class SubProblem : public Problem
79 : {
80 : public:
81 : static InputParameters validParams();
82 :
83 : SubProblem(const InputParameters & parameters);
84 : virtual ~SubProblem();
85 :
86 : virtual libMesh::EquationSystems & es() = 0;
87 : virtual MooseMesh & mesh() = 0;
88 : virtual const MooseMesh & mesh() const = 0;
89 : virtual const MooseMesh & mesh(bool use_displaced) const = 0;
90 :
91 : /**
92 : * @returns whether there will be nonlocal coupling at any point in the simulation, e.g. whether
93 : * there are any active \emph or inactive nonlocal kernels or boundary conditions
94 : */
95 : virtual bool checkNonlocalCouplingRequirement() const = 0;
96 :
97 : /**
98 : * @return whether the given solver system \p sys_num is converged
99 : */
100 0 : virtual bool solverSystemConverged(const unsigned int sys_num) { return converged(sys_num); }
101 :
102 : /**
103 : * @return whether the given nonlinear system \p nl_sys_num is converged.
104 : */
105 : virtual bool nlConverged(const unsigned int nl_sys_num);
106 :
107 : /**
108 : * Eventually we want to convert this virtual over to taking a solver system number argument.
109 : * We will have to first convert apps to use solverSystemConverged, and then once that is done, we
110 : * can change this signature. Then we can go through the apps again and convert back to this
111 : * changed API
112 : */
113 346833 : virtual bool converged(const unsigned int sys_num) { return solverSystemConverged(sys_num); }
114 :
115 : /**
116 : * @return the nonlinear system number corresponding to the provided \p nl_sys_name
117 : */
118 : virtual unsigned int nlSysNum(const NonlinearSystemName & nl_sys_name) const = 0;
119 :
120 : /**
121 : * @return the linear system number corresponding to the provided \p linear_sys_name
122 : */
123 : virtual unsigned int linearSysNum(const LinearSystemName & linear_sys_name) const = 0;
124 :
125 : /**
126 : * @return the solver system number corresponding to the provided \p solver_sys_name
127 : */
128 : virtual unsigned int solverSysNum(const SolverSystemName & solver_sys_name) const = 0;
129 :
130 : virtual void onTimestepBegin() = 0;
131 : virtual void onTimestepEnd() = 0;
132 :
133 : virtual bool isTransient() const = 0;
134 :
135 : /// marks this problem as including/needing finite volume functionality.
136 : virtual void needFV() = 0;
137 :
138 : /// returns true if this problem includes/needs finite volume functionality.
139 : virtual bool haveFV() const = 0;
140 :
141 : /**
142 : * Whether or not the user has requested default ghosting ot be on.
143 : */
144 135412 : bool defaultGhosting() { return _default_ghosting; }
145 :
146 : /**
147 : * Create a Tag. Tags can be associated with Vectors and Matrices and allow objects
148 : * (such as Kernels) to arbitrarily contribute values to any set of vectors/matrics
149 : *
150 : * Note: If the tag is already present then this will simply return the TagID of that Tag, but the
151 : * type must be the same.
152 : *
153 : * @param tag_name The name of the tag to create, the TagID will get automatically generated
154 : * @param type The type of the tag
155 : */
156 : virtual TagID addVectorTag(const TagName & tag_name,
157 : const Moose::VectorTagType type = Moose::VECTOR_TAG_RESIDUAL);
158 :
159 : /**
160 : * Adds a vector tag to the list of vectors that will not be zeroed
161 : * when other tagged vectors are
162 : * @param tag the TagID of the vector that will be manually managed
163 : */
164 : void addNotZeroedVectorTag(const TagID tag);
165 :
166 : /**
167 : * Checks if a vector tag is in the list of vectors that will not be zeroed
168 : * when other tagged vectors are
169 : * @param tag the TagID of the vector that is currently being checked
170 : * @returns false if the tag is not within the set of vectors that are
171 : * intended to not be zero or if the set is empty. returns true otherwise
172 : */
173 : bool vectorTagNotZeroed(const TagID tag) const;
174 :
175 : /**
176 : * Get a VectorTag from a TagID.
177 : */
178 : virtual const VectorTag & getVectorTag(const TagID tag_id) const;
179 : std::vector<VectorTag> getVectorTags(const std::set<TagID> & tag_ids) const;
180 :
181 : /**
182 : * Get a TagID from a TagName.
183 : */
184 : virtual TagID getVectorTagID(const TagName & tag_name) const;
185 :
186 : /**
187 : * Retrieve the name associated with a TagID
188 : */
189 : virtual TagName vectorTagName(const TagID tag) const;
190 :
191 : /**
192 : * Return all vector tags, where a tag is represented by a map from name to ID. Can optionally be
193 : * limited to a vector tag type.
194 : */
195 : virtual const std::vector<VectorTag> &
196 : getVectorTags(const Moose::VectorTagType type = Moose::VECTOR_TAG_ANY) const;
197 :
198 : /**
199 : * Check to see if a particular Tag exists
200 : */
201 3448303706 : virtual bool vectorTagExists(const TagID tag_id) const { return tag_id < _vector_tags.size(); }
202 :
203 : /**
204 : * Check to see if a particular Tag exists by using Tag name
205 : */
206 : virtual bool vectorTagExists(const TagName & tag_name) const;
207 :
208 : /**
209 : * The total number of tags, which can be limited to the tag type
210 : */
211 : virtual unsigned int numVectorTags(const Moose::VectorTagType type = Moose::VECTOR_TAG_ANY) const;
212 :
213 : virtual Moose::VectorTagType vectorTagType(const TagID tag_id) const;
214 :
215 : /**
216 : * Create a Tag. Tags can be associated with Vectors and Matrices and allow objects
217 : * (such as Kernels) to arbitrarily contribute values to any set of vectors/matrics
218 : *
219 : * Note: If the tag is already present then this will simply return the TagID of that Tag
220 : *
221 : * @param tag_name The name of the tag to create, the TagID will get automatically generated
222 : */
223 : virtual TagID addMatrixTag(TagName tag_name);
224 :
225 : /**
226 : * Get a TagID from a TagName.
227 : */
228 : virtual TagID getMatrixTagID(const TagName & tag_name) const;
229 :
230 : /**
231 : * Retrieve the name associated with a TagID
232 : */
233 : virtual TagName matrixTagName(TagID tag);
234 :
235 : /**
236 : * Check to see if a particular Tag exists
237 : */
238 : virtual bool matrixTagExists(const TagName & tag_name) const;
239 :
240 : /**
241 : * Check to see if a particular Tag exists
242 : */
243 : virtual bool matrixTagExists(TagID tag_id) const;
244 :
245 : /**
246 : * The total number of tags
247 : */
248 9803147 : virtual unsigned int numMatrixTags() const { return _matrix_tag_name_to_tag_id.size(); }
249 :
250 : /**
251 : * Return all matrix tags in the system, where a tag is represented by a map from name to ID
252 : */
253 539697 : virtual std::map<TagName, TagID> & getMatrixTags() { return _matrix_tag_name_to_tag_id; }
254 :
255 : /// Whether or not this problem has the variable
256 : virtual bool hasVariable(const std::string & var_name) const = 0;
257 :
258 : /// Whether or not this problem has this linear variable
259 : virtual bool hasLinearVariable(const std::string & var_name) const;
260 :
261 : /// Whether or not this problem has this auxiliary variable
262 : virtual bool hasAuxiliaryVariable(const std::string & var_name) const;
263 :
264 : /**
265 : * Returns the variable reference for requested variable which must
266 : * be of the expected_var_type (Nonlinear vs. Auxiliary) and
267 : * expected_var_field_type (standard, scalar, vector). The default
268 : * values of VAR_ANY and VAR_FIELD_ANY should be used when "any"
269 : * type of variable is acceptable. Throws an error if the variable
270 : * in question is not in the expected System or of the expected
271 : * type.
272 : */
273 : virtual const MooseVariableFieldBase & getVariable(
274 : const THREAD_ID tid,
275 : const std::string & var_name,
276 : Moose::VarKindType expected_var_type = Moose::VarKindType::VAR_ANY,
277 : Moose::VarFieldType expected_var_field_type = Moose::VarFieldType::VAR_FIELD_ANY) const = 0;
278 : virtual MooseVariableFieldBase &
279 4811832 : getVariable(const THREAD_ID tid,
280 : const std::string & var_name,
281 : Moose::VarKindType expected_var_type = Moose::VarKindType::VAR_ANY,
282 : Moose::VarFieldType expected_var_field_type = Moose::VarFieldType::VAR_FIELD_ANY)
283 : {
284 4811832 : return const_cast<MooseVariableFieldBase &>(const_cast<const SubProblem *>(this)->getVariable(
285 4811824 : tid, var_name, expected_var_type, expected_var_field_type));
286 : }
287 :
288 : /// Returns the variable reference for requested MooseVariable which may be in any system
289 : virtual MooseVariable & getStandardVariable(const THREAD_ID tid,
290 : const std::string & var_name) = 0;
291 :
292 : /// Returns the variable reference for requested MooseVariableField which may be in any system
293 : virtual MooseVariableFieldBase & getActualFieldVariable(const THREAD_ID tid,
294 : const std::string & var_name) = 0;
295 :
296 : /// Returns the variable reference for requested VectorMooseVariable which may be in any system
297 : virtual VectorMooseVariable & getVectorVariable(const THREAD_ID tid,
298 : const std::string & var_name) = 0;
299 :
300 : /// Returns the variable reference for requested ArrayMooseVariable which may be in any system
301 : virtual ArrayMooseVariable & getArrayVariable(const THREAD_ID tid,
302 : const std::string & var_name) = 0;
303 :
304 : /// Returns a Boolean indicating whether any system contains a variable with the name provided
305 : virtual bool hasScalarVariable(const std::string & var_name) const = 0;
306 :
307 : /// Returns the scalar variable reference from whichever system contains it
308 : virtual MooseVariableScalar & getScalarVariable(const THREAD_ID tid,
309 : const std::string & var_name) = 0;
310 :
311 : /// Returns the equation system containing the variable provided
312 : virtual libMesh::System & getSystem(const std::string & var_name) = 0;
313 :
314 : /**
315 : * Set the MOOSE variables to be reinited on each element.
316 : * @param moose_vars A set of variables that need to be reinited each time reinit() is called.
317 : *
318 : * @param tid The thread id
319 : */
320 : virtual void
321 : setActiveElementalMooseVariables(const std::set<MooseVariableFieldBase *> & moose_vars,
322 : const THREAD_ID tid);
323 :
324 : /**
325 : * Get the MOOSE variables to be reinited on each element.
326 : *
327 : * @param tid The thread id
328 : */
329 : virtual const std::set<MooseVariableFieldBase *> &
330 : getActiveElementalMooseVariables(const THREAD_ID tid) const;
331 :
332 : /**
333 : * Whether or not a list of active elemental moose variables has been set.
334 : *
335 : * @return True if there has been a list of active elemental moose variables set, False otherwise
336 : */
337 : virtual bool hasActiveElementalMooseVariables(const THREAD_ID tid) const;
338 :
339 : /**
340 : * Clear the active elemental MooseVariableFieldBase. If there are no active variables then they
341 : * will all be reinited. Call this after finishing the computation that was using a restricted set
342 : * of MooseVariableFieldBase
343 : *
344 : * @param tid The thread id
345 : */
346 : virtual void clearActiveElementalMooseVariables(const THREAD_ID tid);
347 :
348 : virtual Assembly & assembly(const THREAD_ID tid, const unsigned int sys_num) = 0;
349 : virtual const Assembly & assembly(const THREAD_ID tid, const unsigned int sys_num) const = 0;
350 :
351 : /**
352 : * Return the nonlinear system object as a base class reference given the system number
353 : */
354 : virtual const SystemBase & systemBaseNonlinear(const unsigned int sys_num) const = 0;
355 : virtual SystemBase & systemBaseNonlinear(const unsigned int sys_num) = 0;
356 : /**
357 : * Return the linear system object as a base class reference given the system number
358 : */
359 : virtual const SystemBase & systemBaseLinear(const unsigned int sys_num) const = 0;
360 : virtual SystemBase & systemBaseLinear(const unsigned int sys_num) = 0;
361 : /**
362 : * Return the solver system object as a base class reference given the system number
363 : */
364 : virtual const SystemBase & systemBaseSolver(const unsigned int sys_num) const = 0;
365 : virtual SystemBase & systemBaseSolver(const unsigned int sys_num) = 0;
366 : /**
367 : * Return the auxiliary system object as a base class reference
368 : */
369 : virtual const SystemBase & systemBaseAuxiliary() const = 0;
370 : virtual SystemBase & systemBaseAuxiliary() = 0;
371 :
372 : virtual void prepareShapes(unsigned int var, const THREAD_ID tid) = 0;
373 : virtual void prepareFaceShapes(unsigned int var, const THREAD_ID tid) = 0;
374 : virtual void prepareNeighborShapes(unsigned int var, const THREAD_ID tid) = 0;
375 : Moose::CoordinateSystemType getCoordSystem(SubdomainID sid) const;
376 :
377 : /**
378 : * Returns the desired radial direction for RZ coordinate transformation
379 : * @return The coordinate direction for the radial direction
380 : */
381 : unsigned int getAxisymmetricRadialCoord() const;
382 :
383 : virtual DiracKernelInfo & diracKernelInfo();
384 : virtual Real finalNonlinearResidual(const unsigned int nl_sys_num) const;
385 : virtual unsigned int nNonlinearIterations(const unsigned int nl_sys_num) const;
386 : virtual unsigned int nLinearIterations(const unsigned int nl_sys_num) const;
387 :
388 : virtual void addResidual(const THREAD_ID tid) = 0;
389 : virtual void addResidualNeighbor(const THREAD_ID tid) = 0;
390 : virtual void addResidualLower(const THREAD_ID tid) = 0;
391 :
392 : virtual void cacheResidual(const THREAD_ID tid);
393 : virtual void cacheResidualNeighbor(const THREAD_ID tid);
394 : virtual void addCachedResidual(const THREAD_ID tid);
395 :
396 : virtual void setResidual(libMesh::NumericVector<libMesh::Number> & residual,
397 : const THREAD_ID tid) = 0;
398 : virtual void setResidualNeighbor(libMesh::NumericVector<libMesh::Number> & residual,
399 : const THREAD_ID tid) = 0;
400 :
401 : virtual void addJacobian(const THREAD_ID tid) = 0;
402 : virtual void addJacobianNeighbor(const THREAD_ID tid) = 0;
403 : virtual void addJacobianNeighborLowerD(const THREAD_ID tid) = 0;
404 : virtual void addJacobianLowerD(const THREAD_ID tid) = 0;
405 : virtual void addJacobianNeighbor(libMesh::SparseMatrix<libMesh::Number> & jacobian,
406 : unsigned int ivar,
407 : unsigned int jvar,
408 : const libMesh::DofMap & dof_map,
409 : std::vector<dof_id_type> & dof_indices,
410 : std::vector<dof_id_type> & neighbor_dof_indices,
411 : const std::set<TagID> & tags,
412 : const THREAD_ID tid) = 0;
413 :
414 : virtual void cacheJacobian(const THREAD_ID tid);
415 : virtual void cacheJacobianNeighbor(const THREAD_ID tid);
416 : virtual void addCachedJacobian(const THREAD_ID tid);
417 :
418 : virtual void prepare(const Elem * elem, const THREAD_ID tid) = 0;
419 : virtual void prepareFace(const Elem * elem, const THREAD_ID tid) = 0;
420 : virtual void prepare(const Elem * elem,
421 : unsigned int ivar,
422 : unsigned int jvar,
423 : const std::vector<dof_id_type> & dof_indices,
424 : const THREAD_ID tid) = 0;
425 : virtual void setCurrentSubdomainID(const Elem * elem, const THREAD_ID tid) = 0;
426 : virtual void
427 : setNeighborSubdomainID(const Elem * elem, unsigned int side, const THREAD_ID tid) = 0;
428 : virtual void prepareAssembly(const THREAD_ID tid) = 0;
429 :
430 : virtual void reinitElem(const Elem * elem, const THREAD_ID tid) = 0;
431 : virtual void reinitElemPhys(const Elem * elem,
432 : const std::vector<Point> & phys_points_in_elem,
433 : const THREAD_ID tid) = 0;
434 : virtual void reinitElemFace(const Elem * elem, unsigned int side, const THREAD_ID tid) = 0;
435 : virtual void reinitLowerDElem(const Elem * lower_d_elem,
436 : const THREAD_ID tid,
437 : const std::vector<Point> * const pts = nullptr,
438 : const std::vector<Real> * const weights = nullptr);
439 : virtual void reinitNode(const Node * node, const THREAD_ID tid) = 0;
440 : virtual void reinitNodeFace(const Node * node, BoundaryID bnd_id, const THREAD_ID tid) = 0;
441 : virtual void reinitNodes(const std::vector<dof_id_type> & nodes, const THREAD_ID tid) = 0;
442 : virtual void reinitNodesNeighbor(const std::vector<dof_id_type> & nodes, const THREAD_ID tid) = 0;
443 : virtual void reinitNeighbor(const Elem * elem, unsigned int side, const THREAD_ID tid) = 0;
444 : virtual void reinitNeighborPhys(const Elem * neighbor,
445 : unsigned int neighbor_side,
446 : const std::vector<Point> & physical_points,
447 : const THREAD_ID tid) = 0;
448 : virtual void reinitNeighborPhys(const Elem * neighbor,
449 : const std::vector<Point> & physical_points,
450 : const THREAD_ID tid) = 0;
451 : virtual void
452 : reinitElemNeighborAndLowerD(const Elem * elem, unsigned int side, const THREAD_ID tid) = 0;
453 : /**
454 : * fills the VariableValue arrays for scalar variables from the solution vector
455 : * @param tid The thread id
456 : * @param reinit_for_derivative_reordering A flag indicating whether we are reinitializing for the
457 : * purpose of re-ordering derivative information for ADNodalBCs
458 : */
459 : virtual void reinitScalars(const THREAD_ID tid,
460 : bool reinit_for_derivative_reordering = false) = 0;
461 : virtual void reinitOffDiagScalars(const THREAD_ID tid) = 0;
462 :
463 : /// sets the current boundary ID in assembly
464 : virtual void setCurrentBoundaryID(BoundaryID bid, const THREAD_ID tid);
465 :
466 : /**
467 : * reinitialize FE objects on a given element on a given side at a given set of reference
468 : * points and then compute variable data. Note that this method makes no assumptions about what's
469 : * been called beforehand, e.g. you don't have to call some prepare method before this one. This
470 : * is an all-in-one reinit
471 : */
472 : virtual void reinitElemFaceRef(const Elem * elem,
473 : unsigned int side,
474 : Real tolerance,
475 : const std::vector<Point> * const pts,
476 : const std::vector<Real> * const weights = nullptr,
477 : const THREAD_ID tid = 0);
478 :
479 : /**
480 : * reinitialize FE objects on a given neighbor element on a given side at a given set of reference
481 : * points and then compute variable data. Note that this method makes no assumptions about what's
482 : * been called beforehand, e.g. you don't have to call some prepare method before this one. This
483 : * is an all-in-one reinit
484 : */
485 : virtual void reinitNeighborFaceRef(const Elem * neighbor_elem,
486 : unsigned int neighbor_side,
487 : Real tolerance,
488 : const std::vector<Point> * const pts,
489 : const std::vector<Real> * const weights = nullptr,
490 : const THREAD_ID tid = 0);
491 :
492 : /**
493 : * reinitialize a neighboring lower dimensional element
494 : */
495 : void reinitNeighborLowerDElem(const Elem * elem, const THREAD_ID tid = 0);
496 :
497 : /**
498 : * Reinit a mortar element to obtain a valid JxW
499 : */
500 : void reinitMortarElem(const Elem * elem, const THREAD_ID tid = 0);
501 :
502 : /**
503 : * Returns true if the Problem has Dirac kernels it needs to compute on elem.
504 : */
505 : virtual bool reinitDirac(const Elem * elem, const THREAD_ID tid) = 0;
506 : /**
507 : * Fills "elems" with the elements that should be looped over for Dirac Kernels
508 : */
509 : virtual void getDiracElements(std::set<const Elem *> & elems) = 0;
510 : /**
511 : * Gets called before Dirac Kernels are asked to add the points they are supposed to be evaluated
512 : * in
513 : */
514 : virtual void clearDiracInfo() = 0;
515 :
516 : // Geom Search
517 : virtual void
518 : updateGeomSearch(GeometricSearchData::GeometricSearchType type = GeometricSearchData::ALL) = 0;
519 :
520 : virtual GeometricSearchData & geomSearchData() = 0;
521 :
522 : /**
523 : * Adds the given material property to a storage map based on block ids
524 : *
525 : * This is method is called from within the Material class when the property
526 : * is first registered.
527 : * @param block_id The block id for the MaterialProperty
528 : * @param name The name of the property
529 : */
530 : virtual void storeSubdomainMatPropName(SubdomainID block_id, const std::string & name);
531 :
532 : /**
533 : * Adds the given material property to a storage map based on boundary ids
534 : *
535 : * This is method is called from within the Material class when the property
536 : * is first registered.
537 : * @param boundary_id The block id for the MaterialProperty
538 : * @param name The name of the property
539 : */
540 : virtual void storeBoundaryMatPropName(BoundaryID boundary_id, const std::string & name);
541 :
542 : /**
543 : * Adds to a map based on block ids of material properties for which a zero
544 : * value can be returned. Thes properties are optional and will not trigger a
545 : * missing material property error.
546 : *
547 : * @param block_id The block id for the MaterialProperty
548 : * @param name The name of the property
549 : */
550 : virtual void storeSubdomainZeroMatProp(SubdomainID block_id, const MaterialPropertyName & name);
551 :
552 : /**
553 : * Adds to a map based on boundary ids of material properties for which a zero
554 : * value can be returned. Thes properties are optional and will not trigger a
555 : * missing material property error.
556 : *
557 : * @param boundary_id The block id for the MaterialProperty
558 : * @param name The name of the property
559 : */
560 : virtual void storeBoundaryZeroMatProp(BoundaryID boundary_id, const MaterialPropertyName & name);
561 :
562 : /**
563 : * Adds to a map based on block ids of material properties to validate
564 : *
565 : * @param block_id The block id for the MaterialProperty
566 : * @param name The name of the property
567 : */
568 : virtual void storeSubdomainDelayedCheckMatProp(const std::string & requestor,
569 : SubdomainID block_id,
570 : const std::string & name);
571 :
572 : /**
573 : * Adds to a map based on boundary ids of material properties to validate
574 : *
575 : * @param requestor The MOOSE object name requesting the material property
576 : * @param boundary_id The block id for the MaterialProperty
577 : * @param name The name of the property
578 : */
579 : virtual void storeBoundaryDelayedCheckMatProp(const std::string & requestor,
580 : BoundaryID boundary_id,
581 : const std::string & name);
582 :
583 : /**
584 : * Checks block material properties integrity
585 : *
586 : * \see FEProblemBase::checkProblemIntegrity
587 : */
588 : virtual void checkBlockMatProps();
589 :
590 : /**
591 : * Checks boundary material properties integrity
592 : *
593 : * \see FEProblemBase::checkProblemIntegrity
594 : */
595 : virtual void checkBoundaryMatProps();
596 :
597 : /**
598 : * Helper method for adding a material property name to the _material_property_requested set
599 : */
600 : virtual void markMatPropRequested(const std::string &);
601 :
602 : /**
603 : * Find out if a material property has been requested by any object
604 : */
605 : virtual bool isMatPropRequested(const std::string & prop_name) const;
606 :
607 : /**
608 : * Helper for tracking the object that is consuming a property for MaterialPropertyDebugOutput
609 : */
610 : void addConsumedPropertyName(const MooseObjectName & obj_name, const std::string & prop_name);
611 :
612 : /**
613 : * Return the map that tracks the object with consumed material properties
614 : */
615 : const std::map<MooseObjectName, std::set<std::string>> & getConsumedPropertyMap() const;
616 :
617 : /**
618 : * Will make sure that all dofs connected to elem_id are ghosted to this processor
619 : */
620 : virtual void addGhostedElem(dof_id_type elem_id) = 0;
621 :
622 : /**
623 : * Will make sure that all necessary elements from boundary_id are ghosted to this processor
624 : */
625 : virtual void addGhostedBoundary(BoundaryID boundary_id) = 0;
626 :
627 : /**
628 : * Causes the boundaries added using addGhostedBoundary to actually be ghosted.
629 : */
630 : virtual void ghostGhostedBoundaries() = 0;
631 :
632 : /**
633 : * Get a vector containing the block ids the material property is defined on.
634 : */
635 : virtual std::set<SubdomainID> getMaterialPropertyBlocks(const std::string & prop_name);
636 :
637 : /**
638 : * Get a vector of block id equivalences that the material property is defined on.
639 : */
640 : virtual std::vector<SubdomainName> getMaterialPropertyBlockNames(const std::string & prop_name);
641 :
642 : /**
643 : * Check if a material property is defined on a block.
644 : */
645 : virtual bool hasBlockMaterialProperty(SubdomainID block_id, const std::string & prop_name);
646 :
647 : /**
648 : * Get a vector containing the block ids the material property is defined on.
649 : */
650 : virtual std::set<BoundaryID> getMaterialPropertyBoundaryIDs(const std::string & prop_name);
651 :
652 : /**
653 : * Get a vector of block id equivalences that the material property is defined on.
654 : */
655 : virtual std::vector<BoundaryName> getMaterialPropertyBoundaryNames(const std::string & prop_name);
656 :
657 : /**
658 : * Check if a material property is defined on a block.
659 : */
660 : virtual bool hasBoundaryMaterialProperty(BoundaryID boundary_id, const std::string & prop_name);
661 :
662 : /**
663 : * Returns true if the problem is in the process of computing it's initial residual.
664 : * @return Whether or not the problem is currently computing the initial residual.
665 : */
666 : virtual bool computingPreSMOResidual(const unsigned int nl_sys_num) const = 0;
667 :
668 : /**
669 : * Return the list of elements that should have their DoFs ghosted to this processor.
670 : * @return The list
671 : */
672 75968 : virtual std::set<dof_id_type> & ghostedElems() { return _ghosted_elems; }
673 :
674 : std::map<std::string, std::vector<dof_id_type>> _var_dof_map;
675 :
676 : /**
677 : * @return the nonlocal coupling matrix for the i'th nonlinear system
678 : */
679 : virtual const libMesh::CouplingMatrix & nonlocalCouplingMatrix(const unsigned i) const = 0;
680 :
681 : /**
682 : * Returns true if the problem is in the process of computing the Jacobian
683 : */
684 453938476 : const bool & currentlyComputingJacobian() const { return _currently_computing_jacobian; };
685 :
686 : /**
687 : * Set whether or not the problem is in the process of computing the Jacobian
688 : */
689 4075358 : void setCurrentlyComputingJacobian(const bool currently_computing_jacobian)
690 : {
691 4075358 : _currently_computing_jacobian = currently_computing_jacobian;
692 4075358 : }
693 :
694 : /**
695 : * Returns true if the problem is in the process of computing the residual and the Jacobian
696 : */
697 : const bool & currentlyComputingResidualAndJacobian() const;
698 :
699 : /**
700 : * Set whether or not the problem is in the process of computing the Jacobian
701 : */
702 : void setCurrentlyComputingResidualAndJacobian(bool currently_computing_residual_and_jacobian);
703 :
704 : /**
705 : * Returns true if the problem is in the process of computing the nonlinear residual
706 : */
707 23415 : bool computingNonlinearResid() const { return _computing_nonlinear_residual; }
708 :
709 : /**
710 : * Set whether or not the problem is in the process of computing the nonlinear residual
711 : */
712 212324 : virtual void computingNonlinearResid(const bool computing_nonlinear_residual)
713 : {
714 212324 : _computing_nonlinear_residual = computing_nonlinear_residual;
715 212324 : }
716 :
717 : /**
718 : * Returns true if the problem is in the process of computing the residual
719 : */
720 73935 : const bool & currentlyComputingResidual() const { return _currently_computing_residual; }
721 :
722 : /**
723 : * Set whether or not the problem is in the process of computing the residual
724 : */
725 595424 : virtual void setCurrentlyComputingResidual(const bool currently_computing_residual)
726 : {
727 595424 : _currently_computing_residual = currently_computing_residual;
728 595424 : }
729 :
730 : /// Is it safe to access the tagged matrices
731 1214449735 : virtual bool safeAccessTaggedMatrices() const { return _safe_access_tagged_matrices; }
732 :
733 : /// Is it safe to access the tagged vectors
734 653435 : virtual bool safeAccessTaggedVectors() const { return _safe_access_tagged_vectors; }
735 :
736 : virtual void clearActiveFEVariableCoupleableMatrixTags(const THREAD_ID tid);
737 :
738 : virtual void clearActiveFEVariableCoupleableVectorTags(const THREAD_ID tid);
739 :
740 : virtual void setActiveFEVariableCoupleableVectorTags(std::set<TagID> & vtags,
741 : const THREAD_ID tid);
742 :
743 : virtual void setActiveFEVariableCoupleableMatrixTags(std::set<TagID> & mtags,
744 : const THREAD_ID tid);
745 :
746 : virtual void clearActiveScalarVariableCoupleableMatrixTags(const THREAD_ID tid);
747 :
748 : virtual void clearActiveScalarVariableCoupleableVectorTags(const THREAD_ID tid);
749 :
750 : virtual void setActiveScalarVariableCoupleableVectorTags(std::set<TagID> & vtags,
751 : const THREAD_ID tid);
752 :
753 : virtual void setActiveScalarVariableCoupleableMatrixTags(std::set<TagID> & mtags,
754 : const THREAD_ID tid);
755 :
756 : const std::set<TagID> & getActiveScalarVariableCoupleableVectorTags(const THREAD_ID tid) const;
757 :
758 : const std::set<TagID> & getActiveScalarVariableCoupleableMatrixTags(const THREAD_ID tid) const;
759 :
760 : const std::set<TagID> & getActiveFEVariableCoupleableVectorTags(const THREAD_ID tid) const;
761 :
762 : const std::set<TagID> & getActiveFEVariableCoupleableMatrixTags(const THREAD_ID tid) const;
763 :
764 : /**
765 : * Method for setting whether we have any ad objects
766 : */
767 560 : virtual void haveADObjects(bool have_ad_objects) { _have_ad_objects = have_ad_objects; }
768 : /**
769 : * Method for reading wehther we have any ad objects
770 : */
771 457060854 : bool haveADObjects() const { return _have_ad_objects; }
772 :
773 : virtual LineSearch * getLineSearch() = 0;
774 :
775 : /**
776 : * The coupling matrix defining what blocks exist in the preconditioning matrix
777 : */
778 : virtual const libMesh::CouplingMatrix * couplingMatrix(const unsigned int nl_sys_num) const = 0;
779 :
780 : private:
781 : /**
782 : * Creates (n_sys - 1) clones of the provided algebraic ghosting functor (corresponding to the
783 : * nonlinear system algebraic ghosting functor), initializes the clone with the appropriate
784 : * DofMap, and then adds the clone to said DofMap
785 : * @param algebraic_gf the (nonlinear system's) algebraic ghosting functor to clone
786 : * @param to_mesh whether the clone should be added to the corresponding DofMap's underlying
787 : * MeshBase (the underlying MeshBase will be the same for every system held by this object's
788 : * EquationSystems object)
789 : */
790 : void cloneAlgebraicGhostingFunctor(libMesh::GhostingFunctor & algebraic_gf, bool to_mesh = true);
791 :
792 : /**
793 : * Creates (n_sys - 1) clones of the provided coupling ghosting functor (corresponding to the
794 : * nonlinear system coupling ghosting functor), initializes the clone with the appropriate
795 : * DofMap, and then adds the clone to said DofMap
796 : * @param coupling_gf the (nonlinear system's) coupling ghosting functor to clone
797 : * @param to_mesh whether the clone should be added to the corresponding DofMap's underlying
798 : * MeshBase (the underlying MeshBase will be the same for every system held by this object's
799 : * EquationSystems object)
800 : */
801 : void cloneCouplingGhostingFunctor(libMesh::GhostingFunctor & coupling_gf, bool to_mesh = true);
802 :
803 : public:
804 : /**
805 : * Add an algebraic ghosting functor to this problem's DofMaps
806 : */
807 : void addAlgebraicGhostingFunctor(libMesh::GhostingFunctor & algebraic_gf, bool to_mesh = true);
808 :
809 : /**
810 : * Add a coupling functor to this problem's DofMaps
811 : */
812 : void addCouplingGhostingFunctor(libMesh::GhostingFunctor & coupling_gf, bool to_mesh = true);
813 :
814 : /**
815 : * Remove an algebraic ghosting functor from this problem's DofMaps
816 : */
817 : void removeAlgebraicGhostingFunctor(libMesh::GhostingFunctor & algebraic_gf);
818 :
819 : /**
820 : * Remove a coupling ghosting functor from this problem's DofMaps
821 : */
822 : void removeCouplingGhostingFunctor(libMesh::GhostingFunctor & coupling_gf);
823 :
824 : /**
825 : * Automatic scaling setter
826 : * @param automatic_scaling A boolean representing whether we are performing automatic scaling
827 : */
828 : virtual void automaticScaling(bool automatic_scaling);
829 :
830 : /**
831 : * Automatic scaling getter
832 : * @return A boolean representing whether we are performing automatic scaling
833 : */
834 : bool automaticScaling() const;
835 :
836 : /**
837 : * Tells this problem that the assembly associated with the given nonlinear system number involves
838 : * a scaling vector
839 : */
840 : void hasScalingVector(const unsigned int nl_sys_num);
841 :
842 : /**
843 : * Whether we have a displaced problem in our simulation
844 : */
845 : virtual bool haveDisplaced() const = 0;
846 :
847 : /**
848 : * Getter for whether we're computing the scaling jacobian
849 : */
850 : virtual bool computingScalingJacobian() const = 0;
851 :
852 : /**
853 : * Getter for whether we're computing the scaling residual
854 : */
855 : virtual bool computingScalingResidual() const = 0;
856 :
857 : /**
858 : * Clear dof indices from variables in nl and aux systems
859 : */
860 : void clearAllDofIndices();
861 :
862 : /**
863 : * @tparam T The type that the functor will return when evaluated, e.g. \p
864 : ADReal or \p Real
865 : * @param name The name of the functor to retrieve
866 : * @param tid The thread ID that we are retrieving the functor property for
867 : * @param requestor_name The name of the object that is requesting this functor property
868 : * @param requestor_is_ad Whether the requesting object is an AD object
869 : * @return a constant reference to the functor
870 : */
871 : template <typename T>
872 : const Moose::Functor<T> & getFunctor(const std::string & name,
873 : const THREAD_ID tid,
874 : const std::string & requestor_name,
875 : bool requestor_is_ad);
876 :
877 : /**
878 : * checks whether we have a functor corresponding to \p name on the thread id \p tid
879 : */
880 : bool hasFunctor(const std::string & name, const THREAD_ID tid) const;
881 :
882 : /**
883 : * checks whether we have a functor of type T corresponding to \p name on the thread id \p tid
884 : */
885 : template <typename T>
886 : bool hasFunctorWithType(const std::string & name, const THREAD_ID tid) const;
887 :
888 : /**
889 : * add a functor to the problem functor container
890 : */
891 : template <typename T>
892 : void
893 : addFunctor(const std::string & name, const Moose::FunctorBase<T> & functor, const THREAD_ID tid);
894 :
895 : /**
896 : * Add a functor that has block-wise lambda definitions, e.g. the evaluations of the functor are
897 : * based on a user-provided lambda expression.
898 : * @param name The name of the functor to add
899 : * @param my_lammy The lambda expression that will be called when the functor is evaluated
900 : * @param clearance_schedule How often to clear functor evaluations. The default value is always,
901 : * which means that the functor will be re-evaluated every time it is called. If it is something
902 : * other than always, than cached values may be returned
903 : * @param mesh The mesh on which this functor operates
904 : * @param block_ids The blocks on which the lambda expression is defined
905 : * @param tid The thread on which the functor we are adding will run
906 : * @return The added functor
907 : */
908 : template <typename T, typename PolymorphicLambda>
909 : const Moose::FunctorBase<T> &
910 : addPiecewiseByBlockLambdaFunctor(const std::string & name,
911 : PolymorphicLambda my_lammy,
912 : const std::set<ExecFlagType> & clearance_schedule,
913 : const MooseMesh & mesh,
914 : const std::set<SubdomainID> & block_ids,
915 : const THREAD_ID tid);
916 :
917 : virtual void initialSetup();
918 : virtual void timestepSetup();
919 : virtual void customSetup(const ExecFlagType & exec_type);
920 : virtual void residualSetup();
921 : virtual void jacobianSetup();
922 :
923 : /// Setter for debug functor output
924 36 : void setFunctorOutput(bool set_output) { _show_functors = set_output; }
925 : /// Setter for debug chain control data output
926 13 : void setChainControlDataOutput(bool set_output) { _show_chain_control_data = set_output; }
927 :
928 : /**
929 : * @return the number of nonlinear systems in the problem
930 : */
931 : virtual std::size_t numNonlinearSystems() const = 0;
932 :
933 : /**
934 : * @return the current nonlinear system number
935 : */
936 : virtual unsigned int currentNlSysNum() const = 0;
937 :
938 : /**
939 : * @return the number of linear systems in the problem
940 : */
941 : virtual std::size_t numLinearSystems() const = 0;
942 :
943 : /**
944 : * @return the number of solver systems in the problem
945 : */
946 : virtual std::size_t numSolverSystems() const = 0;
947 :
948 : /**
949 : * @return the current linear system number
950 : */
951 : virtual unsigned int currentLinearSysNum() const = 0;
952 :
953 : /**
954 : * Register an unfulfilled functor request
955 : */
956 : template <typename T>
957 : void registerUnfilledFunctorRequest(T * functor_interface,
958 : const std::string & functor_name,
959 : const THREAD_ID tid);
960 :
961 : /**
962 : * Return the residual vector tags we are currently computing
963 : */
964 : virtual const std::vector<VectorTag> & currentResidualVectorTags() const = 0;
965 :
966 : /**
967 : * Select the vector tags which belong to a specific system
968 : * @param system Reference to the system
969 : * @param input_vector_tags A vector of vector tags
970 : * @param selected_tags A set which gets populated by the tag-ids that belong to the system
971 : */
972 : static void selectVectorTagsFromSystem(const SystemBase & system,
973 : const std::vector<VectorTag> & input_vector_tags,
974 : std::set<TagID> & selected_tags);
975 :
976 : /**
977 : * Select the matrix tags which belong to a specific system
978 : * @param system Reference to the system
979 : * @param input_matrix_tags A map of matrix tags
980 : * @param selected_tags A set which gets populated by the tag-ids that belong to the system
981 : */
982 : static void selectMatrixTagsFromSystem(const SystemBase & system,
983 : const std::map<TagName, TagID> & input_matrix_tags,
984 : std::set<TagID> & selected_tags);
985 :
986 : /**
987 : * reinitialize the finite volume assembly data for the provided face and thread
988 : */
989 : void reinitFVFace(const THREAD_ID tid, const FaceInfo & fi);
990 :
991 : /**
992 : * Whether the simulation has active nonlocal coupling which should be accounted for in the
993 : * Jacobian. For this to return true, there must be at least one active nonlocal kernel or
994 : * boundary condition
995 : */
996 : virtual bool hasNonlocalCoupling() const = 0;
997 :
998 : /**
999 : * Prepare \p DofMap and \p Assembly classes with our p-refinement information
1000 : */
1001 : void preparePRefinement();
1002 :
1003 : /**
1004 : * @returns whether we're doing p-refinement
1005 : */
1006 : [[nodiscard]] bool doingPRefinement() const;
1007 :
1008 : /**
1009 : * Query whether p-refinement has been requested at any point during the simulation
1010 : */
1011 297525 : [[nodiscard]] bool havePRefinement() const { return _have_p_refinement; }
1012 :
1013 : /**
1014 : * Set the current lower dimensional element. This can be null
1015 : */
1016 : virtual void setCurrentLowerDElem(const Elem * const lower_d_elem, const THREAD_ID tid);
1017 :
1018 : protected:
1019 : /**
1020 : * Helper function called by getVariable that handles the logic for
1021 : * checking whether Variables of the requested type are available.
1022 : */
1023 : template <typename T>
1024 : MooseVariableFieldBase & getVariableHelper(const THREAD_ID tid,
1025 : const std::string & var_name,
1026 : Moose::VarKindType expected_var_type,
1027 : Moose::VarFieldType expected_var_field_type,
1028 : const std::vector<T> & nls,
1029 : const SystemBase & aux) const;
1030 :
1031 : /**
1032 : * Verify the integrity of _vector_tags and _typed_vector_tags
1033 : */
1034 : bool verifyVectorTags() const;
1035 :
1036 : /**
1037 : * Mark a variable family for either disabling or enabling p-refinement with valid parameters of a
1038 : * variable
1039 : */
1040 : void markFamilyPRefinement(const InputParameters & params);
1041 :
1042 : /// The currently declared tags
1043 : std::map<TagName, TagID> _matrix_tag_name_to_tag_id;
1044 :
1045 : /// Reverse map
1046 : std::map<TagID, TagName> _matrix_tag_id_to_tag_name;
1047 :
1048 : /// The Factory for building objects
1049 : Factory & _factory;
1050 :
1051 : DiracKernelInfo _dirac_kernel_info;
1052 :
1053 : /// Map of material properties (block_id -> list of properties)
1054 : std::map<SubdomainID, std::set<std::string>> _map_block_material_props;
1055 :
1056 : /// Map for boundary material properties (boundary_id -> list of properties)
1057 : std::map<BoundaryID, std::set<std::string>> _map_boundary_material_props;
1058 :
1059 : /// Set of properties returned as zero properties
1060 : std::map<SubdomainID, std::set<MaterialPropertyName>> _zero_block_material_props;
1061 : std::map<BoundaryID, std::set<MaterialPropertyName>> _zero_boundary_material_props;
1062 :
1063 : /// set containing all material property names that have been requested by getMaterialProperty*
1064 : std::set<std::string> _material_property_requested;
1065 :
1066 : ///@{
1067 : /**
1068 : * Data structures of the requested material properties. We store them in a map
1069 : * from boundary/block id to multimap. Each of the multimaps is a list of
1070 : * requestor object names to material property names.
1071 : */
1072 : std::map<SubdomainID, std::multimap<std::string, std::string>> _map_block_material_props_check;
1073 : std::map<BoundaryID, std::multimap<std::string, std::string>> _map_boundary_material_props_check;
1074 : ///@}
1075 :
1076 : /// This is the set of MooseVariableFieldBase that will actually get reinited by a call to reinit(elem)
1077 : std::vector<std::set<MooseVariableFieldBase *>> _active_elemental_moose_variables;
1078 :
1079 : /// Whether or not there is currently a list of active elemental moose variables
1080 : /* This needs to remain <unsigned int> for threading purposes */
1081 : std::vector<unsigned int> _has_active_elemental_moose_variables;
1082 :
1083 : std::vector<std::set<TagID>> _active_fe_var_coupleable_matrix_tags;
1084 :
1085 : std::vector<std::set<TagID>> _active_fe_var_coupleable_vector_tags;
1086 :
1087 : std::vector<std::set<TagID>> _active_sc_var_coupleable_matrix_tags;
1088 :
1089 : std::vector<std::set<TagID>> _active_sc_var_coupleable_vector_tags;
1090 :
1091 : /// Whether or not to use default libMesh coupling
1092 : bool _default_ghosting;
1093 :
1094 : /// Elements that should have Dofs ghosted to the local processor
1095 : std::set<dof_id_type> _ghosted_elems;
1096 :
1097 : /// Flag to determine whether the problem is currently computing Jacobian
1098 : bool _currently_computing_jacobian;
1099 :
1100 : /// Flag to determine whether the problem is currently computing the residual and Jacobian
1101 : bool _currently_computing_residual_and_jacobian;
1102 :
1103 : /// Whether the non-linear residual is being evaluated
1104 : bool _computing_nonlinear_residual;
1105 :
1106 : /// Whether the residual is being evaluated
1107 : bool _currently_computing_residual;
1108 :
1109 : /// Is it safe to retrieve data from tagged matrices
1110 : bool _safe_access_tagged_matrices;
1111 :
1112 : /// Is it safe to retrieve data from tagged vectors
1113 : bool _safe_access_tagged_vectors;
1114 :
1115 : /// AD flag indicating whether **any** AD objects have been added
1116 : bool _have_ad_objects;
1117 :
1118 : /// the list of vector tags that will not be zeroed when all other tags are
1119 : std::unordered_set<TagID> _not_zeroed_tagged_vectors;
1120 :
1121 : private:
1122 : /**
1123 : * @return whether a given variable name is in the solver systems (reflected by the first
1124 : * member of the returned pair which is a boolean) and if so, what solver system number it is
1125 : * in (the second member of the returned pair; if the variable is not in the solver systems,
1126 : * then this will be an invalid unsigned integer)
1127 : */
1128 : virtual std::pair<bool, unsigned int>
1129 : determineSolverSystem(const std::string & var_name, bool error_if_not_found = false) const = 0;
1130 :
1131 : enum class TrueFunctorIs
1132 : {
1133 : UNSET,
1134 : NONAD,
1135 : AD
1136 : };
1137 :
1138 : /// A container holding pointers to all the functors in our problem. We hold a tuple where the
1139 : /// zeroth item in the tuple is an enumerator that describes what type of functor the "true"
1140 : /// functor is (either NONAD or AD), the first item in the tuple is the non-AD version of the
1141 : /// functor, and the second item in the tuple is the AD version of the functor
1142 : std::vector<std::multimap<std::string,
1143 : std::tuple<TrueFunctorIs,
1144 : std::unique_ptr<Moose::FunctorEnvelopeBase>,
1145 : std::unique_ptr<Moose::FunctorEnvelopeBase>>>>
1146 : _functors;
1147 :
1148 : /// Container to hold PiecewiseByBlockLambdaFunctors
1149 : std::vector<std::map<std::string, std::unique_ptr<Moose::FunctorAbstract>>> _pbblf_functors;
1150 :
1151 : /// Lists all functors in the problem
1152 : void showFunctors() const;
1153 :
1154 : /// Lists all functors and all the objects that requested them
1155 : void showFunctorRequestors() const;
1156 :
1157 : /// The requestors of functors where the key is the prop name and the value is a set of names of
1158 : /// requestors
1159 : std::map<std::string, std::set<std::string>> _functor_to_requestors;
1160 :
1161 : /// A multimap (for each thread) from unfilled functor requests to whether the requests were for
1162 : /// AD functors and whether the requestor was an AD object
1163 : std::vector<std::multimap<std::string, std::pair<bool, bool>>> _functor_to_request_info;
1164 :
1165 : /// Whether to output a list of the functors used and requested (currently only at initialSetup)
1166 : bool _show_functors;
1167 :
1168 : /// Whether to output a list of all the chain control data
1169 : bool _show_chain_control_data;
1170 :
1171 : /// The declared vector tags
1172 : std::vector<VectorTag> _vector_tags;
1173 :
1174 : /**
1175 : * The vector tags associated with each VectorTagType
1176 : * This is kept separate from _vector_tags for quick access into typed vector tags in places where
1177 : * we don't want to build a new vector every call (like in residual evaluation)
1178 : */
1179 : std::vector<std::vector<VectorTag>> _typed_vector_tags;
1180 :
1181 : /// Map of vector tag TagName to TagID
1182 : std::map<TagName, TagID> _vector_tags_name_map;
1183 :
1184 : ///@{ Helper functions for checking MaterialProperties
1185 : std::string restrictionSubdomainCheckName(SubdomainID check_id);
1186 : std::string restrictionBoundaryCheckName(BoundaryID check_id);
1187 : ///@}
1188 :
1189 : // Contains properties consumed by objects, see addConsumedPropertyName
1190 : std::map<MooseObjectName, std::set<std::string>> _consumed_material_properties;
1191 :
1192 : /// A map from a root algebraic ghosting functor, e.g. the ghosting functor passed into \p
1193 : /// removeAlgebraicGhostingFunctor, to its clones in other systems, e.g. systems other than system
1194 : /// 0
1195 : std::unordered_map<libMesh::GhostingFunctor *,
1196 : std::vector<std::shared_ptr<libMesh::GhostingFunctor>>>
1197 : _root_alg_gf_to_sys_clones;
1198 :
1199 : /// A map from a root coupling ghosting functor, e.g. the ghosting functor passed into \p
1200 : /// removeCouplingGhostingFunctor, to its clones in other systems, e.g. systems other than system
1201 : /// 0
1202 : std::unordered_map<libMesh::GhostingFunctor *,
1203 : std::vector<std::shared_ptr<libMesh::GhostingFunctor>>>
1204 : _root_coupling_gf_to_sys_clones;
1205 :
1206 : /// Whether p-refinement has been requested at any point during the simulation
1207 : bool _have_p_refinement;
1208 :
1209 : /// Indicate whether a family is disabled for p-refinement
1210 : std::unordered_map<FEFamily, bool> _family_for_p_refinement;
1211 : /// The set of variable families by default disable p-refinement
1212 : static const std::unordered_set<FEFamily> _default_families_without_p_refinement;
1213 :
1214 : friend class Restartable;
1215 : };
1216 :
1217 : template <typename T>
1218 : const Moose::Functor<T> &
1219 23768 : SubProblem::getFunctor(const std::string & name,
1220 : const THREAD_ID tid,
1221 : const std::string & requestor_name,
1222 : const bool requestor_is_ad)
1223 : {
1224 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1225 :
1226 : // Log the requestor
1227 23768 : _functor_to_requestors["wraps_" + name].insert(requestor_name);
1228 :
1229 23768 : constexpr bool requested_functor_is_ad =
1230 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1231 :
1232 23768 : auto & functor_to_request_info = _functor_to_request_info[tid];
1233 :
1234 : // Get the requested functor if we already have it
1235 23768 : auto & functors = _functors[tid];
1236 23768 : if (auto find_ret = functors.find("wraps_" + name); find_ret != functors.end())
1237 : {
1238 23621 : if (functors.count("wraps_" + name) > 1)
1239 0 : mooseError("Attempted to get a functor with the name '",
1240 : name,
1241 : "' but multiple (" + std::to_string(functors.count("wraps_" + name)) +
1242 : ") functors match. Make sure that you do not have functor material "
1243 : "properties, functions, postprocessors or variables with the same names.");
1244 :
1245 23621 : auto & [true_functor_is, non_ad_functor, ad_functor] = find_ret->second;
1246 23621 : auto & functor_wrapper = requested_functor_is_ad ? *ad_functor : *non_ad_functor;
1247 :
1248 23621 : auto * const functor = dynamic_cast<Moose::Functor<T> *>(&functor_wrapper);
1249 23621 : if (!functor)
1250 0 : mooseError("A call to SubProblem::getFunctor requested a functor named '",
1251 : name,
1252 : "' that returns the type: '",
1253 : libMesh::demangle(typeid(T).name()),
1254 : "'. However, that functor already exists and returns a different type: '",
1255 0 : functor_wrapper.returnType(),
1256 : "'");
1257 :
1258 23621 : if (functor->template wrapsType<Moose::NullFunctor<T>>())
1259 : // Store for future checking when the actual functor gets added
1260 266 : functor_to_request_info.emplace(name,
1261 532 : std::make_pair(requested_functor_is_ad, requestor_is_ad));
1262 : else
1263 : {
1264 : // We already have the actual functor
1265 23355 : if (true_functor_is == SubProblem::TrueFunctorIs::UNSET)
1266 0 : mooseError("We already have the functor; it should not be unset");
1267 :
1268 : // Check for whether this is a valid request
1269 : // We allow auxiliary variables and linear variables to be retrieved as non AD
1270 20099 : if (!requested_functor_is_ad && requestor_is_ad &&
1271 9413 : true_functor_is == SubProblem::TrueFunctorIs::AD &&
1272 4 : !(hasAuxiliaryVariable(name) || hasLinearVariable(name)))
1273 4 : mooseError("The AD object '",
1274 : requestor_name,
1275 : "' is requesting the functor '",
1276 : name,
1277 : "' as a non-AD functor even though it is truly an AD functor, which is not "
1278 : "allowed, since this may unintentionally drop derivatives.");
1279 : }
1280 :
1281 23617 : return *functor;
1282 : }
1283 :
1284 : // We don't have the functor yet but we could have it in the future. We'll create null functors
1285 : // for now
1286 147 : functor_to_request_info.emplace(name, std::make_pair(requested_functor_is_ad, requestor_is_ad));
1287 : if constexpr (requested_functor_is_ad)
1288 : {
1289 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1290 : typedef T ADType;
1291 :
1292 : auto emplace_ret =
1293 141 : functors.emplace("wraps_" + name,
1294 94 : std::make_tuple(SubProblem::TrueFunctorIs::UNSET,
1295 : std::make_unique<Moose::Functor<NonADType>>(
1296 : std::make_unique<Moose::NullFunctor<NonADType>>()),
1297 : std::make_unique<Moose::Functor<ADType>>(
1298 : std::make_unique<Moose::NullFunctor<ADType>>())));
1299 :
1300 47 : return static_cast<Moose::Functor<T> &>(*(requested_functor_is_ad
1301 47 : ? std::get<2>(emplace_ret->second)
1302 47 : : std::get<1>(emplace_ret->second)));
1303 : }
1304 : else
1305 : {
1306 : typedef T NonADType;
1307 : typedef typename Moose::ADType<T>::type ADType;
1308 :
1309 : auto emplace_ret =
1310 300 : functors.emplace("wraps_" + name,
1311 200 : std::make_tuple(SubProblem::TrueFunctorIs::UNSET,
1312 : std::make_unique<Moose::Functor<NonADType>>(
1313 : std::make_unique<Moose::NullFunctor<NonADType>>()),
1314 : std::make_unique<Moose::Functor<ADType>>(
1315 : std::make_unique<Moose::NullFunctor<ADType>>())));
1316 :
1317 100 : return static_cast<Moose::Functor<T> &>(*(requested_functor_is_ad
1318 : ? std::get<2>(emplace_ret->second)
1319 200 : : std::get<1>(emplace_ret->second)));
1320 : }
1321 : }
1322 :
1323 : template <typename T>
1324 : bool
1325 208 : SubProblem::hasFunctorWithType(const std::string & name, const THREAD_ID tid) const
1326 : {
1327 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1328 208 : auto & functors = _functors[tid];
1329 :
1330 208 : const auto & it = functors.find("wraps_" + name);
1331 208 : constexpr bool requested_functor_is_ad =
1332 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1333 :
1334 208 : if (it == functors.end())
1335 0 : return false;
1336 : else
1337 208 : return dynamic_cast<Moose::Functor<T> *>(
1338 416 : requested_functor_is_ad ? std::get<2>(it->second).get() : std::get<1>(it->second).get());
1339 : }
1340 :
1341 : template <typename T, typename PolymorphicLambda>
1342 : const Moose::FunctorBase<T> &
1343 9175 : SubProblem::addPiecewiseByBlockLambdaFunctor(const std::string & name,
1344 : PolymorphicLambda my_lammy,
1345 : const std::set<ExecFlagType> & clearance_schedule,
1346 : const MooseMesh & mesh,
1347 : const std::set<SubdomainID> & block_ids,
1348 : const THREAD_ID tid)
1349 : {
1350 9175 : auto & pbblf_functors = _pbblf_functors[tid];
1351 :
1352 9175 : auto [it, first_time_added] =
1353 9175 : pbblf_functors.emplace(name,
1354 : std::make_unique<PiecewiseByBlockLambdaFunctor<T>>(
1355 : name, my_lammy, clearance_schedule, mesh, block_ids));
1356 :
1357 9175 : auto * functor = dynamic_cast<PiecewiseByBlockLambdaFunctor<T> *>(it->second.get());
1358 9175 : if (!functor)
1359 : {
1360 0 : if (first_time_added)
1361 0 : mooseError("This should be impossible. If this was the first time we added the functor, then "
1362 : "the dynamic cast absolutely should have succeeded");
1363 : else
1364 0 : mooseError("Attempted to add a lambda functor with the name '",
1365 : name,
1366 : "' but another lambda functor of that name returns a different type");
1367 : }
1368 :
1369 9175 : if (first_time_added)
1370 8136 : addFunctor(name, *functor, tid);
1371 : else
1372 : // The functor already exists
1373 1039 : functor->setFunctor(mesh, block_ids, my_lammy);
1374 :
1375 9171 : return *functor;
1376 : }
1377 :
1378 : template <typename T>
1379 : void
1380 326541 : SubProblem::addFunctor(const std::string & name,
1381 : const Moose::FunctorBase<T> & functor,
1382 : const THREAD_ID tid)
1383 : {
1384 326541 : constexpr bool added_functor_is_ad =
1385 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1386 :
1387 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1388 :
1389 326541 : auto & functor_to_request_info = _functor_to_request_info[tid];
1390 326541 : auto & functors = _functors[tid];
1391 326541 : auto it = functors.find("wraps_" + name);
1392 326541 : if (it != functors.end())
1393 : {
1394 : // We have this functor already. If it's a null functor, we want to replace it with the valid
1395 : // functor we have now. If it's not then we'll add a new entry into the multimap and then we'll
1396 : // error later if a user requests a functor because their request is ambiguous. This is the
1397 : // reason that the functors container is a multimap: for nice error messages
1398 : auto * const existing_wrapper_base =
1399 6938 : added_functor_is_ad ? std::get<2>(it->second).get() : std::get<1>(it->second).get();
1400 6938 : auto * const existing_wrapper = dynamic_cast<Moose::Functor<T> *>(existing_wrapper_base);
1401 6938 : if (existing_wrapper && existing_wrapper->template wrapsType<Moose::NullFunctor<T>>())
1402 : {
1403 : // Sanity check
1404 147 : auto [request_info_it, request_info_end_it] = functor_to_request_info.equal_range(name);
1405 147 : if (request_info_it == request_info_end_it)
1406 0 : mooseError("We are wrapping a NullFunctor but we don't have any unfilled functor request "
1407 : "info. This doesn't make sense.");
1408 :
1409 : // Check for valid requests
1410 560 : while (request_info_it != request_info_end_it)
1411 : {
1412 413 : auto & [requested_functor_is_ad, requestor_is_ad] = request_info_it->second;
1413 5 : if (!requested_functor_is_ad && requestor_is_ad && added_functor_is_ad)
1414 0 : mooseError("We are requesting a non-AD functor '" + name +
1415 : "' from an AD object, but the true functor is AD. This means we could be "
1416 : "dropping important derivatives. We will not allow this");
1417 : // We're going to eventually check whether we've fulfilled all functor requests and our
1418 : // check will be that the multimap is empty. This request is fulfilled, so erase it from the
1419 : // map now
1420 413 : request_info_it = functor_to_request_info.erase(request_info_it);
1421 : }
1422 :
1423 : // Ok we didn't have the functor before, so we will add it now
1424 147 : std::get<0>(it->second) =
1425 : added_functor_is_ad ? SubProblem::TrueFunctorIs::AD : SubProblem::TrueFunctorIs::NONAD;
1426 147 : existing_wrapper->assign(functor);
1427 : // Finally we create the non-AD or AD complement of the just added functor
1428 : if constexpr (added_functor_is_ad)
1429 : {
1430 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1431 5 : auto * const existing_non_ad_wrapper_base = std::get<1>(it->second).get();
1432 5 : auto * const existing_non_ad_wrapper =
1433 5 : dynamic_cast<Moose::Functor<NonADType> *>(existing_non_ad_wrapper_base);
1434 : mooseAssert(existing_non_ad_wrapper->template wrapsType<Moose::NullFunctor<NonADType>>(),
1435 : "Both members of pair should have been wrapping a NullFunctor");
1436 5 : existing_non_ad_wrapper->assign(
1437 : std::make_unique<Moose::RawValueFunctor<NonADType>>(functor));
1438 : }
1439 : else
1440 : {
1441 : typedef typename Moose::ADType<T>::type ADType;
1442 142 : auto * const existing_ad_wrapper_base = std::get<2>(it->second).get();
1443 142 : auto * const existing_ad_wrapper =
1444 142 : dynamic_cast<Moose::Functor<ADType> *>(existing_ad_wrapper_base);
1445 : mooseAssert(existing_ad_wrapper->template wrapsType<Moose::NullFunctor<ADType>>(),
1446 : "Both members of pair should have been wrapping a NullFunctor");
1447 142 : existing_ad_wrapper->assign(std::make_unique<Moose::ADWrapperFunctor<ADType>>(functor));
1448 : }
1449 147 : return;
1450 : }
1451 6791 : else if (!existing_wrapper)
1452 : {
1453 : // Functor was emplaced but the cast failed. This could be a double definition with
1454 : // different types, or it could be a request with one type then a definition with another
1455 : // type. Either way it is going to error later, but it is cleaner to catch it now
1456 0 : mooseError("Functor '",
1457 : name,
1458 : "' is being added with return type '",
1459 : MooseUtils::prettyCppType<T>(),
1460 : "' but it has already been defined or requested with return type '",
1461 0 : existing_wrapper_base->returnType(),
1462 : "'.");
1463 : }
1464 : }
1465 :
1466 : // We are a new functor, create the opposite ADType one and store it with other functors
1467 : if constexpr (added_functor_is_ad)
1468 : {
1469 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1470 196951 : auto new_non_ad_wrapper = std::make_unique<Moose::Functor<NonADType>>(
1471 : std::make_unique<Moose::RawValueFunctor<NonADType>>(functor));
1472 196951 : auto new_ad_wrapper = std::make_unique<Moose::Functor<T>>(functor);
1473 393902 : _functors[tid].emplace("wraps_" + name,
1474 196951 : std::make_tuple(SubProblem::TrueFunctorIs::AD,
1475 196951 : std::move(new_non_ad_wrapper),
1476 196951 : std::move(new_ad_wrapper)));
1477 196951 : }
1478 : else
1479 : {
1480 : typedef typename Moose::ADType<T>::type ADType;
1481 129443 : auto new_non_ad_wrapper = std::make_unique<Moose::Functor<T>>((functor));
1482 129443 : auto new_ad_wrapper = std::make_unique<Moose::Functor<ADType>>(
1483 : std::make_unique<Moose::ADWrapperFunctor<ADType>>(functor));
1484 258886 : _functors[tid].emplace("wraps_" + name,
1485 129443 : std::make_tuple(SubProblem::TrueFunctorIs::NONAD,
1486 129443 : std::move(new_non_ad_wrapper),
1487 129443 : std::move(new_ad_wrapper)));
1488 129443 : }
1489 : }
1490 :
1491 : inline const bool &
1492 298732 : SubProblem::currentlyComputingResidualAndJacobian() const
1493 : {
1494 298732 : return _currently_computing_residual_and_jacobian;
1495 : }
1496 :
1497 : inline void
1498 4052010 : SubProblem::setCurrentlyComputingResidualAndJacobian(
1499 : const bool currently_computing_residual_and_jacobian)
1500 : {
1501 4052010 : _currently_computing_residual_and_jacobian = currently_computing_residual_and_jacobian;
1502 4052010 : }
1503 :
1504 : namespace Moose
1505 : {
1506 : void initial_condition(libMesh::EquationSystems & es, const std::string & system_name);
1507 : } // namespace Moose
|
