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 <vector>
13 :
14 : #include "DataIO.h"
15 : #include "MooseTypes.h"
16 : #include "VariableWarehouse.h"
17 : #include "InputParameters.h"
18 : #include "MooseVariableBase.h"
19 : #include "ConsoleStreamInterface.h"
20 :
21 : // libMesh
22 : #include "libmesh/exodusII_io.h"
23 : #include "libmesh/parallel_object.h"
24 : #include "libmesh/numeric_vector.h"
25 : #include "libmesh/sparse_matrix.h"
26 :
27 : // Forward declarations
28 : class Factory;
29 : class MooseApp;
30 : class MooseVariableFieldBase;
31 : template <typename>
32 : class MooseVariableFE;
33 : typedef MooseVariableFE<Real> MooseVariable;
34 : typedef MooseVariableFE<VectorValue<Real>> VectorMooseVariable;
35 : class MooseMesh;
36 : class SubProblem;
37 : class SystemBase;
38 : class TimeIntegrator;
39 : class InputParameters;
40 : class FEProblemBase;
41 :
42 : // libMesh forward declarations
43 : namespace libMesh
44 : {
45 : class System;
46 : class DofMap;
47 : class FEType;
48 : }
49 :
50 : /**
51 : * ///< Type of coordinate system
52 : */
53 : void extraSendList(std::vector<dof_id_type> & send_list, void * context);
54 :
55 : /**
56 : * Free function used for a libMesh callback
57 : */
58 : void extraSparsity(libMesh::SparsityPattern::Graph & sparsity,
59 : std::vector<dof_id_type> & n_nz,
60 : std::vector<dof_id_type> & n_oz,
61 : void * context);
62 :
63 : /**
64 : * Information about variables that will be copied
65 : */
66 : struct VarCopyInfo
67 : {
68 430 : VarCopyInfo(const std::string & dest_name,
69 : const std::string & source_name,
70 : const std::string & timestep)
71 430 : : _dest_name(dest_name), _source_name(source_name), _timestep(timestep)
72 : {
73 430 : }
74 :
75 : std::string _dest_name;
76 : std::string _source_name;
77 : std::string _timestep;
78 : };
79 :
80 : /**
81 : * Base class for a system (of equations)
82 : *
83 : */
84 : class SystemBase : public libMesh::ParallelObject, public ConsoleStreamInterface
85 :
86 : {
87 : public:
88 : SystemBase(SubProblem & subproblem,
89 : FEProblemBase & fe_problem,
90 : const std::string & name,
91 : Moose::VarKindType var_kind);
92 122512 : virtual ~SystemBase() {}
93 :
94 : /**
95 : * Gets the number of this system
96 : * @return The number of this system
97 : */
98 : unsigned int number() const;
99 44467066 : MooseMesh & mesh() { return _mesh; }
100 : const MooseMesh & mesh() const { return _mesh; }
101 168558128 : SubProblem & subproblem() { return _subproblem; }
102 1711660 : const SubProblem & subproblem() const { return _subproblem; }
103 1200 : FEProblemBase & feProblem() { return _fe_problem; }
104 : const FEProblemBase & feProblem() const { return _fe_problem; }
105 :
106 : /**
107 : * Applies scaling factors to the system's variables
108 : * @param inverse_scaling_factors A vector containing the inverse of each variable's scaling
109 : * factor, e.g. 1 / scaling_factor
110 : */
111 : void applyScalingFactors(const std::vector<Real> & inverse_scaling_factors);
112 :
113 : /**
114 : * Whether we are computing an initial Jacobian for automatic variable scaling
115 : */
116 : bool computingScalingJacobian() const;
117 :
118 : /**
119 : * Getter for whether we are performing automatic scaling
120 : * @return whether we are performing automatic scaling
121 : */
122 14955 : bool automaticScaling() const { return _automatic_scaling; }
123 :
124 : /**
125 : * Setter for whether we are performing automatic scaling
126 : * @param automatic_scaling A boolean representing whether we are performing automatic scaling
127 : */
128 66022 : void automaticScaling(bool automatic_scaling) { _automatic_scaling = automatic_scaling; }
129 :
130 : /**
131 : * Sets the verbose flag
132 : * @param[in] verbose Verbose flag
133 : */
134 61845 : void setVerboseFlag(const bool & verbose) { _verbose = verbose; }
135 :
136 : /**
137 : * Gets writeable reference to the dof map
138 : */
139 : virtual libMesh::DofMap & dofMap();
140 :
141 : /**
142 : * Gets const reference to the dof map
143 : */
144 : virtual const libMesh::DofMap & dofMap() const;
145 :
146 : /**
147 : * Get the reference to the libMesh system
148 : */
149 : virtual libMesh::System & system() = 0;
150 : virtual const libMesh::System & system() const = 0;
151 :
152 : /**
153 : * This is called prior to the libMesh system has been init'd. MOOSE system wrappers can use this
154 : * method to add vectors and matrices to the libMesh system
155 : */
156 127465 : virtual void preInit() {}
157 :
158 : /*
159 : * This is called after the libMesh system has been init'd. This can be used to initialize MOOSE
160 : * system data that relies on the libMesh system data already being initialized
161 : */
162 127465 : virtual void postInit() {}
163 :
164 : /**
165 : * Reinitialize the system when the degrees of freedom in this system have changed. This is called
166 : * after the libMesh system has been reinit'd
167 : */
168 8890 : virtual void reinit() {}
169 :
170 : /**
171 : * Called only once, just before the solve begins so objects can do some precalculations
172 : */
173 0 : virtual void initializeObjects() {}
174 :
175 : /**
176 : * Update the system (doing libMesh magic)
177 : */
178 : void update();
179 :
180 : /**
181 : * Solve the system (using libMesh magic)
182 : */
183 : virtual void solve();
184 :
185 : virtual void copyOldSolutions();
186 : virtual void copyPreviousNonlinearSolutions();
187 : virtual void copyPreviousFixedPointSolutions();
188 : virtual void restoreSolutions();
189 :
190 : /**
191 : * The solution vector that is currently being operated on.
192 : * This is typically a ghosted vector that comes in from the Nonlinear solver.
193 : */
194 : virtual const NumericVector<Number> * const & currentSolution() const = 0;
195 :
196 60764678 : NumericVector<Number> & solution() { return solutionState(0); }
197 183748 : NumericVector<Number> & solutionOld() { return solutionState(1); }
198 13032 : NumericVector<Number> & solutionOlder() { return solutionState(2); }
199 100430 : const NumericVector<Number> & solution() const { return solutionState(0); }
200 109207 : const NumericVector<Number> & solutionOld() const { return solutionState(1); }
201 2319 : const NumericVector<Number> & solutionOlder() const { return solutionState(2); }
202 :
203 : virtual const NumericVector<Number> * solutionPreviousNewton() const;
204 : virtual NumericVector<Number> * solutionPreviousNewton();
205 :
206 : /**
207 : * Initializes the solution state.
208 : */
209 : virtual void initSolutionState();
210 :
211 : /**
212 : * Get a state of the solution (0 = current, 1 = old, 2 = older, etc).
213 : *
214 : * If the state does not exist, it will be initialized in addition to any newer
215 : * states before it that have not been initialized.
216 : */
217 : virtual NumericVector<Number> &
218 : solutionState(const unsigned int state,
219 : Moose::SolutionIterationType iteration_type = Moose::SolutionIterationType::Time);
220 :
221 : /**
222 : * Get a state of the solution (0 = current, 1 = old, 2 = older, etc).
223 : */
224 : virtual const NumericVector<Number> & solutionState(
225 : const unsigned int state,
226 : Moose::SolutionIterationType iteration_type = Moose::SolutionIterationType::Time) const;
227 :
228 : /**
229 : * Returns the parallel type of the given solution state
230 : */
231 : libMesh::ParallelType
232 : solutionStateParallelType(const unsigned int state,
233 : const Moose::SolutionIterationType iteration_type) const;
234 :
235 : /**
236 : * Registers that the solution state \p state is needed.
237 : */
238 : virtual void needSolutionState(
239 : const unsigned int state,
240 : Moose::SolutionIterationType iteration_type = Moose::SolutionIterationType::Time,
241 : libMesh::ParallelType parallel_type = GHOSTED);
242 :
243 : /**
244 : * Whether or not the system has the solution state (0 = current, 1 = old, 2 = older, etc).
245 : */
246 : virtual bool hasSolutionState(
247 : const unsigned int state,
248 : Moose::SolutionIterationType iteration_type = Moose::SolutionIterationType::Time) const;
249 :
250 : /**
251 : * Add u_dot, u_dotdot, u_dot_old and u_dotdot_old
252 : * vectors if requested by the time integrator
253 : */
254 : virtual void addDotVectors();
255 :
256 62628 : virtual std::vector<Number> & duDotDus() { return _du_dot_du; }
257 55806 : virtual Number & duDotDotDu() { return _du_dotdot_du; }
258 : virtual const Number & duDotDu(unsigned int var_num = 0) const;
259 569490 : virtual const Number & duDotDotDu() const { return _du_dotdot_du; }
260 :
261 367289702 : virtual NumericVector<Number> * solutionUDot() { return _u_dot; }
262 86631 : virtual NumericVector<Number> * solutionUDotDot() { return _u_dotdot; }
263 1372272 : virtual NumericVector<Number> * solutionUDotOld() { return _u_dot_old; }
264 1372272 : virtual NumericVector<Number> * solutionUDotDotOld() { return _u_dotdot_old; }
265 569490 : virtual const NumericVector<Number> * solutionUDot() const { return _u_dot; }
266 569490 : virtual const NumericVector<Number> * solutionUDotDot() const { return _u_dotdot; }
267 569490 : virtual const NumericVector<Number> * solutionUDotOld() const { return _u_dot_old; }
268 569490 : virtual const NumericVector<Number> * solutionUDotDotOld() const { return _u_dotdot_old; }
269 :
270 : virtual void saveOldSolutions();
271 : virtual void restoreOldSolutions();
272 :
273 : /**
274 : * Check if the named vector exists in the system.
275 : */
276 : bool hasVector(const std::string & tag_name) const;
277 :
278 : /**
279 : * Check if the tagged vector exists in the system.
280 : */
281 6837286309 : virtual bool hasVector(TagID tag_id) const
282 : {
283 6837286309 : return tag_id < _tagged_vectors.size() && _tagged_vectors[tag_id];
284 : }
285 :
286 : /**
287 : * Ideally, we should not need this API.
288 : * There exists a really bad API "addCachedResidualDirectly " in FEProblem and DisplacedProblem
289 : * This API should go away once addCachedResidualDirectly is removed in the future
290 : * Return Tag ID for Time
291 : */
292 0 : virtual TagID timeVectorTag() const { mooseError("Not implemented yet"); }
293 :
294 : /**
295 : * Return the Matrix Tag ID for System
296 : */
297 0 : virtual TagID systemMatrixTag() const { mooseError("Not implemented yet"); }
298 :
299 : /*
300 : * Return TagID for nontime
301 : */
302 0 : virtual TagID nonTimeVectorTag() const { mooseError("Not implemented yet"); }
303 :
304 : /*
305 : * Return TagID for nontime
306 : */
307 0 : virtual TagID residualVectorTag() const { mooseError("Not implemented yet"); }
308 :
309 : /**
310 : * Get the default vector tags associated with this system
311 : */
312 31846 : virtual std::set<TagID> defaultVectorTags() const
313 : {
314 95538 : return {timeVectorTag(), nonTimeVectorTag(), residualVectorTag()};
315 : }
316 : /**
317 : * Get the default matrix tags associted with this system
318 : */
319 12300 : virtual std::set<TagID> defaultMatrixTags() const { return {systemMatrixTag()}; }
320 :
321 : /**
322 : * Get a raw NumericVector by name
323 : */
324 : ///@{
325 : virtual NumericVector<Number> & getVector(const std::string & name);
326 : virtual const NumericVector<Number> & getVector(const std::string & name) const;
327 : ///@}
328 :
329 : /**
330 : * Get a raw NumericVector by tag
331 : */
332 : ///@{
333 : virtual NumericVector<Number> & getVector(TagID tag);
334 : virtual const NumericVector<Number> & getVector(TagID tag) const;
335 : ///@}
336 :
337 : /**
338 : * Associate a vector for a given tag
339 : */
340 : virtual void associateVectorToTag(NumericVector<Number> & vec, TagID tag);
341 :
342 : /**
343 : * Disassociate a given vector from a given tag
344 : */
345 : virtual void disassociateVectorFromTag(NumericVector<Number> & vec, TagID tag);
346 :
347 : /**
348 : * Disassociate any vector that is associated with a given tag
349 : */
350 : virtual void disassociateVectorFromTag(TagID tag);
351 :
352 : /**
353 : * Disassociate the vectors associated with the default vector tags of this system
354 : */
355 : virtual void disassociateDefaultVectorTags();
356 :
357 : /**
358 : * Check if the tagged matrix exists in the system.
359 : */
360 2504440987 : virtual bool hasMatrix(TagID tag) const
361 : {
362 2504440987 : return tag < _tagged_matrices.size() && _tagged_matrices[tag];
363 : }
364 :
365 : /**
366 : * Get a raw SparseMatrix
367 : */
368 : virtual libMesh::SparseMatrix<Number> & getMatrix(TagID tag);
369 :
370 : /**
371 : * Get a raw SparseMatrix
372 : */
373 : virtual const libMesh::SparseMatrix<Number> & getMatrix(TagID tag) const;
374 :
375 : /**
376 : * Make all existing matrices active
377 : */
378 : virtual void activateAllMatrixTags();
379 :
380 : /**
381 : * If or not a matrix tag is active
382 : */
383 : virtual bool matrixTagActive(TagID tag) const;
384 :
385 : /**
386 : * Make matrices inactive
387 : */
388 : virtual void deactivateAllMatrixTags();
389 :
390 : /**
391 : * Close all matrices associated the tags
392 : */
393 : void closeTaggedMatrices(const std::set<TagID> & tags);
394 :
395 : /**
396 : * flushes all matrices associated to tags. Flush assembles the matrix but doesn't shrink memory
397 : * allocation
398 : */
399 : void flushTaggedMatrices(const std::set<TagID> & tags);
400 :
401 : /**
402 : * Associate a matrix to a tag
403 : */
404 : virtual void associateMatrixToTag(libMesh::SparseMatrix<Number> & matrix, TagID tag);
405 :
406 : /**
407 : * Disassociate a matrix from a tag
408 : */
409 : virtual void disassociateMatrixFromTag(libMesh::SparseMatrix<Number> & matrix, TagID tag);
410 :
411 : /**
412 : * Disassociate any matrix that is associated with a given tag
413 : */
414 : virtual void disassociateMatrixFromTag(TagID tag);
415 :
416 : /**
417 : * Disassociate the matrices associated with the default matrix tags of this system
418 : */
419 : virtual void disassociateDefaultMatrixTags();
420 :
421 : /**
422 : * Returns a reference to a serialized version of the solution vector for this subproblem
423 : */
424 : virtual NumericVector<Number> & serializedSolution();
425 :
426 0 : virtual NumericVector<Number> & residualCopy()
427 : {
428 0 : mooseError("This system does not support getting a copy of the residual");
429 : }
430 0 : virtual NumericVector<Number> & residualGhosted()
431 : {
432 0 : mooseError("This system does not support getting a ghosted copy of the residual");
433 : }
434 :
435 : /**
436 : * Will modify the send_list to add all of the extra ghosted dofs for this system
437 : */
438 : virtual void augmentSendList(std::vector<dof_id_type> & send_list);
439 :
440 : /**
441 : * Will modify the sparsity pattern to add logical geometric connections
442 : */
443 : virtual void augmentSparsity(libMesh::SparsityPattern::Graph & sparsity,
444 : std::vector<dof_id_type> & n_nz,
445 : std::vector<dof_id_type> & n_oz) = 0;
446 :
447 : /**
448 : * Canonical method for adding a variable
449 : * @param var_type the type of the variable, e.g. MooseVariableScalar
450 : * @param var_name the variable name, e.g. 'u'
451 : * @param params the InputParameters from which to construct the variable
452 : */
453 : virtual void addVariable(const std::string & var_type,
454 : const std::string & var_name,
455 : InputParameters & parameters);
456 :
457 : /**
458 : * If a variable is an array variable
459 : */
460 : virtual bool isArrayVariable(const std::string & var_name) const;
461 :
462 : ///@{
463 : /**
464 : * Query a system for a variable
465 : *
466 : * @param var_name name of the variable
467 : * @return true if the variable exists
468 : */
469 : virtual bool hasVariable(const std::string & var_name) const;
470 : virtual bool hasScalarVariable(const std::string & var_name) const;
471 : ///@}
472 :
473 : virtual bool isScalarVariable(unsigned int var_name) const;
474 :
475 : /**
476 : * Gets a reference to a variable of with specified name
477 : *
478 : * @param tid Thread id
479 : * @param var_name variable name
480 : * @return reference the variable (class)
481 : */
482 : MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string & var_name) const;
483 :
484 : /**
485 : * Gets a reference to a variable with specified number
486 : *
487 : * @param tid Thread id
488 : * @param var_number libMesh variable number
489 : * @return reference the variable (class)
490 : */
491 : MooseVariableFieldBase & getVariable(THREAD_ID tid, unsigned int var_number) const;
492 :
493 : /**
494 : * Gets a reference to a variable of with specified name
495 : *
496 : * This excludes and cannot return finite volume variables.
497 : *
498 : * @param tid Thread id
499 : * @param var_name variable name
500 : * @return reference the variable (class)
501 : */
502 : template <typename T>
503 : MooseVariableFE<T> & getFieldVariable(THREAD_ID tid, const std::string & var_name);
504 :
505 : /**
506 : * Returns a field variable pointer - this includes finite volume variables.
507 : */
508 : template <typename T>
509 : MooseVariableField<T> & getActualFieldVariable(THREAD_ID tid, const std::string & var_name);
510 :
511 : /**
512 : * Gets a reference to a variable with specified number
513 : *
514 : * This excludes and cannot return finite volume variables.
515 : *
516 : * @param tid Thread id
517 : * @param var_number libMesh variable number
518 : * @return reference the variable (class)
519 : */
520 : template <typename T>
521 : MooseVariableFE<T> & getFieldVariable(THREAD_ID tid, unsigned int var_number);
522 :
523 : /**
524 : * Returns a field variable pointer - this includes finite volume variables.
525 : */
526 : template <typename T>
527 : MooseVariableField<T> & getActualFieldVariable(THREAD_ID tid, unsigned int var_number);
528 :
529 : /**
530 : * Return a finite volume variable
531 : */
532 : template <typename T>
533 : MooseVariableFV<T> & getFVVariable(THREAD_ID tid, const std::string & var_name);
534 :
535 : /**
536 : * Gets a reference to a scalar variable with specified number
537 : *
538 : * @param tid Thread id
539 : * @param var_name A string which is the name of the variable to get.
540 : * @return reference the variable (class)
541 : */
542 : virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid,
543 : const std::string & var_name) const;
544 :
545 : /**
546 : * Gets a reference to a variable with specified number
547 : *
548 : * @param tid Thread id
549 : * @param var_number libMesh variable number
550 : * @return reference the variable (class)
551 : */
552 : virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, unsigned int var_number) const;
553 :
554 : /**
555 : * Get the block where a variable of this system is defined
556 : *
557 : * @param var_number The number of the variable
558 : * @return the set of subdomain ids where the variable is active (defined)
559 : */
560 : virtual const std::set<SubdomainID> * getVariableBlocks(unsigned int var_number);
561 :
562 : /**
563 : * Get the number of variables in this system
564 : * @return the number of variables
565 : */
566 : virtual unsigned int nVariables() const;
567 :
568 : /**
569 : * Get the number of field variables in this system
570 : * @return the number of field variables
571 : */
572 : unsigned int nFieldVariables() const;
573 :
574 : /**
575 : * Get the number of finite volume variables in this system
576 : * @return the number of finite volume variables
577 : */
578 : unsigned int nFVVariables() const;
579 :
580 : /**
581 : * Gets the maximum number of dofs used by any one variable on any one element
582 : *
583 : * @return The max
584 : */
585 84 : std::size_t getMaxVarNDofsPerElem() const { return _max_var_n_dofs_per_elem; }
586 :
587 : /**
588 : * Gets the maximum number of dofs used by any one variable on any one node
589 : *
590 : * @return The max
591 : */
592 : std::size_t getMaxVarNDofsPerNode() const { return _max_var_n_dofs_per_node; }
593 :
594 : /**
595 : * assign the maximum element dofs
596 : */
597 63271 : void assignMaxVarNDofsPerElem(std::size_t max_dofs) { _max_var_n_dofs_per_elem = max_dofs; }
598 :
599 : /**
600 : * assign the maximum node dofs
601 : */
602 63271 : void assignMaxVarNDofsPerNode(std::size_t max_dofs) { _max_var_n_dofs_per_node = max_dofs; }
603 :
604 : /**
605 : * Adds this variable to the list of variables to be zeroed during each residual evaluation.
606 : * @param var_name The name of the variable to be zeroed.
607 : */
608 : virtual void addVariableToZeroOnResidual(std::string var_name);
609 :
610 : /**
611 : * Adds this variable to the list of variables to be zeroed during each Jacobian evaluation.
612 : * @param var_name The name of the variable to be zeroed.
613 : */
614 : virtual void addVariableToZeroOnJacobian(std::string var_name);
615 :
616 : /**
617 : * Zero out the solution for the list of variables passed in.
618 : *
619 : * @ param vars_to_be_zeroed The variable names in this vector will have their solutions set to
620 : * zero after this call
621 : */
622 : virtual void zeroVariables(std::vector<std::string> & vars_to_be_zeroed);
623 :
624 : /**
625 : * Zero out the solution for the variables that were registered as needing to have their solutions
626 : * zeroed on out on residual evaluation by a call to addVariableToZeroOnResidual()
627 : */
628 : virtual void zeroVariablesForResidual();
629 :
630 : /**
631 : * Zero out the solution for the variables that were registered as needing to have their solutions
632 : * zeroed on out on Jacobian evaluation by a call to addVariableToZeroOnResidual()
633 : */
634 : virtual void zeroVariablesForJacobian();
635 :
636 : /**
637 : * Get minimal quadrature order needed for integrating variables in this system
638 : * @return The minimal order of quadrature
639 : */
640 : virtual libMesh::Order getMinQuadratureOrder();
641 :
642 : /**
643 : * Prepare the system for use
644 : * @param tid ID of the thread
645 : */
646 : virtual void prepare(THREAD_ID tid);
647 :
648 : /**
649 : * Prepare the system for use on sides
650 : *
651 : * This will try to reuse the preparation done on the element.
652 : *
653 : * @param tid ID of the thread
654 : * @param resize_data Pass True if this system needs to resize residual and jacobian
655 : * datastructures based on preparing this face
656 : */
657 : virtual void prepareFace(THREAD_ID tid, bool resize_data);
658 :
659 : /**
660 : * Prepare the system for use
661 : * @param tid ID of the thread
662 : */
663 : virtual void prepareNeighbor(THREAD_ID tid);
664 :
665 : /**
666 : * Prepare the system for use for lower dimensional elements
667 : * @param tid ID of the thread
668 : */
669 : virtual void prepareLowerD(THREAD_ID tid);
670 :
671 : /**
672 : * Reinit an element assembly info
673 : * @param elem Which element we are reinitializing for
674 : * @param tid ID of the thread
675 : */
676 : virtual void reinitElem(const Elem * elem, THREAD_ID tid);
677 :
678 : /**
679 : * Reinit assembly info for a side of an element
680 : * @param elem The element
681 : * @param side Side of of the element
682 : * @param tid Thread ID
683 : */
684 : virtual void reinitElemFace(const Elem * elem, unsigned int side, THREAD_ID tid);
685 :
686 : /**
687 : * Compute the values of the variables at all the current points.
688 : */
689 : virtual void reinitNeighborFace(const Elem * elem, unsigned int side, THREAD_ID tid);
690 :
691 : /**
692 : * Compute the values of the variables at all the current points.
693 : */
694 : virtual void reinitNeighbor(const Elem * elem, THREAD_ID tid);
695 :
696 : /**
697 : * Compute the values of the variables on the lower dimensional element
698 : */
699 : virtual void reinitLowerD(THREAD_ID tid);
700 :
701 : /**
702 : * Reinit nodal assembly info
703 : * @param node Node to reinit for
704 : * @param tid Thread ID
705 : */
706 : virtual void reinitNode(const Node * node, THREAD_ID tid);
707 :
708 : /**
709 : * Reinit nodal assembly info on a face
710 : * @param node Node to reinit
711 : * @param bnd_id Boundary ID
712 : * @param tid Thread ID
713 : */
714 : virtual void reinitNodeFace(const Node * node, BoundaryID bnd_id, THREAD_ID tid);
715 :
716 : /**
717 : * Reinit variables at a set of nodes
718 : * @param nodes List of node ids to reinit
719 : * @param tid Thread ID
720 : */
721 : virtual void reinitNodes(const std::vector<dof_id_type> & nodes, THREAD_ID tid);
722 :
723 : /**
724 : * Reinit variables at a set of neighbor nodes
725 : * @param nodes List of node ids to reinit
726 : * @param tid Thread ID
727 : */
728 : virtual void reinitNodesNeighbor(const std::vector<dof_id_type> & nodes, THREAD_ID tid);
729 :
730 : /**
731 : * Reinit scalar varaibles
732 : * @param tid Thread ID
733 : * @param reinit_for_derivative_reordering A flag indicating whether we are reinitializing for the
734 : * purpose of re-ordering derivative information for ADNodalBCs
735 : */
736 : virtual void reinitScalars(THREAD_ID tid, bool reinit_for_derivative_reordering = false);
737 :
738 : /**
739 : * Add info about variable that will be copied
740 : *
741 : * @param dest_name Name of the nodal variable being used for copying into (name is from the
742 : * exodusII file)
743 : * @param source_name Name of the nodal variable being used for copying from (name is from the
744 : * exodusII file)
745 : * @param timestep Timestep in the file being used
746 : */
747 : virtual void addVariableToCopy(const std::string & dest_name,
748 : const std::string & source_name,
749 : const std::string & timestep);
750 :
751 930831129 : const std::vector<MooseVariableFieldBase *> & getVariables(THREAD_ID tid)
752 : {
753 930831129 : return _vars[tid].fieldVariables();
754 : }
755 :
756 83164412 : const std::vector<MooseVariableScalar *> & getScalarVariables(THREAD_ID tid)
757 : {
758 83164412 : return _vars[tid].scalars();
759 : }
760 :
761 363 : const std::set<SubdomainID> & getSubdomainsForVar(unsigned int var_number) const
762 : {
763 363 : return _var_map.at(var_number);
764 : }
765 :
766 : /**
767 : * Get the block where a variable of this system is defined
768 : *
769 : * @param var_name The name of the variable
770 : * @return the set of subdomain ids where the variable is active (defined)
771 : */
772 : const std::set<SubdomainID> & getSubdomainsForVar(const std::string & var_name) const;
773 :
774 : /**
775 : * Remove a vector from the system with the given name.
776 : */
777 : void removeVector(const std::string & name);
778 :
779 : /**
780 : * Adds a solution length vector to the system.
781 : *
782 : * @param vector_name The name of the vector.
783 : * @param project Whether or not to project this vector when doing mesh refinement.
784 : * If the vector is just going to be recomputed then there is no need to project
785 : * it.
786 : * @param type What type of parallel vector. This is usually either PARALLEL or GHOSTED.
787 : * GHOSTED is needed if you are going to be accessing
788 : * off-processor entries.
789 : * The ghosting pattern is the same as the solution
790 : * vector.
791 : */
792 : NumericVector<Number> &
793 : addVector(const std::string & vector_name, const bool project, const libMesh::ParallelType type);
794 :
795 : /**
796 : * Adds a solution length vector to the system with the specified TagID
797 : *
798 : * @param tag_name The name of the tag
799 : * @param project Whether or not to project this vector when doing mesh refinement.
800 : * If the vector is just going to be recomputed then there is no need to project
801 : * it.
802 : * @param type What type of parallel vector. This is usually either PARALLEL or GHOSTED.
803 : * GHOSTED is needed if you are going to be accessing
804 : * off-processor entries.
805 : * The ghosting pattern is the same as the solution
806 : * vector.
807 : */
808 : NumericVector<Number> &
809 : addVector(TagID tag, const bool project, const libMesh::ParallelType type);
810 :
811 : /**
812 : * Close vector with the given tag
813 : */
814 : void closeTaggedVector(const TagID tag);
815 : /**
816 : * Close all vectors for given tags
817 : */
818 : void closeTaggedVectors(const std::set<TagID> & tags);
819 :
820 : /**
821 : * Zero vector with the given tag
822 : */
823 : void zeroTaggedVector(const TagID tag);
824 : /**
825 : * Zero all vectors for given tags
826 : */
827 : void zeroTaggedVectors(const std::set<TagID> & tags);
828 :
829 : /**
830 : * Remove a solution length vector from the system with the specified TagID
831 : *
832 : * @param tag_id Tag ID
833 : */
834 : void removeVector(TagID tag_id);
835 :
836 : /// set all the global dof indices for a variable
837 : /// @param var_name The name of the variable
838 : void setVariableGlobalDoFs(const std::string & var_name);
839 :
840 : /// Get the global dof indices of a variable, this needs to be called
841 : /// after the indices have been set by `setVariableGlobalDoFs`
842 116 : const std::vector<dof_id_type> & getVariableGlobalDoFs() { return _var_all_dof_indices; }
843 :
844 : /**
845 : * Adds a matrix with a given tag
846 : *
847 : * @param tag_name The name of the tag
848 : */
849 : libMesh::SparseMatrix<Number> & addMatrix(TagID tag);
850 :
851 : /**
852 : * Removes a matrix with a given tag
853 : *
854 : * @param tag_name The name of the tag
855 : */
856 : void removeMatrix(TagID tag);
857 :
858 : virtual const std::string & name() const;
859 :
860 7283438 : const std::vector<VariableName> & getVariableNames() const { return _vars[0].names(); }
861 :
862 : void getStandardFieldVariableNames(std::vector<VariableName> & std_field_variables) const;
863 :
864 : /**
865 : * Returns the maximum number of all variables on the system
866 : */
867 : unsigned int getMaxVariableNumber() const { return _max_var_number; }
868 :
869 0 : virtual void computeVariables(const NumericVector<Number> & /*soln*/) {}
870 :
871 : void copyVars(libMesh::ExodusII_IO & io);
872 :
873 : /**
874 : * Copy current solution into old and older
875 : */
876 : virtual void copySolutionsBackwards();
877 :
878 : void addTimeIntegrator(const std::string & type,
879 : const std::string & name,
880 : InputParameters & parameters);
881 :
882 : /// Whether or not there are variables to be restarted from an Exodus mesh file
883 112128 : bool hasVarCopy() const { return _var_to_copy.size() > 0; }
884 :
885 : /**
886 : * Add the scaling factor vector to the system
887 : */
888 : void addScalingVector();
889 :
890 : /**
891 : * Whether or not the solution states have been initialized via initSolutionState()
892 : *
893 : * After the solution states have been initialized, additional solution
894 : * states cannot be added.
895 : */
896 154 : bool solutionStatesInitialized() const { return _solution_states_initialized; }
897 :
898 : /// Setup Functions
899 : virtual void initialSetup();
900 : virtual void timestepSetup();
901 : virtual void customSetup(const ExecFlagType & exec_type);
902 : virtual void subdomainSetup();
903 : virtual void residualSetup();
904 : virtual void jacobianSetup();
905 :
906 : /**
907 : * Clear all dof indices from moose variables
908 : */
909 : void clearAllDofIndices();
910 :
911 : /**
912 : * Set the active vector tags for the variables
913 : */
914 : void setActiveVariableCoupleableVectorTags(const std::set<TagID> & vtags, THREAD_ID tid);
915 :
916 : /**
917 : * Set the active vector tags for the scalar variables
918 : */
919 : void setActiveScalarVariableCoupleableVectorTags(const std::set<TagID> & vtags, THREAD_ID tid);
920 :
921 : /**
922 : * @return the type of variables this system holds, e.g. nonlinear or auxiliary
923 : */
924 5177 : Moose::VarKindType varKind() const { return _var_kind; }
925 :
926 : /**
927 : * Reference to the container vector which hold gradients at dofs (if it can be interpreted).
928 : * Mainly used for finite volume systems.
929 : */
930 1559 : const std::vector<std::unique_ptr<NumericVector<Number>>> & gradientContainer() const
931 : {
932 1559 : return _raw_grad_container;
933 : }
934 :
935 : /**
936 : * Compute time derivatives, auxiliary variables, etc.
937 : * @param type Our current execution stage
938 : */
939 : virtual void compute(ExecFlagType type) = 0;
940 :
941 : /**
942 : * Copy time integrators from another system
943 : */
944 : void copyTimeIntegrators(const SystemBase & other_sys);
945 :
946 : /**
947 : * Retrieve the time integrator that integrates the given variable's equation
948 : */
949 : const TimeIntegrator & getTimeIntegrator(const unsigned int var_num) const;
950 :
951 : /**
952 : * Retrieve the time integrator that integrates the given variable's equation. If no suitable time
953 : * integrator is found (this could happen for instance if we're solving a non-transient problem),
954 : * then a nullptr will be returned
955 : */
956 : const TimeIntegrator * queryTimeIntegrator(const unsigned int var_num) const;
957 :
958 : /**
959 : * @returns All the time integrators owned by this system
960 : */
961 : const std::vector<std::shared_ptr<TimeIntegrator>> & getTimeIntegrators();
962 :
963 : /**
964 : * @returns The prefix used for this system for solver settings for PETSc. This prefix is used to
965 : * prevent collision of solver settings for different systems. Note that this prefix does not have
966 : * a leading dash so it's appropriate for passage straight to PETSc APIs
967 : */
968 : std::string prefix() const;
969 :
970 : /**
971 : * size the matrix data for each variable for the number of matrix tags we have
972 : */
973 : void sizeVariableMatrixData();
974 :
975 : protected:
976 : /**
977 : * Internal getter for solution owned by libMesh.
978 : */
979 : virtual NumericVector<Number> & solutionInternal() const = 0;
980 :
981 : /// The subproblem for whom this class holds variable data, etc; this can either be the governing
982 : /// finite element/volume problem or a subjugate displaced problem
983 : SubProblem & _subproblem;
984 :
985 : /// the governing finite element/volume problem
986 : FEProblemBase & _fe_problem;
987 :
988 : MooseApp & _app;
989 : Factory & _factory;
990 :
991 : MooseMesh & _mesh;
992 : /// The name of this system
993 : std::string _name;
994 :
995 : /// Variable warehouses (one for each thread)
996 : std::vector<VariableWarehouse> _vars;
997 : /// Map of variables (variable id -> array of subdomains where it lives)
998 : std::map<unsigned int, std::set<SubdomainID>> _var_map;
999 : /// Maximum variable number
1000 : unsigned int _max_var_number;
1001 :
1002 : std::vector<std::string> _vars_to_be_zeroed_on_residual;
1003 : std::vector<std::string> _vars_to_be_zeroed_on_jacobian;
1004 :
1005 : /// solution vector for u^dot
1006 : NumericVector<Number> * _u_dot;
1007 : /// solution vector for u^dotdot
1008 : NumericVector<Number> * _u_dotdot;
1009 :
1010 : /// old solution vector for u^dot
1011 : NumericVector<Number> * _u_dot_old;
1012 : /// old solution vector for u^dotdot
1013 : NumericVector<Number> * _u_dotdot_old;
1014 :
1015 : /// Derivative of time derivative of u with respect to uj. This depends on the time integration
1016 : /// scheme
1017 : std::vector<Real> _du_dot_du;
1018 : Real _du_dotdot_du;
1019 :
1020 : /// Tagged vectors (pointer)
1021 : std::vector<NumericVector<Number> *> _tagged_vectors;
1022 : /// Tagged matrices (pointer)
1023 : std::vector<libMesh::SparseMatrix<Number> *> _tagged_matrices;
1024 : /// Active tagged matrices. A matrix is active if its tag-matrix pair is present in the map. We use a map instead of a vector so that users can easily add and remove to this container with calls to (de)activateMatrixTag
1025 : std::unordered_map<TagID, libMesh::SparseMatrix<Number> *> _active_tagged_matrices;
1026 : /// Active flags for tagged matrices
1027 : std::vector<bool> _matrix_tag_active_flags;
1028 :
1029 : // Used for saving old solutions so that they wont be accidentally changed
1030 : NumericVector<Real> * _saved_old;
1031 : NumericVector<Real> * _saved_older;
1032 :
1033 : // Used for saving old u_dot and u_dotdot so that they wont be accidentally changed
1034 : NumericVector<Real> * _saved_dot_old;
1035 : NumericVector<Real> * _saved_dotdot_old;
1036 :
1037 : /// default kind of variables in this system
1038 : Moose::VarKindType _var_kind;
1039 :
1040 : std::vector<VarCopyInfo> _var_to_copy;
1041 :
1042 : /// Maximum number of dofs for any one variable on any one element
1043 : size_t _max_var_n_dofs_per_elem;
1044 :
1045 : /// Maximum number of dofs for any one variable on any one node
1046 : size_t _max_var_n_dofs_per_node;
1047 :
1048 : /// Time integrator
1049 : std::vector<std::shared_ptr<TimeIntegrator>> _time_integrators;
1050 :
1051 : /// Map variable number to its pointer
1052 : std::vector<std::vector<MooseVariableFieldBase *>> _numbered_vars;
1053 :
1054 : /// Whether to automatically scale the variables
1055 : bool _automatic_scaling;
1056 :
1057 : /// True if printing out additional information
1058 : bool _verbose;
1059 :
1060 : /// Whether or not the solution states have been initialized
1061 : bool _solution_states_initialized;
1062 :
1063 : /// Container for the dof indices of a given variable
1064 : std::vector<dof_id_type> _var_all_dof_indices;
1065 :
1066 : /// Serialized version of the solution vector, or nullptr if a
1067 : /// serialized solution is not needed
1068 : std::unique_ptr<NumericVector<Number>> _serialized_solution;
1069 :
1070 : /// A cache for storing gradients at dof locations. We store it on the system
1071 : /// because we create copies of variables on each thread and that would
1072 : /// lead to increased data duplication when using threading-based parallelism.
1073 : std::vector<std::unique_ptr<NumericVector<Number>>> _raw_grad_container;
1074 :
1075 : private:
1076 : /**
1077 : * Gets the vector name used for an old (not current) solution state.
1078 : */
1079 : TagName oldSolutionStateVectorName(const unsigned int,
1080 : Moose::SolutionIterationType iteration_type) const;
1081 :
1082 : /// 2D array of solution state vector pointers; first index corresponds to
1083 : /// SolutionIterationType, second index corresponds to state index (0=current, 1=old, 2=older)
1084 : std::array<std::vector<NumericVector<Number> *>, 3> _solution_states;
1085 : /// The saved solution states (0 = current, 1 = old, 2 = older, etc)
1086 : std::vector<NumericVector<Number> *> _saved_solution_states;
1087 : };
1088 :
1089 : inline bool
1090 229279710 : SystemBase::hasSolutionState(const unsigned int state,
1091 : const Moose::SolutionIterationType iteration_type) const
1092 : {
1093 229279710 : return _solution_states[static_cast<unsigned short>(iteration_type)].size() > state;
1094 : }
1095 :
1096 : #define PARALLEL_TRY
1097 :
1098 : #define PARALLEL_CATCH _fe_problem.checkExceptionAndStopSolve();
|
