www.mooseframework.org
MultiAppNearestNodeTransfer.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 "DisplacedProblem.h"
19 #include "FEProblem.h"
20 #include "MooseMesh.h"
21 #include "MooseTypes.h"
22 #include "MooseVariable.h"
23 
24 #include "libmesh/system.h"
25 #include "libmesh/mesh_tools.h"
26 #include "libmesh/id_types.h"
27 #include "libmesh/parallel_algebra.h"
28 
29 template <>
32 {
34 
35  params.addRequiredParam<AuxVariableName>(
36  "variable", "The auxiliary variable to store the transferred values in.");
37  params.addRequiredParam<VariableName>("source_variable", "The variable to transfer from.");
38  params.addParam<BoundaryName>(
39  "source_boundary",
40  "The boundary we are transferring from (if not specified, whole domain is used).");
41  params.addParam<BoundaryName>(
42  "target_boundary",
43  "The boundary we are transferring to (if not specified, whole domain is used).");
44  params.addParam<bool>("displaced_source_mesh",
45  false,
46  "Whether or not to use the displaced mesh for the source mesh.");
47  params.addParam<bool>("displaced_target_mesh",
48  false,
49  "Whether or not to use the displaced mesh for the target mesh.");
50  params.addParam<bool>("fixed_meshes",
51  false,
52  "Set to true when the meshes are not changing (ie, "
53  "no movement or adaptivity). This will cache "
54  "nearest node neighbors to greatly speed up the "
55  "transfer.");
56 
57  return params;
58 }
59 
61  : MultiAppTransfer(parameters),
62  _to_var_name(getParam<AuxVariableName>("variable")),
63  _from_var_name(getParam<VariableName>("source_variable")),
64  _fixed_meshes(getParam<bool>("fixed_meshes")),
65  _node_map(declareRestartableData<std::map<dof_id_type, Node *>>("node_map")),
66  _distance_map(declareRestartableData<std::map<dof_id_type, Real>>("distance_map")),
67  _neighbors_cached(declareRestartableData<bool>("neighbors_cached", false)),
68  _cached_froms(declareRestartableData<std::vector<std::vector<unsigned int>>>("cached_froms")),
69  _cached_dof_ids(
70  declareRestartableData<std::vector<std::vector<dof_id_type>>>("cached_dof_ids")),
71  _cached_from_inds(
72  declareRestartableData<std::map<dof_id_type, unsigned int>>("cached_from_ids")),
73  _cached_qp_inds(declareRestartableData<std::map<dof_id_type, unsigned int>>("cached_qp_inds"))
74 {
75  // This transfer does not work with DistributedMesh
76  _displaced_source_mesh = getParam<bool>("displaced_source_mesh");
77  _displaced_target_mesh = getParam<bool>("displaced_target_mesh");
78 }
79 
80 void
82 {
83  if (_direction == TO_MULTIAPP)
85  else
87 }
88 
89 void
91 {
92  _console << "Beginning NearestNodeTransfer " << name() << std::endl;
93 
94  getAppInfo();
95 
96  // Get the bounding boxes for the "from" domains.
97  std::vector<BoundingBox> bboxes = getFromBoundingBoxes();
98 
99  // Figure out how many "from" domains each processor owns.
100  std::vector<unsigned int> froms_per_proc = getFromsPerProc();
101 
103  // For every point in the local "to" domain, figure out which "from" domains
104  // might contain it's nearest neighbor, and send that point to the processors
105  // that own those "from" domains.
106  //
107  // How do we know which "from" domains might contain the nearest neighbor, you
108  // ask? Well, consider two "from" domains, A and B. If every point in A is
109  // closer than every point in B, then we know that B cannot possibly contain
110  // the nearest neighbor. Hence, we'll only check A for the nearest neighbor.
111  // We'll use the functions bboxMaxDistance and bboxMinDistance to figure out
112  // if every point in A is closer than every point in B.
114 
115  // outgoing_qps = nodes/centroids we'll send to other processors.
116  std::vector<std::vector<Point>> outgoing_qps(n_processors());
117  // When we get results back, node_index_map will tell us which results go with
118  // which points
119  std::vector<std::map<std::pair<unsigned int, unsigned int>, unsigned int>> node_index_map(
120  n_processors());
121 
122  if (!_neighbors_cached)
123  {
124  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
125  {
126  System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
127  unsigned int sys_num = to_sys->number();
128  unsigned int var_num = to_sys->variable_number(_to_var_name);
129  MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
130  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
131 
132  if (is_nodal)
133  {
134  std::vector<Node *> target_local_nodes;
135 
136  if (isParamValid("target_boundary"))
137  {
138  BoundaryID target_bnd_id =
139  _to_meshes[i_to]->getBoundaryID(getParam<BoundaryName>("target_boundary"));
140 
141  ConstBndNodeRange & bnd_nodes = *(_to_meshes[i_to])->getBoundaryNodeRange();
142  for (const auto & bnode : bnd_nodes)
143  if (bnode->_bnd_id == target_bnd_id && bnode->_node->processor_id() == processor_id())
144  target_local_nodes.push_back(bnode->_node);
145  }
146  else
147  {
148  target_local_nodes.resize(to_mesh->n_local_nodes());
149  MeshBase::const_node_iterator nodes_begin = to_mesh->local_nodes_begin();
150  MeshBase::const_node_iterator nodes_end = to_mesh->local_nodes_end();
151 
152  unsigned int i = 0;
153  for (MeshBase::const_node_iterator nodes_it = nodes_begin; nodes_it != nodes_end;
154  ++nodes_it, ++i)
155  target_local_nodes[i] = *nodes_it;
156  }
157 
158  for (const auto & node : target_local_nodes)
159  {
160  // Skip this node if the variable has no dofs at it.
161  if (node->n_dofs(sys_num, var_num) < 1)
162  continue;
163 
164  // Find which bboxes might have the nearest node to this point.
165  Real nearest_max_distance = std::numeric_limits<Real>::max();
166  for (const auto & bbox : bboxes)
167  {
168  Real distance = bboxMaxDistance(*node, bbox);
169  if (distance < nearest_max_distance)
170  nearest_max_distance = distance;
171  }
172 
173  unsigned int from0 = 0;
174  for (processor_id_type i_proc = 0; i_proc < n_processors();
175  from0 += froms_per_proc[i_proc], i_proc++)
176  {
177  bool qp_found = false;
178  for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc] && !qp_found;
179  i_from++)
180  {
181  Real distance = bboxMinDistance(*node, bboxes[i_from]);
182  if (distance < nearest_max_distance || bboxes[i_from].contains_point(*node))
183  {
184  std::pair<unsigned int, unsigned int> key(i_to, node->id());
185  node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
186  outgoing_qps[i_proc].push_back(*node + _to_positions[i_to]);
187  qp_found = true;
188  }
189  }
190  }
191  }
192  }
193  else // Elemental
194  {
195  MeshBase::const_element_iterator elem_it = to_mesh->local_elements_begin();
196  MeshBase::const_element_iterator elem_end = to_mesh->local_elements_end();
197 
198  for (; elem_it != elem_end; ++elem_it)
199  {
200  Elem * elem = *elem_it;
201 
202  Point centroid = elem->centroid();
203 
204  // Skip this element if the variable has no dofs at it.
205  if (elem->n_dofs(sys_num, var_num) < 1)
206  continue;
207 
208  // Find which bboxes might have the nearest node to this point.
209  Real nearest_max_distance = std::numeric_limits<Real>::max();
210  for (const auto & bbox : bboxes)
211  {
212  Real distance = bboxMaxDistance(centroid, bbox);
213  if (distance < nearest_max_distance)
214  nearest_max_distance = distance;
215  }
216 
217  unsigned int from0 = 0;
218  for (processor_id_type i_proc = 0; i_proc < n_processors();
219  from0 += froms_per_proc[i_proc], i_proc++)
220  {
221  bool qp_found = false;
222  for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc] && !qp_found;
223  i_from++)
224  {
225  Real distance = bboxMinDistance(centroid, bboxes[i_from]);
226  if (distance < nearest_max_distance || bboxes[i_from].contains_point(centroid))
227  {
228  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
229  node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
230  outgoing_qps[i_proc].push_back(centroid + _to_positions[i_to]);
231  qp_found = true;
232  }
233  }
234  }
235  }
236  }
237  }
238  }
239 
241  // Send local node/centroid positions off to the other processors and take
242  // care of points sent to this processor. We'll need to check the points
243  // against all of the "from" domains that this processor owns. For each
244  // point, we'll find the nearest node, then we'll send the value at that node
245  // and the distance between the node and the point back to the processor that
246  // requested that point.
248 
249  std::vector<std::vector<Real>> incoming_evals(n_processors());
250  std::vector<Parallel::Request> send_qps(n_processors());
251  std::vector<Parallel::Request> send_evals(n_processors());
252 
253  // Create these here so that they live the entire life of this function
254  // and are NOT reused per processor.
255  std::vector<std::vector<Real>> processor_outgoing_evals(n_processors());
256 
257  if (!_neighbors_cached)
258  {
259  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
260  {
261  if (i_proc == processor_id())
262  continue;
263  _communicator.send(i_proc, outgoing_qps[i_proc], send_qps[i_proc]);
264  }
265 
266  // Build an array of pointers to all of this processor's local nodes. We
267  // need to do this to avoid the expense of using LibMesh iterators. This
268  // step also takes care of limiting the search to boundary nodes, if
269  // applicable.
270  std::vector<std::vector<Node *>> local_nodes(froms_per_proc[processor_id()]);
271  for (unsigned int i = 0; i < froms_per_proc[processor_id()]; i++)
272  {
273  getLocalNodes(_from_meshes[i], local_nodes[i]);
274  }
275 
276  if (_fixed_meshes)
277  {
278  _cached_froms.resize(n_processors());
279  _cached_dof_ids.resize(n_processors());
280  }
281 
282  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
283  {
284  std::vector<Point> incoming_qps;
285  if (i_proc == processor_id())
286  incoming_qps = outgoing_qps[i_proc];
287  else
288  _communicator.receive(i_proc, incoming_qps);
289 
290  if (_fixed_meshes)
291  {
292  _cached_froms[i_proc].resize(incoming_qps.size());
293  _cached_dof_ids[i_proc].resize(incoming_qps.size());
294  }
295 
296  std::vector<Real> & outgoing_evals = processor_outgoing_evals[i_proc];
297  outgoing_evals.resize(2 * incoming_qps.size());
298 
299  for (unsigned int qp = 0; qp < incoming_qps.size(); qp++)
300  {
301  Point qpt = incoming_qps[qp];
302  outgoing_evals[2 * qp] = std::numeric_limits<Real>::max();
303  for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
304  i_local_from++)
305  {
306  MooseVariable & from_var = _from_problems[i_local_from]->getVariable(0, _from_var_name);
307  System & from_sys = from_var.sys().system();
308  unsigned int from_sys_num = from_sys.number();
309  unsigned int from_var_num = from_sys.variable_number(from_var.name());
310 
311  for (unsigned int i_node = 0; i_node < local_nodes[i_local_from].size(); i_node++)
312  {
313  Real current_distance =
314  (qpt - *(local_nodes[i_local_from][i_node]) - _from_positions[i_local_from]).norm();
315  if (current_distance < outgoing_evals[2 * qp])
316  {
317  // Assuming LAGRANGE!
318  if (local_nodes[i_local_from][i_node]->n_dofs(from_sys_num, from_var_num) > 0)
319  {
320  dof_id_type from_dof =
321  local_nodes[i_local_from][i_node]->dof_number(from_sys_num, from_var_num, 0);
322 
323  outgoing_evals[2 * qp] = current_distance;
324  outgoing_evals[2 * qp + 1] = (*from_sys.solution)(from_dof);
325 
326  if (_fixed_meshes)
327  {
328  // Cache the nearest nodes.
329  _cached_froms[i_proc][qp] = i_local_from;
330  _cached_dof_ids[i_proc][qp] = from_dof;
331  }
332  }
333  }
334  }
335  }
336  }
337 
338  if (i_proc == processor_id())
339  incoming_evals[i_proc] = outgoing_evals;
340  else
341  _communicator.send(i_proc, outgoing_evals, send_evals[i_proc]);
342  }
343  }
344 
345  else // We've cached the nearest nodes.
346  {
347  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
348  {
349  std::vector<Real> & outgoing_evals = processor_outgoing_evals[i_proc];
350  outgoing_evals.resize(_cached_froms[i_proc].size());
351 
352  for (unsigned int qp = 0; qp < outgoing_evals.size(); qp++)
353  {
354  MooseVariable & from_var =
355  _from_problems[_cached_froms[i_proc][qp]]->getVariable(0, _from_var_name);
356  System & from_sys = from_var.sys().system();
357  dof_id_type from_dof = _cached_dof_ids[i_proc][qp];
358  // outgoing_evals[qp] = (*from_sys.solution)(_cached_dof_ids[i_proc][qp]);
359  outgoing_evals[qp] = (*from_sys.solution)(from_dof);
360  }
361 
362  if (i_proc == processor_id())
363  incoming_evals[i_proc] = outgoing_evals;
364  else
365  _communicator.send(i_proc, outgoing_evals, send_evals[i_proc]);
366  }
367  }
368 
370  // Gather all of the evaluations, find the nearest one for each node/element,
371  // and apply the values.
373 
374  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
375  {
376  if (i_proc == processor_id())
377  continue;
378 
379  _communicator.receive(i_proc, incoming_evals[i_proc]);
380  }
381 
382  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
383  {
384  // Loop over the master nodes and set the value of the variable
385  System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
386 
387  unsigned int sys_num = to_sys->number();
388  unsigned int var_num = to_sys->variable_number(_to_var_name);
389 
390  NumericVector<Real> * solution = nullptr;
391  switch (_direction)
392  {
393  case TO_MULTIAPP:
394  solution = &getTransferVector(i_to, _to_var_name);
395  break;
396  case FROM_MULTIAPP:
397  solution = to_sys->solution.get();
398  break;
399  default:
400  mooseError("Unknown direction");
401  }
402 
403  MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
404 
405  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
406 
407  if (is_nodal)
408  {
409  std::vector<Node *> target_local_nodes;
410 
411  if (isParamValid("target_boundary"))
412  {
413  BoundaryID target_bnd_id =
414  _to_meshes[i_to]->getBoundaryID(getParam<BoundaryName>("target_boundary"));
415 
416  ConstBndNodeRange & bnd_nodes = *(_to_meshes[i_to])->getBoundaryNodeRange();
417  for (const auto & bnode : bnd_nodes)
418  if (bnode->_bnd_id == target_bnd_id && bnode->_node->processor_id() == processor_id())
419  target_local_nodes.push_back(bnode->_node);
420  }
421  else
422  {
423  target_local_nodes.resize(to_mesh->n_local_nodes());
424  MeshBase::const_node_iterator nodes_begin = to_mesh->local_nodes_begin();
425  MeshBase::const_node_iterator nodes_end = to_mesh->local_nodes_end();
426 
427  unsigned int i = 0;
428  for (MeshBase::const_node_iterator nodes_it = nodes_begin; nodes_it != nodes_end;
429  ++nodes_it, ++i)
430  target_local_nodes[i] = *nodes_it;
431  }
432 
433  for (const auto & node : target_local_nodes)
434  {
435  // Skip this node if the variable has no dofs at it.
436  if (node->n_dofs(sys_num, var_num) < 1)
437  continue;
438 
439  Real best_val = 0;
440  if (!_neighbors_cached)
441  {
442  Real min_dist = std::numeric_limits<Real>::max();
443  for (unsigned int i_from = 0; i_from < incoming_evals.size(); i_from++)
444  {
445  std::pair<unsigned int, unsigned int> key(i_to, node->id());
446  if (node_index_map[i_from].find(key) == node_index_map[i_from].end())
447  continue;
448  unsigned int qp_ind = node_index_map[i_from][key];
449  if (incoming_evals[i_from][2 * qp_ind] >= min_dist)
450  continue;
451  min_dist = incoming_evals[i_from][2 * qp_ind];
452  best_val = incoming_evals[i_from][2 * qp_ind + 1];
453 
454  if (_fixed_meshes)
455  {
456  // Cache these indices.
457  _cached_from_inds[node->id()] = i_from;
458  _cached_qp_inds[node->id()] = qp_ind;
459  }
460  }
461  }
462 
463  else
464  {
465  best_val = incoming_evals[_cached_from_inds[node->id()]][_cached_qp_inds[node->id()]];
466  }
467 
468  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
469  solution->set(dof, best_val);
470  }
471  }
472  else // Elemental
473  {
474  MeshBase::const_element_iterator elem_it = to_mesh->local_elements_begin();
475  MeshBase::const_element_iterator elem_end = to_mesh->local_elements_end();
476 
477  for (; elem_it != elem_end; ++elem_it)
478  {
479  Elem * elem = *elem_it;
480 
481  // Skip this element if the variable has no dofs at it.
482  if (elem->n_dofs(sys_num, var_num) < 1)
483  continue;
484 
485  Real best_val = 0;
486  if (!_neighbors_cached)
487  {
488  Real min_dist = std::numeric_limits<Real>::max();
489  for (unsigned int i_from = 0; i_from < incoming_evals.size(); i_from++)
490  {
491  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
492  if (node_index_map[i_from].find(key) == node_index_map[i_from].end())
493  continue;
494  unsigned int qp_ind = node_index_map[i_from][key];
495  if (incoming_evals[i_from][2 * qp_ind] >= min_dist)
496  continue;
497  min_dist = incoming_evals[i_from][2 * qp_ind];
498  best_val = incoming_evals[i_from][2 * qp_ind + 1];
499 
500  if (_fixed_meshes)
501  {
502  // Cache these indices.
503  _cached_from_inds[elem->id()] = i_from;
504  _cached_qp_inds[elem->id()] = qp_ind;
505  }
506  }
507  }
508 
509  else
510  {
511  best_val = incoming_evals[_cached_from_inds[elem->id()]][_cached_qp_inds[elem->id()]];
512  }
513 
514  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
515  solution->set(dof, best_val);
516  }
517  }
518  solution->close();
519  to_sys->update();
520  }
521 
522  if (_fixed_meshes)
523  _neighbors_cached = true;
524 
525  // Make sure all our sends succeeded.
526  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
527  {
528  if (i_proc == processor_id())
529  continue;
530  send_qps[i_proc].wait();
531  send_evals[i_proc].wait();
532  }
533 
534  _console << "Finished NearestNodeTransfer " << name() << std::endl;
535 }
536 
537 Node *
539  Real & distance,
540  MooseMesh * mesh,
541  bool local)
542 {
543  distance = std::numeric_limits<Real>::max();
544  Node * nearest = NULL;
545 
546  if (isParamValid("source_boundary"))
547  {
548  BoundaryID src_bnd_id = mesh->getBoundaryID(getParam<BoundaryName>("source_boundary"));
549 
550  ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
551  for (const auto & bnode : bnd_nodes)
552  {
553  if (bnode->_bnd_id == src_bnd_id)
554  {
555  Node * node = bnode->_node;
556  Real current_distance = (p - *node).norm();
557 
558  if (current_distance < distance)
559  {
560  distance = current_distance;
561  nearest = node;
562  }
563  }
564  }
565  }
566  else
567  {
568  SimpleRange<MeshBase::const_node_iterator> range(
569  local ? mesh->localNodesBegin() : mesh->getMesh().nodes_begin(),
570  local ? mesh->localNodesEnd() : mesh->getMesh().nodes_end());
571 
572  for (auto & node : range)
573  {
574  Real current_distance = (p - *node).norm();
575 
576  if (current_distance < distance)
577  {
578  distance = current_distance;
579  nearest = node;
580  }
581  }
582  }
583 
584  return nearest;
585 }
586 
587 Real
589 {
590  std::vector<Point> source_points = {bbox.first, bbox.second};
591 
592  std::vector<Point> all_points(8);
593  for (unsigned int x = 0; x < 2; x++)
594  for (unsigned int y = 0; y < 2; y++)
595  for (unsigned int z = 0; z < 2; z++)
596  all_points[x + 2 * y + 4 * z] =
597  Point(source_points[x](0), source_points[y](1), source_points[z](2));
598 
599  Real max_distance = 0.;
600 
601  for (unsigned int i = 0; i < 8; i++)
602  {
603  Real distance = (p - all_points[i]).norm();
604  if (distance > max_distance)
605  max_distance = distance;
606  }
607  return max_distance;
608 }
609 
610 Real
612 {
613  std::vector<Point> source_points = {bbox.first, bbox.second};
614 
615  std::vector<Point> all_points(8);
616  for (unsigned int x = 0; x < 2; x++)
617  for (unsigned int y = 0; y < 2; y++)
618  for (unsigned int z = 0; z < 2; z++)
619  all_points[x + 2 * y + 4 * z] =
620  Point(source_points[x](0), source_points[y](1), source_points[z](2));
621 
622  Real min_distance = 0.;
623 
624  for (unsigned int i = 0; i < 8; i++)
625  {
626  Real distance = (p - all_points[i]).norm();
627  if (distance < min_distance)
628  min_distance = distance;
629  }
630  return min_distance;
631 }
632 
633 void
634 MultiAppNearestNodeTransfer::getLocalNodes(MooseMesh * mesh, std::vector<Node *> & local_nodes)
635 {
636  if (isParamValid("source_boundary"))
637  {
638  BoundaryID src_bnd_id = mesh->getBoundaryID(getParam<BoundaryName>("source_boundary"));
639 
640  ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
641  for (const auto & bnode : bnd_nodes)
642  if (bnode->_bnd_id == src_bnd_id && bnode->_node->processor_id() == processor_id())
643  local_nodes.push_back(bnode->_node);
644  }
645  else
646  {
647  local_nodes.resize(mesh->getMesh().n_local_nodes());
648 
649  MeshBase::const_node_iterator nodes_begin = mesh->localNodesBegin();
650  MeshBase::const_node_iterator nodes_end = mesh->localNodesEnd();
651 
652  unsigned int i = 0;
653  for (MeshBase::const_node_iterator node_it = nodes_begin; node_it != nodes_end; ++node_it, ++i)
654  local_nodes[i] = *node_it;
655  }
656 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
NumericVector< Real > & getTransferVector(unsigned int i_local, std::string var_name)
If we are transferring to a multiapp, return the appropriate solution vector.
void variableIntegrityCheck(const AuxVariableName &var_name) const
Utility to verify that the vEariable in the destination system exists.
Class for stuff related to variables.
Definition: MooseVariable.h:43
MultiAppNearestNodeTransfer(const InputParameters &parameters)
std::vector< EquationSystems * > _to_es
MeshBase::const_node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2053
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
static PetscErrorCode Vec x
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.
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
Definition: MooseMesh.C:968
std::vector< std::vector< dof_id_type > > & _cached_dof_ids
virtual void execute() override
Execute the transfer.
MeshBase::const_node_iterator localNodesEnd()
Definition: MooseMesh.C:2059
std::map< dof_id_type, unsigned int > & _cached_from_inds
std::vector< MooseMesh * > _from_meshes
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2408
MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
virtual void initialSetup() override
Method called at the beginning of the simulation for checking integrity or doing one-time setup...
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
std::vector< unsigned int > getFromsPerProc()
Return the number of "from" domains that each processor owns.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:67
virtual System & system()=0
Get the reference to the libMesh system.
std::map< dof_id_type, unsigned int > & _cached_qp_inds
Real bboxMinDistance(Point p, BoundingBox bbox)
Return the distance between the given point and the nearest corner of the given bounding box...
std::vector< Point > _to_positions
std::vector< std::vector< unsigned int > > & _cached_froms
InputParameters validParams< MultiAppNearestNodeTransfer >()
InputParameters validParams< MultiAppTransfer >()
Base class for all MultiAppTransfer objects.
Node * getNearestNode(const Point &p, Real &distance, MooseMesh *mesh, bool local)
Return the nearest node to the point p.
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 mooseError(Args &&...args) const
Definition: MooseObject.h:80
Real bboxMaxDistance(Point p, BoundingBox bbox)
Return the distance between the given point and the farthest corner of the given bounding box...
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
static System * find_sys(EquationSystems &es, const std::string &var_name)
Small helper function for finding the system containing the variable.
Definition: Transfer.C:70
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > ConstBndNodeRange
Some useful StoredRange typedefs.
Definition: MooseMesh.h:1213
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.
void getLocalNodes(MooseMesh *mesh, std::vector< Node * > &local_nodes)
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
Definition: MooseMesh.C:763
std::vector< FEProblemBase * > _from_problems
bool _fixed_meshes
If true then node connections will be cached.
std::vector< MooseMesh * > _to_meshes
boundary_id_type BoundaryID
Definition: MooseTypes.h:75