Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
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 :
10 : #include "MultiAppNearestNodeTransfer.h"
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 : #include "MooseAppCoordTransform.h"
19 :
20 : #include "libmesh/system.h"
21 : #include "libmesh/mesh_tools.h"
22 : #include "libmesh/id_types.h"
23 : #include "libmesh/parallel_algebra.h"
24 : #include "libmesh/dof_object.h"
25 :
26 : // TIMPI includes
27 : #include "timpi/parallel_sync.h"
28 :
29 : registerMooseObjectDeprecated("MooseApp", MultiAppNearestNodeTransfer, "12/31/2024 24:00");
30 :
31 : InputParameters
32 15529 : MultiAppNearestNodeTransfer::validParams()
33 : {
34 15529 : InputParameters params = MultiAppConservativeTransfer::validParams();
35 15529 : params.addClassDescription(
36 : "Transfer the value to the target domain from the nearest node in the source domain.");
37 :
38 15529 : params.addParam<BoundaryName>(
39 : "source_boundary",
40 : "The boundary we are transferring from (if not specified, whole domain is used).");
41 15529 : params.addParam<std::vector<BoundaryName>>(
42 : "target_boundary",
43 : "The boundary we are transferring to (if not specified, whole domain is used).");
44 46587 : params.addParam<bool>("fixed_meshes",
45 31058 : false,
46 : "Set to true when the meshes are not changing (ie, "
47 : "no movement or adaptivity). This will cache "
48 : "nearest node neighbors to greatly speed up the "
49 : "transfer.");
50 :
51 15529 : MultiAppTransfer::addBBoxFactorParam(params);
52 15529 : return params;
53 0 : }
54 :
55 628 : MultiAppNearestNodeTransfer::MultiAppNearestNodeTransfer(const InputParameters & parameters)
56 : : MultiAppConservativeTransfer(parameters),
57 628 : _fixed_meshes(getParam<bool>("fixed_meshes")),
58 628 : _node_map(declareRestartableData<std::map<dof_id_type, Node *>>("node_map")),
59 628 : _distance_map(declareRestartableData<std::map<dof_id_type, Real>>("distance_map")),
60 628 : _neighbors_cached(declareRestartableData<bool>("neighbors_cached", false)),
61 628 : _cached_froms(declareRestartableData<std::map<processor_id_type, std::vector<unsigned int>>>(
62 : "cached_froms")),
63 628 : _cached_dof_ids(declareRestartableData<std::map<processor_id_type, std::vector<dof_id_type>>>(
64 : "cached_dof_ids")),
65 628 : _cached_from_inds(
66 628 : declareRestartableData<std::map<std::pair<unsigned int, dof_id_type>, unsigned int>>(
67 : "cached_from_ids")),
68 628 : _cached_qp_inds(
69 628 : declareRestartableData<std::map<std::pair<unsigned int, dof_id_type>, unsigned int>>(
70 1256 : "cached_qp_inds"))
71 : {
72 628 : mooseDeprecated("MultiAppNearestNodeTransfer is deprecated. Use "
73 : "MultiAppGeneralFieldNearestNodeTransfer instead and adapt the parameters");
74 :
75 628 : if (_to_var_names.size() != 1)
76 0 : paramError("variable", " Support single to-variable only");
77 :
78 628 : if (_from_var_names.size() != 1)
79 0 : paramError("source_variable", " Support single from-variable only");
80 628 : }
81 :
82 : void
83 1308 : MultiAppNearestNodeTransfer::execute()
84 : {
85 1308 : TIME_SECTION(
86 : "MultiAppNearestNodeTransfer::execute()", 5, "Transferring variables based on nearest nodes");
87 :
88 : // Get the bounding boxes for the "from" domains.
89 1308 : std::vector<BoundingBox> bboxes;
90 1308 : if (isParamValid("source_boundary"))
91 : {
92 93 : if (_from_meshes.size())
93 : {
94 93 : const auto & sb = getParam<BoundaryName>("source_boundary");
95 93 : if (!MooseMeshUtils::hasBoundaryName(_from_meshes[0]->getMesh(), sb))
96 8 : paramError("source_boundary", "The boundary '", sb, "' was not found in the mesh");
97 :
98 85 : bboxes = getFromBoundingBoxes(_from_meshes[0]->getBoundaryID(sb));
99 : }
100 : else
101 0 : bboxes = getFromBoundingBoxes(Moose::INVALID_BOUNDARY_ID);
102 : }
103 : else
104 1215 : bboxes = getFromBoundingBoxes();
105 :
106 : // Figure out how many "from" domains each processor owns.
107 1300 : std::vector<unsigned int> froms_per_proc = getFromsPerProc();
108 :
109 : ////////////////////
110 : // For every point in the local "to" domain, figure out which "from" domains
111 : // might contain its nearest neighbor, and send that point to the processors
112 : // that own those "from" domains.
113 : //
114 : // How do we know which "from" domains might contain the nearest neighbor, you
115 : // ask? Well, consider two "from" domains, A and B. If every point in A is
116 : // closer than every point in B, then we know that B cannot possibly contain
117 : // the nearest neighbor. Hence, we'll only check A for the nearest neighbor.
118 : // We'll use the functions bboxMaxDistance and bboxMinDistance to figure out
119 : // if every point in A is closer than every point in B.
120 : ////////////////////
121 :
122 : // outgoing_qps = nodes/centroids we'll send to other processors.
123 : // Map processor to quadrature points. We will send these points to remote processors
124 1300 : std::map<processor_id_type, std::vector<Point>> outgoing_qps;
125 : // When we get results back, node_index_map will tell us which results go with
126 : // which points
127 : // <processor, <system_id, node_i>> --> point_id
128 : std::map<processor_id_type, std::map<std::pair<unsigned int, dof_id_type>, dof_id_type>>
129 1300 : node_index_map;
130 :
131 1300 : if (!_neighbors_cached)
132 : {
133 2811 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
134 : {
135 1596 : System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
136 1596 : unsigned int sys_num = to_sys->number();
137 1596 : unsigned int var_num = to_sys->variable_number(_to_var_name);
138 1596 : MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
139 : const auto to_global_num =
140 1596 : _current_direction == FROM_MULTIAPP ? 0 : _to_local2global_map[i_to];
141 1596 : const auto & to_transform = *_to_transforms[to_global_num];
142 1596 : auto & fe_type = to_sys->variable_type(var_num);
143 1596 : bool is_constant = fe_type.order == CONSTANT;
144 1596 : bool is_to_nodal = fe_type.family == LAGRANGE;
145 :
146 : // We support L2_LAGRANGE elemental variable with the first order
147 1596 : if (fe_type.order > FIRST && !is_to_nodal)
148 0 : mooseError("We don't currently support second order or higher elemental variable ");
149 :
150 1596 : if (!is_to_nodal && isParamValid("target_boundary"))
151 0 : paramWarning("target_boundary",
152 : "Setting a target boundary is only valid for receiving "
153 : "variables of the LAGRANGE basis");
154 :
155 1596 : if (is_to_nodal)
156 : {
157 1059 : const std::vector<Node *> & target_local_nodes = getTargetLocalNodes(i_to);
158 :
159 : // For error checking: keep track of all target_local_nodes
160 : // which are successfully mapped to at least one domain where
161 : // the nearest neighbor might be found.
162 1051 : std::set<Node *> local_nodes_found;
163 :
164 174853 : for (const auto & node : target_local_nodes)
165 : {
166 : // Skip this node if the variable has no dofs at it.
167 173802 : if (node->n_dofs(sys_num, var_num) < 1)
168 7040 : continue;
169 :
170 166762 : const auto transformed_node = to_transform(*node);
171 :
172 : // Find which bboxes might have the nearest node to this point.
173 166762 : Real nearest_max_distance = std::numeric_limits<Real>::max();
174 470950 : for (const auto & bbox : bboxes)
175 : {
176 304188 : Real distance = bboxMaxDistance(transformed_node, bbox);
177 304188 : if (distance < nearest_max_distance)
178 226756 : nearest_max_distance = distance;
179 : }
180 :
181 166762 : unsigned int from0 = 0;
182 448426 : for (processor_id_type i_proc = 0; i_proc < n_processors();
183 281664 : from0 += froms_per_proc[i_proc], i_proc++)
184 : {
185 585852 : for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc]; i_from++)
186 : {
187 304188 : Real distance = bboxMinDistance(transformed_node, bboxes[i_from]);
188 :
189 315085 : if (distance <= nearest_max_distance ||
190 10897 : bboxes[i_from].contains_point(transformed_node))
191 : {
192 293291 : std::pair<unsigned int, dof_id_type> key(i_to, node->id());
193 : // Record a local ID for each quadrature point
194 293291 : node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
195 293291 : outgoing_qps[i_proc].push_back(transformed_node);
196 293291 : local_nodes_found.insert(node);
197 : }
198 : }
199 : }
200 : }
201 :
202 : // By the time we get to here, we should have found at least
203 : // one candidate BoundingBox for every node in the
204 : // target_local_nodes array that has dofs for the current
205 : // variable in the current System.
206 174853 : for (const auto & node : target_local_nodes)
207 173802 : if (node->n_dofs(sys_num, var_num) && !local_nodes_found.count(node))
208 0 : mooseError("In ",
209 0 : name(),
210 : ": No candidate BoundingBoxes found for node ",
211 0 : node->id(),
212 : " at position ",
213 0 : to_transform(*node));
214 1051 : }
215 : else // Elemental
216 : {
217 : // For error checking: keep track of all local elements
218 : // which are successfully mapped to at least one domain where
219 : // the nearest neighbor might be found.
220 537 : std::set<Elem *> local_elems_found;
221 537 : std::vector<Point> points;
222 537 : std::vector<dof_id_type> point_ids;
223 28935 : for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
224 : {
225 : // Skip this element if the variable has no dofs at it.
226 14199 : if (elem->n_dofs(sys_num, var_num) < 1)
227 0 : continue;
228 :
229 14199 : points.clear();
230 14199 : point_ids.clear();
231 : // For constant monomial, we take the centroid of element
232 14199 : if (is_constant)
233 : {
234 14199 : points.push_back(to_transform(elem->vertex_average()));
235 14199 : point_ids.push_back(elem->id());
236 : }
237 :
238 : // For L2_LAGRANGE, we take all the nodes of element
239 : else
240 0 : for (auto & node : elem->node_ref_range())
241 : {
242 0 : points.push_back(to_transform(node));
243 0 : point_ids.push_back(node.id());
244 : }
245 :
246 14199 : unsigned int offset = 0;
247 28398 : for (auto & point : points)
248 : {
249 : // Find which bboxes might have the nearest node to this point.
250 14199 : Real nearest_max_distance = std::numeric_limits<Real>::max();
251 42387 : for (const auto & bbox : bboxes)
252 : {
253 28188 : Real distance = bboxMaxDistance(point, bbox);
254 28188 : if (distance < nearest_max_distance)
255 18295 : nearest_max_distance = distance;
256 : }
257 :
258 14199 : unsigned int from0 = 0;
259 33737 : for (processor_id_type i_proc = 0; i_proc < n_processors();
260 19538 : from0 += froms_per_proc[i_proc], i_proc++)
261 : {
262 47726 : for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc]; i_from++)
263 : {
264 28188 : Real distance = bboxMinDistance(point, bboxes[i_from]);
265 28188 : if (distance <= nearest_max_distance || bboxes[i_from].contains_point(point))
266 : {
267 : std::pair<unsigned int, dof_id_type> key(
268 : i_to,
269 25260 : point_ids[offset]); // Create an unique ID
270 : // If this point already exist, we skip it
271 25260 : if (node_index_map[i_proc].find(key) != node_index_map[i_proc].end())
272 6052 : continue;
273 19208 : node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
274 19208 : outgoing_qps[i_proc].push_back(point);
275 19208 : local_elems_found.insert(elem);
276 : } // if distance
277 : } // for i_from
278 : } // for i_proc
279 14199 : offset++;
280 : } // point
281 537 : } // for elem
282 :
283 : // Verify that we found at least one candidate bounding
284 : // box for each local element with dofs for the current
285 : // variable in the current System.
286 28935 : for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
287 14199 : if (elem->n_dofs(sys_num, var_num) && !local_elems_found.count(elem))
288 0 : mooseError("In ",
289 0 : name(),
290 : ": No candidate BoundingBoxes found for Elem ",
291 0 : elem->id(),
292 : ", centroid = ",
293 537 : to_transform(elem->vertex_average()));
294 537 : }
295 : }
296 : }
297 :
298 : ////////////////////
299 : // Send local node/centroid positions off to the other processors and take
300 : // care of points sent to this processor. We'll need to check the points
301 : // against all of the "from" domains that this processor owns. For each
302 : // point, we'll find the nearest node, then we'll send the value at that node
303 : // and the distance between the node and the point back to the processor that
304 : // requested that point.
305 : ////////////////////
306 :
307 1292 : std::map<processor_id_type, std::vector<Real>> incoming_evals;
308 :
309 : // Create these here so that they live the entire life of this function
310 : // and are NOT reused per processor.
311 1292 : std::map<processor_id_type, std::vector<Real>> processor_outgoing_evals;
312 :
313 1292 : if (!_neighbors_cached)
314 : {
315 : // Build an array of pointers to all of this processor's local entities (nodes or
316 : // elements). We need to do this to avoid the expense of using LibMesh iterators.
317 : // This step also takes care of limiting the search to boundary nodes, if
318 : // applicable.
319 : std::vector<std::vector<std::pair<Point, DofObject *>>> local_entities(
320 1215 : froms_per_proc[processor_id()]);
321 :
322 1215 : std::vector<std::vector<unsigned int>> local_comps(froms_per_proc[processor_id()]);
323 :
324 : // Local array of all from Variable references
325 1215 : std::vector<std::reference_wrapper<MooseVariableFEBase>> _from_vars;
326 :
327 2868 : for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
328 : i_local_from++)
329 : {
330 1653 : MooseVariableFEBase & from_var = _from_problems[i_local_from]->getVariable(
331 : 0, _from_var_name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD);
332 1653 : auto & from_fe_type = from_var.feType();
333 1653 : bool is_constant = from_fe_type.order == CONSTANT;
334 1653 : bool is_to_nodal = from_fe_type.family == LAGRANGE;
335 :
336 : // We support L2_LAGRANGE elemental variable with the first order
337 1653 : if (from_fe_type.order > FIRST && !is_to_nodal)
338 0 : mooseError("We don't currently support second order or higher elemental variable ");
339 :
340 1653 : _from_vars.emplace_back(from_var);
341 1653 : getLocalEntitiesAndComponents(_from_meshes[i_local_from],
342 1653 : local_entities[i_local_from],
343 1653 : local_comps[i_local_from],
344 : is_to_nodal,
345 : is_constant);
346 : }
347 :
348 : // Quadrature points I will receive from remote processors
349 1215 : std::map<processor_id_type, std::vector<Point>> incoming_qps;
350 1860 : auto qps_action_functor = [&incoming_qps](processor_id_type pid, const std::vector<Point> & qps)
351 : {
352 : // Quadrature points from processor 'pid'
353 1860 : auto & incoming_qps_from_pid = incoming_qps[pid];
354 : // Store data for late use
355 1860 : incoming_qps_from_pid.reserve(incoming_qps_from_pid.size() + qps.size());
356 1860 : std::copy(qps.begin(), qps.end(), std::back_inserter(incoming_qps_from_pid));
357 3075 : };
358 :
359 1215 : Parallel::push_parallel_vector_data(comm(), outgoing_qps, qps_action_functor);
360 :
361 3075 : for (auto & qps : incoming_qps)
362 : {
363 1860 : const processor_id_type pid = qps.first;
364 :
365 1860 : if (_fixed_meshes)
366 : {
367 51 : auto & froms = _cached_froms[pid];
368 51 : froms.resize(qps.second.size());
369 51 : std::fill(froms.begin(), froms.end(), libMesh::invalid_uint);
370 :
371 51 : auto & dof_ids = _cached_dof_ids[pid];
372 51 : dof_ids.resize(qps.second.size());
373 51 : std::fill(dof_ids.begin(), dof_ids.end(), DofObject::invalid_id);
374 : }
375 :
376 1860 : std::vector<Real> & outgoing_evals = processor_outgoing_evals[pid];
377 : // Resize this vector to two times the size of the incoming_qps
378 : // vector because we are going to store both the value from the nearest
379 : // local node *and* the distance between the incoming_qp and that node
380 : // for later comparison purposes.
381 1860 : outgoing_evals.resize(2 * qps.second.size());
382 :
383 314359 : for (std::size_t qp = 0; qp < qps.second.size(); qp++)
384 : {
385 312499 : const Point & qpt = qps.second[qp];
386 312499 : outgoing_evals[2 * qp] = std::numeric_limits<Real>::max();
387 768575 : for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
388 : i_local_from++)
389 : {
390 456076 : MooseVariableFEBase & from_var = _from_vars[i_local_from];
391 456076 : System & from_sys = from_var.sys().system();
392 456076 : unsigned int from_sys_num = from_sys.number();
393 456076 : unsigned int from_var_num = from_sys.variable_number(from_var.name());
394 : const auto from_global_num =
395 456076 : _current_direction == TO_MULTIAPP ? 0 : _from_local2global_map[i_local_from];
396 456076 : const auto & from_transform = *_from_transforms[from_global_num];
397 :
398 32138490 : for (unsigned int i_node = 0; i_node < local_entities[i_local_from].size(); i_node++)
399 : {
400 : // Compute distance between the current incoming_qp to local node i_node. No need to
401 : // transform the 'to' because we already did it
402 : Real current_distance =
403 31682414 : (qpt - from_transform(local_entities[i_local_from][i_node].first)).norm();
404 :
405 : // If an incoming_qp is equally close to two or more local nodes, then
406 : // the first one we test will "win", even though any of the others could
407 : // also potentially be chosen instead... there's no way to decide among
408 : // the set of all equidistant points.
409 : //
410 : // outgoing_evals[2 * qp] is the current closest distance between a local point and
411 : // the incoming_qp.
412 31682414 : if (current_distance < outgoing_evals[2 * qp])
413 : {
414 : // Assuming LAGRANGE!
415 9696267 : if (local_entities[i_local_from][i_node].second->n_dofs(from_sys_num, from_var_num) >
416 : 0)
417 : {
418 19342024 : dof_id_type from_dof = local_entities[i_local_from][i_node].second->dof_number(
419 9671012 : from_sys_num, from_var_num, local_comps[i_local_from][i_node]);
420 :
421 : // The indexing of the outgoing_evals vector looks
422 : // like [(distance, value), (distance, value), ...]
423 : // for each incoming_qp. We only keep the value from
424 : // the node with the smallest distance to the
425 : // incoming_qp, and then we compare across all
426 : // processors later and pick the closest one.
427 9671012 : outgoing_evals[2 * qp] = current_distance;
428 9671012 : outgoing_evals[2 * qp + 1] = (*from_sys.solution)(from_dof);
429 :
430 9671012 : if (_fixed_meshes)
431 : {
432 : // Cache the nearest nodes.
433 50000 : _cached_froms[pid][qp] = i_local_from;
434 50000 : _cached_dof_ids[pid][qp] = from_dof;
435 : }
436 : }
437 : }
438 : }
439 : }
440 : }
441 : }
442 1215 : }
443 :
444 : else // We've cached the nearest nodes.
445 : {
446 196 : for (auto & problem_from : _cached_froms)
447 : {
448 119 : const processor_id_type pid = problem_from.first;
449 119 : std::vector<Real> & outgoing_evals = processor_outgoing_evals[pid];
450 119 : outgoing_evals.resize(problem_from.second.size());
451 :
452 8864 : for (unsigned int qp = 0; qp < outgoing_evals.size(); qp++)
453 : {
454 8745 : const auto from_problem = problem_from.second[qp];
455 8745 : if (from_problem == libMesh::invalid_uint)
456 : {
457 : mooseAssert(_cached_dof_ids[pid][qp] == DofObject::invalid_id,
458 : "The state of the from problem and dof id should match.");
459 0 : continue;
460 : }
461 :
462 : MooseVariableFEBase & from_var =
463 8745 : _from_problems[from_problem]->getVariable(0,
464 : _from_var_name,
465 : Moose::VarKindType::VAR_ANY,
466 : Moose::VarFieldType::VAR_FIELD_STANDARD);
467 8745 : System & from_sys = from_var.sys().system();
468 8745 : dof_id_type from_dof = _cached_dof_ids[pid][qp];
469 8745 : outgoing_evals[qp] = (*from_sys.solution)(from_dof);
470 : }
471 : }
472 : }
473 :
474 : auto evals_action_functor =
475 1979 : [&incoming_evals](processor_id_type pid, const std::vector<Real> & evals)
476 : {
477 : // evals for processor 'pid'
478 1979 : auto & incoming_evals_for_pid = incoming_evals[pid];
479 : // Copy evals for late use
480 1979 : incoming_evals_for_pid.reserve(incoming_evals_for_pid.size() + evals.size());
481 1979 : std::copy(evals.begin(), evals.end(), std::back_inserter(incoming_evals_for_pid));
482 3271 : };
483 :
484 1292 : Parallel::push_parallel_vector_data(comm(), processor_outgoing_evals, evals_action_functor);
485 :
486 : ////////////////////
487 : // Gather all of the evaluations, find the nearest one for each node/element,
488 : // and apply the values.
489 : ////////////////////
490 :
491 3042 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
492 : {
493 : // Loop over the master nodes and set the value of the variable
494 1750 : System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
495 :
496 1750 : unsigned int sys_num = to_sys->number();
497 1750 : unsigned int var_num = to_sys->variable_number(_to_var_name);
498 :
499 1750 : NumericVector<Real> * solution = nullptr;
500 1750 : switch (_current_direction)
501 : {
502 1020 : case TO_MULTIAPP:
503 1020 : solution = &getTransferVector(i_to, _to_var_name);
504 1020 : break;
505 730 : case FROM_MULTIAPP:
506 730 : solution = to_sys->solution.get();
507 730 : break;
508 0 : default:
509 0 : mooseError("Unknown direction");
510 : }
511 :
512 1750 : const MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
513 :
514 1750 : auto & fe_type = to_sys->variable_type(var_num);
515 1750 : bool is_constant = fe_type.order == CONSTANT;
516 1750 : bool is_to_nodal = fe_type.family == LAGRANGE;
517 :
518 : // We support L2_LAGRANGE elemental variable with the first order
519 1750 : if (fe_type.order > FIRST && !is_to_nodal)
520 0 : mooseError("We don't currently support second order or higher elemental variable ");
521 :
522 1750 : if (is_to_nodal)
523 : {
524 1180 : const std::vector<Node *> & target_local_nodes = getTargetLocalNodes(i_to);
525 :
526 178942 : for (const auto & node : target_local_nodes)
527 : {
528 : // Skip this node if the variable has no dofs at it.
529 177762 : if (node->n_dofs(sys_num, var_num) < 1)
530 7040 : continue;
531 :
532 : // If nothing is in the node_index_map for a given local node,
533 : // it will get the value 0.
534 170722 : Real best_val = 0;
535 170722 : if (!_neighbors_cached)
536 : {
537 : // Search through all the incoming evaluation points from
538 : // different processors for the one with the closest
539 : // point. If there are multiple values from other processors
540 : // which are equidistant, the first one we check will "win".
541 166762 : Real min_dist = std::numeric_limits<Real>::max();
542 448282 : for (auto & evals : incoming_evals)
543 : {
544 : // processor Id
545 281520 : const processor_id_type pid = evals.first;
546 281520 : std::pair<unsigned int, dof_id_type> key(i_to, node->id());
547 281520 : if (node_index_map[pid].find(key) == node_index_map[pid].end())
548 61119 : continue;
549 278469 : unsigned int qp_ind = node_index_map[pid][key];
550 : // Distances
551 278469 : if (evals.second[2 * qp_ind] >= min_dist)
552 58068 : continue;
553 :
554 : // If we made it here, we are going to set a new value and
555 : // distance because we found one that was closer.
556 220401 : min_dist = evals.second[2 * qp_ind];
557 220401 : best_val = evals.second[2 * qp_ind + 1];
558 :
559 220401 : if (_fixed_meshes)
560 : {
561 : // Cache these indices.
562 2365 : _cached_from_inds[std::make_pair(i_to, node->id())] = pid;
563 2365 : _cached_qp_inds[std::make_pair(i_to, node->id())] = qp_ind;
564 : }
565 : }
566 : }
567 :
568 : else
569 : {
570 3960 : best_val = incoming_evals[_cached_from_inds[std::make_pair(i_to, node->id())]]
571 3960 : [_cached_qp_inds[std::make_pair(i_to, node->id())]];
572 : }
573 :
574 170722 : dof_id_type dof = node->dof_number(sys_num, var_num, 0);
575 170722 : solution->set(dof, best_val);
576 : }
577 : }
578 : else // Elemental
579 : {
580 570 : std::vector<Point> points;
581 570 : std::vector<dof_id_type> point_ids;
582 33768 : for (auto & elem : to_mesh.active_local_element_ptr_range())
583 : {
584 : // Skip this element if the variable has no dofs at it.
585 16599 : if (elem->n_dofs(sys_num, var_num) < 1)
586 0 : continue;
587 :
588 16599 : points.clear();
589 16599 : point_ids.clear();
590 : // grap sample points
591 : // for constant shape function, we take the element centroid
592 16599 : if (is_constant)
593 : {
594 16599 : points.push_back(elem->vertex_average());
595 16599 : point_ids.push_back(elem->id());
596 : }
597 : // for higher order method, we take all nodes of element
598 : // this works for the first order L2 Lagrange. Might not work
599 : // with something higher than the first order
600 : else
601 0 : for (auto & node : elem->node_ref_range())
602 : {
603 0 : points.push_back(node);
604 0 : point_ids.push_back(node.id());
605 : }
606 :
607 16599 : unsigned int n_comp = elem->n_comp(sys_num, var_num);
608 : // We assume each point corresponds to one component of elemental variable
609 16599 : if (points.size() != n_comp)
610 0 : mooseError(" Number of points ",
611 0 : points.size(),
612 : " does not equal to number of variable components ",
613 : n_comp);
614 :
615 33198 : for (MooseIndex(points) offset = 0; offset < points.size(); offset++)
616 : {
617 16599 : dof_id_type point_id = point_ids[offset];
618 16599 : Real best_val = 0;
619 16599 : if (!_neighbors_cached)
620 : {
621 14199 : Real min_dist = std::numeric_limits<Real>::max();
622 33737 : for (auto & evals : incoming_evals)
623 : {
624 19538 : const processor_id_type pid = evals.first;
625 :
626 19538 : std::pair<unsigned int, dof_id_type> key(i_to, point_id);
627 19538 : if (node_index_map[pid].find(key) == node_index_map[pid].end())
628 3321 : continue;
629 :
630 19208 : unsigned int qp_ind = node_index_map[pid][key];
631 19208 : if (evals.second[2 * qp_ind] >= min_dist)
632 2991 : continue;
633 :
634 16217 : min_dist = evals.second[2 * qp_ind];
635 16217 : best_val = evals.second[2 * qp_ind + 1];
636 :
637 16217 : if (_fixed_meshes)
638 : {
639 : // Cache these indices.
640 920 : _cached_from_inds[std::make_pair(i_to, point_id)] = pid;
641 920 : _cached_qp_inds[std::make_pair(i_to, point_id)] = qp_ind;
642 : } // if _fixed_meshes
643 : } // i_from
644 : } //
645 : else
646 : {
647 2400 : best_val = incoming_evals[_cached_from_inds[std::make_pair(i_to, point_id)]]
648 2400 : [_cached_qp_inds[std::make_pair(i_to, point_id)]];
649 : }
650 16599 : dof_id_type dof = elem->dof_number(sys_num, var_num, offset);
651 16599 : solution->set(dof, best_val);
652 : } // for offset
653 570 : }
654 570 : }
655 1750 : solution->close();
656 1750 : to_sys->update();
657 : }
658 :
659 1292 : if (_fixed_meshes)
660 110 : _neighbors_cached = true;
661 :
662 1292 : postExecute();
663 1292 : }
664 :
665 : Real
666 332376 : MultiAppNearestNodeTransfer::bboxMaxDistance(const Point & p, const BoundingBox & bbox)
667 : {
668 332376 : std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
669 :
670 332376 : std::array<Point, 8> all_points;
671 997128 : for (unsigned int x = 0; x < 2; x++)
672 1994256 : for (unsigned int y = 0; y < 2; y++)
673 3988512 : for (unsigned int z = 0; z < 2; z++)
674 2659008 : all_points[x + 2 * y + 4 * z] =
675 5318016 : Point(source_points[x](0), source_points[y](1), source_points[z](2));
676 :
677 332376 : Real max_distance = 0.;
678 :
679 2991384 : for (unsigned int i = 0; i < 8; i++)
680 : {
681 2659008 : Real distance = (p - all_points[i]).norm();
682 2659008 : if (distance > max_distance)
683 578088 : max_distance = distance;
684 : }
685 :
686 332376 : return max_distance;
687 : }
688 :
689 : Real
690 332376 : MultiAppNearestNodeTransfer::bboxMinDistance(const Point & p, const BoundingBox & bbox)
691 : {
692 332376 : std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
693 :
694 332376 : std::array<Point, 8> all_points;
695 997128 : for (unsigned int x = 0; x < 2; x++)
696 1994256 : for (unsigned int y = 0; y < 2; y++)
697 3988512 : for (unsigned int z = 0; z < 2; z++)
698 2659008 : all_points[x + 2 * y + 4 * z] =
699 5318016 : Point(source_points[x](0), source_points[y](1), source_points[z](2));
700 :
701 332376 : Real min_distance = std::numeric_limits<Real>::max();
702 :
703 2991384 : for (unsigned int i = 0; i < 8; i++)
704 : {
705 2659008 : Real distance = (p - all_points[i]).norm();
706 2659008 : if (distance < min_distance)
707 572339 : min_distance = distance;
708 : }
709 :
710 332376 : return min_distance;
711 : }
712 :
713 : void
714 1653 : MultiAppNearestNodeTransfer::getLocalEntitiesAndComponents(
715 : MooseMesh * mesh,
716 : std::vector<std::pair<Point, DofObject *>> & local_entities,
717 : std::vector<unsigned int> & local_comps,
718 : bool is_nodal,
719 : bool is_constant)
720 : {
721 : mooseAssert(mesh, "mesh should not be a nullptr");
722 : mooseAssert(local_entities.empty(), "local_entities should be empty");
723 1653 : MeshBase & mesh_base = mesh->getMesh();
724 :
725 1653 : if (isParamValid("source_boundary"))
726 : {
727 87 : const auto & sb = getParam<BoundaryName>("source_boundary");
728 87 : BoundaryID src_bnd_id = mesh->getBoundaryID(sb);
729 87 : if (!MooseMeshUtils::hasBoundaryName(mesh_base, sb))
730 0 : paramError("source_boundary", "The boundary '", sb, "' was not found in the mesh");
731 :
732 87 : if (is_nodal)
733 : {
734 87 : const ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
735 4849 : for (const auto & bnode : bnd_nodes)
736 : {
737 4762 : unsigned int comp = 0;
738 5950 : if (bnode->_bnd_id == src_bnd_id &&
739 1188 : bnode->_node->processor_id() == mesh_base.processor_id())
740 : {
741 1144 : local_entities.emplace_back(*bnode->_node, bnode->_node);
742 1144 : local_comps.push_back(comp++);
743 : }
744 : }
745 : }
746 : else
747 : {
748 0 : const ConstBndElemRange & bnd_elems = *mesh->getBoundaryElementRange();
749 0 : for (const auto & belem : bnd_elems)
750 : {
751 0 : unsigned int comp = 0;
752 0 : if (belem->_bnd_id == src_bnd_id &&
753 0 : belem->_elem->processor_id() == mesh_base.processor_id())
754 : {
755 : // CONSTANT Monomial
756 0 : if (is_constant)
757 : {
758 0 : local_entities.emplace_back(belem->_elem->vertex_average(), belem->_elem);
759 0 : local_comps.push_back(comp++);
760 : }
761 : // L2_LAGRANGE
762 : else
763 : {
764 0 : for (auto & node : belem->_elem->node_ref_range())
765 : {
766 0 : local_entities.emplace_back(node, belem->_elem);
767 0 : local_comps.push_back(comp++);
768 : }
769 : }
770 : }
771 : }
772 : }
773 : }
774 : else
775 : {
776 1566 : if (is_nodal)
777 : {
778 1276 : local_entities.reserve(mesh_base.n_local_nodes());
779 114328 : for (auto & node : mesh_base.local_node_ptr_range())
780 : {
781 113052 : unsigned int comp = 0;
782 113052 : local_entities.emplace_back(*node, node);
783 113052 : local_comps.push_back(comp++);
784 1276 : }
785 : }
786 : else
787 : {
788 290 : local_entities.reserve(mesh_base.n_local_elem());
789 17730 : for (auto & elem : mesh_base.active_local_element_ptr_range())
790 : {
791 8720 : unsigned int comp = 0;
792 : // CONSTANT Monomial
793 8720 : if (is_constant)
794 : {
795 8720 : local_entities.emplace_back(elem->vertex_average(), elem);
796 8720 : local_comps.push_back(comp++);
797 : }
798 : // L2_LAGRANGE
799 : else
800 : {
801 0 : for (auto & node : elem->node_ref_range())
802 : {
803 0 : local_entities.emplace_back(node, elem);
804 0 : local_comps.push_back(comp++);
805 : }
806 : }
807 290 : }
808 : }
809 : }
810 1653 : }
811 :
812 : const std::vector<Node *> &
813 2239 : MultiAppNearestNodeTransfer::getTargetLocalNodes(const unsigned int to_problem_id)
814 : {
815 2239 : _target_local_nodes.clear();
816 2239 : MeshBase & to_mesh = _to_meshes[to_problem_id]->getMesh();
817 :
818 2239 : if (isParamValid("target_boundary"))
819 : {
820 : const std::vector<BoundaryName> & target_boundaries =
821 128 : getParam<std::vector<BoundaryName>>("target_boundary");
822 292 : for (const auto & b : target_boundaries)
823 172 : if (!MooseMeshUtils::hasBoundaryName(to_mesh, b))
824 8 : paramError("target_boundary", "The boundary '", b, "' was not found in the mesh");
825 :
826 120 : ConstBndNodeRange & bnd_nodes = *(_to_meshes[to_problem_id])->getBoundaryNodeRange();
827 :
828 284 : for (const auto & t : target_boundaries)
829 : {
830 164 : BoundaryID target_bnd_id = _to_meshes[to_problem_id]->getBoundaryID(t);
831 :
832 8564 : for (const auto & bnode : bnd_nodes)
833 10510 : if (bnode->_bnd_id == target_bnd_id &&
834 2110 : bnode->_node->processor_id() == _to_meshes[to_problem_id]->processor_id())
835 1952 : _target_local_nodes.push_back(bnode->_node);
836 : }
837 : }
838 : else
839 : {
840 2111 : _target_local_nodes.resize(to_mesh.n_local_nodes());
841 2111 : unsigned int i = 0;
842 351723 : for (auto & node : to_mesh.local_node_ptr_range())
843 351723 : _target_local_nodes[i++] = node;
844 : }
845 :
846 2231 : return _target_local_nodes;
847 : }
|
