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