www.mooseframework.org
MultiAppNearestNodeTransfer.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 "DisplacedProblem.h"
14 #include "FEProblem.h"
15 #include "MooseMesh.h"
16 #include "MooseTypes.h"
17 #include "MooseVariableFE.h"
18 
19 #include "libmesh/system.h"
20 #include "libmesh/mesh_tools.h"
21 #include "libmesh/id_types.h"
22 #include "libmesh/parallel_algebra.h"
23 #include "libmesh/dof_object.h"
24 
26 
27 template <>
30 {
32 
33  params.addParam<BoundaryName>(
34  "source_boundary",
35  "The boundary we are transferring from (if not specified, whole domain is used).");
36  params.addParam<BoundaryName>(
37  "target_boundary",
38  "The boundary we are transferring to (if not specified, whole domain is used).");
39  params.addParam<bool>("fixed_meshes",
40  false,
41  "Set to true when the meshes are not changing (ie, "
42  "no movement or adaptivity). This will cache "
43  "nearest node neighbors to greatly speed up the "
44  "transfer.");
45 
46  return params;
47 }
48 
50  : MultiAppFieldTransfer(parameters),
51  _fixed_meshes(getParam<bool>("fixed_meshes")),
52  _node_map(declareRestartableData<std::map<dof_id_type, Node *>>("node_map")),
53  _distance_map(declareRestartableData<std::map<dof_id_type, Real>>("distance_map")),
54  _neighbors_cached(declareRestartableData<bool>("neighbors_cached", false)),
55  _cached_froms(declareRestartableData<std::vector<std::vector<unsigned int>>>("cached_froms")),
56  _cached_dof_ids(
57  declareRestartableData<std::vector<std::vector<dof_id_type>>>("cached_dof_ids")),
58  _cached_from_inds(
59  declareRestartableData<std::map<dof_id_type, unsigned int>>("cached_from_ids")),
60  _cached_qp_inds(declareRestartableData<std::map<dof_id_type, unsigned int>>("cached_qp_inds"))
61 {
62  if (_to_var_names.size() != 1 && _from_var_names.size() != 1)
63  mooseError(" Support single variable only");
64 }
65 
66 void
68 {
69  _console << "Beginning NearestNodeTransfer " << name() << std::endl;
70 
71  getAppInfo();
72 
73  // Get the bounding boxes for the "from" domains.
74  std::vector<BoundingBox> bboxes;
75  if (isParamValid("source_boundary"))
76  bboxes = getFromBoundingBoxes(
77  _from_meshes[0]->getBoundaryID(getParam<BoundaryName>("source_boundary")));
78  else
79  bboxes = getFromBoundingBoxes();
80 
81  // Figure out how many "from" domains each processor owns.
82  std::vector<unsigned int> froms_per_proc = getFromsPerProc();
83 
85  // For every point in the local "to" domain, figure out which "from" domains
86  // might contain it's nearest neighbor, and send that point to the processors
87  // that own those "from" domains.
88  //
89  // How do we know which "from" domains might contain the nearest neighbor, you
90  // ask? Well, consider two "from" domains, A and B. If every point in A is
91  // closer than every point in B, then we know that B cannot possibly contain
92  // the nearest neighbor. Hence, we'll only check A for the nearest neighbor.
93  // We'll use the functions bboxMaxDistance and bboxMinDistance to figure out
94  // if every point in A is closer than every point in B.
96 
97  // outgoing_qps = nodes/centroids we'll send to other processors.
98  std::vector<std::vector<Point>> outgoing_qps(n_processors());
99  // When we get results back, node_index_map will tell us which results go with
100  // which points
101  std::vector<std::map<std::pair<unsigned int, unsigned int>, unsigned int>> node_index_map(
102  n_processors());
103 
104  if (!_neighbors_cached)
105  {
106  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
107  {
108  System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
109  unsigned int sys_num = to_sys->number();
110  unsigned int var_num = to_sys->variable_number(_to_var_name);
111  MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
112  bool is_to_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
113 
114  if (is_to_nodal)
115  {
116  std::vector<Node *> target_local_nodes;
117 
118  if (isParamValid("target_boundary"))
119  {
120  BoundaryID target_bnd_id =
121  _to_meshes[i_to]->getBoundaryID(getParam<BoundaryName>("target_boundary"));
122 
123  ConstBndNodeRange & bnd_nodes = *(_to_meshes[i_to])->getBoundaryNodeRange();
124  for (const auto & bnode : bnd_nodes)
125  if (bnode->_bnd_id == target_bnd_id && bnode->_node->processor_id() == processor_id())
126  target_local_nodes.push_back(bnode->_node);
127  }
128  else
129  {
130  target_local_nodes.resize(to_mesh->n_local_nodes());
131  unsigned int i = 0;
132  for (auto & node : to_mesh->local_node_ptr_range())
133  target_local_nodes[i++] = node;
134  }
135 
136  // For error checking: keep track of all target_local_nodes
137  // which are successfully mapped to at least one domain where
138  // the nearest neighbor might be found.
139  std::set<Node *> local_nodes_found;
140 
141  for (const auto & node : target_local_nodes)
142  {
143  // Skip this node if the variable has no dofs at it.
144  if (node->n_dofs(sys_num, var_num) < 1)
145  continue;
146 
147  // Find which bboxes might have the nearest node to this point.
148  Real nearest_max_distance = std::numeric_limits<Real>::max();
149  for (const auto & bbox : bboxes)
150  {
151  Real distance = bboxMaxDistance(*node, bbox);
152  if (distance < nearest_max_distance)
153  nearest_max_distance = distance;
154  }
155 
156  unsigned int from0 = 0;
157  for (processor_id_type i_proc = 0; i_proc < n_processors();
158  from0 += froms_per_proc[i_proc], i_proc++)
159  {
160  bool qp_found = false;
161 
162  for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc] && !qp_found;
163  i_from++)
164  {
165 
166  Real distance = bboxMinDistance(*node, bboxes[i_from]);
167 
168  if (distance <= nearest_max_distance || bboxes[i_from].contains_point(*node))
169  {
170  std::pair<unsigned int, unsigned int> key(i_to, node->id());
171  node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
172  outgoing_qps[i_proc].push_back(*node + _to_positions[i_to]);
173  qp_found = true;
174  local_nodes_found.insert(node);
175  }
176  }
177  }
178  }
179 
180  // By the time we get to here, we should have found at least
181  // one candidate BoundingBox for every node in the
182  // target_local_nodes array that has dofs for the current
183  // variable in the current System.
184  for (const auto & node : target_local_nodes)
185  if (node->n_dofs(sys_num, var_num) && !local_nodes_found.count(node))
186  mooseError("In ",
187  name(),
188  ": No candidate BoundingBoxes found for node ",
189  node->id(),
190  " at position ",
191  *node);
192  }
193  else // Elemental
194  {
195  // For error checking: keep track of all local elements
196  // which are successfully mapped to at least one domain where
197  // the nearest neighbor might be found.
198  std::set<Elem *> local_elems_found;
199 
200  for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
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  local_elems_found.insert(elem);
233  }
234  }
235  }
236  }
237 
238  // Verify that we found at least one candidate bounding
239  // box for each local element with dofs for the current
240  // variable in the current System.
241  for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
242  if (elem->n_dofs(sys_num, var_num) && !local_elems_found.count(elem))
243  mooseError("In ",
244  name(),
245  ": No candidate BoundingBoxes found for Elem ",
246  elem->id(),
247  ", centroid = ",
248  elem->centroid());
249  }
250  }
251  }
252 
254  // Send local node/centroid positions off to the other processors and take
255  // care of points sent to this processor. We'll need to check the points
256  // against all of the "from" domains that this processor owns. For each
257  // point, we'll find the nearest node, then we'll send the value at that node
258  // and the distance between the node and the point back to the processor that
259  // requested that point.
261 
262  std::vector<std::vector<Real>> incoming_evals(n_processors());
263  std::vector<Parallel::Request> send_qps(n_processors());
264  std::vector<Parallel::Request> send_evals(n_processors());
265 
266  // Create these here so that they live the entire life of this function
267  // and are NOT reused per processor.
268  std::vector<std::vector<Real>> processor_outgoing_evals(n_processors());
269 
270  if (!_neighbors_cached)
271  {
272  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
273  {
274  if (i_proc == processor_id())
275  continue;
276  _communicator.send(i_proc, outgoing_qps[i_proc], send_qps[i_proc]);
277  }
278 
279  // Build an array of pointers to all of this processor's local entities (nodes or
280  // elements). We need to do this to avoid the expense of using LibMesh iterators.
281  // This step also takes care of limiting the search to boundary nodes, if
282  // applicable.
283  std::vector<std::vector<std::pair<Point, DofObject *>>> local_entities(
284  froms_per_proc[processor_id()]);
285 
286  // Local array of all from Variable references
287  std::vector<std::reference_wrapper<MooseVariableFEBase>> _from_vars;
288 
289  for (unsigned int i = 0; i < froms_per_proc[processor_id()]; i++)
290  {
291  MooseVariableFEBase & from_var = _from_problems[i]->getVariable(
293  bool is_to_nodal = from_var.feType().family == LAGRANGE;
294 
295  _from_vars.emplace_back(from_var);
296  getLocalEntities(_from_meshes[i], local_entities[i], is_to_nodal);
297  }
298 
299  if (_fixed_meshes)
300  {
301  _cached_froms.resize(n_processors());
302  _cached_dof_ids.resize(n_processors());
303  }
304 
305  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
306  {
307  // We either use our own outgoing_qps or receive them from
308  // another processor.
309  std::vector<Point> incoming_qps;
310  if (i_proc == processor_id())
311  incoming_qps = outgoing_qps[i_proc];
312  else
313  _communicator.receive(i_proc, incoming_qps);
314 
315  if (_fixed_meshes)
316  {
317  _cached_froms[i_proc].resize(incoming_qps.size());
318  _cached_dof_ids[i_proc].resize(incoming_qps.size());
319  }
320 
321  std::vector<Real> & outgoing_evals = processor_outgoing_evals[i_proc];
322  // Resize this vector to two times the size of the incoming_qps
323  // vector because we are going to store both the value from the nearest
324  // local node *and* the distance between the incoming_qp and that node
325  // for later comparison purposes.
326  outgoing_evals.resize(2 * incoming_qps.size());
327 
328  for (unsigned int qp = 0; qp < incoming_qps.size(); qp++)
329  {
330  const Point & qpt = incoming_qps[qp];
331  outgoing_evals[2 * qp] = std::numeric_limits<Real>::max();
332  for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
333  i_local_from++)
334  {
335  MooseVariableFEBase & from_var = _from_vars[i_local_from];
336  System & from_sys = from_var.sys().system();
337  unsigned int from_sys_num = from_sys.number();
338  unsigned int from_var_num = from_sys.variable_number(from_var.name());
339 
340  for (unsigned int i_node = 0; i_node < local_entities[i_local_from].size(); i_node++)
341  {
342  // Compute distance between the current incoming_qp to local node i_node.
343  Real current_distance =
344  (qpt - local_entities[i_local_from][i_node].first - _from_positions[i_local_from])
345  .norm();
346 
347  // If an incoming_qp is equally close to two or more local nodes, then
348  // the first one we test will "win", even though any of the others could
349  // also potentially be chosen instead... there's no way to decide among
350  // the set of all equidistant points.
351  //
352  // outgoing_evals[2 * qp] is the current closest distance between a local point and
353  // the incoming_qp.
354  if (current_distance < outgoing_evals[2 * qp])
355  {
356  // Assuming LAGRANGE!
357  if (local_entities[i_local_from][i_node].second->n_dofs(from_sys_num, from_var_num) >
358  0)
359  {
360  dof_id_type from_dof = local_entities[i_local_from][i_node].second->dof_number(
361  from_sys_num, from_var_num, 0);
362 
363  // The indexing of the outgoing_evals vector looks
364  // like [(distance, value), (distance, value), ...]
365  // for each incoming_qp. We only keep the value from
366  // the node with the smallest distance to the
367  // incoming_qp, and then we compare across all
368  // processors later and pick the closest one.
369  outgoing_evals[2 * qp] = current_distance;
370  outgoing_evals[2 * qp + 1] = (*from_sys.solution)(from_dof);
371 
372  if (_fixed_meshes)
373  {
374  // Cache the nearest nodes.
375  _cached_froms[i_proc][qp] = i_local_from;
376  _cached_dof_ids[i_proc][qp] = from_dof;
377  }
378  }
379  }
380  }
381  }
382  }
383 
384  if (i_proc == processor_id())
385  incoming_evals[i_proc] = outgoing_evals;
386  else
387  _communicator.send(i_proc, outgoing_evals, send_evals[i_proc]);
388  }
389  }
390 
391  else // We've cached the nearest nodes.
392  {
393  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
394  {
395  std::vector<Real> & outgoing_evals = processor_outgoing_evals[i_proc];
396  outgoing_evals.resize(_cached_froms[i_proc].size());
397 
398  for (unsigned int qp = 0; qp < outgoing_evals.size(); qp++)
399  {
400  MooseVariableFEBase & from_var = _from_problems[_cached_froms[i_proc][qp]]->getVariable(
401  0,
405  System & from_sys = from_var.sys().system();
406  dof_id_type from_dof = _cached_dof_ids[i_proc][qp];
407  // outgoing_evals[qp] = (*from_sys.solution)(_cached_dof_ids[i_proc][qp]);
408  outgoing_evals[qp] = (*from_sys.solution)(from_dof);
409  }
410 
411  if (i_proc == processor_id())
412  incoming_evals[i_proc] = outgoing_evals;
413  else
414  _communicator.send(i_proc, outgoing_evals, send_evals[i_proc]);
415  }
416  }
417 
419  // Gather all of the evaluations, find the nearest one for each node/element,
420  // and apply the values.
422 
423  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
424  {
425  if (i_proc == processor_id())
426  continue;
427 
428  _communicator.receive(i_proc, incoming_evals[i_proc]);
429  }
430 
431  for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
432  {
433  // Loop over the master nodes and set the value of the variable
434  System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
435 
436  unsigned int sys_num = to_sys->number();
437  unsigned int var_num = to_sys->variable_number(_to_var_name);
438 
439  NumericVector<Real> * solution = nullptr;
440  switch (_direction)
441  {
442  case TO_MULTIAPP:
443  solution = &getTransferVector(i_to, _to_var_name);
444  break;
445  case FROM_MULTIAPP:
446  solution = to_sys->solution.get();
447  break;
448  default:
449  mooseError("Unknown direction");
450  }
451 
452  const MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
453 
454  bool is_to_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
455 
456  if (is_to_nodal)
457  {
458  std::vector<Node *> target_local_nodes;
459 
460  if (isParamValid("target_boundary"))
461  {
462  BoundaryID target_bnd_id =
463  _to_meshes[i_to]->getBoundaryID(getParam<BoundaryName>("target_boundary"));
464 
465  ConstBndNodeRange & bnd_nodes = *(_to_meshes[i_to])->getBoundaryNodeRange();
466  for (const auto & bnode : bnd_nodes)
467  if (bnode->_bnd_id == target_bnd_id && bnode->_node->processor_id() == processor_id())
468  target_local_nodes.push_back(bnode->_node);
469  }
470  else
471  {
472  target_local_nodes.resize(to_mesh.n_local_nodes());
473  unsigned int i = 0;
474  for (auto & node : to_mesh.local_node_ptr_range())
475  target_local_nodes[i++] = node;
476  }
477 
478  for (const auto & node : target_local_nodes)
479  {
480  // Skip this node if the variable has no dofs at it.
481  if (node->n_dofs(sys_num, var_num) < 1)
482  continue;
483 
484  // If nothing is in the node_index_map for a given local node,
485  // it will get the value 0.
486  Real best_val = 0;
487  if (!_neighbors_cached)
488  {
489  // Search through all the incoming evaluation points from
490  // different processors for the one with the closest
491  // point. If there are multiple values from other processors
492  // which are equidistant, the first one we check will "win".
493  Real min_dist = std::numeric_limits<Real>::max();
494  for (unsigned int i_from = 0; i_from < incoming_evals.size(); i_from++)
495  {
496  std::pair<unsigned int, unsigned int> key(i_to, node->id());
497  if (node_index_map[i_from].find(key) == node_index_map[i_from].end())
498  continue;
499  unsigned int qp_ind = node_index_map[i_from][key];
500  if (incoming_evals[i_from][2 * qp_ind] >= min_dist)
501  continue;
502 
503  // If we made it here, we are going set a new value and
504  // distance because we found one that was closer.
505  min_dist = incoming_evals[i_from][2 * qp_ind];
506  best_val = incoming_evals[i_from][2 * qp_ind + 1];
507 
508  if (_fixed_meshes)
509  {
510  // Cache these indices.
511  _cached_from_inds[node->id()] = i_from;
512  _cached_qp_inds[node->id()] = qp_ind;
513  }
514  }
515  }
516 
517  else
518  {
519  best_val = incoming_evals[_cached_from_inds[node->id()]][_cached_qp_inds[node->id()]];
520  }
521 
522  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
523  solution->set(dof, best_val);
524  }
525  }
526  else // Elemental
527  {
528  for (auto & elem : to_mesh.active_local_element_ptr_range())
529  {
530  // Skip this element if the variable has no dofs at it.
531  if (elem->n_dofs(sys_num, var_num) < 1)
532  continue;
533 
534  Real best_val = 0;
535  if (!_neighbors_cached)
536  {
537  Real min_dist = std::numeric_limits<Real>::max();
538  for (unsigned int i_from = 0; i_from < incoming_evals.size(); i_from++)
539  {
540  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
541  if (node_index_map[i_from].find(key) == node_index_map[i_from].end())
542  continue;
543  unsigned int qp_ind = node_index_map[i_from][key];
544  if (incoming_evals[i_from][2 * qp_ind] >= min_dist)
545  continue;
546  min_dist = incoming_evals[i_from][2 * qp_ind];
547  best_val = incoming_evals[i_from][2 * qp_ind + 1];
548 
549  if (_fixed_meshes)
550  {
551  // Cache these indices.
552  _cached_from_inds[elem->id()] = i_from;
553  _cached_qp_inds[elem->id()] = qp_ind;
554  }
555  }
556  }
557 
558  else
559  {
560  best_val = incoming_evals[_cached_from_inds[elem->id()]][_cached_qp_inds[elem->id()]];
561  }
562 
563  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
564  solution->set(dof, best_val);
565  }
566  }
567  solution->close();
568  to_sys->update();
569  }
570 
571  if (_fixed_meshes)
572  _neighbors_cached = true;
573 
574  // Make sure all our sends succeeded.
575  for (processor_id_type i_proc = 0; i_proc < n_processors(); i_proc++)
576  {
577  if (i_proc == processor_id())
578  continue;
579  send_qps[i_proc].wait();
580  send_evals[i_proc].wait();
581  }
582 
583  _console << "Finished NearestNodeTransfer " << name() << std::endl;
584 
585  postExecute();
586 }
587 
588 Node *
590  Real & distance,
591  MooseMesh * mesh,
592  bool local)
593 {
594  distance = std::numeric_limits<Real>::max();
595  Node * nearest = nullptr;
596 
597  if (isParamValid("source_boundary"))
598  {
599  BoundaryID src_bnd_id = mesh->getBoundaryID(getParam<BoundaryName>("source_boundary"));
600 
601  const ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
602  for (const auto & bnode : bnd_nodes)
603  {
604  if (bnode->_bnd_id == src_bnd_id)
605  {
606  Node * node = bnode->_node;
607  Real current_distance = (p - *node).norm();
608 
609  if (current_distance < distance)
610  {
611  distance = current_distance;
612  nearest = node;
613  }
614  }
615  }
616  }
617  else
618  {
619  for (auto & node : as_range(local ? mesh->localNodesBegin() : mesh->getMesh().nodes_begin(),
620  local ? mesh->localNodesEnd() : mesh->getMesh().nodes_end()))
621  {
622  Real current_distance = (p - *node).norm();
623 
624  if (current_distance < distance)
625  {
626  distance = current_distance;
627  nearest = node;
628  }
629  }
630  }
631 
632  return nearest;
633 }
634 
635 Real
636 MultiAppNearestNodeTransfer::bboxMaxDistance(const Point & p, const BoundingBox & bbox)
637 {
638  std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
639 
640  std::array<Point, 8> all_points;
641  for (unsigned int x = 0; x < 2; x++)
642  for (unsigned int y = 0; y < 2; y++)
643  for (unsigned int z = 0; z < 2; z++)
644  all_points[x + 2 * y + 4 * z] =
645  Point(source_points[x](0), source_points[y](1), source_points[z](2));
646 
647  Real max_distance = 0.;
648 
649  for (unsigned int i = 0; i < 8; i++)
650  {
651  Real distance = (p - all_points[i]).norm();
652  if (distance > max_distance)
653  max_distance = distance;
654  }
655 
656  return max_distance;
657 }
658 
659 Real
660 MultiAppNearestNodeTransfer::bboxMinDistance(const Point & p, const BoundingBox & bbox)
661 {
662  std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
663 
664  std::array<Point, 8> all_points;
665  for (unsigned int x = 0; x < 2; x++)
666  for (unsigned int y = 0; y < 2; y++)
667  for (unsigned int z = 0; z < 2; z++)
668  all_points[x + 2 * y + 4 * z] =
669  Point(source_points[x](0), source_points[y](1), source_points[z](2));
670 
671  Real min_distance = std::numeric_limits<Real>::max();
672 
673  for (unsigned int i = 0; i < 8; i++)
674  {
675  Real distance = (p - all_points[i]).norm();
676  if (distance < min_distance)
677  min_distance = distance;
678  }
679 
680  return min_distance;
681 }
682 
683 void
685  MooseMesh * mesh, std::vector<std::pair<Point, DofObject *>> & local_entities, bool is_nodal)
686 {
687  mooseAssert(local_entities.empty(), "local_entities should be empty");
688  const MeshBase & mesh_base = mesh->getMesh();
689 
690  if (isParamValid("source_boundary"))
691  {
692  BoundaryID src_bnd_id = mesh->getBoundaryID(getParam<BoundaryName>("source_boundary"));
693  auto proc_id = processor_id();
694  if (is_nodal)
695  {
696  const ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
697  for (const auto & bnode : bnd_nodes)
698  if (bnode->_bnd_id == src_bnd_id && bnode->_node->processor_id() == proc_id)
699  local_entities.emplace_back(*bnode->_node, bnode->_node);
700  }
701  else
702  {
703  const ConstBndElemRange & bnd_elems = *mesh->getBoundaryElementRange();
704  for (const auto & belem : bnd_elems)
705  if (belem->_bnd_id == src_bnd_id && belem->_elem->processor_id() == proc_id)
706  local_entities.emplace_back(belem->_elem->centroid(), belem->_elem);
707  }
708  }
709  else
710  {
711  if (is_nodal)
712  {
713  local_entities.reserve(mesh_base.n_local_nodes());
714  for (auto & node : mesh_base.local_node_ptr_range())
715  local_entities.emplace_back(*node, node);
716  }
717  else
718  {
719  local_entities.reserve(mesh_base.n_local_elem());
720  for (auto & elem : mesh_base.active_local_element_ptr_range())
721  local_entities.emplace_back(elem->centroid(), elem);
722  }
723  }
724 }
NumericVector< Real > & getTransferVector(unsigned int i_local, std::string var_name)
If we are transferring to a multiapp, return the appropriate solution vector.
MultiAppNearestNodeTransfer(const InputParameters &parameters)
Real bboxMaxDistance(const Point &p, const BoundingBox &bbox)
Return the distance between the given point and the farthest corner of the given bounding box...
virtual void postExecute()
Add some extra work if necessary after execute().
StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > ConstBndElemRange
Definition: MooseMesh.h:1259
std::vector< EquationSystems * > _to_es
MeshBase::const_node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2205
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
Real bboxMinDistance(const Point &p, const BoundingBox &bbox)
Return the distance between the given point and the nearest corner of the given bounding box...
static PetscErrorCode Vec x
const FEType & feType() const
Get the type of finite element object.
Copy the value to the target domain from the nearest node in the source domain.
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:2211
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:2567
registerMooseObject("MooseApp", MultiAppNearestNodeTransfer)
This serves an interface for MultiAppInterpolationTransfer, MultiAppNearestNodeTransfer and so on...
boundary_id_type BoundaryID
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.
virtual System & system()=0
Get the reference to the libMesh system.
std::map< dof_id_type, unsigned int > & _cached_qp_inds
std::vector< Point > _to_positions
std::vector< std::vector< unsigned int > > & _cached_froms
InputParameters validParams< MultiAppNearestNodeTransfer >()
const std::string & name() const
Get the variable name.
Node * getNearestNode(const Point &p, Real &distance, MooseMesh *mesh, bool local)
Return the nearest node to the point p.
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...
StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > * getBoundaryElementRange()
Definition: MooseMesh.C:801
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:65
InputParameters validParams< MultiAppFieldTransfer >()
void getLocalEntities(MooseMesh *mesh, std::vector< std::pair< Point, DofObject *>> &local_entities, bool nodal)
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > ConstBndNodeRange
Some useful StoredRange typedefs.
Definition: MooseMesh.h:1258
const MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
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.
std::vector< VariableName > _from_var_names
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
Definition: MooseMesh.C:788
std::vector< FEProblemBase * > _from_problems
bool _fixed_meshes
If true then node connections will be cached.
std::vector< MooseMesh * > _to_meshes
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
Definition: MooseMesh.C:1007