www.mooseframework.org
MultiAppProjectionTransfer.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 
16 
17 // MOOSE includes
18 #include "AddVariableAction.h"
19 #include "FEProblem.h"
20 #include "MooseMesh.h"
21 #include "MooseVariable.h"
22 #include "SystemBase.h"
23 
24 #include "libmesh/dof_map.h"
25 #include "libmesh/linear_implicit_system.h"
26 #include "libmesh/mesh_function.h"
27 #include "libmesh/mesh_tools.h"
28 #include "libmesh/numeric_vector.h"
29 #include "libmesh/parallel_algebra.h"
30 #include "libmesh/quadrature_gauss.h"
31 #include "libmesh/sparse_matrix.h"
32 #include "libmesh/string_to_enum.h"
33 
34 void
35 assemble_l2(EquationSystems & es, const std::string & system_name)
36 {
37  MultiAppProjectionTransfer * transfer =
38  es.parameters.get<MultiAppProjectionTransfer *>("transfer");
39  transfer->assembleL2(es, system_name);
40 }
41 
42 template <>
45 {
47  params.addRequiredParam<AuxVariableName>(
48  "variable", "The auxiliary variable to store the transferred values in.");
49  params.addRequiredParam<VariableName>("source_variable", "The variable to transfer from.");
50 
51  MooseEnum proj_type("l2", "l2");
52  params.addParam<MooseEnum>("proj_type", proj_type, "The type of the projection.");
53 
54  params.addParam<bool>("fixed_meshes",
55  false,
56  "Set to true when the meshes are not changing (ie, "
57  "no movement or adaptivity). This will cache some "
58  "information to speed up the transfer.");
59 
60  return params;
61 }
62 
64  : MultiAppTransfer(parameters),
65  _to_var_name(getParam<AuxVariableName>("variable")),
66  _from_var_name(getParam<VariableName>("source_variable")),
67  _proj_type(getParam<MooseEnum>("proj_type")),
68  _compute_matrix(true),
69  _fixed_meshes(getParam<bool>("fixed_meshes")),
70  _qps_cached(false)
71 {
72 }
73 
74 void
76 {
77  getAppInfo();
78 
79  _proj_sys.resize(_to_problems.size(), NULL);
80 
81  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
82  {
83  FEProblemBase & to_problem = *_to_problems[i_to];
84  EquationSystems & to_es = to_problem.es();
85 
86  // Add the projection system.
87  FEType fe_type = to_problem.getVariable(0, _to_var_name).feType();
88  LinearImplicitSystem & proj_sys = to_es.add_system<LinearImplicitSystem>("proj-sys-" + name());
89  _proj_var_num = proj_sys.add_variable("var", fe_type);
90  proj_sys.attach_assemble_function(assemble_l2);
91  _proj_sys[i_to] = &proj_sys;
92 
93  // Prevent the projection system from being written to checkpoint
94  // files. In the event of a recover or restart, we'll read the checkpoint
95  // before this initialSetup method is called. As a result, we'll find
96  // systems in the checkpoint (the projection systems) that we don't know
97  // what to do with, and there will be a crash. We could fix this by making
98  // the systems in the constructor, except we don't know how many sub apps
99  // there are at the time of construction. So instead, we'll just nuke the
100  // projection system and rebuild it from scratch every recover/restart.
101  proj_sys.hide_output() = true;
102 
103  // Reinitialize EquationSystems since we added a system.
104  to_es.reinit();
105  }
106 
107  if (_fixed_meshes)
108  {
109  _cached_qps.resize(n_processors());
110  _cached_index_map.resize(n_processors());
111  }
112 }
113 
114 void
115 MultiAppProjectionTransfer::assembleL2(EquationSystems & es, const std::string & system_name)
116 {
117  // Get the system and mesh from the input arguments.
118  LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>(system_name);
119  MeshBase & to_mesh = es.get_mesh();
120 
121  // Get the meshfunction evaluations and the map that was stashed in the es.
122  std::vector<Real> & final_evals = *es.parameters.get<std::vector<Real> *>("final_evals");
123  std::map<dof_id_type, unsigned int> & element_map =
124  *es.parameters.get<std::map<dof_id_type, unsigned int> *>("element_map");
125 
126  // Setup system vectors and matrices.
127  FEType fe_type = system.variable_type(0);
128  std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
129  QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
130  fe->attach_quadrature_rule(&qrule);
131  const DofMap & dof_map = system.get_dof_map();
132  DenseMatrix<Number> Ke;
133  DenseVector<Number> Fe;
134  std::vector<dof_id_type> dof_indices;
135  const std::vector<Real> & JxW = fe->get_JxW();
136  const std::vector<std::vector<Real>> & phi = fe->get_phi();
137 
138  for (const auto & elem : to_mesh.active_local_element_ptr_range())
139  {
140  fe->reinit(elem);
141 
142  dof_map.dof_indices(elem, dof_indices);
143  Ke.resize(dof_indices.size(), dof_indices.size());
144  Fe.resize(dof_indices.size());
145 
146  for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
147  {
148  Real meshfun_eval = 0.;
149  if (element_map.find(elem->id()) != element_map.end())
150  {
151  // We have evaluations for this element.
152  meshfun_eval = final_evals[element_map[elem->id()] + qp];
153  }
154 
155  // Now compute the element matrix and RHS contributions.
156  for (unsigned int i = 0; i < phi.size(); i++)
157  {
158  // RHS
159  Fe(i) += JxW[qp] * (meshfun_eval * phi[i][qp]);
160 
161  if (_compute_matrix)
162  for (unsigned int j = 0; j < phi.size(); j++)
163  {
164  // The matrix contribution
165  Ke(i, j) += JxW[qp] * (phi[i][qp] * phi[j][qp]);
166  }
167  }
168  dof_map.constrain_element_matrix_and_vector(Ke, Fe, dof_indices);
169 
170  if (_compute_matrix)
171  system.matrix->add_matrix(Ke, dof_indices);
172  system.rhs->add_vector(Fe, dof_indices);
173  }
174  }
175 }
176 
177 void
179 {
180  _console << "Beginning projection transfer " << name() << std::endl;
181 
182  getAppInfo();
183 
185  // We are going to project the solutions by solving some linear systems. In
186  // order to assemble the systems, we need to evaluate the "from" domain
187  // solutions at quadrature points in the "to" domain. Some parallel
188  // communication is necessary because each processor doesn't necessarily have
189  // all the "from" information it needs to set its "to" values. We don't want
190  // to use a bunch of big all-to-all broadcasts, so we'll use bounding boxes to
191  // figure out which processors have the information we need and only
192  // communicate with those processors.
193  //
194  // Each processor will
195  // 1. Check its local quadrature points in the "to" domains to see which
196  // "from" domains they might be in.
197  // 2. Send quadrature points to the processors with "from" domains that might
198  // contain those points.
199  // 3. Recieve quadrature points from other processors, evaluate its mesh
200  // functions at those points, and send the values back to the proper
201  // processor
202  // 4. Recieve mesh function evaluations from all relevant processors and
203  // decide which one to use at every quadrature point (the lowest global app
204  // index always wins)
205  // 5. And use the mesh function evaluations to assemble and solve an L2
206  // projection system on its local elements.
208 
210  // For every combination of global "from" problem and local "to" problem, find
211  // which "from" bounding boxes overlap with which "to" elements. Keep track
212  // of which processors own bounding boxes that overlap with which elements.
213  // Build vectors of quadrature points to send to other processors for mesh
214  // function evaluations.
216 
217  // Get the bounding boxes for the "from" domains.
218  std::vector<BoundingBox> bboxes = getFromBoundingBoxes();
219 
220  // Figure out how many "from" domains each processor owns.
221  std::vector<unsigned int> froms_per_proc = getFromsPerProc();
222 
223  std::vector<std::vector<Point>> outgoing_qps(n_processors());
224  std::vector<std::map<std::pair<unsigned int, unsigned int>, unsigned int>> element_index_map(
225  n_processors());
226  // element_index_map[i_to, element_id] = index
227  // outgoing_qps[index] is the first quadrature point in element
228 
229  if (!_qps_cached)
230  {
231  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
232  {
233  MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
234 
235  LinearImplicitSystem & system = *_proj_sys[i_to];
236 
237  FEType fe_type = system.variable_type(0);
238  std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
239  QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
240  fe->attach_quadrature_rule(&qrule);
241  const std::vector<Point> & xyz = fe->get_xyz();
242 
243  MeshBase::const_element_iterator el = to_mesh.local_elements_begin();
244  const MeshBase::const_element_iterator end_el = to_mesh.local_elements_end();
245 
246  unsigned int from0 = 0;
247  for (processor_id_type i_proc = 0; i_proc < n_processors();
248  from0 += froms_per_proc[i_proc], i_proc++)
249  {
250  for (el = to_mesh.local_elements_begin(); el != end_el; el++)
251  {
252  const Elem * elem = *el;
253  fe->reinit(elem);
254 
255  bool qp_hit = false;
256  for (unsigned int i_from = 0; i_from < froms_per_proc[i_proc] && !qp_hit; i_from++)
257  {
258  for (unsigned int qp = 0; qp < qrule.n_points() && !qp_hit; qp++)
259  {
260  Point qpt = xyz[qp];
261  if (bboxes[from0 + i_from].contains_point(qpt + _to_positions[i_to]))
262  qp_hit = true;
263  }
264  }
265 
266  if (qp_hit)
267  {
268  // The selected processor's bounding box contains at least one
269  // quadrature point from this element. Send all qps from this element
270  // and remember where they are in the array using the map.
271  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
272  element_index_map[i_proc][key] = outgoing_qps[i_proc].size();
273  for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
274  {
275  Point qpt = xyz[qp];
276  outgoing_qps[i_proc].push_back(qpt + _to_positions[i_to]);
277  }
278  }
279  }
280  }
281  }
282 
283  if (_fixed_meshes)
284  _cached_index_map = element_index_map;
285  }
286  else
287  {
288  element_index_map = _cached_index_map;
289  }
290 
292  // Request quadrature point evaluations from other processors and handle
293  // requests sent to this processor.
295 
296  // Non-blocking send quadrature points to other processors.
297  std::vector<Parallel::Request> send_qps(n_processors());
298  if (!_qps_cached)
299  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
300  if (i_proc != processor_id())
301  _communicator.send(i_proc, outgoing_qps[i_proc], send_qps[i_proc]);
302 
303  // Get the local bounding boxes.
304  std::vector<BoundingBox> local_bboxes(froms_per_proc[processor_id()]);
305  {
306  // Find the index to the first of this processor's local bounding boxes.
307  unsigned int local_start = 0;
308  for (processor_id_type i_proc = 0; i_proc < n_processors() && i_proc != processor_id();
309  i_proc++)
310  local_start += froms_per_proc[i_proc];
311 
312  // Extract the local bounding boxes.
313  for (unsigned int i_from = 0; i_from < froms_per_proc[processor_id()]; i_from++)
314  local_bboxes[i_from] = bboxes[local_start + i_from];
315  }
316 
317  // Setup the local mesh functions.
318  std::vector<MeshFunction *> local_meshfuns(froms_per_proc[processor_id()], NULL);
319  for (unsigned int i_from = 0; i_from < _from_problems.size(); i_from++)
320  {
321  FEProblemBase & from_problem = *_from_problems[i_from];
322  MooseVariable & from_var = from_problem.getVariable(0, _from_var_name);
323  System & from_sys = from_var.sys().system();
324  unsigned int from_var_num = from_sys.variable_number(from_var.name());
325 
326  MeshFunction * from_func = new MeshFunction(
327  from_problem.es(), *from_sys.current_local_solution, from_sys.get_dof_map(), from_var_num);
328  from_func->init(Trees::ELEMENTS);
329  from_func->enable_out_of_mesh_mode(OutOfMeshValue);
330  local_meshfuns[i_from] = from_func;
331  }
332 
333  // Recieve quadrature points from other processors, evaluate mesh frunctions
334  // at those points, and send the values back.
335  std::vector<Parallel::Request> send_evals(n_processors());
336  std::vector<Parallel::Request> send_ids(n_processors());
337  std::vector<std::vector<Real>> outgoing_evals(n_processors());
338  std::vector<std::vector<unsigned int>> outgoing_ids(n_processors());
339  std::vector<std::vector<Real>> incoming_evals(n_processors());
340  std::vector<std::vector<unsigned int>> incoming_app_ids(n_processors());
341  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
342  {
343  // Use the cached qps if they're available.
344  std::vector<Point> incoming_qps;
345  if (!_qps_cached)
346  {
347  if (i_proc == processor_id())
348  incoming_qps = outgoing_qps[i_proc];
349  else
350  _communicator.receive(i_proc, incoming_qps);
351  // Cache these qps for later if _fixed_meshes
352  if (_fixed_meshes)
353  _cached_qps[i_proc] = incoming_qps;
354  }
355  else
356  {
357  incoming_qps = _cached_qps[i_proc];
358  }
359 
360  outgoing_evals[i_proc].resize(incoming_qps.size(), OutOfMeshValue);
361  if (_direction == FROM_MULTIAPP)
362  outgoing_ids[i_proc].resize(incoming_qps.size(), libMesh::invalid_uint);
363  for (unsigned int qp = 0; qp < incoming_qps.size(); qp++)
364  {
365  Point qpt = incoming_qps[qp];
366 
367  // Loop until we've found the lowest-ranked app that actually contains
368  // the quadrature point.
369  for (unsigned int i_from = 0; i_from < _from_problems.size(); i_from++)
370  {
371  if (local_bboxes[i_from].contains_point(qpt))
372  {
373  outgoing_evals[i_proc][qp] = (*local_meshfuns[i_from])(qpt - _from_positions[i_from]);
374  if (_direction == FROM_MULTIAPP)
375  outgoing_ids[i_proc][qp] = _local2global_map[i_from];
376  }
377  }
378  }
379 
380  if (i_proc == processor_id())
381  {
382  incoming_evals[i_proc] = outgoing_evals[i_proc];
383  if (_direction == FROM_MULTIAPP)
384  incoming_app_ids[i_proc] = outgoing_ids[i_proc];
385  }
386  else
387  {
388  _communicator.send(i_proc, outgoing_evals[i_proc], send_evals[i_proc]);
389  if (_direction == FROM_MULTIAPP)
390  _communicator.send(i_proc, outgoing_ids[i_proc], send_ids[i_proc]);
391  }
392  }
393 
395  // Gather all of the qp evaluations and pick out the best ones for each qp.
397  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
398  {
399  if (i_proc == processor_id())
400  continue;
401  _communicator.receive(i_proc, incoming_evals[i_proc]);
402  if (_direction == FROM_MULTIAPP)
403  _communicator.receive(i_proc, incoming_app_ids[i_proc]);
404  }
405 
406  std::vector<std::vector<Real>> final_evals(_to_problems.size());
407  std::vector<std::map<dof_id_type, unsigned int>> trimmed_element_maps(_to_problems.size());
408 
409  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
410  {
411  MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
412  LinearImplicitSystem & system = *_proj_sys[i_to];
413 
414  FEType fe_type = system.variable_type(0);
415  std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
416  QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
417  fe->attach_quadrature_rule(&qrule);
418  const std::vector<Point> & xyz = fe->get_xyz();
419 
420  for (const auto & elem : to_mesh.active_local_element_ptr_range())
421  {
422  fe->reinit(elem);
423 
424  bool element_is_evaled = false;
425  std::vector<Real> evals(qrule.n_points(), 0.);
426 
427  for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
428  {
429  Point qpt = xyz[qp];
430 
431  unsigned int lowest_app_rank = libMesh::invalid_uint;
432  for (unsigned int i_proc = 0; i_proc < n_processors(); i_proc++)
433  {
434  // Ignore the selected processor if the element wasn't found in it's
435  // bounding box.
436  std::map<std::pair<unsigned int, unsigned int>, unsigned int> & map =
437  element_index_map[i_proc];
438  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
439  if (map.find(key) == map.end())
440  continue;
441  unsigned int qp0 = map[key];
442 
443  // Ignore the selected processor if it's app has a higher rank than the
444  // previously found lowest app rank.
445  if (_direction == FROM_MULTIAPP)
446  if (incoming_app_ids[i_proc][qp0 + qp] >= lowest_app_rank)
447  continue;
448 
449  // Ignore the selected processor if the qp was actually outside the
450  // processor's subapp's mesh.
451  if (incoming_evals[i_proc][qp0 + qp] == OutOfMeshValue)
452  continue;
453 
454  // This is the best meshfunction evaluation so far, save it.
455  element_is_evaled = true;
456  evals[qp] = incoming_evals[i_proc][qp0 + qp];
457  }
458  }
459 
460  // If we found good evaluations for any of the qps in this element, save
461  // those evaluations for later.
462  if (element_is_evaled)
463  {
464  trimmed_element_maps[i_to][elem->id()] = final_evals[i_to].size();
465  for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
466  final_evals[i_to].push_back(evals[qp]);
467  }
468  }
469  }
470 
472  // We now have just one or zero mesh function values at all of our local
473  // quadrature points. Stash those values (and a map linking them to element
474  // ids) in the equation systems parameters and project the solution.
476 
477  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
478  {
479  _to_es[i_to]->parameters.set<std::vector<Real> *>("final_evals") = &final_evals[i_to];
480  _to_es[i_to]->parameters.set<std::map<dof_id_type, unsigned int> *>("element_map") =
481  &trimmed_element_maps[i_to];
482  projectSolution(i_to);
483  _to_es[i_to]->parameters.set<std::vector<Real> *>("final_evals") = NULL;
484  _to_es[i_to]->parameters.set<std::map<dof_id_type, unsigned int> *>("element_map") = NULL;
485  }
486 
487  for (unsigned int i = 0; i < _from_problems.size(); i++)
488  delete local_meshfuns[i];
489 
490  // Make sure all our sends succeeded.
491  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
492  {
493  if (i_proc == processor_id())
494  continue;
495  if (!_qps_cached)
496  send_qps[i_proc].wait();
497  send_evals[i_proc].wait();
498  if (_direction == FROM_MULTIAPP)
499  send_ids[i_proc].wait();
500  }
501 
502  if (_fixed_meshes)
503  _qps_cached = true;
504 
505  _console << "Finished projection transfer " << name() << std::endl;
506 }
507 
508 void
510 {
511  FEProblemBase & to_problem = *_to_problems[i_to];
512  EquationSystems & proj_es = to_problem.es();
513  LinearImplicitSystem & ls = *_proj_sys[i_to];
514  // activate the current transfer
515  proj_es.parameters.set<MultiAppProjectionTransfer *>("transfer") = this;
516 
517  // TODO: specify solver params in an input file
518  // solver tolerance
519  Real tol = proj_es.parameters.get<Real>("linear solver tolerance");
520  proj_es.parameters.set<Real>("linear solver tolerance") = 1e-10; // set our tolerance
521  // solve it
522  ls.solve();
523  proj_es.parameters.set<Real>("linear solver tolerance") = tol; // restore the original tolerance
524 
525  // copy projected solution into target es
526  MeshBase & to_mesh = proj_es.get_mesh();
527 
528  MooseVariable & to_var = to_problem.getVariable(0, _to_var_name);
529  System & to_sys = to_var.sys().system();
530  NumericVector<Number> * to_solution = to_sys.solution.get();
531 
532  {
533  MeshBase::const_node_iterator it = to_mesh.local_nodes_begin();
534  const MeshBase::const_node_iterator end_it = to_mesh.local_nodes_end();
535  for (; it != end_it; ++it)
536  {
537  const Node * node = *it;
538  for (unsigned int comp = 0; comp < node->n_comp(to_sys.number(), to_var.number()); comp++)
539  {
540  const dof_id_type proj_index = node->dof_number(ls.number(), _proj_var_num, comp);
541  const dof_id_type to_index = node->dof_number(to_sys.number(), to_var.number(), comp);
542  to_solution->set(to_index, (*ls.solution)(proj_index));
543  }
544  }
545  }
546  for (const auto & elem : to_mesh.active_local_element_ptr_range())
547  for (unsigned int comp = 0; comp < elem->n_comp(to_sys.number(), to_var.number()); comp++)
548  {
549  const dof_id_type proj_index = elem->dof_number(ls.number(), _proj_var_num, comp);
550  const dof_id_type to_index = elem->dof_number(to_sys.number(), to_var.number(), comp);
551  to_solution->set(to_index, (*ls.solution)(proj_index));
552  }
553 
554  to_solution->close();
555  to_sys.update();
556 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
void assemble_l2(EquationSystems &es, const std::string &system_name)
const FEType & feType() const
Get the type of finite element object.
Class for stuff related to variables.
Definition: MooseVariable.h:43
std::vector< std::map< std::pair< unsigned int, unsigned int >, unsigned int > > _cached_index_map
std::vector< EquationSystems * > _to_es
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
std::vector< BoundingBox > getFromBoundingBoxes()
Return the bounding boxes of all the "from" domains, including all the domains not local to this proc...
std::vector< FEProblemBase * > _to_problems
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::vector< Point > _from_positions
InputParameters validParams< MultiAppProjectionTransfer >()
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
void assembleL2(EquationSystems &es, const std::string &system_name)
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
const std::string & name() const
Get the variable name.
virtual EquationSystems & es() override
nl system()
MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
Project values from one domain to another.
bool _compute_matrix
True, if we need to recompute the projection matrix.
DofMap & dof_map
std::vector< unsigned int > getFromsPerProc()
Return the number of "from" domains that each processor owns.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
virtual void execute() override
Execute the transfer.
virtual System & system()=0
Get the reference to the libMesh system.
std::vector< unsigned int > _local2global_map
static const Number OutOfMeshValue
Definition: Transfer.h:81
MultiAppProjectionTransfer(const InputParameters &parameters)
unsigned int number() const
Get variable number coming from libMesh.
std::vector< Point > _to_positions
InputParameters validParams< MultiAppTransfer >()
Base class for all MultiAppTransfer objects.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:53
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void projectSolution(unsigned int to_problem)
std::vector< LinearImplicitSystem * > _proj_sys
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
virtual void initialSetup() override
Method called at the beginning of the simulation for checking integrity or doing one-time setup...
unsigned int _proj_var_num
Having one projection variable number seems weird, but there is always one variable in every system b...
void getAppInfo()
This method will fill information into the convenience member variables (_to_problems, _from_meshes, etc.)
SystemBase & sys()
Get the system this variable is part of.
friend void assemble_l2(EquationSystems &es, const std::string &system_name)
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name) override
Returns the variable reference for requested variable which may be in any system. ...
std::vector< FEProblemBase * > _from_problems
std::vector< std::vector< Point > > _cached_qps
std::vector< MooseMesh * > _to_meshes