www.mooseframework.org
SystemBase.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "MooseApp.h"
16 #include "SystemBase.h"
17 #include "Factory.h"
18 #include "SubProblem.h"
19 #include "MooseVariable.h"
20 #include "MooseVariableScalar.h"
22 #include "Conversion.h"
23 #include "Parser.h"
25 #include "MooseTypes.h"
26 #include "InitialCondition.h"
27 #include "ScalarInitialCondition.h"
28 #include "Assembly.h"
29 #include "MooseMesh.h"
30 
32 void
33 extraSendList(std::vector<dof_id_type> & send_list, void * context)
34 {
35  SystemBase * sys = static_cast<SystemBase *>(context);
36  sys->augmentSendList(send_list);
37 }
38 
40 void
41 extraSparsity(SparsityPattern::Graph & sparsity,
42  std::vector<dof_id_type> & n_nz,
43  std::vector<dof_id_type> & n_oz,
44  void * context)
45 {
46  SystemBase * sys = static_cast<SystemBase *>(context);
47  sys->augmentSparsity(sparsity, n_nz, n_oz);
48 }
49 
50 template <>
51 void
52 dataStore(std::ostream & stream, SystemBase & system_base, void * context)
53 {
54  System & libmesh_system = system_base.system();
55 
56  NumericVector<Real> & solution = *(libmesh_system.solution.get());
57 
58  dataStore(stream, solution, context);
59 
60  for (System::vectors_iterator it = libmesh_system.vectors_begin();
61  it != libmesh_system.vectors_end();
62  it++)
63  dataStore(stream, *(it->second), context);
64 }
65 
66 template <>
67 void
68 dataLoad(std::istream & stream, SystemBase & system_base, void * context)
69 {
70  System & libmesh_system = system_base.system();
71 
72  NumericVector<Real> & solution = *(libmesh_system.solution.get());
73 
74  dataLoad(stream, solution, context);
75 
76  for (System::vectors_iterator it = libmesh_system.vectors_begin();
77  it != libmesh_system.vectors_end();
78  it++)
79  dataLoad(stream, *(it->second), context);
80 
81  system_base.update();
82 }
83 
85  const std::string & name,
86  Moose::VarKindType var_kind)
87  : libMesh::ParallelObject(subproblem),
88  _subproblem(subproblem),
89  _app(subproblem.getMooseApp()),
90  _factory(_app.getFactory()),
91  _mesh(subproblem.mesh()),
92  _name(name),
93  _vars(libMesh::n_threads()),
94  _var_map(),
95  _dummy_vec(NULL),
96  _saved_old(NULL),
97  _saved_older(NULL),
98  _var_kind(var_kind)
99 {
100 }
101 
103 SystemBase::getVariable(THREAD_ID tid, const std::string & var_name)
104 {
105  MooseVariable * var = dynamic_cast<MooseVariable *>(_vars[tid].getVariable(var_name));
106  if (var == NULL)
107  mooseError("Variable '" + var_name + "' does not exist in this system");
108  return *var;
109 }
110 
112 SystemBase::getVariable(THREAD_ID tid, unsigned int var_number)
113 {
114  MooseVariable * var = dynamic_cast<MooseVariable *>(_vars[tid].getVariable(var_number));
115  if (var == NULL)
116  mooseError("variable #" + Moose::stringify(var_number) + " does not exist in this system");
117  return *var;
118 }
119 
121 SystemBase::getScalarVariable(THREAD_ID tid, const std::string & var_name)
122 {
123  MooseVariableScalar * var = dynamic_cast<MooseVariableScalar *>(_vars[tid].getVariable(var_name));
124  if (var == NULL)
125  mooseError("Scalar variable '" + var_name + "' does not exist in this system");
126  return *var;
127 }
128 
130 SystemBase::getScalarVariable(THREAD_ID tid, unsigned int var_number)
131 {
132  MooseVariableScalar * var =
133  dynamic_cast<MooseVariableScalar *>(_vars[tid].getVariable(var_number));
134  if (var == NULL)
135  mooseError("variable #" + Moose::stringify(var_number) + " does not exist in this system");
136  return *var;
137 }
138 
139 const std::set<SubdomainID> *
140 SystemBase::getVariableBlocks(unsigned int var_number)
141 {
142  mooseAssert(_var_map.find(var_number) != _var_map.end(), "Variable does not exist.");
143  if (_var_map[var_number].empty())
144  return NULL;
145  else
146  return &_var_map[var_number];
147 }
148 
149 void
151 {
152  _vars_to_be_zeroed_on_residual.push_back(var_name);
153 }
154 
155 void
157 {
158  _vars_to_be_zeroed_on_jacobian.push_back(var_name);
159 }
160 
161 void
162 SystemBase::zeroVariables(std::vector<std::string> & vars_to_be_zeroed)
163 {
164  if (vars_to_be_zeroed.size() > 0)
165  {
166  NumericVector<Number> & solution = this->solution();
167 
168  AllLocalDofIndicesThread aldit(system(), vars_to_be_zeroed);
169  ConstElemRange & elem_range = *_mesh.getActiveLocalElementRange();
170  Threads::parallel_reduce(elem_range, aldit);
171 
172  std::set<dof_id_type> dof_indices_to_zero = aldit._all_dof_indices;
173 
174  solution.close();
175 
176  for (const auto & dof : dof_indices_to_zero)
177  solution.set(dof, 0);
178 
179  solution.close();
180 
181  // Call update to update the current_local_solution for this system
182  system().update();
183  }
184 }
185 
186 void
188 {
190 }
191 
192 void
194 {
196 }
197 
198 Order
200 {
201  Order order = CONSTANT;
202  std::vector<MooseVariable *> vars = _vars[0].variables();
203  for (const auto & var : vars)
204  {
205  FEType fe_type = var->feType();
206  if (fe_type.default_quadrature_order() > order)
207  order = fe_type.default_quadrature_order();
208  }
209 
210  return order;
211 }
212 
213 void
215 {
217  {
218  const std::set<MooseVariable *> & active_elemental_moose_variables =
220  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
221  for (const auto & var : vars)
222  var->clearDofIndices();
223 
224  for (const auto & var : active_elemental_moose_variables)
225  if (&(var->sys()) == this)
226  var->prepare();
227  }
228  else
229  {
230  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
231  for (const auto & var : vars)
232  var->prepare();
233  }
234 }
235 
236 void
237 SystemBase::prepareFace(THREAD_ID tid, bool resize_data)
238 {
240  tid)) // We only need to do something if the element prepare was restricted
241  {
242  const std::set<MooseVariable *> & active_elemental_moose_variables =
244 
245  std::vector<MooseVariable *> newly_prepared_vars;
246 
247  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
248  for (const auto & var : vars)
249  {
250  if (&(var->sys()) == this &&
251  !active_elemental_moose_variables.count(
252  var)) // If it wasnt in the active list we need to prepare it
253  {
254  var->prepare();
255  newly_prepared_vars.push_back(var);
256  }
257  }
258 
259  // Make sure to resize the residual and jacobian datastructures for all the new variables
260  if (resize_data)
261  for (unsigned int i = 0; i < newly_prepared_vars.size(); i++)
262  {
263  _subproblem.assembly(tid).prepareVariable(newly_prepared_vars[i]);
265  _subproblem.assembly(tid).prepareVariableNonlocal(newly_prepared_vars[i]);
266  }
267  }
268 }
269 
270 void
272 {
273  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
274  for (const auto & var : vars)
275  var->prepareNeighbor();
276 }
277 
278 void
279 SystemBase::reinitElem(const Elem * /*elem*/, THREAD_ID tid)
280 {
281 
283  {
284  const std::set<MooseVariable *> & active_elemental_moose_variables =
286  for (const auto & var : active_elemental_moose_variables)
287  if (&(var->sys()) == this)
288  var->computeElemValues();
289  }
290  else
291  {
292  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
293  for (const auto & var : vars)
294  var->computeElemValues();
295  }
296 }
297 
298 void
299 SystemBase::reinitElemFace(const Elem * /*elem*/,
300  unsigned int /*side*/,
301  BoundaryID /*bnd_id*/,
302  THREAD_ID tid)
303 {
304  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
305  for (const auto & var : vars)
306  var->computeElemValuesFace();
307 }
308 
309 void
310 SystemBase::reinitNeighborFace(const Elem * /*elem*/,
311  unsigned int /*side*/,
312  BoundaryID /*bnd_id*/,
313  THREAD_ID tid)
314 {
315  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
316  for (const auto & var : vars)
317  var->computeNeighborValuesFace();
318 }
319 
320 void
321 SystemBase::reinitNeighbor(const Elem * /*elem*/, THREAD_ID tid)
322 {
323  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
324  for (const auto & var : vars)
325  var->computeNeighborValues();
326 }
327 
328 void
329 SystemBase::reinitNode(const Node * /*node*/, THREAD_ID tid)
330 {
331  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
332  for (const auto & var : vars)
333  {
334  if (var->isNodal())
335  {
336  var->reinitNode();
337  var->computeNodalValues();
338  }
339  }
340 }
341 
342 void
343 SystemBase::reinitNodeFace(const Node * /*node*/, BoundaryID /*bnd_id*/, THREAD_ID tid)
344 {
345  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
346  for (const auto & var : vars)
347  {
348  if (var->isNodal())
349  {
350  var->reinitNode();
351  var->computeNodalValues();
352  }
353  }
354 }
355 
356 void
357 SystemBase::reinitNodeNeighbor(const Node * /*node*/, THREAD_ID tid)
358 {
359  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
360  for (const auto & var : vars)
361  {
362  if (var->isNodal())
363  {
364  var->reinitNodeNeighbor();
365  var->computeNodalNeighborValues();
366  }
367  }
368 }
369 
370 void
371 SystemBase::reinitNodes(const std::vector<dof_id_type> & nodes, THREAD_ID tid)
372 {
373  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
374  for (const auto & var : vars)
375  {
376  var->reinitNodes(nodes);
377  var->computeNodalValues();
378  }
379 }
380 
381 void
382 SystemBase::reinitNodesNeighbor(const std::vector<dof_id_type> & nodes, THREAD_ID tid)
383 {
384  const std::vector<MooseVariable *> & vars = _vars[tid].variables();
385  for (const auto & var : vars)
386  {
387  var->reinitNodesNeighbor(nodes);
388  var->computeNodalNeighborValues();
389  }
390 }
391 
392 void
394 {
395  const std::vector<MooseVariableScalar *> & vars = _vars[tid].scalars();
396  for (const auto & var : vars)
397  var->reinit();
398 }
399 
400 void
401 SystemBase::augmentSendList(std::vector<dof_id_type> & send_list)
402 {
403  std::set<dof_id_type> & ghosted_elems = _subproblem.ghostedElems();
404 
405  DofMap & dof_map = dofMap();
406 
407  std::vector<dof_id_type> dof_indices;
408 
409  System & sys = system();
410 
411  unsigned int sys_num = sys.number();
412 
413  unsigned int n_vars = sys.n_vars();
414 
415  for (const auto & elem_id : ghosted_elems)
416  {
417  Elem * elem = _mesh.elemPtr(elem_id);
418 
419  if (elem->active())
420  {
421  dof_map.dof_indices(elem, dof_indices);
422 
423  // Only need to ghost it if it's actually not on this processor
424  for (const auto & dof : dof_indices)
425  if (dof < dof_map.first_dof() || dof >= dof_map.end_dof())
426  send_list.push_back(dof);
427 
428  // Now add the DoFs from all of the nodes. This is necessary because of block
429  // restricted variables. A variable might not live _on_ this element but it
430  // might live on nodes connected to this element.
431  for (unsigned int n = 0; n < elem->n_nodes(); n++)
432  {
433  Node * node = elem->node_ptr(n);
434 
435  // Have to get each variable's dofs
436  for (unsigned int v = 0; v < n_vars; v++)
437  {
438  const Variable & var = sys.variable(v);
439  unsigned int var_num = var.number();
440  unsigned int n_comp = var.n_components();
441 
442  // See if this variable has any dofs at this node
443  if (node->n_dofs(sys_num, var_num) > 0)
444  {
445  // Loop over components of the variable
446  for (unsigned int c = 0; c < n_comp; c++)
447  send_list.push_back(node->dof_number(sys_num, var_num, c));
448  }
449  }
450  }
451  }
452  }
453 }
454 
458 void
460 {
461  if (!_saved_old)
462  _saved_old = &addVector("save_solution_old", false, PARALLEL);
463  if (!_saved_older)
464  _saved_older = &addVector("save_solution_older", false, PARALLEL);
465  *_saved_old = solutionOld();
467 }
468 
472 void
474 {
475  if (_saved_old)
476  {
477  solutionOld() = *_saved_old;
478  removeVector("save_solution_old");
479  _saved_old = NULL;
480  }
481  if (_saved_older)
482  {
484  removeVector("save_solution_older");
485  _saved_older = NULL;
486  }
487 }
488 
489 NumericVector<Number> &
490 SystemBase::addVector(const std::string & vector_name, const bool project, const ParallelType type)
491 {
492  if (hasVector(vector_name))
493  return getVector(vector_name);
494 
495  NumericVector<Number> * vec = &system().add_vector(vector_name, project, type);
496  return *vec;
497 }
498 
499 void
500 SystemBase::addVariable(const std::string & var_name,
501  const FEType & type,
502  Real scale_factor,
503  const std::set<SubdomainID> * const active_subdomains)
504 {
505  unsigned int var_num = system().add_variable(var_name, type, active_subdomains);
506  for (THREAD_ID tid = 0; tid < libMesh::n_threads(); tid++)
507  {
508  // FIXME: we cannot refer fetype in libMesh at this point, so we will just make a copy in
509  // MooseVariableBase.
510  MooseVariableBase * var;
511  if (type == FEType(0, MONOMIAL))
512  var = new MooseVariableConstMonomial(
513  var_num, type, *this, _subproblem.assembly(tid), _var_kind);
514  else
515  var = new MooseVariable(var_num, type, *this, _subproblem.assembly(tid), _var_kind);
516 
517  var->scalingFactor(scale_factor);
518  _vars[tid].add(var_name, var);
519  }
520  if (active_subdomains == NULL)
521  _var_map[var_num] = std::set<SubdomainID>();
522  else
523  for (std::set<SubdomainID>::iterator it = active_subdomains->begin();
524  it != active_subdomains->end();
525  ++it)
526  _var_map[var_num].insert(*it);
527 }
528 
529 void
530 SystemBase::addScalarVariable(const std::string & var_name,
531  Order order,
532  Real scale_factor,
533  const std::set<SubdomainID> * const active_subdomains)
534 {
535  FEType type(order, SCALAR);
536  unsigned int var_num = system().add_variable(var_name, type, active_subdomains);
537  for (THREAD_ID tid = 0; tid < libMesh::n_threads(); tid++)
538  {
539  // FIXME: we cannot refer fetype in libMesh at this point, so we will just make a copy in
540  // MooseVariableBase.
541  MooseVariableScalar * var =
542  new MooseVariableScalar(var_num, type, *this, _subproblem.assembly(tid), _var_kind);
543  var->scalingFactor(scale_factor);
544  _vars[tid].add(var_name, var);
545  }
546  if (active_subdomains == NULL)
547  _var_map[var_num] = std::set<SubdomainID>();
548  else
549  for (std::set<SubdomainID>::iterator it = active_subdomains->begin();
550  it != active_subdomains->end();
551  ++it)
552  _var_map[var_num].insert(*it);
553 }
554 
555 bool
556 SystemBase::hasVariable(const std::string & var_name)
557 {
558  if (system().has_variable(var_name))
559  return system().variable_type(var_name).family != SCALAR;
560  else
561  return false;
562 }
563 
564 bool
565 SystemBase::hasScalarVariable(const std::string & var_name)
566 {
567  if (system().has_variable(var_name))
568  return system().variable_type(var_name).family == SCALAR;
569  else
570  return false;
571 }
572 
573 bool
574 SystemBase::isScalarVariable(unsigned int var_num)
575 {
576  return (system().variable(var_num).type().family == SCALAR);
577 }
578 
579 unsigned int
581 {
582  return _vars[0].names().size();
583 }
584 
588 bool
589 SystemBase::hasVector(const std::string & name)
590 {
591  return system().have_vector(name);
592 }
593 
597 NumericVector<Number> &
598 SystemBase::getVector(const std::string & name)
599 {
600  return system().get_vector(name);
601 }
602 
603 unsigned int
605 {
606  return system().number();
607 }
608 
609 DofMap &
611 {
612  return system().get_dof_map();
613 }
614 
615 void
616 SystemBase::addVariableToCopy(const std::string & dest_name,
617  const std::string & source_name,
618  const std::string & timestep)
619 {
620  _var_to_copy.push_back(VarCopyInfo(dest_name, source_name, timestep));
621 }
622 
623 void
624 SystemBase::copyVars(ExodusII_IO & io)
625 {
626  int n_steps = io.get_num_time_steps();
627 
628  bool did_copy = false;
629  for (std::vector<VarCopyInfo>::iterator it = _var_to_copy.begin(); it != _var_to_copy.end(); ++it)
630  {
631  VarCopyInfo & vci = *it;
632  int timestep = -1;
633 
634  if (vci._timestep == "LATEST")
635  // Use the last time step in the file from which to retrieve the solution
636  timestep = n_steps;
637  else
638  {
639  std::istringstream ss(vci._timestep);
640  if (!((ss >> timestep) && ss.eof()) || timestep > n_steps)
641  mooseError("Invalid value passed as \"initial_from_file_timestep\". Expected \"LATEST\" or "
642  "a valid integer between 1 and ",
643  n_steps,
644  " inclusive, received ",
645  vci._timestep);
646  }
647 
648  did_copy = true;
649  if (getVariable(0, vci._dest_name).isNodal())
650  io.copy_nodal_solution(system(), vci._dest_name, vci._source_name, timestep);
651  else
652  io.copy_elemental_solution(system(), vci._dest_name, vci._source_name, timestep);
653  }
654 
655  if (did_copy)
656  solution().close();
657 }
658 
659 void
661 {
662 }
663 
664 void
666 {
667  system().update();
668 }
669 
670 void
672 {
673  system().solve();
674 }
675 
679 void
681 {
682  system().update();
687 }
688 
692 void
694 {
699 }
700 
704 void
706 {
707  *(const_cast<NumericVector<Number> *&>(currentSolution())) = solutionOld();
708  solution() = solutionOld();
711  system().update();
712 }
virtual const std::set< SubdomainID > * getVariableBlocks(unsigned int var_number)
Get the block where a variable of this system is defined.
Definition: SystemBase.C:140
virtual void reinitNode(const Node *node, THREAD_ID tid)
Reinit nodal assembly info.
Definition: SystemBase.C:329
ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
Definition: MooseMesh.C:685
virtual void reinitScalars(THREAD_ID tid)
Reinit scalar varaibles.
Definition: SystemBase.C:393
virtual void reinitNeighborFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid)
Compute the values of the variables at all the current points.
Definition: SystemBase.C:310
std::string _source_name
Definition: SystemBase.h:83
virtual void copyOldSolutions()
Shifts the solutions backwards in time.
Definition: SystemBase.C:693
Class for stuff related to variables.
Definition: MooseVariable.h:43
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:2062
virtual void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid)
Reinit assembly info for a side of an element.
Definition: SystemBase.C:299
virtual NumericVector< Number > & addVector(const std::string &vector_name, const bool project, const ParallelType type)
Adds a solution length vector to the system.
Definition: SystemBase.C:490
virtual bool hasScalarVariable(const std::string &var_name)
Definition: SystemBase.C:565
std::vector< std::string > _vars_to_be_zeroed_on_residual
Definition: SystemBase.h:492
virtual void addVariable(const std::string &var_name, const FEType &type, Real scale_factor, const std::set< SubdomainID > *const active_subdomains=NULL)
Adds a variable to the system.
Definition: SystemBase.C:500
void dataStore(std::ostream &stream, SystemBase &system_base, void *context)
IO Methods for restart, backup and restore.
Definition: SystemBase.C:52
virtual NumericVector< Number > & solutionOld()=0
virtual const std::string & name()
Definition: SystemBase.h:453
virtual void reinitElem(const Elem *elem, THREAD_ID tid)
Reinit an element assembly info.
Definition: SystemBase.C:279
void mooseError(Args &&...args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:182
virtual Assembly & assembly(THREAD_ID tid)=0
virtual Order getMinQuadratureOrder()
Get minimal quadrature order needed for integrating variables in this system.
Definition: SystemBase.C:199
virtual void addVariableToCopy(const std::string &dest_name, const std::string &source_name, const std::string &timestep)
Add info about variable that will be copied.
Definition: SystemBase.C:616
virtual bool hasActiveElementalMooseVariables(THREAD_ID tid)
Whether or not a list of active elemental moose variables has been set.
Definition: SubProblem.C:76
virtual void prepareFace(THREAD_ID tid, bool resize_data)
Prepare the system for use on sides.
Definition: SystemBase.C:237
virtual void removeVector(const std::string &name)
Remove a vector from the system with the given name.
Definition: SystemBase.h:435
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
Base class for a system (of equations)
Definition: SystemBase.h:91
virtual void reinitNodes(const std::vector< dof_id_type > &nodes, THREAD_ID tid)
Reinit variables at a set of nodes.
Definition: SystemBase.C:371
virtual void reinitNodesNeighbor(const std::vector< dof_id_type > &nodes, THREAD_ID tid)
Reinit variables at a set of neighbor nodes.
Definition: SystemBase.C:382
NumericVector< Real > * _saved_older
Definition: SystemBase.h:501
virtual NumericVector< Number > * solutionPreviousNewton()=0
void prepareVariable(MooseVariable *var)
Used for preparing the dense residual and jacobian blocks for one particular variable.
Definition: Assembly.C:1076
virtual void update()
Update the system (doing libMesh magic)
Definition: SystemBase.C:665
Grab all the local dof indices for the variables passed in, in the system passed in.
virtual bool isNodal() const override
Is this variable nodal.
std::vector< VarCopyInfo > _var_to_copy
Definition: SystemBase.h:506
virtual DofMap & dofMap()
Gets the dof map.
Definition: SystemBase.C:610
virtual bool hasVector(const std::string &name)
Check if the named vector exists in the system.
Definition: SystemBase.C:589
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
virtual void reinitNodeNeighbor(const Node *node, THREAD_ID tid)
Reinit nodal assembly info for neighbor node.
Definition: SystemBase.C:357
virtual bool isScalarVariable(unsigned int var_name)
Definition: SystemBase.C:574
virtual void zeroVariablesForResidual()
Zero out the solution for the variables that were registered as needing to have their solutions zeroe...
Definition: SystemBase.C:187
std::set< dof_id_type > _all_dof_indices
std::vector< std::string > _vars_to_be_zeroed_on_jacobian
Definition: SystemBase.h:493
Information about variables that will be copied.
Definition: SystemBase.h:73
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:156
virtual void prepareNeighbor(THREAD_ID tid)
Prepare the system for use.
Definition: SystemBase.C:271
DofMap & dof_map
const std::vector< numeric_index_type > & n_nz
std::string _dest_name
Definition: SystemBase.h:82
SubProblem & _subproblem
Definition: SystemBase.h:478
virtual void augmentSparsity(SparsityPattern::Graph &sparsity, std::vector< dof_id_type > &n_nz, std::vector< dof_id_type > &n_oz)=0
Will modify the sparsity pattern to add logical geometric connections.
Moose::VarKindType _var_kind
default kind of variables in this system
Definition: SystemBase.h:504
virtual bool checkNonlocalCouplingRequirement()
Definition: SubProblem.h:62
virtual std::set< dof_id_type > & ghostedElems()
Return the list of elements that should have their DoFs ghosted to this processor.
Definition: SubProblem.h:400
virtual System & system()=0
Get the reference to the libMesh system.
std::string _timestep
Definition: SystemBase.h:84
virtual bool hasVariable(const std::string &var_name)
Query a system for a variable.
Definition: SystemBase.C:556
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:66
void copyVars(ExodusII_IO &io)
Definition: SystemBase.C:624
virtual void reinitNeighbor(const Elem *elem, THREAD_ID tid)
Compute the values of the variables at all the current points.
Definition: SystemBase.C:321
void extraSendList(std::vector< dof_id_type > &send_list, void *context)
Free function used for a libMesh callback.
Definition: SystemBase.C:33
std::map< unsigned int, std::set< SubdomainID > > _var_map
Map of variables (variable id -> array of subdomains where it lives)
Definition: SystemBase.h:490
virtual void copySolutionsBackwards()
Copy current solution into old and older.
Definition: SystemBase.C:680
virtual void solve()
Solve the system (using libMesh magic)
Definition: SystemBase.C:671
virtual unsigned int number()
Gets the number of this system.
Definition: SystemBase.C:604
virtual void saveOldSolutions()
Save the old and older solutions.
Definition: SystemBase.C:459
std::vector< VariableWarehouse > _vars
Variable warehouses (one for each thread)
Definition: SystemBase.h:488
MatType type
virtual void addExtraVectors()
Method called during initialSetup to add extra system vector if they are required by the simulation...
Definition: SystemBase.C:660
SystemBase(SubProblem &subproblem, const std::string &name, Moose::VarKindType var_kind)
Definition: SystemBase.C:84
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
void prepareVariableNonlocal(MooseVariable *var)
Definition: Assembly.C:1098
MooseMesh & _mesh
Definition: SystemBase.h:483
virtual void addVariableToZeroOnResidual(std::string var_name)
Adds this variable to the list of variables to be zeroed during each residual evaluation.
Definition: SystemBase.C:150
PetscInt n
virtual void reinitNodeFace(const Node *node, BoundaryID bnd_id, THREAD_ID tid)
Reinit nodal assembly info on a face.
Definition: SystemBase.C:343
virtual void zeroVariablesForJacobian()
Zero out the solution for the variables that were registered as needing to have their solutions zeroe...
Definition: SystemBase.C:193
NumericVector< Real > * _saved_old
Definition: SystemBase.h:500
Class for scalar variables (they are different).
virtual NumericVector< Number > & solutionOlder()=0
virtual NumericVector< Number > & solution()=0
virtual const std::set< MooseVariable * > & getActiveElementalMooseVariables(THREAD_ID tid)
Get the MOOSE variables to be reinited on each element.
Definition: SubProblem.C:70
virtual void prepare(THREAD_ID tid)
Prepare the system for use.
Definition: SystemBase.C:214
virtual void addVariableToZeroOnJacobian(std::string var_name)
Adds this variable to the list of variables to be zeroed during each Jacobian evaluation.
Definition: SystemBase.C:156
virtual void augmentSendList(std::vector< dof_id_type > &send_list)
Will modify the send_list to add all of the extra ghosted dofs for this system.
Definition: SystemBase.C:401
virtual void restoreSolutions()
Restore current solutions (call after your solve failed)
Definition: SystemBase.C:705
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector.
Definition: SystemBase.C:598
void extraSparsity(SparsityPattern::Graph &sparsity, std::vector< dof_id_type > &n_nz, std::vector< dof_id_type > &n_oz, void *context)
Free function used for a libMesh callback.
Definition: SystemBase.C:41
void dataLoad(std::istream &stream, SystemBase &system_base, void *context)
IO Methods for restart, backup and restore.
Definition: SystemBase.C:68
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:121
virtual void zeroVariables(std::vector< std::string > &vars_to_be_zeroed)
Zero out the solution for the list of variables passed in.
Definition: SystemBase.C:162
virtual void restoreOldSolutions()
Restore the old and older solutions when the saved solutions present.
Definition: SystemBase.C:473
virtual unsigned int nVariables()
Get the number of variables in this system.
Definition: SystemBase.C:580
const std::vector< numeric_index_type > & n_oz
virtual void addScalarVariable(const std::string &var_name, Order order, Real scale_factor, const std::set< SubdomainID > *const active_subdomains=NULL)
Adds a scalar variable.
Definition: SystemBase.C:530
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
unsigned int THREAD_ID
Definition: MooseTypes.h:79
void scalingFactor(Real factor)
Set the scaling factor for this variable.
virtual const NumericVector< Number > *& currentSolution()=0
The solution vector that is currently being operated on.