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 310072 : 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 121478 : 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 3146819739 : 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 9032361 : 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 487556 : 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 4322272 : 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 4322272 : return const_cast<MooseVariableFieldBase &>(const_cast<const SubProblem *>(this)->getVariable(
285 4322264 : 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 71184 : 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 402381613 : 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 3704828 : void setCurrentlyComputingJacobian(const bool currently_computing_jacobian)
690 : {
691 3704828 : _currently_computing_jacobian = currently_computing_jacobian;
692 3704828 : }
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 21781 : 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 197082 : virtual void computingNonlinearResid(const bool computing_nonlinear_residual)
713 : {
714 197082 : _computing_nonlinear_residual = computing_nonlinear_residual;
715 197082 : }
716 :
717 : /**
718 : * Returns true if the problem is in the process of computing the residual
719 : */
720 66739 : 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 548912 : virtual void setCurrentlyComputingResidual(const bool currently_computing_residual)
726 : {
727 548912 : _currently_computing_residual = currently_computing_residual;
728 548912 : }
729 :
730 : /// Is it safe to access the tagged matrices
731 1111507935 : virtual bool safeAccessTaggedMatrices() const { return _safe_access_tagged_matrices; }
732 :
733 : /// Is it safe to access the tagged vectors
734 588273 : 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 522 : 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 407560380 : 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 34 : void setFunctorOutput(bool set_output) { _output_functors = set_output; }
925 :
926 : /**
927 : * @return the number of nonlinear systems in the problem
928 : */
929 : virtual std::size_t numNonlinearSystems() const = 0;
930 :
931 : /**
932 : * @return the current nonlinear system number
933 : */
934 : virtual unsigned int currentNlSysNum() const = 0;
935 :
936 : /**
937 : * @return the number of linear systems in the problem
938 : */
939 : virtual std::size_t numLinearSystems() const = 0;
940 :
941 : /**
942 : * @return the number of solver systems in the problem
943 : */
944 : virtual std::size_t numSolverSystems() const = 0;
945 :
946 : /**
947 : * @return the current linear system number
948 : */
949 : virtual unsigned int currentLinearSysNum() const = 0;
950 :
951 : /**
952 : * Register an unfulfilled functor request
953 : */
954 : template <typename T>
955 : void registerUnfilledFunctorRequest(T * functor_interface,
956 : const std::string & functor_name,
957 : const THREAD_ID tid);
958 :
959 : /**
960 : * Return the residual vector tags we are currently computing
961 : */
962 : virtual const std::vector<VectorTag> & currentResidualVectorTags() const = 0;
963 :
964 : /**
965 : * Select the vector tags which belong to a specific system
966 : * @param system Reference to the system
967 : * @param input_vector_tags A vector of vector tags
968 : * @param selected_tags A set which gets populated by the tag-ids that belong to the system
969 : */
970 : static void selectVectorTagsFromSystem(const SystemBase & system,
971 : const std::vector<VectorTag> & input_vector_tags,
972 : std::set<TagID> & selected_tags);
973 :
974 : /**
975 : * Select the matrix tags which belong to a specific system
976 : * @param system Reference to the system
977 : * @param input_matrix_tags A map of matrix tags
978 : * @param selected_tags A set which gets populated by the tag-ids that belong to the system
979 : */
980 : static void selectMatrixTagsFromSystem(const SystemBase & system,
981 : const std::map<TagName, TagID> & input_matrix_tags,
982 : std::set<TagID> & selected_tags);
983 :
984 : /**
985 : * reinitialize the finite volume assembly data for the provided face and thread
986 : */
987 : void reinitFVFace(const THREAD_ID tid, const FaceInfo & fi);
988 :
989 : /**
990 : * Whether the simulation has active nonlocal coupling which should be accounted for in the
991 : * Jacobian. For this to return true, there must be at least one active nonlocal kernel or
992 : * boundary condition
993 : */
994 : virtual bool hasNonlocalCoupling() const = 0;
995 :
996 : /**
997 : * Prepare \p DofMap and \p Assembly classes with our p-refinement information
998 : */
999 : void preparePRefinement();
1000 :
1001 : /**
1002 : * @returns whether we're doing p-refinement
1003 : */
1004 : [[nodiscard]] bool doingPRefinement() const;
1005 :
1006 : /**
1007 : * Query whether p-refinement has been requested at any point during the simulation
1008 : */
1009 272744 : [[nodiscard]] bool havePRefinement() const { return _have_p_refinement; }
1010 :
1011 : /**
1012 : * Set the current lower dimensional element. This can be null
1013 : */
1014 : virtual void setCurrentLowerDElem(const Elem * const lower_d_elem, const THREAD_ID tid);
1015 :
1016 : protected:
1017 : /**
1018 : * Helper function called by getVariable that handles the logic for
1019 : * checking whether Variables of the requested type are available.
1020 : */
1021 : template <typename T>
1022 : MooseVariableFieldBase & getVariableHelper(const THREAD_ID tid,
1023 : const std::string & var_name,
1024 : Moose::VarKindType expected_var_type,
1025 : Moose::VarFieldType expected_var_field_type,
1026 : const std::vector<T> & nls,
1027 : const SystemBase & aux) const;
1028 :
1029 : /**
1030 : * Verify the integrity of _vector_tags and _typed_vector_tags
1031 : */
1032 : bool verifyVectorTags() const;
1033 :
1034 : /**
1035 : * Mark a variable family for either disabling or enabling p-refinement with valid parameters of a
1036 : * variable
1037 : */
1038 : void markFamilyPRefinement(const InputParameters & params);
1039 :
1040 : /// The currently declared tags
1041 : std::map<TagName, TagID> _matrix_tag_name_to_tag_id;
1042 :
1043 : /// Reverse map
1044 : std::map<TagID, TagName> _matrix_tag_id_to_tag_name;
1045 :
1046 : /// The Factory for building objects
1047 : Factory & _factory;
1048 :
1049 : DiracKernelInfo _dirac_kernel_info;
1050 :
1051 : /// Map of material properties (block_id -> list of properties)
1052 : std::map<SubdomainID, std::set<std::string>> _map_block_material_props;
1053 :
1054 : /// Map for boundary material properties (boundary_id -> list of properties)
1055 : std::map<BoundaryID, std::set<std::string>> _map_boundary_material_props;
1056 :
1057 : /// Set of properties returned as zero properties
1058 : std::map<SubdomainID, std::set<MaterialPropertyName>> _zero_block_material_props;
1059 : std::map<BoundaryID, std::set<MaterialPropertyName>> _zero_boundary_material_props;
1060 :
1061 : /// set containing all material property names that have been requested by getMaterialProperty*
1062 : std::set<std::string> _material_property_requested;
1063 :
1064 : ///@{
1065 : /**
1066 : * Data structures of the requested material properties. We store them in a map
1067 : * from boundary/block id to multimap. Each of the multimaps is a list of
1068 : * requestor object names to material property names.
1069 : */
1070 : std::map<SubdomainID, std::multimap<std::string, std::string>> _map_block_material_props_check;
1071 : std::map<BoundaryID, std::multimap<std::string, std::string>> _map_boundary_material_props_check;
1072 : ///@}
1073 :
1074 : /// This is the set of MooseVariableFieldBase that will actually get reinited by a call to reinit(elem)
1075 : std::vector<std::set<MooseVariableFieldBase *>> _active_elemental_moose_variables;
1076 :
1077 : /// Whether or not there is currently a list of active elemental moose variables
1078 : /* This needs to remain <unsigned int> for threading purposes */
1079 : std::vector<unsigned int> _has_active_elemental_moose_variables;
1080 :
1081 : std::vector<std::set<TagID>> _active_fe_var_coupleable_matrix_tags;
1082 :
1083 : std::vector<std::set<TagID>> _active_fe_var_coupleable_vector_tags;
1084 :
1085 : std::vector<std::set<TagID>> _active_sc_var_coupleable_matrix_tags;
1086 :
1087 : std::vector<std::set<TagID>> _active_sc_var_coupleable_vector_tags;
1088 :
1089 : /// Whether or not to use default libMesh coupling
1090 : bool _default_ghosting;
1091 :
1092 : /// Elements that should have Dofs ghosted to the local processor
1093 : std::set<dof_id_type> _ghosted_elems;
1094 :
1095 : /// Flag to determine whether the problem is currently computing Jacobian
1096 : bool _currently_computing_jacobian;
1097 :
1098 : /// Flag to determine whether the problem is currently computing the residual and Jacobian
1099 : bool _currently_computing_residual_and_jacobian;
1100 :
1101 : /// Whether the non-linear residual is being evaluated
1102 : bool _computing_nonlinear_residual;
1103 :
1104 : /// Whether the residual is being evaluated
1105 : bool _currently_computing_residual;
1106 :
1107 : /// Is it safe to retrieve data from tagged matrices
1108 : bool _safe_access_tagged_matrices;
1109 :
1110 : /// Is it safe to retrieve data from tagged vectors
1111 : bool _safe_access_tagged_vectors;
1112 :
1113 : /// AD flag indicating whether **any** AD objects have been added
1114 : bool _have_ad_objects;
1115 :
1116 : /// the list of vector tags that will not be zeroed when all other tags are
1117 : std::unordered_set<TagID> _not_zeroed_tagged_vectors;
1118 :
1119 : private:
1120 : /**
1121 : * @return whether a given variable name is in the solver systems (reflected by the first
1122 : * member of the returned pair which is a boolean) and if so, what solver system number it is
1123 : * in (the second member of the returned pair; if the variable is not in the solver systems,
1124 : * then this will be an invalid unsigned integer)
1125 : */
1126 : virtual std::pair<bool, unsigned int>
1127 : determineSolverSystem(const std::string & var_name, bool error_if_not_found = false) const = 0;
1128 :
1129 : enum class TrueFunctorIs
1130 : {
1131 : UNSET,
1132 : NONAD,
1133 : AD
1134 : };
1135 :
1136 : /// A container holding pointers to all the functors in our problem. We hold a tuple where the
1137 : /// zeroth item in the tuple is an enumerator that describes what type of functor the "true"
1138 : /// functor is (either NONAD or AD), the first item in the tuple is the non-AD version of the
1139 : /// functor, and the second item in the tuple is the AD version of the functor
1140 : std::vector<std::multimap<std::string,
1141 : std::tuple<TrueFunctorIs,
1142 : std::unique_ptr<Moose::FunctorEnvelopeBase>,
1143 : std::unique_ptr<Moose::FunctorEnvelopeBase>>>>
1144 : _functors;
1145 :
1146 : /// Container to hold PiecewiseByBlockLambdaFunctors
1147 : std::vector<std::map<std::string, std::unique_ptr<Moose::FunctorAbstract>>> _pbblf_functors;
1148 :
1149 : /// Lists all functors in the problem
1150 : void showFunctors() const;
1151 :
1152 : /// Lists all functors and all the objects that requested them
1153 : void showFunctorRequestors() const;
1154 :
1155 : /// The requestors of functors where the key is the prop name and the value is a set of names of
1156 : /// requestors
1157 : std::map<std::string, std::set<std::string>> _functor_to_requestors;
1158 :
1159 : /// A multimap (for each thread) from unfilled functor requests to whether the requests were for
1160 : /// AD functors and whether the requestor was an AD object
1161 : std::vector<std::multimap<std::string, std::pair<bool, bool>>> _functor_to_request_info;
1162 :
1163 : /// Whether to output a list of the functors used and requested (currently only at initialSetup)
1164 : bool _output_functors;
1165 :
1166 : /// The declared vector tags
1167 : std::vector<VectorTag> _vector_tags;
1168 :
1169 : /**
1170 : * The vector tags associated with each VectorTagType
1171 : * This is kept separate from _vector_tags for quick access into typed vector tags in places where
1172 : * we don't want to build a new vector every call (like in residual evaluation)
1173 : */
1174 : std::vector<std::vector<VectorTag>> _typed_vector_tags;
1175 :
1176 : /// Map of vector tag TagName to TagID
1177 : std::map<TagName, TagID> _vector_tags_name_map;
1178 :
1179 : ///@{ Helper functions for checking MaterialProperties
1180 : std::string restrictionSubdomainCheckName(SubdomainID check_id);
1181 : std::string restrictionBoundaryCheckName(BoundaryID check_id);
1182 : ///@}
1183 :
1184 : // Contains properties consumed by objects, see addConsumedPropertyName
1185 : std::map<MooseObjectName, std::set<std::string>> _consumed_material_properties;
1186 :
1187 : /// A map from a root algebraic ghosting functor, e.g. the ghosting functor passed into \p
1188 : /// removeAlgebraicGhostingFunctor, to its clones in other systems, e.g. systems other than system
1189 : /// 0
1190 : std::unordered_map<libMesh::GhostingFunctor *,
1191 : std::vector<std::shared_ptr<libMesh::GhostingFunctor>>>
1192 : _root_alg_gf_to_sys_clones;
1193 :
1194 : /// A map from a root coupling ghosting functor, e.g. the ghosting functor passed into \p
1195 : /// removeCouplingGhostingFunctor, to its clones in other systems, e.g. systems other than system
1196 : /// 0
1197 : std::unordered_map<libMesh::GhostingFunctor *,
1198 : std::vector<std::shared_ptr<libMesh::GhostingFunctor>>>
1199 : _root_coupling_gf_to_sys_clones;
1200 :
1201 : /// Whether p-refinement has been requested at any point during the simulation
1202 : bool _have_p_refinement;
1203 :
1204 : /// Indicate whether a family is disabled for p-refinement
1205 : std::unordered_map<FEFamily, bool> _family_for_p_refinement;
1206 : /// The set of variable families by default disable p-refinement
1207 : static const std::unordered_set<FEFamily> _default_families_without_p_refinement;
1208 :
1209 : friend class Restartable;
1210 : };
1211 :
1212 : template <typename T>
1213 : const Moose::Functor<T> &
1214 21844 : SubProblem::getFunctor(const std::string & name,
1215 : const THREAD_ID tid,
1216 : const std::string & requestor_name,
1217 : const bool requestor_is_ad)
1218 : {
1219 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1220 :
1221 : // Log the requestor
1222 21844 : _functor_to_requestors["wraps_" + name].insert(requestor_name);
1223 :
1224 21844 : constexpr bool requested_functor_is_ad =
1225 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1226 :
1227 21844 : auto & functor_to_request_info = _functor_to_request_info[tid];
1228 :
1229 : // Get the requested functor if we already have it
1230 21844 : auto & functors = _functors[tid];
1231 21844 : if (auto find_ret = functors.find("wraps_" + name); find_ret != functors.end())
1232 : {
1233 21706 : if (functors.count("wraps_" + name) > 1)
1234 0 : mooseError("Attempted to get a functor with the name '",
1235 : name,
1236 : "' but multiple (" + std::to_string(functors.count("wraps_" + name)) +
1237 : ") functors match. Make sure that you do not have functor material "
1238 : "properties, functions, postprocessors or variables with the same names.");
1239 :
1240 21706 : auto & [true_functor_is, non_ad_functor, ad_functor] = find_ret->second;
1241 21706 : auto & functor_wrapper = requested_functor_is_ad ? *ad_functor : *non_ad_functor;
1242 :
1243 21706 : auto * const functor = dynamic_cast<Moose::Functor<T> *>(&functor_wrapper);
1244 21706 : if (!functor)
1245 0 : mooseError("A call to SubProblem::getFunctor requested a functor named '",
1246 : name,
1247 : "' that returns the type: '",
1248 : libMesh::demangle(typeid(T).name()),
1249 : "'. However, that functor already exists and returns a different type: '",
1250 0 : functor_wrapper.returnType(),
1251 : "'");
1252 :
1253 21706 : if (functor->template wrapsType<Moose::NullFunctor<T>>())
1254 : // Store for future checking when the actual functor gets added
1255 247 : functor_to_request_info.emplace(name,
1256 494 : std::make_pair(requested_functor_is_ad, requestor_is_ad));
1257 : else
1258 : {
1259 : // We already have the actual functor
1260 21459 : if (true_functor_is == SubProblem::TrueFunctorIs::UNSET)
1261 0 : mooseError("We already have the functor; it should not be unset");
1262 :
1263 : // Check for whether this is a valid request
1264 : // We allow auxiliary variables and linear variables to be retrieved as non AD
1265 18958 : if (!requested_functor_is_ad && requestor_is_ad &&
1266 8861 : true_functor_is == SubProblem::TrueFunctorIs::AD &&
1267 4 : !(hasAuxiliaryVariable(name) || hasLinearVariable(name)))
1268 4 : mooseError("The AD object '",
1269 : requestor_name,
1270 : "' is requesting the functor '",
1271 : name,
1272 : "' as a non-AD functor even though it is truly an AD functor, which is not "
1273 : "allowed, since this may unintentionally drop derivatives.");
1274 : }
1275 :
1276 21702 : return *functor;
1277 : }
1278 :
1279 : // We don't have the functor yet but we could have it in the future. We'll create null functors
1280 : // for now
1281 138 : functor_to_request_info.emplace(name, std::make_pair(requested_functor_is_ad, requestor_is_ad));
1282 : if constexpr (requested_functor_is_ad)
1283 : {
1284 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1285 : typedef T ADType;
1286 :
1287 : auto emplace_ret =
1288 132 : functors.emplace("wraps_" + name,
1289 88 : std::make_tuple(SubProblem::TrueFunctorIs::UNSET,
1290 : std::make_unique<Moose::Functor<NonADType>>(
1291 : std::make_unique<Moose::NullFunctor<NonADType>>()),
1292 : std::make_unique<Moose::Functor<ADType>>(
1293 : std::make_unique<Moose::NullFunctor<ADType>>())));
1294 :
1295 44 : return static_cast<Moose::Functor<T> &>(*(requested_functor_is_ad
1296 44 : ? std::get<2>(emplace_ret->second)
1297 44 : : std::get<1>(emplace_ret->second)));
1298 : }
1299 : else
1300 : {
1301 : typedef T NonADType;
1302 : typedef typename Moose::ADType<T>::type ADType;
1303 :
1304 : auto emplace_ret =
1305 282 : functors.emplace("wraps_" + name,
1306 188 : std::make_tuple(SubProblem::TrueFunctorIs::UNSET,
1307 : std::make_unique<Moose::Functor<NonADType>>(
1308 : std::make_unique<Moose::NullFunctor<NonADType>>()),
1309 : std::make_unique<Moose::Functor<ADType>>(
1310 : std::make_unique<Moose::NullFunctor<ADType>>())));
1311 :
1312 94 : return static_cast<Moose::Functor<T> &>(*(requested_functor_is_ad
1313 : ? std::get<2>(emplace_ret->second)
1314 188 : : std::get<1>(emplace_ret->second)));
1315 : }
1316 : }
1317 :
1318 : template <typename T>
1319 : bool
1320 432 : SubProblem::hasFunctorWithType(const std::string & name, const THREAD_ID tid) const
1321 : {
1322 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1323 432 : auto & functors = _functors[tid];
1324 :
1325 432 : const auto & it = functors.find("wraps_" + name);
1326 432 : constexpr bool requested_functor_is_ad =
1327 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1328 :
1329 432 : if (it == functors.end())
1330 240 : return false;
1331 : else
1332 192 : return dynamic_cast<Moose::Functor<T> *>(
1333 384 : requested_functor_is_ad ? std::get<2>(it->second).get() : std::get<1>(it->second).get());
1334 : }
1335 :
1336 : template <typename T, typename PolymorphicLambda>
1337 : const Moose::FunctorBase<T> &
1338 8432 : SubProblem::addPiecewiseByBlockLambdaFunctor(const std::string & name,
1339 : PolymorphicLambda my_lammy,
1340 : const std::set<ExecFlagType> & clearance_schedule,
1341 : const MooseMesh & mesh,
1342 : const std::set<SubdomainID> & block_ids,
1343 : const THREAD_ID tid)
1344 : {
1345 8432 : auto & pbblf_functors = _pbblf_functors[tid];
1346 :
1347 8432 : auto [it, first_time_added] =
1348 8432 : pbblf_functors.emplace(name,
1349 : std::make_unique<PiecewiseByBlockLambdaFunctor<T>>(
1350 : name, my_lammy, clearance_schedule, mesh, block_ids));
1351 :
1352 8432 : auto * functor = dynamic_cast<PiecewiseByBlockLambdaFunctor<T> *>(it->second.get());
1353 8432 : if (!functor)
1354 : {
1355 0 : if (first_time_added)
1356 0 : mooseError("This should be impossible. If this was the first time we added the functor, then "
1357 : "the dynamic cast absolutely should have succeeded");
1358 : else
1359 0 : mooseError("Attempted to add a lambda functor with the name '",
1360 : name,
1361 : "' but another lambda functor of that name returns a different type");
1362 : }
1363 :
1364 8432 : if (first_time_added)
1365 7443 : addFunctor(name, *functor, tid);
1366 : else
1367 : // The functor already exists
1368 989 : functor->setFunctor(mesh, block_ids, my_lammy);
1369 :
1370 8428 : return *functor;
1371 : }
1372 :
1373 : template <typename T>
1374 : void
1375 295276 : SubProblem::addFunctor(const std::string & name,
1376 : const Moose::FunctorBase<T> & functor,
1377 : const THREAD_ID tid)
1378 : {
1379 295276 : constexpr bool added_functor_is_ad =
1380 : !std::is_same<T, typename MetaPhysicL::RawType<T>::value_type>::value;
1381 :
1382 : mooseAssert(tid < _functors.size(), "Too large a thread ID");
1383 :
1384 295276 : auto & functor_to_request_info = _functor_to_request_info[tid];
1385 295276 : auto & functors = _functors[tid];
1386 295276 : auto it = functors.find("wraps_" + name);
1387 295276 : if (it != functors.end())
1388 : {
1389 : // We have this functor already. If it's a null functor, we want to replace it with the valid
1390 : // functor we have now. If it's not then we'll add a new entry into the multimap and then we'll
1391 : // error later if a user requests a functor because their request is ambiguous. This is the
1392 : // reason that the functors container is a multimap: for nice error messages
1393 : auto * const existing_wrapper_base =
1394 6265 : added_functor_is_ad ? std::get<2>(it->second).get() : std::get<1>(it->second).get();
1395 6265 : auto * const existing_wrapper = dynamic_cast<Moose::Functor<T> *>(existing_wrapper_base);
1396 6265 : if (existing_wrapper && existing_wrapper->template wrapsType<Moose::NullFunctor<T>>())
1397 : {
1398 : // Sanity check
1399 138 : auto [request_info_it, request_info_end_it] = functor_to_request_info.equal_range(name);
1400 138 : if (request_info_it == request_info_end_it)
1401 0 : mooseError("We are wrapping a NullFunctor but we don't have any unfilled functor request "
1402 : "info. This doesn't make sense.");
1403 :
1404 : // Check for valid requests
1405 523 : while (request_info_it != request_info_end_it)
1406 : {
1407 385 : auto & [requested_functor_is_ad, requestor_is_ad] = request_info_it->second;
1408 5 : if (!requested_functor_is_ad && requestor_is_ad && added_functor_is_ad)
1409 0 : mooseError("We are requesting a non-AD functor '" + name +
1410 : "' from an AD object, but the true functor is AD. This means we could be "
1411 : "dropping important derivatives. We will not allow this");
1412 : // We're going to eventually check whether we've fulfilled all functor requests and our
1413 : // check will be that the multimap is empty. This request is fulfilled, so erase it from the
1414 : // map now
1415 385 : request_info_it = functor_to_request_info.erase(request_info_it);
1416 : }
1417 :
1418 : // Ok we didn't have the functor before, so we will add it now
1419 138 : std::get<0>(it->second) =
1420 : added_functor_is_ad ? SubProblem::TrueFunctorIs::AD : SubProblem::TrueFunctorIs::NONAD;
1421 138 : existing_wrapper->assign(functor);
1422 : // Finally we create the non-AD or AD complement of the just added functor
1423 : if constexpr (added_functor_is_ad)
1424 : {
1425 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1426 5 : auto * const existing_non_ad_wrapper_base = std::get<1>(it->second).get();
1427 5 : auto * const existing_non_ad_wrapper =
1428 5 : dynamic_cast<Moose::Functor<NonADType> *>(existing_non_ad_wrapper_base);
1429 : mooseAssert(existing_non_ad_wrapper->template wrapsType<Moose::NullFunctor<NonADType>>(),
1430 : "Both members of pair should have been wrapping a NullFunctor");
1431 5 : existing_non_ad_wrapper->assign(
1432 : std::make_unique<Moose::RawValueFunctor<NonADType>>(functor));
1433 : }
1434 : else
1435 : {
1436 : typedef typename Moose::ADType<T>::type ADType;
1437 133 : auto * const existing_ad_wrapper_base = std::get<2>(it->second).get();
1438 133 : auto * const existing_ad_wrapper =
1439 133 : dynamic_cast<Moose::Functor<ADType> *>(existing_ad_wrapper_base);
1440 : mooseAssert(existing_ad_wrapper->template wrapsType<Moose::NullFunctor<ADType>>(),
1441 : "Both members of pair should have been wrapping a NullFunctor");
1442 133 : existing_ad_wrapper->assign(std::make_unique<Moose::ADWrapperFunctor<ADType>>(functor));
1443 : }
1444 138 : return;
1445 : }
1446 6127 : else if (!existing_wrapper)
1447 : {
1448 : // Functor was emplaced but the cast failed. This could be a double definition with
1449 : // different types, or it could be a request with one type then a definition with another
1450 : // type. Either way it is going to error later, but it is cleaner to catch it now
1451 0 : mooseError("Functor '",
1452 : name,
1453 : "' is being added with return type '",
1454 : MooseUtils::prettyCppType<T>(),
1455 : "' but it has already been defined or requested with return type '",
1456 0 : existing_wrapper_base->returnType(),
1457 : "'.");
1458 : }
1459 : }
1460 :
1461 : // We are a new functor, create the opposite ADType one and store it with other functors
1462 : if constexpr (added_functor_is_ad)
1463 : {
1464 : typedef typename MetaPhysicL::RawType<T>::value_type NonADType;
1465 176474 : auto new_non_ad_wrapper = std::make_unique<Moose::Functor<NonADType>>(
1466 : std::make_unique<Moose::RawValueFunctor<NonADType>>(functor));
1467 176474 : auto new_ad_wrapper = std::make_unique<Moose::Functor<T>>(functor);
1468 352948 : _functors[tid].emplace("wraps_" + name,
1469 176474 : std::make_tuple(SubProblem::TrueFunctorIs::AD,
1470 176474 : std::move(new_non_ad_wrapper),
1471 176474 : std::move(new_ad_wrapper)));
1472 176474 : }
1473 : else
1474 : {
1475 : typedef typename Moose::ADType<T>::type ADType;
1476 118664 : auto new_non_ad_wrapper = std::make_unique<Moose::Functor<T>>((functor));
1477 118664 : auto new_ad_wrapper = std::make_unique<Moose::Functor<ADType>>(
1478 : std::make_unique<Moose::ADWrapperFunctor<ADType>>(functor));
1479 237328 : _functors[tid].emplace("wraps_" + name,
1480 118664 : std::make_tuple(SubProblem::TrueFunctorIs::NONAD,
1481 118664 : std::move(new_non_ad_wrapper),
1482 118664 : std::move(new_ad_wrapper)));
1483 118664 : }
1484 : }
1485 :
1486 : inline const bool &
1487 294699 : SubProblem::currentlyComputingResidualAndJacobian() const
1488 : {
1489 294699 : return _currently_computing_residual_and_jacobian;
1490 : }
1491 :
1492 : inline void
1493 3683624 : SubProblem::setCurrentlyComputingResidualAndJacobian(
1494 : const bool currently_computing_residual_and_jacobian)
1495 : {
1496 3683624 : _currently_computing_residual_and_jacobian = currently_computing_residual_and_jacobian;
1497 3683624 : }
1498 :
1499 : namespace Moose
1500 : {
1501 : void initial_condition(libMesh::EquationSystems & es, const std::string & system_name);
1502 : } // namespace Moose
|
