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 4299 : MultiAppNearestNodeTransfer::validParams()
33 : {
34 4299 : InputParameters params = MultiAppConservativeTransfer::validParams();
35 8598 : params.addClassDescription(
36 : "Transfer the value to the target domain from the nearest node in the source domain.");
37 :
38 17196 : params.addParam<BoundaryName>(
39 : "source_boundary",
40 : "The boundary we are transferring from (if not specified, whole domain is used).");
41 17196 : params.addParam<std::vector<BoundaryName>>(
42 : "target_boundary",
43 : "The boundary we are transferring to (if not specified, whole domain is used).");
44 8598 : params.addParam<bool>("fixed_meshes",
45 8598 : 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 4299 : MultiAppTransfer::addBBoxFactorParam(params);
52 4299 : return params;
53 0 : }
54 :
55 616 : MultiAppNearestNodeTransfer::MultiAppNearestNodeTransfer(const InputParameters & parameters)
56 : : MultiAppConservativeTransfer(parameters),
57 616 : _fixed_meshes(getParam<bool>("fixed_meshes")),
58 1232 : _node_map(declareRestartableData<std::map<dof_id_type, Node *>>("node_map")),
59 1232 : _distance_map(declareRestartableData<std::map<dof_id_type, Real>>("distance_map")),
60 1232 : _neighbors_cached(declareRestartableData<bool>("neighbors_cached", false)),
61 1232 : _cached_froms(declareRestartableData<std::map<processor_id_type, std::vector<unsigned int>>>(
62 : "cached_froms")),
63 1232 : _cached_dof_ids(declareRestartableData<std::map<processor_id_type, std::vector<dof_id_type>>>(
64 : "cached_dof_ids")),
65 616 : _cached_from_inds(
66 1232 : declareRestartableData<std::map<std::pair<unsigned int, dof_id_type>, unsigned int>>(
67 : "cached_from_ids")),
68 616 : _cached_qp_inds(
69 1232 : declareRestartableData<std::map<std::pair<unsigned int, dof_id_type>, unsigned int>>(
70 1232 : "cached_qp_inds"))
71 : {
72 616 : mooseDeprecated("MultiAppNearestNodeTransfer is deprecated. Use "
73 : "MultiAppGeneralFieldNearestNodeTransfer instead and adapt the parameters");
74 :
75 616 : if (_to_var_names.size() != 1)
76 0 : paramError("variable", " Support single to-variable only");
77 :
78 616 : if (_from_var_names.size() != 1)
79 0 : paramError("source_variable", " Support single from-variable only");
80 616 : }
81 :
82 : void
83 1347 : MultiAppNearestNodeTransfer::execute()
84 : {
85 6735 : TIME_SECTION(
86 : "MultiAppNearestNodeTransfer::execute()", 5, "Transferring variables based on nearest nodes");
87 :
88 : // Get the bounding boxes for the "from" domains.
89 1347 : std::vector<BoundingBox> bboxes;
90 4041 : if (isParamValid("source_boundary"))
91 : {
92 94 : if (_from_meshes.size())
93 : {
94 188 : const auto & sb = getParam<BoundaryName>("source_boundary");
95 94 : if (!MooseMeshUtils::hasBoundaryName(_from_meshes[0]->getMesh(), sb))
96 12 : paramError("source_boundary", "The boundary '", sb, "' was not found in the mesh");
97 :
98 88 : bboxes = getFromBoundingBoxes(_from_meshes[0]->getBoundaryID(sb));
99 : }
100 : else
101 0 : bboxes = getFromBoundingBoxes(Moose::INVALID_BOUNDARY_ID);
102 : }
103 : else
104 1253 : bboxes = getFromBoundingBoxes();
105 :
106 : // Figure out how many "from" domains each processor owns.
107 1341 : 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 1341 : 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 1341 : node_index_map;
130 :
131 1341 : if (!_neighbors_cached)
132 : {
133 2909 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
134 : {
135 1654 : System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
136 1654 : unsigned int sys_num = to_sys->number();
137 1654 : unsigned int var_num = to_sys->variable_number(_to_var_name);
138 1654 : MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
139 : const auto to_global_num =
140 1654 : _current_direction == FROM_MULTIAPP ? 0 : _to_local2global_map[i_to];
141 1654 : const auto & to_transform = *_to_transforms[to_global_num];
142 1654 : auto & fe_type = to_sys->variable_type(var_num);
143 1654 : bool is_constant = fe_type.order == CONSTANT;
144 1654 : bool is_to_nodal = fe_type.family == LAGRANGE;
145 :
146 : // We support L2_LAGRANGE elemental variable with the first order
147 1654 : if (fe_type.order > FIRST && !is_to_nodal)
148 0 : mooseError("We don't currently support second order or higher elemental variable ");
149 :
150 2774 : 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 1654 : if (is_to_nodal)
156 : {
157 1094 : 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 1088 : std::set<Node *> local_nodes_found;
163 :
164 189351 : for (const auto & node : target_local_nodes)
165 : {
166 : // Skip this node if the variable has no dofs at it.
167 188263 : if (node->n_dofs(sys_num, var_num) < 1)
168 7920 : continue;
169 :
170 180343 : const auto transformed_node = to_transform(*node);
171 :
172 : // Find which bboxes might have the nearest node to this point.
173 180343 : Real nearest_max_distance = std::numeric_limits<Real>::max();
174 509952 : for (const auto & bbox : bboxes)
175 : {
176 329609 : Real distance = bboxMaxDistance(transformed_node, bbox);
177 329609 : if (distance < nearest_max_distance)
178 245782 : nearest_max_distance = distance;
179 : }
180 :
181 180343 : unsigned int from0 = 0;
182 485789 : for (processor_id_type i_proc = 0; i_proc < n_processors();
183 305446 : from0 += froms_per_proc[i_proc], i_proc++)
184 : {
185 635055 : for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc]; i_from++)
186 : {
187 329609 : Real distance = bboxMinDistance(transformed_node, bboxes[i_from]);
188 :
189 340953 : if (distance <= nearest_max_distance ||
190 11344 : bboxes[i_from].contains_point(transformed_node))
191 : {
192 318265 : std::pair<unsigned int, dof_id_type> key(i_to, node->id());
193 : // Record a local ID for each quadrature point
194 318265 : node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
195 318265 : outgoing_qps[i_proc].push_back(transformed_node);
196 318265 : 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 189351 : for (const auto & node : target_local_nodes)
207 188263 : 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 1088 : }
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 560 : std::set<Elem *> local_elems_found;
221 560 : std::vector<Point> points;
222 560 : std::vector<dof_id_type> point_ids;
223 15931 : 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 15371 : if (elem->n_dofs(sys_num, var_num) < 1)
227 0 : continue;
228 :
229 15371 : points.clear();
230 15371 : point_ids.clear();
231 : // For constant monomial, we take the centroid of element
232 15371 : if (is_constant)
233 : {
234 15371 : points.push_back(to_transform(elem->vertex_average()));
235 15371 : 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 15371 : unsigned int offset = 0;
247 30742 : for (auto & point : points)
248 : {
249 : // Find which bboxes might have the nearest node to this point.
250 15371 : Real nearest_max_distance = std::numeric_limits<Real>::max();
251 45706 : for (const auto & bbox : bboxes)
252 : {
253 30335 : Real distance = bboxMaxDistance(point, bbox);
254 30335 : if (distance < nearest_max_distance)
255 19701 : nearest_max_distance = distance;
256 : }
257 :
258 15371 : unsigned int from0 = 0;
259 36081 : for (processor_id_type i_proc = 0; i_proc < n_processors();
260 20710 : from0 += froms_per_proc[i_proc], i_proc++)
261 : {
262 51045 : for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc]; i_from++)
263 : {
264 30335 : Real distance = bboxMinDistance(point, bboxes[i_from]);
265 30335 : 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 27041 : point_ids[offset]); // Create an unique ID
270 : // If this point already exist, we skip it
271 27041 : if (node_index_map[i_proc].find(key) != node_index_map[i_proc].end())
272 6661 : continue;
273 20380 : node_index_map[i_proc][key] = outgoing_qps[i_proc].size();
274 20380 : outgoing_qps[i_proc].push_back(point);
275 20380 : local_elems_found.insert(elem);
276 : } // if distance
277 : } // for i_from
278 : } // for i_proc
279 15371 : offset++;
280 : } // point
281 560 : } // 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 15931 : for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
287 15371 : 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 560 : to_transform(elem->vertex_average()));
294 560 : }
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 1335 : 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 1335 : std::map<processor_id_type, std::vector<Real>> processor_outgoing_evals;
312 :
313 1335 : 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 2510 : froms_per_proc[processor_id()]);
321 :
322 1255 : std::vector<std::vector<unsigned int>> local_comps(froms_per_proc[processor_id()]);
323 :
324 : // Local array of all from Variable references
325 1255 : std::vector<std::reference_wrapper<MooseVariableFEBase>> _from_vars;
326 :
327 2984 : for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
328 : i_local_from++)
329 : {
330 1729 : MooseVariableFEBase & from_var = _from_problems[i_local_from]->getVariable(
331 : 0, _from_var_name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD);
332 1729 : auto & from_fe_type = from_var.feType();
333 1729 : bool is_constant = from_fe_type.order == CONSTANT;
334 1729 : bool is_to_nodal = from_fe_type.family == LAGRANGE;
335 :
336 : // We support L2_LAGRANGE elemental variable with the first order
337 1729 : 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 1729 : _from_vars.emplace_back(from_var);
341 1729 : getLocalEntitiesAndComponents(_from_meshes[i_local_from],
342 1729 : local_entities[i_local_from],
343 1729 : local_comps[i_local_from],
344 : is_to_nodal,
345 : is_constant);
346 : }
347 :
348 : // Quadrature points I will receive from remote processors
349 1255 : std::map<processor_id_type, std::vector<Point>> incoming_qps;
350 1902 : auto qps_action_functor = [&incoming_qps](processor_id_type pid, const std::vector<Point> & qps)
351 : {
352 : // Quadrature points from processor 'pid'
353 1902 : auto & incoming_qps_from_pid = incoming_qps[pid];
354 : // Store data for late use
355 1902 : incoming_qps_from_pid.reserve(incoming_qps_from_pid.size() + qps.size());
356 1902 : std::copy(qps.begin(), qps.end(), std::back_inserter(incoming_qps_from_pid));
357 3157 : };
358 :
359 1255 : Parallel::push_parallel_vector_data(comm(), outgoing_qps, qps_action_functor);
360 :
361 3157 : for (auto & qps : incoming_qps)
362 : {
363 1902 : const processor_id_type pid = qps.first;
364 :
365 1902 : 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 1902 : 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 1902 : outgoing_evals.resize(2 * qps.second.size());
382 :
383 340547 : for (std::size_t qp = 0; qp < qps.second.size(); qp++)
384 : {
385 338645 : const Point & qpt = qps.second[qp];
386 338645 : outgoing_evals[2 * qp] = std::numeric_limits<Real>::max();
387 832881 : for (unsigned int i_local_from = 0; i_local_from < froms_per_proc[processor_id()];
388 : i_local_from++)
389 : {
390 494236 : MooseVariableFEBase & from_var = _from_vars[i_local_from];
391 494236 : System & from_sys = from_var.sys().system();
392 494236 : unsigned int from_sys_num = from_sys.number();
393 494236 : unsigned int from_var_num = from_sys.variable_number(from_var.name());
394 : const auto from_global_num =
395 494236 : _current_direction == TO_MULTIAPP ? 0 : _from_local2global_map[i_local_from];
396 494236 : const auto & from_transform = *_from_transforms[from_global_num];
397 :
398 35156430 : 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 34662194 : (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 34662194 : if (current_distance < outgoing_evals[2 * qp])
413 : {
414 : // Assuming LAGRANGE!
415 10729897 : if (local_entities[i_local_from][i_node].second->n_dofs(from_sys_num, from_var_num) >
416 : 0)
417 : {
418 21407656 : dof_id_type from_dof = local_entities[i_local_from][i_node].second->dof_number(
419 10703828 : 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 10703828 : outgoing_evals[2 * qp] = current_distance;
428 10703828 : outgoing_evals[2 * qp + 1] = (*from_sys.solution)(from_dof);
429 :
430 10703828 : 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 1255 : }
443 :
444 : else // We've cached the nearest nodes.
445 : {
446 202 : for (auto & problem_from : _cached_froms)
447 : {
448 122 : const processor_id_type pid = problem_from.first;
449 122 : std::vector<Real> & outgoing_evals = processor_outgoing_evals[pid];
450 122 : outgoing_evals.resize(problem_from.second.size());
451 :
452 9220 : for (unsigned int qp = 0; qp < outgoing_evals.size(); qp++)
453 : {
454 9098 : const auto from_problem = problem_from.second[qp];
455 9098 : 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 9098 : _from_problems[from_problem]->getVariable(0,
464 : _from_var_name,
465 : Moose::VarKindType::VAR_ANY,
466 : Moose::VarFieldType::VAR_FIELD_STANDARD);
467 9098 : System & from_sys = from_var.sys().system();
468 9098 : dof_id_type from_dof = _cached_dof_ids[pid][qp];
469 9098 : outgoing_evals[qp] = (*from_sys.solution)(from_dof);
470 : }
471 : }
472 : }
473 :
474 : auto evals_action_functor =
475 2024 : [&incoming_evals](processor_id_type pid, const std::vector<Real> & evals)
476 : {
477 : // evals for processor 'pid'
478 2024 : auto & incoming_evals_for_pid = incoming_evals[pid];
479 : // Copy evals for late use
480 2024 : incoming_evals_for_pid.reserve(incoming_evals_for_pid.size() + evals.size());
481 2024 : std::copy(evals.begin(), evals.end(), std::back_inserter(incoming_evals_for_pid));
482 3359 : };
483 :
484 1335 : 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 3159 : 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 1824 : System * to_sys = find_sys(*_to_es[i_to], _to_var_name);
495 :
496 1824 : unsigned int sys_num = to_sys->number();
497 1824 : unsigned int var_num = to_sys->variable_number(_to_var_name);
498 :
499 1824 : NumericVector<Real> * solution = nullptr;
500 1824 : switch (_current_direction)
501 : {
502 1073 : case TO_MULTIAPP:
503 1073 : solution = &getTransferVector(i_to, _to_var_name);
504 1073 : break;
505 751 : case FROM_MULTIAPP:
506 751 : solution = to_sys->solution.get();
507 751 : break;
508 0 : default:
509 0 : mooseError("Unknown direction");
510 : }
511 :
512 1824 : const MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
513 :
514 1824 : auto & fe_type = to_sys->variable_type(var_num);
515 1824 : bool is_constant = fe_type.order == CONSTANT;
516 1824 : bool is_to_nodal = fe_type.family == LAGRANGE;
517 :
518 : // We support L2_LAGRANGE elemental variable with the first order
519 1824 : if (fe_type.order > FIRST && !is_to_nodal)
520 0 : mooseError("We don't currently support second order or higher elemental variable ");
521 :
522 1824 : if (is_to_nodal)
523 : {
524 1230 : const std::vector<Node *> & target_local_nodes = getTargetLocalNodes(i_to);
525 :
526 193706 : for (const auto & node : target_local_nodes)
527 : {
528 : // Skip this node if the variable has no dofs at it.
529 192476 : if (node->n_dofs(sys_num, var_num) < 1)
530 7920 : continue;
531 :
532 : // If nothing is in the node_index_map for a given local node,
533 : // it will get the value 0.
534 184556 : Real best_val = 0;
535 184556 : 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 180343 : Real min_dist = std::numeric_limits<Real>::max();
542 485645 : for (auto & evals : incoming_evals)
543 : {
544 : // processor Id
545 305302 : const processor_id_type pid = evals.first;
546 305302 : std::pair<unsigned int, dof_id_type> key(i_to, node->id());
547 305302 : if (node_index_map[pid].find(key) == node_index_map[pid].end())
548 66270 : continue;
549 302251 : unsigned int qp_ind = node_index_map[pid][key];
550 : // Distances
551 302251 : if (evals.second[2 * qp_ind] >= min_dist)
552 63219 : 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 239032 : min_dist = evals.second[2 * qp_ind];
557 239032 : best_val = evals.second[2 * qp_ind + 1];
558 :
559 239032 : 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 4213 : best_val = incoming_evals[_cached_from_inds[std::make_pair(i_to, node->id())]]
571 4213 : [_cached_qp_inds[std::make_pair(i_to, node->id())]];
572 : }
573 :
574 184556 : dof_id_type dof = node->dof_number(sys_num, var_num, 0);
575 184556 : solution->set(dof, best_val);
576 : }
577 : }
578 : else // Elemental
579 : {
580 594 : std::vector<Point> points;
581 594 : std::vector<dof_id_type> point_ids;
582 18465 : for (auto & elem : to_mesh.active_local_element_ptr_range())
583 : {
584 : // Skip this element if the variable has no dofs at it.
585 17871 : if (elem->n_dofs(sys_num, var_num) < 1)
586 0 : continue;
587 :
588 17871 : points.clear();
589 17871 : point_ids.clear();
590 : // grap sample points
591 : // for constant shape function, we take the element centroid
592 17871 : if (is_constant)
593 : {
594 17871 : points.push_back(elem->vertex_average());
595 17871 : 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 17871 : unsigned int n_comp = elem->n_comp(sys_num, var_num);
608 : // We assume each point corresponds to one component of elemental variable
609 17871 : 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 35742 : for (MooseIndex(points) offset = 0; offset < points.size(); offset++)
616 : {
617 17871 : dof_id_type point_id = point_ids[offset];
618 17871 : Real best_val = 0;
619 17871 : if (!_neighbors_cached)
620 : {
621 15371 : Real min_dist = std::numeric_limits<Real>::max();
622 36081 : for (auto & evals : incoming_evals)
623 : {
624 20710 : const processor_id_type pid = evals.first;
625 :
626 20710 : std::pair<unsigned int, dof_id_type> key(i_to, point_id);
627 20710 : if (node_index_map[pid].find(key) == node_index_map[pid].end())
628 3321 : continue;
629 :
630 20380 : unsigned int qp_ind = node_index_map[pid][key];
631 20380 : if (evals.second[2 * qp_ind] >= min_dist)
632 2991 : continue;
633 :
634 17389 : min_dist = evals.second[2 * qp_ind];
635 17389 : best_val = evals.second[2 * qp_ind + 1];
636 :
637 17389 : 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 2500 : best_val = incoming_evals[_cached_from_inds[std::make_pair(i_to, point_id)]]
648 2500 : [_cached_qp_inds[std::make_pair(i_to, point_id)]];
649 : }
650 17871 : dof_id_type dof = elem->dof_number(sys_num, var_num, offset);
651 17871 : solution->set(dof, best_val);
652 : } // for offset
653 594 : }
654 594 : }
655 1824 : solution->close();
656 1824 : to_sys->update();
657 : }
658 :
659 1335 : if (_fixed_meshes)
660 113 : _neighbors_cached = true;
661 :
662 1335 : postExecute();
663 1335 : }
664 :
665 : Real
666 359944 : MultiAppNearestNodeTransfer::bboxMaxDistance(const Point & p, const BoundingBox & bbox)
667 : {
668 359944 : std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
669 :
670 359944 : std::array<Point, 8> all_points;
671 1079832 : for (unsigned int x = 0; x < 2; x++)
672 2159664 : for (unsigned int y = 0; y < 2; y++)
673 4319328 : for (unsigned int z = 0; z < 2; z++)
674 2879552 : all_points[x + 2 * y + 4 * z] =
675 5759104 : Point(source_points[x](0), source_points[y](1), source_points[z](2));
676 :
677 359944 : Real max_distance = 0.;
678 :
679 3239496 : for (unsigned int i = 0; i < 8; i++)
680 : {
681 2879552 : Real distance = (p - all_points[i]).norm();
682 2879552 : if (distance > max_distance)
683 623018 : max_distance = distance;
684 : }
685 :
686 359944 : return max_distance;
687 : }
688 :
689 : Real
690 359944 : MultiAppNearestNodeTransfer::bboxMinDistance(const Point & p, const BoundingBox & bbox)
691 : {
692 359944 : std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
693 :
694 359944 : std::array<Point, 8> all_points;
695 1079832 : for (unsigned int x = 0; x < 2; x++)
696 2159664 : for (unsigned int y = 0; y < 2; y++)
697 4319328 : for (unsigned int z = 0; z < 2; z++)
698 2879552 : all_points[x + 2 * y + 4 * z] =
699 5759104 : Point(source_points[x](0), source_points[y](1), source_points[z](2));
700 :
701 359944 : Real min_distance = std::numeric_limits<Real>::max();
702 :
703 3239496 : for (unsigned int i = 0; i < 8; i++)
704 : {
705 2879552 : Real distance = (p - all_points[i]).norm();
706 2879552 : if (distance < min_distance)
707 616658 : min_distance = distance;
708 : }
709 :
710 359944 : return min_distance;
711 : }
712 :
713 : void
714 1729 : 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 1729 : MeshBase & mesh_base = mesh->getMesh();
724 :
725 5187 : if (isParamValid("source_boundary"))
726 : {
727 184 : const auto & sb = getParam<BoundaryName>("source_boundary");
728 92 : BoundaryID src_bnd_id = mesh->getBoundaryID(sb);
729 92 : if (!MooseMeshUtils::hasBoundaryName(mesh_base, sb))
730 0 : paramError("source_boundary", "The boundary '", sb, "' was not found in the mesh");
731 :
732 92 : if (is_nodal)
733 : {
734 92 : const ConstBndNodeRange & bnd_nodes = *mesh->getBoundaryNodeRange();
735 4930 : for (const auto & bnode : bnd_nodes)
736 : {
737 4838 : unsigned int comp = 0;
738 6045 : if (bnode->_bnd_id == src_bnd_id &&
739 1207 : bnode->_node->processor_id() == mesh_base.processor_id())
740 : {
741 1163 : local_entities.emplace_back(*bnode->_node, bnode->_node);
742 1163 : 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 1637 : if (is_nodal)
777 : {
778 1328 : local_entities.reserve(mesh_base.n_local_nodes());
779 120323 : for (auto & node : mesh_base.local_node_ptr_range())
780 : {
781 118995 : unsigned int comp = 0;
782 118995 : local_entities.emplace_back(*node, node);
783 118995 : local_comps.push_back(comp++);
784 1328 : }
785 : }
786 : else
787 : {
788 309 : local_entities.reserve(mesh_base.n_local_elem());
789 9774 : for (auto & elem : mesh_base.active_local_element_ptr_range())
790 : {
791 9465 : unsigned int comp = 0;
792 : // CONSTANT Monomial
793 9465 : if (is_constant)
794 : {
795 9465 : local_entities.emplace_back(elem->vertex_average(), elem);
796 9465 : 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 309 : }
808 : }
809 : }
810 1729 : }
811 :
812 : const std::vector<Node *> &
813 2324 : MultiAppNearestNodeTransfer::getTargetLocalNodes(const unsigned int to_problem_id)
814 : {
815 2324 : _target_local_nodes.clear();
816 2324 : MeshBase & to_mesh = _to_meshes[to_problem_id]->getMesh();
817 :
818 6972 : if (isParamValid("target_boundary"))
819 : {
820 : const std::vector<BoundaryName> & target_boundaries =
821 264 : getParam<std::vector<BoundaryName>>("target_boundary");
822 306 : for (const auto & b : target_boundaries)
823 180 : if (!MooseMeshUtils::hasBoundaryName(to_mesh, b))
824 12 : paramError("target_boundary", "The boundary '", b, "' was not found in the mesh");
825 :
826 126 : ConstBndNodeRange & bnd_nodes = *(_to_meshes[to_problem_id])->getBoundaryNodeRange();
827 :
828 300 : for (const auto & t : target_boundaries)
829 : {
830 174 : BoundaryID target_bnd_id = _to_meshes[to_problem_id]->getBoundaryID(t);
831 :
832 8806 : for (const auto & bnode : bnd_nodes)
833 10800 : if (bnode->_bnd_id == target_bnd_id &&
834 2168 : bnode->_node->processor_id() == _to_meshes[to_problem_id]->processor_id())
835 2010 : _target_local_nodes.push_back(bnode->_node);
836 : }
837 : }
838 : else
839 : {
840 2192 : _target_local_nodes.resize(to_mesh.n_local_nodes());
841 2192 : unsigned int i = 0;
842 380921 : for (auto & node : to_mesh.local_node_ptr_range())
843 380921 : _target_local_nodes[i++] = node;
844 : }
845 :
846 2318 : return _target_local_nodes;
847 : }
|
