https://mooseframework.inl.gov
NearestNodeLocator.C
Go to the documentation of this file.
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"
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 
27  MooseMesh & mesh,
28  BoundaryID boundary1,
29  BoundaryID boundary2)
30  : Restartable(subproblem.getMooseApp(),
31  Moose::stringify(boundary1) + Moose::stringify(boundary2),
32  "NearestNodeLocator",
33  0),
34  PerfGraphInterface(subproblem.getMooseApp().perfGraph(),
35  "NearestNodeLocator_" + Moose::stringify(boundary1) + "_" +
36  Moose::stringify(boundary2)),
37  _subproblem(subproblem),
38  _mesh(mesh),
39  _boundary1(boundary1),
40  _boundary2(boundary2),
41  _first(true),
42  _reinit_iteration(true),
43  _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 }
59 
61 
62 void
64 {
65  TIME_SECTION("findNodes", 3, "Finding Nearest Nodes");
66 
71  const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map = _mesh.nodeToElemMap();
72 
74  {
75  _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  std::vector<dof_id_type> trial_secondary_nodes;
82  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  std::unique_ptr<BoundingBox> my_inflated_box = nullptr;
86 
87  const std::vector<Real> & inflation = _mesh.getGhostedBoundaryInflation();
88 
89  // This means there was a user specified inflation... so we can build a BB
90  if (inflation.size() > 0)
91  {
92  BoundingBox my_box = MeshTools::create_local_bounding_box(_mesh);
93 
94  Point distance;
95  for (unsigned int i = 0; i < inflation.size(); ++i)
96  distance(i) = inflation[i];
97 
98  my_inflated_box =
99  std::make_unique<BoundingBox>(my_box.first - distance, my_box.second + distance);
100  }
101 
102  // Data structures to hold the boundary nodes
104  for (const auto & bnode : bnd_nodes)
105  {
106  BoundaryID boundary_id = bnode->_bnd_id;
107  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  if (!my_inflated_box || (my_inflated_box->contains_point(*bnode->_node)))
111  {
112  if (boundary_id == _boundary1)
113  trial_primary_nodes.push_back(node_id);
114  else if (boundary_id == _boundary2)
115  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  std::vector<Point> primary_points(trial_primary_nodes.size());
122  for (unsigned int i = 0; i < trial_primary_nodes.size(); ++i)
123  {
124  const Node & node = _mesh.nodeRef(trial_primary_nodes[i]);
125  primary_points[i] = node;
126  }
127 
128  // Create object kd_tree of class KDTree using the coordinates of trial
129  // primary nodes.
130  KDTree kd_tree(primary_points, _mesh.getMaxLeafSize());
131 
132  NodeIdRange trial_secondary_node_range(
133  trial_secondary_nodes.begin(), trial_secondary_nodes.end(), 1);
134 
136  _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getPatchSize(), kd_tree);
137 
138  Threads::parallel_reduce(trial_secondary_node_range, snt);
139 
140  _secondary_nodes = snt._secondary_nodes;
141  _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.
149  {
150  SecondaryNeighborhoodThread snt_ghosting(
151  _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getGhostingPatchSize(), kd_tree);
152 
153  Threads::parallel_reduce(trial_secondary_node_range, snt_ghosting);
154 
155  for (const auto & dof : snt_ghosting._ghosted_elems)
157  }
158  else
159  {
160  for (const auto & dof : snt._ghosted_elems)
162  }
163 
164  // Cache the secondary_node_range so we don't have to build it each time
166  std::make_unique<NodeIdRange>(_secondary_nodes.begin(), _secondary_nodes.end(), 1);
167  }
168 
169  _nearest_node_info.clear();
170 
172 
173  Threads::parallel_reduce(*_secondary_node_range, nnt);
174 
176 
178 
180  {
181  // Get the set of elements that are currently being ghosted
182  std::set<dof_id_type> ghost = _subproblem.ghostedElems();
183 
184  for (const auto & node_id : *_secondary_node_range)
185  {
186  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  auto node_to_elem_pair = node_to_elem_map.find(nearest_node->id());
191 
192  if (node_to_elem_pair != node_to_elem_map.end())
193  {
194  const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
195  for (const auto & dof : elems_connected_to_node)
196  if (std::find(ghost.begin(), ghost.end(), dof) == ghost.end() &&
197  _mesh.elemPtr(dof)->processor_id() != _mesh.processor_id())
198  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  }
204 }
205 
206 void
208 {
209  TIME_SECTION("reinit", 3, "Reinitializing Nearest Node Search");
210 
211  // Reset all data
212  _secondary_node_range.reset();
213  _nearest_node_info.clear();
214 
215  _first = true;
216 
217  _secondary_nodes.clear();
218  _neighbor_nodes.clear();
219 
220  _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  _reinit_iteration = true;
225 
226  // Redo the search
227  findNodes();
228 
229  _reinit_iteration = false;
230 }
231 
232 Real
234 {
235  return _nearest_node_info[node_id]._distance;
236 }
237 
238 const Node *
240 {
241  const Node * returnval = _nearest_node_info[node_id]._nearest_node;
242  libmesh_assert(_mesh.getMesh().get_boundary_info().has_boundary_id(returnval, _boundary1));
243  return returnval;
244 }
245 
246 void
247 NearestNodeLocator::updatePatch(std::vector<dof_id_type> & secondary_nodes)
248 {
249  TIME_SECTION("updatePatch", 3, "Updating Nearest Node Search Patch");
250 
251  std::vector<dof_id_type> trial_primary_nodes;
252 
253  // Build a bounding box. No reason to consider nodes outside of our inflated BB
254  std::unique_ptr<BoundingBox> my_inflated_box = nullptr;
255 
256  const std::vector<Real> & inflation = _mesh.getGhostedBoundaryInflation();
257 
258  // This means there was a user specified inflation... so we can build a BB
259  if (inflation.size() > 0)
260  {
261  BoundingBox my_box = MeshTools::create_local_bounding_box(_mesh);
262 
263  Point distance;
264  for (unsigned int i = 0; i < inflation.size(); ++i)
265  distance(i) = inflation[i];
266 
267  my_inflated_box =
268  std::make_unique<BoundingBox>(my_box.first - distance, my_box.second + distance);
269  }
270 
271  // Data structures to hold the boundary nodes
273  for (const auto & bnode : bnd_nodes)
274  {
275  BoundaryID boundary_id = bnode->_bnd_id;
276  dof_id_type node_id = bnode->_node->id();
277 
278  // If we have a BB only consider saving this node if it's in our inflated BB
279  if (!my_inflated_box || (my_inflated_box->contains_point(*bnode->_node)))
280  {
281  if (boundary_id == _boundary1)
282  trial_primary_nodes.push_back(node_id);
283  }
284  }
285 
286  // Convert trial primary nodes to a vector of Points. This will be used to construct the KDTree.
287  std::vector<Point> primary_points(trial_primary_nodes.size());
288  for (unsigned int i = 0; i < trial_primary_nodes.size(); ++i)
289  {
290  const Node & node = _mesh.nodeRef(trial_primary_nodes[i]);
291  primary_points[i] = node;
292  }
293 
294  const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map = _mesh.nodeToElemMap();
295 
296  // Create object kd_tree of class KDTree using the coordinates of trial
297  // primary nodes.
298  KDTree kd_tree(primary_points, _mesh.getMaxLeafSize());
299 
300  NodeIdRange secondary_node_range(secondary_nodes.begin(), secondary_nodes.end(), 1);
301 
303  _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getPatchSize(), kd_tree);
304 
305  Threads::parallel_reduce(secondary_node_range, snt);
306 
307  // Calculate new ghosting patch for the secondary_node_range
308  SecondaryNeighborhoodThread snt_ghosting(
309  _mesh, trial_primary_nodes, node_to_elem_map, _mesh.getGhostingPatchSize(), kd_tree);
310 
311  Threads::parallel_reduce(secondary_node_range, snt_ghosting);
312 
313  // Add the new set of elements that need to be ghosted into _new_ghosted_elems
314  for (const auto & dof : snt_ghosting._ghosted_elems)
315  _new_ghosted_elems.push_back(dof);
316 
317  std::vector<dof_id_type> tracked_secondary_nodes = snt._secondary_nodes;
318 
319  // Update the neighbor nodes (patch) for these tracked secondary nodes
320  for (const auto & node_id : tracked_secondary_nodes)
321  _neighbor_nodes[node_id] = snt._neighbor_nodes[node_id];
322 
323  NodeIdRange tracked_secondary_node_range(
324  tracked_secondary_nodes.begin(), tracked_secondary_nodes.end(), 1);
325 
326  NearestNodeThread nnt(_mesh, snt._neighbor_nodes);
327 
328  Threads::parallel_reduce(tracked_secondary_node_range, nnt);
329 
330  _max_patch_percentage = nnt._max_patch_percentage;
331 
332  // Get the set of elements that are currently being ghosted
333  std::set<dof_id_type> ghost = _subproblem.ghostedElems();
334 
335  // Update the nearest node information corresponding to these tracked secondary nodes
336  for (const auto & node_id : tracked_secondary_node_range)
337  {
338  _nearest_node_info[node_id] = nnt._nearest_node_info[node_id];
339 
340  // Check if the elements attached to the nearest node are within the ghosted
341  // set of elements. If not produce an error.
342  const Node * nearest_node = nnt._nearest_node_info[node_id]._nearest_node;
343 
344  auto node_to_elem_pair = node_to_elem_map.find(nearest_node->id());
345 
346  if (node_to_elem_pair != node_to_elem_map.end())
347  {
348  const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
349  for (const auto & dof : elems_connected_to_node)
350  if (std::find(ghost.begin(), ghost.end(), dof) == ghost.end() &&
351  _mesh.elemPtr(dof)->processor_id() != _mesh.processor_id())
352  mooseError("Error in NearestNodeLocator: The nearest neighbor lies outside the ghosted "
353  "set of elements. Increase the ghosting_patch_size parameter in the mesh "
354  "block and try again.");
355  }
356  }
357 }
358 
359 void
361 {
362  TIME_SECTION("updateGhostedElems", 5, "Updating Nearest Node Search Because of Ghosting");
363 
364  // When 'iteration' patch update strategy is used, add the elements in
365  // _new_ghosted_elems, which were accumulated in the nonlinear iterations
366  // during the previous time step, to the list of ghosted elements. Also clear
367  // the _new_ghosted_elems array for storing the ghosted elements from the
368  // nonlinear iterations in the current time step.
369 
370  for (const auto & dof : _new_ghosted_elems)
372 
373  _new_ghosted_elems.clear();
374 }
375 //===================================================================
377  : _nearest_node(nullptr), _distance(std::numeric_limits<Real>::max())
378 {
379 }
std::map< dof_id_type, std::vector< dof_id_type > > _neighbor_nodes
void findNodes()
This is the main method that is going to start the search.
KOKKOS_INLINE_FUNCTION const T * find(const T &target, const T *const begin, const T *const end)
Find a value in an array.
Definition: KokkosUtils.h:30
A class for creating restricted objects.
Definition: Restartable.h:28
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:3153
Definition: KDTree.h:28
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
void updatePatch(std::vector< dof_id_type > &secondary_nodes)
Reconstructs the KDtree, updates the patch for the nodes in secondary_nodes, and updates the closest ...
unsigned int getGhostingPatchSize() const
Getter for the ghosting_patch_size parameter.
Definition: MooseMesh.h:630
MeshBase & mesh
std::map< dof_id_type, NearestNodeLocator::NearestNodeInfo > _nearest_node_info
Real distance(dof_id_type node_id)
Valid to call this after findNodes() has been called to get the distance to the nearest node...
const std::vector< Real > & getGhostedBoundaryInflation() const
Return a writable reference to the _ghosted_boundaries_inflation vector.
Definition: MooseMesh.C:3287
void reinit()
Completely redo the search from scratch.
virtual const Node & nodeRef(const dof_id_type i) const
Definition: MooseMesh.C:849
NearestNodeLocator(SubProblem &subproblem, MooseMesh &mesh, BoundaryID boundary1, BoundaryID boundary2)
auto max(const L &left, const R &right)
std::map< dof_id_type, NearestNodeInfo > _nearest_node_info
std::unique_ptr< NodeIdRange > _secondary_node_range
std::vector< dof_id_type > _secondary_nodes
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:3488
boundary_id_type BoundaryID
unsigned int getMaxLeafSize() const
Getter for the maximum leaf size parameter.
Definition: MooseMesh.h:635
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:92
libmesh_assert(ctx)
const Moose::PatchUpdateType _patch_update_strategy
void updateGhostedElems()
Updates the ghosted elements at the start of the time step for iteration patch update strategy...
virtual std::set< dof_id_type > & ghostedElems()
Return the list of elements that should have their DoFs ghosted to this processor.
Definition: SubProblem.h:672
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
unsigned int getPatchSize() const
Getter for the patch_size parameter.
Definition: MooseMesh.C:3441
Interface for objects interacting with the PerfGraph.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:78
const Node * nearestNode(dof_id_type node_id)
Valid to call this after findNodes() has been called to get a pointer to the nearest node...
virtual void addGhostedElem(dof_id_type elem_id)=0
Will make sure that all dofs connected to elem_id are ghosted to this processor.
SubProblem & _subproblem
std::vector< dof_id_type > _new_ghosted_elems
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
processor_id_type processor_id() const
libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
Definition: MooseMesh.C:1327
uint8_t dof_id_type
const std::map< dof_id_type, std::vector< dof_id_type > > & nodeToElemMap()
If not already created, creates a map from every node to all elements to which they are connected...
Definition: MooseMesh.C:1216