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 "NearestNodeLocator.h"
11 : #include "MooseMesh.h"
12 : #include "SubProblem.h"
13 : #include "SecondaryNeighborhoodThread.h"
14 : #include "NearestNodeThread.h"
15 : #include "Moose.h"
16 : #include "KDTree.h"
17 : #include "Conversion.h"
18 : #include "MooseApp.h"
19 :
20 : // libMesh
21 : #include "libmesh/boundary_info.h"
22 : #include "libmesh/elem.h"
23 : #include "libmesh/plane.h"
24 : #include "libmesh/mesh_tools.h"
25 :
26 1768 : NearestNodeLocator::NearestNodeLocator(SubProblem & subproblem,
27 : MooseMesh & mesh,
28 : BoundaryID boundary1,
29 1768 : BoundaryID boundary2)
30 : : Restartable(subproblem.getMooseApp(),
31 3536 : Moose::stringify(boundary1) + Moose::stringify(boundary2),
32 : "NearestNodeLocator",
33 : 0),
34 1768 : PerfGraphInterface(subproblem.getMooseApp().perfGraph(),
35 3536 : "NearestNodeLocator_" + Moose::stringify(boundary1) + "_" +
36 1768 : Moose::stringify(boundary2)),
37 1768 : _subproblem(subproblem),
38 1768 : _mesh(mesh),
39 1768 : _boundary1(boundary1),
40 1768 : _boundary2(boundary2),
41 1768 : _first(true),
42 1768 : _reinit_iteration(true),
43 8840 : _patch_update_strategy(_mesh.getPatchUpdateStrategy())
44 : {
45 : /*
46 : //sanity check on boundary ids
47 : const std::set<BoundaryID>& bids=_mesh.getBoundaryIDs();
48 : std::set<BoundaryID>::const_iterator sit;
49 : sit=bids.find(_boundary1);
50 : if (sit == bids.end())
51 : mooseError("NearestNodeLocator being created for boundaries ", _boundary1, " and ", _boundary2,
52 : ", but boundary ", _boundary1, " does not exist");
53 : sit=bids.find(_boundary2);
54 : if (sit == bids.end())
55 : mooseError("NearestNodeLocator being created for boundaries ", _boundary1, " and ", _boundary2,
56 : ", but boundary ", _boundary2, " does not exist");
57 : */
58 1768 : }
59 :
60 3520 : NearestNodeLocator::~NearestNodeLocator() = default;
61 :
62 : void
63 173096 : NearestNodeLocator::findNodes()
64 : {
65 173096 : TIME_SECTION("findNodes", 3, "Finding Nearest Nodes");
66 :
67 : /**
68 : * If this is the first time through we're going to build up a "neighborhood" of nodes
69 : * surrounding each of the secondary nodes. This will speed searching later.
70 : */
71 173096 : const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map = _mesh.nodeToElemMap();
72 :
73 173096 : if (_first || (_reinit_iteration && _patch_update_strategy == Moose::Iteration))
74 : {
75 11364 : _first = false;
76 :
77 : // Trial secondary nodes are all the nodes on the secondary side
78 : // We only keep the ones that are either on this processor or are likely
79 : // to interact with elements on this processor (ie nodes owned by this processor
80 : // are in the "neighborhood" of the secondary node
81 11364 : std::vector<dof_id_type> trial_secondary_nodes;
82 11364 : std::vector<dof_id_type> trial_primary_nodes;
83 :
84 : // Build a bounding box. No reason to consider nodes outside of our inflated BB
85 11364 : std::unique_ptr<BoundingBox> my_inflated_box = nullptr;
86 :
87 11364 : const std::vector<Real> & inflation = _mesh.getGhostedBoundaryInflation();
88 :
89 : // This means there was a user specified inflation... so we can build a BB
90 11364 : if (inflation.size() > 0)
91 : {
92 0 : BoundingBox my_box = MeshTools::create_local_bounding_box(_mesh);
93 :
94 0 : Point distance;
95 0 : for (unsigned int i = 0; i < inflation.size(); ++i)
96 0 : distance(i) = inflation[i];
97 :
98 : my_inflated_box =
99 0 : std::make_unique<BoundingBox>(my_box.first - distance, my_box.second + distance);
100 : }
101 :
102 : // Data structures to hold the boundary nodes
103 11364 : ConstBndNodeRange & bnd_nodes = *_mesh.getBoundaryNodeRange();
104 2196962 : for (const auto & bnode : bnd_nodes)
105 : {
106 2185598 : BoundaryID boundary_id = bnode->_bnd_id;
107 2185598 : dof_id_type node_id = bnode->_node->id();
108 :
109 : // If we have a BB only consider saving this node if it's in our inflated BB
110 2185598 : if (!my_inflated_box || (my_inflated_box->contains_point(*bnode->_node)))
111 : {
112 2185598 : if (boundary_id == _boundary1)
113 1004525 : trial_primary_nodes.push_back(node_id);
114 1181073 : else if (boundary_id == _boundary2)
115 61895 : trial_secondary_nodes.push_back(node_id);
116 : }
117 : }
118 :
119 : // Convert trial primary nodes to a vector of Points. This will be used to
120 : // construct the Kdtree.
121 11364 : std::vector<Point> primary_points(trial_primary_nodes.size());
122 1015889 : for (unsigned int i = 0; i < trial_primary_nodes.size(); ++i)
123 : {
124 1004525 : const Node & node = _mesh.nodeRef(trial_primary_nodes[i]);
125 1004525 : primary_points[i] = node;
126 : }
127 :
128 : // Create object kd_tree of class KDTree using the coordinates of trial
129 : // primary nodes.
130 11364 : KDTree kd_tree(primary_points, _mesh.getMaxLeafSize());
131 :
132 : NodeIdRange trial_secondary_node_range(
133 11364 : trial_secondary_nodes.begin(), trial_secondary_nodes.end(), 1);
134 :
135 : SecondaryNeighborhoodThread snt(
136 11364 : _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getPatchSize(), kd_tree);
137 :
138 11364 : Threads::parallel_reduce(trial_secondary_node_range, snt);
139 :
140 11364 : _secondary_nodes = snt._secondary_nodes;
141 11364 : _neighbor_nodes = snt._neighbor_nodes;
142 :
143 : // If 'iteration' patch update strategy is used, a second neighborhood
144 : // search using the ghosting_patch_size, which is larger than the regular
145 : // patch_size used for contact search, is conducted. The ghosted element set
146 : // given by this search is used for ghosting the elements connected to the
147 : // secondary and neighboring primary nodes.
148 11364 : if (_patch_update_strategy == Moose::Iteration)
149 : {
150 : SecondaryNeighborhoodThread snt_ghosting(
151 8630 : _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getGhostingPatchSize(), kd_tree);
152 :
153 8630 : Threads::parallel_reduce(trial_secondary_node_range, snt_ghosting);
154 :
155 768938 : for (const auto & dof : snt_ghosting._ghosted_elems)
156 760308 : _subproblem.addGhostedElem(dof);
157 8630 : }
158 : else
159 : {
160 100655 : for (const auto & dof : snt._ghosted_elems)
161 97921 : _subproblem.addGhostedElem(dof);
162 : }
163 :
164 : // Cache the secondary_node_range so we don't have to build it each time
165 : _secondary_node_range =
166 11364 : std::make_unique<NodeIdRange>(_secondary_nodes.begin(), _secondary_nodes.end(), 1);
167 11364 : }
168 :
169 173096 : _nearest_node_info.clear();
170 :
171 173096 : NearestNodeThread nnt(_mesh, _neighbor_nodes);
172 :
173 173096 : Threads::parallel_reduce(*_secondary_node_range, nnt);
174 :
175 173096 : _max_patch_percentage = nnt._max_patch_percentage;
176 :
177 173096 : _nearest_node_info = nnt._nearest_node_info;
178 :
179 173096 : if (_patch_update_strategy == Moose::Iteration)
180 : {
181 : // Get the set of elements that are currently being ghosted
182 8630 : std::set<dof_id_type> ghost = _subproblem.ghostedElems();
183 :
184 34131 : for (const auto & node_id : *_secondary_node_range)
185 : {
186 25501 : const Node * nearest_node = _nearest_node_info[node_id]._nearest_node;
187 :
188 : // Check if the elements attached to the nearest node are within the ghosted
189 : // set of elements. If not produce an error.
190 25501 : auto node_to_elem_pair = node_to_elem_map.find(nearest_node->id());
191 :
192 25501 : if (node_to_elem_pair != node_to_elem_map.end())
193 : {
194 25501 : const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
195 76503 : for (const auto & dof : elems_connected_to_node)
196 85540 : if (std::find(ghost.begin(), ghost.end(), dof) == ghost.end() &&
197 34538 : _mesh.elemPtr(dof)->processor_id() != _mesh.processor_id())
198 0 : mooseError("Error in NearestNodeLocator: The nearest neighbor lies outside the "
199 : "ghosted set of elements. Increase the ghosting_patch_size parameter in the "
200 : "mesh block and try again.");
201 : }
202 : }
203 8630 : }
204 173096 : }
205 :
206 : void
207 990 : NearestNodeLocator::reinit()
208 : {
209 990 : TIME_SECTION("reinit", 3, "Reinitializing Nearest Node Search");
210 :
211 : // Reset all data
212 990 : _secondary_node_range.reset();
213 990 : _nearest_node_info.clear();
214 :
215 990 : _first = true;
216 :
217 990 : _secondary_nodes.clear();
218 990 : _neighbor_nodes.clear();
219 :
220 990 : _new_ghosted_elems.clear();
221 :
222 : // After a call from system reinit, mesh has been updated with initial adaptivity.
223 : // Moose::Iteration relies on data generated for ghosting (i.e. trial_primary_nodes)
224 990 : _reinit_iteration = true;
225 :
226 : // Redo the search
227 990 : findNodes();
228 :
229 990 : _reinit_iteration = false;
230 990 : }
231 :
232 : Real
233 515820 : NearestNodeLocator::distance(dof_id_type node_id)
234 : {
235 515820 : return _nearest_node_info[node_id]._distance;
236 : }
237 :
238 : const Node *
239 561953 : NearestNodeLocator::nearestNode(dof_id_type node_id)
240 : {
241 561953 : return _nearest_node_info[node_id]._nearest_node;
242 : }
243 :
244 : void
245 0 : NearestNodeLocator::updatePatch(std::vector<dof_id_type> & secondary_nodes)
246 : {
247 0 : TIME_SECTION("updatePatch", 3, "Updating Nearest Node Search Patch");
248 :
249 0 : std::vector<dof_id_type> trial_primary_nodes;
250 :
251 : // Build a bounding box. No reason to consider nodes outside of our inflated BB
252 0 : std::unique_ptr<BoundingBox> my_inflated_box = nullptr;
253 :
254 0 : const std::vector<Real> & inflation = _mesh.getGhostedBoundaryInflation();
255 :
256 : // This means there was a user specified inflation... so we can build a BB
257 0 : if (inflation.size() > 0)
258 : {
259 0 : BoundingBox my_box = MeshTools::create_local_bounding_box(_mesh);
260 :
261 0 : Point distance;
262 0 : for (unsigned int i = 0; i < inflation.size(); ++i)
263 0 : distance(i) = inflation[i];
264 :
265 : my_inflated_box =
266 0 : std::make_unique<BoundingBox>(my_box.first - distance, my_box.second + distance);
267 : }
268 :
269 : // Data structures to hold the boundary nodes
270 0 : ConstBndNodeRange & bnd_nodes = *_mesh.getBoundaryNodeRange();
271 0 : for (const auto & bnode : bnd_nodes)
272 : {
273 0 : BoundaryID boundary_id = bnode->_bnd_id;
274 0 : dof_id_type node_id = bnode->_node->id();
275 :
276 : // If we have a BB only consider saving this node if it's in our inflated BB
277 0 : if (!my_inflated_box || (my_inflated_box->contains_point(*bnode->_node)))
278 : {
279 0 : if (boundary_id == _boundary1)
280 0 : trial_primary_nodes.push_back(node_id);
281 : }
282 : }
283 :
284 : // Convert trial primary nodes to a vector of Points. This will be used to construct the KDTree.
285 0 : std::vector<Point> primary_points(trial_primary_nodes.size());
286 0 : for (unsigned int i = 0; i < trial_primary_nodes.size(); ++i)
287 : {
288 0 : const Node & node = _mesh.nodeRef(trial_primary_nodes[i]);
289 0 : primary_points[i] = node;
290 : }
291 :
292 0 : const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map = _mesh.nodeToElemMap();
293 :
294 : // Create object kd_tree of class KDTree using the coordinates of trial
295 : // primary nodes.
296 0 : KDTree kd_tree(primary_points, _mesh.getMaxLeafSize());
297 :
298 0 : NodeIdRange secondary_node_range(secondary_nodes.begin(), secondary_nodes.end(), 1);
299 :
300 : SecondaryNeighborhoodThread snt(
301 0 : _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getPatchSize(), kd_tree);
302 :
303 0 : Threads::parallel_reduce(secondary_node_range, snt);
304 :
305 : // Calculate new ghosting patch for the secondary_node_range
306 : SecondaryNeighborhoodThread snt_ghosting(
307 0 : _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getGhostingPatchSize(), kd_tree);
308 :
309 0 : Threads::parallel_reduce(secondary_node_range, snt_ghosting);
310 :
311 : // Add the new set of elements that need to be ghosted into _new_ghosted_elems
312 0 : for (const auto & dof : snt_ghosting._ghosted_elems)
313 0 : _new_ghosted_elems.push_back(dof);
314 :
315 0 : std::vector<dof_id_type> tracked_secondary_nodes = snt._secondary_nodes;
316 :
317 : // Update the neighbor nodes (patch) for these tracked secondary nodes
318 0 : for (const auto & node_id : tracked_secondary_nodes)
319 0 : _neighbor_nodes[node_id] = snt._neighbor_nodes[node_id];
320 :
321 : NodeIdRange tracked_secondary_node_range(
322 0 : tracked_secondary_nodes.begin(), tracked_secondary_nodes.end(), 1);
323 :
324 0 : NearestNodeThread nnt(_mesh, snt._neighbor_nodes);
325 :
326 0 : Threads::parallel_reduce(tracked_secondary_node_range, nnt);
327 :
328 0 : _max_patch_percentage = nnt._max_patch_percentage;
329 :
330 : // Get the set of elements that are currently being ghosted
331 0 : std::set<dof_id_type> ghost = _subproblem.ghostedElems();
332 :
333 : // Update the nearest node information corresponding to these tracked secondary nodes
334 0 : for (const auto & node_id : tracked_secondary_node_range)
335 : {
336 0 : _nearest_node_info[node_id] = nnt._nearest_node_info[node_id];
337 :
338 : // Check if the elements attached to the nearest node are within the ghosted
339 : // set of elements. If not produce an error.
340 0 : const Node * nearest_node = nnt._nearest_node_info[node_id]._nearest_node;
341 :
342 0 : auto node_to_elem_pair = node_to_elem_map.find(nearest_node->id());
343 :
344 0 : if (node_to_elem_pair != node_to_elem_map.end())
345 : {
346 0 : const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
347 0 : for (const auto & dof : elems_connected_to_node)
348 0 : if (std::find(ghost.begin(), ghost.end(), dof) == ghost.end() &&
349 0 : _mesh.elemPtr(dof)->processor_id() != _mesh.processor_id())
350 0 : mooseError("Error in NearestNodeLocator: The nearest neighbor lies outside the ghosted "
351 : "set of elements. Increase the ghosting_patch_size parameter in the mesh "
352 : "block and try again.");
353 : }
354 : }
355 0 : }
356 :
357 : void
358 660 : NearestNodeLocator::updateGhostedElems()
359 : {
360 660 : TIME_SECTION("updateGhostedElems", 5, "Updating Nearest Node Search Because of Ghosting");
361 :
362 : // When 'iteration' patch update strategy is used, add the elements in
363 : // _new_ghosted_elems, which were accumulated in the nonlinear iterations
364 : // during the previous time step, to the list of ghosted elements. Also clear
365 : // the _new_ghosted_elems array for storing the ghosted elements from the
366 : // nonlinear iterations in the current time step.
367 :
368 660 : for (const auto & dof : _new_ghosted_elems)
369 0 : _subproblem.addGhostedElem(dof);
370 :
371 660 : _new_ghosted_elems.clear();
372 660 : }
373 : //===================================================================
374 1194360 : NearestNodeLocator::NearestNodeInfo::NearestNodeInfo()
375 1194360 : : _nearest_node(nullptr), _distance(std::numeric_limits<Real>::max())
376 : {
377 1194360 : }
|
