libMesh
boundary_info.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // C++ includes
21 #include <iterator> // std::distance
22 
23 // Local includes
24 #include "libmesh/libmesh_config.h"
25 
26 #include "libmesh/boundary_info.h"
27 #include "libmesh/distributed_mesh.h"
28 #include "libmesh/elem.h"
29 #include "libmesh/mesh_communication.h"
30 #include "libmesh/mesh_serializer.h"
31 #include "libmesh/parallel.h"
32 #include "libmesh/partitioner.h"
33 #include "libmesh/remote_elem.h"
34 #include "libmesh/unstructured_mesh.h"
35 
36 namespace libMesh
37 {
38 
39 
40 
41 //------------------------------------------------------
42 // BoundaryInfo static member initializations
44 
45 
46 
47 //------------------------------------------------------
48 // BoundaryInfo functions
50  ParallelObject(m.comm()),
51  _mesh (m)
52 {
53 }
54 
55 BoundaryInfo & BoundaryInfo::operator=(const BoundaryInfo & other_boundary_info)
56 {
57  // Overwrite any preexisting boundary info
58  this->clear();
59 
70  // Copy node boundary info
71  {
72  boundary_node_iter it = other_boundary_info._boundary_node_id.begin();
73  const boundary_node_iter end = other_boundary_info._boundary_node_id.end();
74 
75  for (; it != end; ++it)
76  {
77  const Node * other_node = it->first;
78  _boundary_node_id.insert(std::make_pair(_mesh.node_ptr(other_node->id()),
79  it->second));
80  }
81  }
82 
83  // Copy edge boundary info
84  {
85  boundary_edge_iter it = other_boundary_info._boundary_edge_id.begin();
86  const boundary_edge_iter end = other_boundary_info._boundary_edge_id.end();
87 
88  for (; it != end; ++it)
89  {
90  const Elem * other_elem = it->first;
91  _boundary_edge_id.insert(std::make_pair(_mesh.elem_ptr(other_elem->id()),
92  it->second));
93  }
94  }
95 
96  // Copy shellface boundary info
97  {
98  boundary_shellface_iter it = other_boundary_info._boundary_shellface_id.begin();
99  const boundary_shellface_iter end = other_boundary_info._boundary_shellface_id.end();
100 
101  for (; it != end; ++it)
102  {
103  const Elem * other_elem = it->first;
104  _boundary_shellface_id.insert(std::make_pair(_mesh.elem_ptr(other_elem->id()),
105  it->second));
106  }
107  }
108 
109  // Copy side boundary info
110  {
111  boundary_side_iter it = other_boundary_info._boundary_side_id.begin();
112  const boundary_side_iter end = other_boundary_info._boundary_side_id.end();
113 
114  for (; it != end; ++it)
115  {
116  const Elem * other_elem = it->first;
117  _boundary_side_id.insert(std::make_pair(_mesh.elem_ptr(other_elem->id()),
118  it->second));
119  }
120  }
121 
122  _boundary_ids = other_boundary_info._boundary_ids;
123  _side_boundary_ids = other_boundary_info._side_boundary_ids;
124  _node_boundary_ids = other_boundary_info._node_boundary_ids;
125  _edge_boundary_ids = other_boundary_info._edge_boundary_ids;
126  _shellface_boundary_ids = other_boundary_info._shellface_boundary_ids;
127 
128  return *this;
129 }
130 
131 
133 {
134  this->clear();
135 }
136 
137 
138 
140 {
141  _boundary_node_id.clear();
142  _boundary_side_id.clear();
143  _boundary_edge_id.clear();
144  _boundary_shellface_id.clear();
145  _boundary_ids.clear();
146  _side_boundary_ids.clear();
147  _node_boundary_ids.clear();
148  _edge_boundary_ids.clear();
149  _shellface_boundary_ids.clear();
150 }
151 
152 
153 
155 {
156  // Clear the old caches
157  _boundary_ids.clear();
158  _side_boundary_ids.clear();
159  _node_boundary_ids.clear();
160  _edge_boundary_ids.clear();
161  _shellface_boundary_ids.clear();
162 
163  // Loop over id maps to regenerate each set.
164  for (boundary_node_iter it = _boundary_node_id.begin(),
165  end = _boundary_node_id.end();
166  it != end; ++it)
167  {
168  const boundary_id_type id = it->second;
169  _boundary_ids.insert(id);
170  _node_boundary_ids.insert(id);
171  }
172 
173  for (boundary_edge_iter it = _boundary_edge_id.begin(),
174  end = _boundary_edge_id.end();
175  it != end; ++it)
176  {
177  const boundary_id_type id = it->second.second;
178  _boundary_ids.insert(id);
179  _edge_boundary_ids.insert(id);
180  }
181 
182  for (boundary_side_iter it = _boundary_side_id.begin(),
183  end = _boundary_side_id.end();
184  it != end; ++it)
185  {
186  const boundary_id_type id = it->second.second;
187  _boundary_ids.insert(id);
188  _side_boundary_ids.insert(id);
189  }
190 
192  end = _boundary_shellface_id.end();
193  it != end; ++it)
194  {
195  const boundary_id_type id = it->second.second;
196  _boundary_ids.insert(id);
197  _shellface_boundary_ids.insert(id);
198  }
199 }
200 
201 
202 
203 void BoundaryInfo::sync (UnstructuredMesh & boundary_mesh)
204 {
205  std::set<boundary_id_type> request_boundary_ids(_boundary_ids);
206  request_boundary_ids.insert(invalid_id);
207  if (!_mesh.is_serial())
208  this->comm().set_union(request_boundary_ids);
209 
210  this->sync(request_boundary_ids,
211  boundary_mesh);
212 }
213 
214 
215 
216 void BoundaryInfo::sync (const std::set<boundary_id_type> & requested_boundary_ids,
217  UnstructuredMesh & boundary_mesh)
218 {
219  // Call the 3 argument version of this function with a dummy value for the third set.
220  std::set<subdomain_id_type> subdomains_relative_to;
221  subdomains_relative_to.insert(Elem::invalid_subdomain_id);
222 
223  this->sync(requested_boundary_ids,
224  boundary_mesh,
225  subdomains_relative_to);
226 }
227 
228 
229 
230 void BoundaryInfo::sync (const std::set<boundary_id_type> & requested_boundary_ids,
231  UnstructuredMesh & boundary_mesh,
232  const std::set<subdomain_id_type> & subdomains_relative_to)
233 {
234  LOG_SCOPE("sync()", "BoundaryInfo");
235 
236  boundary_mesh.clear();
237 
243  if (!_mesh.is_serial())
244  boundary_mesh.delete_remote_elements();
245 
252  MeshSerializer serializer
253  (const_cast<MeshBase &>(_mesh), boundary_mesh.is_serial());
254 
259  boundary_mesh.set_n_partitions() = _mesh.n_partitions();
260 
261  std::map<dof_id_type, dof_id_type> node_id_map;
262  std::map<std::pair<dof_id_type, unsigned char>, dof_id_type> side_id_map;
263 
264  this->_find_id_maps(requested_boundary_ids, 0, &node_id_map, 0, &side_id_map, subdomains_relative_to);
265 
266  // Let's add all the boundary nodes we found to the boundary mesh
267  for (const auto & node : _mesh.node_ptr_range())
268  {
269  dof_id_type node_id = node->id();
270  if (node_id_map.count(node_id))
271  {
272  boundary_mesh.add_point(*node, node_id_map[node_id], node->processor_id());
273 
274  // Copy over all the node's boundary IDs to boundary_mesh
275  std::vector<boundary_id_type> node_boundary_ids;
276  this->boundary_ids(node, node_boundary_ids);
277  for (std::size_t index=0; index<node_boundary_ids.size(); index++)
278  {
279  boundary_mesh.get_boundary_info().add_node(node_id_map[node_id],
280  node_boundary_ids[index]);
281  }
282  }
283  }
284 
285  // Let's add the elements
286  this->add_elements (requested_boundary_ids, boundary_mesh, subdomains_relative_to);
287 
288  // The new elements are currently using the interior mesh's nodes;
289  // we want them to use the boundary mesh's nodes instead.
290 
291  // This side's Node pointers still point to the nodes of the original mesh.
292  // We need to re-point them to the boundary mesh's nodes! Since we copied *ALL* of
293  // the original mesh's nodes over, we should be guaranteed to have the same ordering.
294  for (auto & new_elem : boundary_mesh.element_ptr_range())
295  {
296  for (auto nn : new_elem->node_index_range())
297  {
298  // Get the correct node pointer, based on the id()
299  Node * new_node =
300  boundary_mesh.node_ptr(node_id_map[new_elem->node_id(nn)]);
301 
302  // sanity check: be sure that the new Node exists and its
303  // global id really matches
304  libmesh_assert (new_node);
305  libmesh_assert_equal_to (new_node->id(),
306  node_id_map[new_elem->node_id(nn)]);
307 
308  // Assign the new node pointer
309  new_elem->set_node(nn) = new_node;
310  }
311  }
312 
313  // Don't repartition this mesh; we want it to stay in sync with the
314  // interior partitioning.
315  boundary_mesh.partitioner().reset(libmesh_nullptr);
316 
317  // Make boundary_mesh nodes and elements contiguous
318  boundary_mesh.prepare_for_use(/*skip_renumber =*/ false);
319 
320  // and finally distribute element partitioning to the nodes
322 }
323 
324 
326  std::map<dof_id_type, dof_id_type> & node_id_map,
327  std::map<dof_id_type, unsigned char> & side_id_map,
328  Real tolerance)
329 {
330  LOG_SCOPE("get_side_and_node_maps()", "BoundaryInfo");
331 
332  node_id_map.clear();
333  side_id_map.clear();
334 
335  for (const auto & boundary_elem : boundary_mesh.active_element_ptr_range())
336  {
337  const Elem * interior_parent = boundary_elem->interior_parent();
338 
339  // Find out which side of interior_parent boundary_elem corresponds to.
340  // Use centroid comparison as a way to check.
341  unsigned char interior_parent_side_index = 0;
342  bool found_matching_sides = false;
343  for (auto side : interior_parent->side_index_range())
344  {
345  UniquePtr<const Elem> interior_parent_side = interior_parent->build_side_ptr(side);
346  Real centroid_distance = (boundary_elem->centroid() - interior_parent_side->centroid()).norm();
347 
348  if (centroid_distance < (tolerance * boundary_elem->hmin()))
349  {
350  interior_parent_side_index = side;
351  found_matching_sides = true;
352  break;
353  }
354  }
355 
356  if (!found_matching_sides)
357  libmesh_error_msg("No matching side found within the specified tolerance");
358 
359  side_id_map[boundary_elem->id()] = interior_parent_side_index;
360 
361  UniquePtr<const Elem> interior_parent_side = interior_parent->build_side_ptr(interior_parent_side_index);
362  for (auto local_node_index : boundary_elem->node_index_range())
363  {
364  dof_id_type boundary_node_id = boundary_elem->node_id(local_node_index);
365  dof_id_type interior_node_id = interior_parent_side->node_id(local_node_index);
366 
367  node_id_map[interior_node_id] = boundary_node_id;
368  }
369  }
370 }
371 
372 
373 
374 void BoundaryInfo::add_elements(const std::set<boundary_id_type> & requested_boundary_ids,
375  UnstructuredMesh & boundary_mesh)
376 {
377  // Call the 3 argument version of this function with a dummy value for the third arg.
378  std::set<subdomain_id_type> subdomains_relative_to;
379  subdomains_relative_to.insert(Elem::invalid_subdomain_id);
380 
381  this->add_elements(requested_boundary_ids,
382  boundary_mesh,
383  subdomains_relative_to);
384 }
385 
386 
387 
388 void BoundaryInfo::add_elements(const std::set<boundary_id_type> & requested_boundary_ids,
389  UnstructuredMesh & boundary_mesh,
390  const std::set<subdomain_id_type> & subdomains_relative_to)
391 {
392  LOG_SCOPE("add_elements()", "BoundaryInfo");
393 
394  // We're not prepared to mix serial and distributed meshes in this
395  // method, so make sure they match from the start.
396  libmesh_assert_equal_to(_mesh.is_serial(),
397  boundary_mesh.is_serial());
398 
399  std::map<std::pair<dof_id_type, unsigned char>, dof_id_type> side_id_map;
400  this->_find_id_maps(requested_boundary_ids,
401  0,
403  boundary_mesh.max_elem_id(),
404  &side_id_map,
405  subdomains_relative_to);
406 
407  // We have to add sides *outside* any element loop, because if
408  // boundary_mesh and _mesh are the same then those additions can
409  // invalidate our element iterators. So we just use the element
410  // loop to make a list of sides to add.
411  typedef std::vector<std::pair<dof_id_type, unsigned char>>
412  side_container;
413  side_container sides_to_add;
414 
415  for (const auto & elem : _mesh.element_ptr_range())
416  {
417  // If the subdomains_relative_to container has the
418  // invalid_subdomain_id, we fall back on the "old" behavior of
419  // adding sides regardless of this Elem's subdomain. Otherwise,
420  // if the subdomains_relative_to container doesn't contain the
421  // current Elem's subdomain_id(), we won't add any sides from
422  // it.
423  if (!subdomains_relative_to.count(Elem::invalid_subdomain_id) &&
424  !subdomains_relative_to.count(elem->subdomain_id()))
425  continue;
426 
427  // Get the top-level parent for this element
428  const Elem * top_parent = elem->top_parent();
429 
430  // Find all the boundary side ids for this Elem.
431  const std::pair<boundary_side_iter, boundary_side_iter>
432  bounds = _boundary_side_id.equal_range(top_parent);
433 
434  for (auto s : elem->side_index_range())
435  {
436  bool add_this_side = false;
437  boundary_id_type this_bcid = invalid_id;
438 
439  for (const auto & pr : as_range(bounds))
440  {
441  this_bcid = pr.second.second;
442 
443  // if this side is flagged with a boundary condition
444  // and the user wants this id
445  if ((pr.second.first == s) &&
446  (requested_boundary_ids.count(this_bcid)))
447  {
448  add_this_side = true;
449  break;
450  }
451  }
452 
453  // We may still want to add this side if the user called
454  // sync() with no requested_boundary_ids. This corresponds
455  // to the "old" style of calling sync() in which the entire
456  // boundary was copied to the BoundaryMesh, and handles the
457  // case where elements on the geometric boundary are not in
458  // any sidesets.
459  if (bounds.first == bounds.second &&
460  requested_boundary_ids.count(invalid_id) &&
461  elem->neighbor_ptr(s) == libmesh_nullptr)
462  add_this_side = true;
463 
464  if (add_this_side)
465  sides_to_add.push_back(std::make_pair(elem->id(), s));
466  }
467  }
468 
469 #ifdef LIBMESH_ENABLE_UNIQUE_ID
470  unique_id_type old_max_unique_id = boundary_mesh.parallel_max_unique_id();
471 #endif
472 
473  for (side_container::const_iterator it = sides_to_add.begin();
474  it != sides_to_add.end(); ++it)
475  {
476  const dof_id_type elem_id = it->first;
477  const unsigned char s = it->second;
478  Elem * elem = _mesh.elem_ptr(elem_id);
479 
480  // Build the side - do not use a "proxy" element here:
481  // This will be going into the boundary_mesh and needs to
482  // stand on its own.
483  UniquePtr<Elem> side (elem->build_side_ptr(s, false));
484 
485  side->processor_id() = elem->processor_id();
486 
487  const std::pair<dof_id_type, unsigned char> side_pair(elem_id, s);
488 
489  libmesh_assert(side_id_map.count(side_pair));
490 
491  const dof_id_type new_side_id = side_id_map[side_pair];
492 
493  side->set_id(new_side_id);
494 
495 #ifdef LIBMESH_ENABLE_UNIQUE_ID
496  side->set_unique_id() = old_max_unique_id + new_side_id;
497 #endif
498 
499  // Add the side
500  Elem * new_elem = boundary_mesh.add_elem(side.release());
501 
502 #ifdef LIBMESH_ENABLE_AMR
503  // Set parent links
504  if (elem->parent())
505  {
506  const std::pair<dof_id_type, unsigned char> parent_side_pair(elem->parent()->id(), s);
507 
508  libmesh_assert(side_id_map.count(parent_side_pair));
509 
510  Elem * side_parent = boundary_mesh.elem_ptr(side_id_map[parent_side_pair]);
511 
512  libmesh_assert(side_parent);
513 
514  new_elem->set_parent(side_parent);
515 
516  side_parent->set_refinement_flag(Elem::INACTIVE);
517 
518  // Figuring out which child we are of our parent
519  // is a trick. Due to libMesh child numbering
520  // conventions, if we are an element on a vertex,
521  // then we share that vertex with our parent, with
522  // the same local index.
523  bool found_child = false;
524  for (unsigned int v=0; v != new_elem->n_vertices(); ++v)
525  if (new_elem->node_ptr(v) == side_parent->node_ptr(v))
526  {
527  side_parent->add_child(new_elem, v);
528  found_child = true;
529  }
530 
531  // If we don't share any vertex with our parent,
532  // then we're the fourth child (index 3) of a
533  // triangle.
534  if (!found_child)
535  {
536  libmesh_assert_equal_to (new_elem->n_vertices(), 3);
537  side_parent->add_child(new_elem, 3);
538  }
539  }
540 #endif
541 
542  new_elem->set_interior_parent (const_cast<Elem *>(elem));
543 
544  // On non-local elements on DistributedMesh we might have
545  // RemoteElem neighbor links to construct
546  if (!_mesh.is_serial() &&
547  (elem->processor_id() != this->processor_id()))
548  {
549  const unsigned short n_nodes = elem->n_nodes();
550 
551  const unsigned short bdy_n_sides = new_elem->n_sides();
552  const unsigned short bdy_n_nodes = new_elem->n_nodes();
553 
554  // Check every interior side for a RemoteElem
555  for (auto interior_side : elem->side_index_range())
556  {
557  // Might this interior side have a RemoteElem that
558  // needs a corresponding Remote on a boundary side?
559  if (elem->neighbor_ptr(interior_side) != remote_elem)
560  continue;
561 
562  // Which boundary side?
563  for (unsigned short boundary_side = 0;
564  boundary_side != bdy_n_sides; ++boundary_side)
565  {
566  // Look for matching node points. This is safe in
567  // *this* context.
568  bool found_all_nodes = true;
569  for (unsigned short boundary_node = 0;
570  boundary_node != bdy_n_nodes; ++boundary_node)
571  {
572  if (!new_elem->is_node_on_side(boundary_node,
573  boundary_side))
574  continue;
575 
576  bool found_this_node = false;
577  for (unsigned short interior_node = 0;
578  interior_node != n_nodes; ++interior_node)
579  {
580  if (!elem->is_node_on_side(interior_node,
581  interior_side))
582  continue;
583 
584  if (new_elem->point(boundary_node) ==
585  elem->point(interior_node))
586  {
587  found_this_node = true;
588  break;
589  }
590  }
591  if (!found_this_node)
592  {
593  found_all_nodes = false;
594  break;
595  }
596  }
597 
598  if (found_all_nodes)
599  {
600  new_elem->set_neighbor
601  (boundary_side,
602  const_cast<RemoteElem *>(remote_elem));
603  break;
604  }
605  }
606  }
607  }
608  }
609 
610  // We haven't been bothering to keep unique ids consistent on ghost
611  // elements
612  if (!boundary_mesh.is_serial())
614 
615  // Make sure we didn't add ids inconsistently
616 #ifdef DEBUG
617 # ifdef LIBMESH_HAVE_RTTI
618  DistributedMesh * parmesh = dynamic_cast<DistributedMesh *>(&boundary_mesh);
619  if (parmesh)
621 # endif
622 #endif
623 }
624 
625 
626 
628  const boundary_id_type id)
629 {
630  this->add_node (_mesh.node_ptr(node), id);
631 }
632 
633 
634 
635 void BoundaryInfo::add_node(const Node * node,
636  const boundary_id_type id)
637 {
638  if (id == invalid_id)
639  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
640  << invalid_id \
641  << "\n That is reserved for internal use.");
642 
643  // Don't add the same ID twice
644  for (const auto & pr : as_range(_boundary_node_id.equal_range(node)))
645  if (pr.second == id)
646  return;
647 
648  _boundary_node_id.insert(std::make_pair(node, id));
649  _boundary_ids.insert(id);
650  _node_boundary_ids.insert(id); // Also add this ID to the set of node boundary IDs
651 }
652 
653 
654 
655 void BoundaryInfo::add_node(const Node * node,
656  const std::vector<boundary_id_type> & ids)
657 {
658  if (ids.empty())
659  return;
660 
661  libmesh_assert(node);
662 
663  // Don't add the same ID twice
664  std::pair<boundary_node_iter, boundary_node_iter> bounds = _boundary_node_id.equal_range(node);
665 
666  // The entries in the ids vector may be non-unique. If we expected
667  // *lots* of ids, it might be fastest to construct a std::set from
668  // the entries, but for a small number of entries, which is more
669  // typical, it is probably faster to copy the vector and do sort+unique.
670  // http://stackoverflow.com/questions/1041620/whats-the-most-efficient-way-to-erase-duplicates-and-sort-a-vector
671  std::vector<boundary_id_type> unique_ids(ids.begin(), ids.end());
672  std::sort(unique_ids.begin(), unique_ids.end());
673  std::vector<boundary_id_type>::iterator new_end =
674  std::unique(unique_ids.begin(), unique_ids.end());
675 
676  std::vector<boundary_id_type>::iterator it = unique_ids.begin();
677  for (; it != new_end; ++it)
678  {
679  boundary_id_type id = *it;
680 
681  if (id == invalid_id)
682  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
683  << invalid_id \
684  << "\n That is reserved for internal use.");
685 
686  bool already_inserted = false;
687  for (const auto & pr : as_range(bounds))
688  if (pr.second == id)
689  {
690  already_inserted = true;
691  break;
692  }
693  if (already_inserted)
694  continue;
695 
696  _boundary_node_id.insert(std::make_pair(node,id));
697  _boundary_ids.insert(id);
698  _node_boundary_ids.insert(id); // Also add this ID to the set of node boundary IDs
699  }
700 }
701 
702 
703 
705 {
706  _boundary_node_id.clear();
707 }
708 
710  const unsigned short int edge,
711  const boundary_id_type id)
712 {
713  this->add_edge (_mesh.elem_ptr(e), edge, id);
714 }
715 
716 
717 
718 void BoundaryInfo::add_edge(const Elem * elem,
719  const unsigned short int edge,
720  const boundary_id_type id)
721 {
722  libmesh_assert(elem);
723 
724  // Only add BCs for level-0 elements.
725  libmesh_assert_equal_to (elem->level(), 0);
726 
727  if (id == invalid_id)
728  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
729  << invalid_id \
730  << "\n That is reserved for internal use.");
731 
732  // Don't add the same ID twice
733  for (const auto & pr : as_range(_boundary_edge_id.equal_range(elem)))
734  if (pr.second.first == edge &&
735  pr.second.second == id)
736  return;
737 
738  _boundary_edge_id.insert(std::make_pair(elem, std::make_pair(edge, id)));
739  _boundary_ids.insert(id);
740  _edge_boundary_ids.insert(id); // Also add this ID to the set of edge boundary IDs
741 }
742 
743 
744 
745 void BoundaryInfo::add_edge(const Elem * elem,
746  const unsigned short int edge,
747  const std::vector<boundary_id_type> & ids)
748 {
749  if (ids.empty())
750  return;
751 
752  libmesh_assert(elem);
753 
754  // Only add BCs for level-0 elements.
755  libmesh_assert_equal_to (elem->level(), 0);
756 
757  // Don't add the same ID twice
758  std::pair<boundary_edge_iter, boundary_edge_iter> bounds = _boundary_edge_id.equal_range(elem);
759 
760  // The entries in the ids vector may be non-unique. If we expected
761  // *lots* of ids, it might be fastest to construct a std::set from
762  // the entries, but for a small number of entries, which is more
763  // typical, it is probably faster to copy the vector and do sort+unique.
764  // http://stackoverflow.com/questions/1041620/whats-the-most-efficient-way-to-erase-duplicates-and-sort-a-vector
765  std::vector<boundary_id_type> unique_ids(ids.begin(), ids.end());
766  std::sort(unique_ids.begin(), unique_ids.end());
767  std::vector<boundary_id_type>::iterator new_end =
768  std::unique(unique_ids.begin(), unique_ids.end());
769 
770  std::vector<boundary_id_type>::iterator it = unique_ids.begin();
771  for (; it != new_end; ++it)
772  {
773  boundary_id_type id = *it;
774 
775  if (id == invalid_id)
776  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
777  << invalid_id \
778  << "\n That is reserved for internal use.");
779 
780  bool already_inserted = false;
781  for (const auto & pr : as_range(bounds))
782  if (pr.second.first == edge &&
783  pr.second.second == id)
784  {
785  already_inserted = true;
786  break;
787  }
788  if (already_inserted)
789  continue;
790 
791  _boundary_edge_id.insert(std::make_pair(elem, std::make_pair(edge, id)));
792  _boundary_ids.insert(id);
793  _edge_boundary_ids.insert(id); // Also add this ID to the set of edge boundary IDs
794  }
795 }
796 
797 
798 
800  const unsigned short int shellface,
801  const boundary_id_type id)
802 {
803  this->add_shellface (_mesh.elem_ptr(e), shellface, id);
804 }
805 
806 
807 
809  const unsigned short int shellface,
810  const boundary_id_type id)
811 {
812  libmesh_assert(elem);
813 
814  // Only add BCs for level-0 elements.
815  libmesh_assert_equal_to (elem->level(), 0);
816 
817  // Shells only have 2 faces
818  libmesh_assert_less(shellface, 2);
819 
820  if (id == invalid_id)
821  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
822  << invalid_id \
823  << "\n That is reserved for internal use.");
824 
825  // Don't add the same ID twice
826  for (const auto & pr : as_range(_boundary_shellface_id.equal_range(elem)))
827  if (pr.second.first == shellface &&
828  pr.second.second == id)
829  return;
830 
831  _boundary_shellface_id.insert(std::make_pair(elem, std::make_pair(shellface, id)));
832  _boundary_ids.insert(id);
833  _shellface_boundary_ids.insert(id); // Also add this ID to the set of shellface boundary IDs
834 }
835 
836 
837 
839  const unsigned short int shellface,
840  const std::vector<boundary_id_type> & ids)
841 {
842  if (ids.empty())
843  return;
844 
845  libmesh_assert(elem);
846 
847  // Only add BCs for level-0 elements.
848  libmesh_assert_equal_to (elem->level(), 0);
849 
850  // Shells only have 2 faces
851  libmesh_assert_less(shellface, 2);
852 
853  // Don't add the same ID twice
854  std::pair<boundary_shellface_iter, boundary_shellface_iter> bounds = _boundary_shellface_id.equal_range(elem);
855 
856  // The entries in the ids vector may be non-unique. If we expected
857  // *lots* of ids, it might be fastest to construct a std::set from
858  // the entries, but for a small number of entries, which is more
859  // typical, it is probably faster to copy the vector and do sort+unique.
860  // http://stackoverflow.com/questions/1041620/whats-the-most-efficient-way-to-erase-duplicates-and-sort-a-vector
861  std::vector<boundary_id_type> unique_ids(ids.begin(), ids.end());
862  std::sort(unique_ids.begin(), unique_ids.end());
863  std::vector<boundary_id_type>::iterator new_end =
864  std::unique(unique_ids.begin(), unique_ids.end());
865 
866  std::vector<boundary_id_type>::iterator it = unique_ids.begin();
867  for (; it != new_end; ++it)
868  {
869  boundary_id_type id = *it;
870 
871  if (id == invalid_id)
872  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
873  << invalid_id \
874  << "\n That is reserved for internal use.");
875 
876  bool already_inserted = false;
877  for (const auto & pr : as_range(bounds))
878  if (pr.second.first == shellface &&
879  pr.second.second == id)
880  {
881  already_inserted = true;
882  break;
883  }
884  if (already_inserted)
885  continue;
886 
887  _boundary_shellface_id.insert(std::make_pair(elem, std::make_pair(shellface, id)));
888  _boundary_ids.insert(id);
889  _shellface_boundary_ids.insert(id); // Also add this ID to the set of shellface boundary IDs
890  }
891 }
892 
893 
895  const unsigned short int side,
896  const boundary_id_type id)
897 {
898  this->add_side (_mesh.elem_ptr(e), side, id);
899 }
900 
901 
902 
903 void BoundaryInfo::add_side(const Elem * elem,
904  const unsigned short int side,
905  const boundary_id_type id)
906 {
907  libmesh_assert(elem);
908 
909  // Only add BCs for level-0 elements.
910  libmesh_assert_equal_to (elem->level(), 0);
911 
912  if (id == invalid_id)
913  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
914  << invalid_id \
915  << "\n That is reserved for internal use.");
916 
917  // Don't add the same ID twice
918  for (const auto & pr : as_range(_boundary_side_id.equal_range(elem)))
919  if (pr.second.first == side &&
920  pr.second.second == id)
921  return;
922 
923  _boundary_side_id.insert(std::make_pair(elem, std::make_pair(side, id)));
924  _boundary_ids.insert(id);
925  _side_boundary_ids.insert(id); // Also add this ID to the set of side boundary IDs
926 }
927 
928 
929 
930 void BoundaryInfo::add_side(const Elem * elem,
931  const unsigned short int side,
932  const std::vector<boundary_id_type> & ids)
933 {
934  if (ids.empty())
935  return;
936 
937  libmesh_assert(elem);
938 
939  // Only add BCs for level-0 elements.
940  libmesh_assert_equal_to (elem->level(), 0);
941 
942  // Don't add the same ID twice
943  std::pair<boundary_side_iter, boundary_side_iter> bounds = _boundary_side_id.equal_range(elem);
944 
945  // The entries in the ids vector may be non-unique. If we expected
946  // *lots* of ids, it might be fastest to construct a std::set from
947  // the entries, but for a small number of entries, which is more
948  // typical, it is probably faster to copy the vector and do sort+unique.
949  // http://stackoverflow.com/questions/1041620/whats-the-most-efficient-way-to-erase-duplicates-and-sort-a-vector
950  std::vector<boundary_id_type> unique_ids(ids.begin(), ids.end());
951  std::sort(unique_ids.begin(), unique_ids.end());
952  std::vector<boundary_id_type>::iterator new_end =
953  std::unique(unique_ids.begin(), unique_ids.end());
954 
955  std::vector<boundary_id_type>::const_iterator it = unique_ids.begin();
956  for (; it != new_end; ++it)
957  {
958  boundary_id_type id = *it;
959 
960  if (id == invalid_id)
961  libmesh_error_msg("ERROR: You may not set a boundary ID of " \
962  << invalid_id \
963  << "\n That is reserved for internal use.");
964 
965  bool already_inserted = false;
966  for (const auto & pr : as_range(bounds))
967  if (pr.second.first == side && pr.second.second == id)
968  {
969  already_inserted = true;
970  break;
971  }
972  if (already_inserted)
973  continue;
974 
975  _boundary_side_id.insert(std::make_pair(elem, std::make_pair(side, id)));
976  _boundary_ids.insert(id);
977  _side_boundary_ids.insert(id); // Also add this ID to the set of side boundary IDs
978  }
979 }
980 
981 
982 
983 bool BoundaryInfo::has_boundary_id(const Node * const node,
984  const boundary_id_type id) const
985 {
986  for (const auto & pr : as_range(_boundary_node_id.equal_range(node)))
987  if (pr.second == id)
988  return true;
989 
990  return false;
991 }
992 
993 
994 
995 #ifdef LIBMESH_ENABLE_DEPRECATED
996 std::vector<boundary_id_type> BoundaryInfo::boundary_ids(const Node * node) const
997 {
998  libmesh_deprecated();
999 
1000  std::vector<boundary_id_type> ids;
1001  this->boundary_ids(node, ids);
1002  return ids;
1003 }
1004 #endif
1005 
1006 
1007 
1009  std::vector<boundary_id_type> & vec_to_fill) const
1010 {
1011  // Clear out any previous contents
1012  vec_to_fill.clear();
1013 
1014  for (const auto & pr : as_range(_boundary_node_id.equal_range(node)))
1015  vec_to_fill.push_back(pr.second);
1016 }
1017 
1018 
1019 
1020 unsigned int BoundaryInfo::n_boundary_ids(const Node * node) const
1021 {
1022  std::pair<boundary_node_iter, boundary_node_iter> pos = _boundary_node_id.equal_range(node);
1023 
1024  return cast_int<unsigned int>(std::distance(pos.first, pos.second));
1025 }
1026 
1027 
1028 
1029 #ifdef LIBMESH_ENABLE_DEPRECATED
1030 std::vector<boundary_id_type> BoundaryInfo::edge_boundary_ids (const Elem * const elem,
1031  const unsigned short int edge) const
1032 {
1033  libmesh_deprecated();
1034 
1035  std::vector<boundary_id_type> ids;
1036  this->edge_boundary_ids(elem, edge, ids);
1037  return ids;
1038 }
1039 #endif
1040 
1041 
1042 
1043 void BoundaryInfo::edge_boundary_ids (const Elem * const elem,
1044  const unsigned short int edge,
1045  std::vector<boundary_id_type> & vec_to_fill) const
1046 {
1047  libmesh_assert(elem);
1048 
1049  // Clear out any previous contents
1050  vec_to_fill.clear();
1051 
1052  // Only level-0 elements store BCs. If this is not a level-0
1053  // element get its level-0 parent and infer the BCs.
1054  const Elem * searched_elem = elem;
1055 #ifdef LIBMESH_ENABLE_AMR
1056  if (elem->level() != 0)
1057  {
1058  // Find all the sides that contain edge. If one of those is a boundary
1059  // side, then this must be a boundary edge. In that case, we just use the
1060  // top-level parent.
1061  bool found_boundary_edge = false;
1062  for (auto side : elem->side_index_range())
1063  {
1064  if (elem->is_edge_on_side(edge,side))
1065  {
1066  if (elem->neighbor_ptr(side) == libmesh_nullptr)
1067  {
1068  searched_elem = elem->top_parent ();
1069  found_boundary_edge = true;
1070  break;
1071  }
1072  }
1073  }
1074 
1075  if (!found_boundary_edge)
1076  {
1077  // Child element is not on external edge, but it may have internal
1078  // "boundary" IDs. We will walk up the tree, at each level checking that
1079  // the current child is actually on the same edge of the parent that is
1080  // currently being searched for (i.e. that was passed in as "edge").
1081  while (searched_elem->parent() != libmesh_nullptr)
1082  {
1083  const Elem * parent = searched_elem->parent();
1084  if (parent->is_child_on_edge(parent->which_child_am_i(searched_elem), edge) == false)
1085  return;
1086  searched_elem = parent;
1087  }
1088  }
1089  }
1090 #endif
1091 
1092  // Check each element in the range to see if its edge matches the requested edge.
1093  for (const auto & pr : as_range(_boundary_edge_id.equal_range(searched_elem)))
1094  if (pr.second.first == edge)
1095  vec_to_fill.push_back(pr.second.second);
1096 }
1097 
1098 
1099 
1100 unsigned int BoundaryInfo::n_edge_boundary_ids (const Elem * const elem,
1101  const unsigned short int edge) const
1102 {
1103  std::vector<boundary_id_type> ids;
1104  this->edge_boundary_ids(elem, edge, ids);
1105  return ids.size();
1106 }
1107 
1108 
1109 
1110 #ifdef LIBMESH_ENABLE_DEPRECATED
1111 std::vector<boundary_id_type> BoundaryInfo::raw_edge_boundary_ids (const Elem * const elem,
1112  const unsigned short int edge) const
1113 {
1114  libmesh_deprecated();
1115 
1116  std::vector<boundary_id_type> ids;
1117  this->raw_edge_boundary_ids(elem, edge, ids);
1118  return ids;
1119 }
1120 #endif
1121 
1122 
1123 
1125  const unsigned short int edge,
1126  std::vector<boundary_id_type> & vec_to_fill) const
1127 {
1128  libmesh_assert(elem);
1129 
1130  // Clear out any previous contents
1131  vec_to_fill.clear();
1132 
1133  // Only level-0 elements store BCs.
1134  if (elem->parent())
1135  return;
1136 
1137  // Check each element in the range to see if its edge matches the requested edge.
1138  for (const auto & pr : as_range(_boundary_edge_id.equal_range(elem)))
1139  if (pr.second.first == edge)
1140  vec_to_fill.push_back(pr.second.second);
1141 }
1142 
1143 
1144 
1146  const unsigned short int shellface,
1147  std::vector<boundary_id_type> & vec_to_fill) const
1148 {
1149  libmesh_assert(elem);
1150 
1151  // Shells only have 2 faces
1152  libmesh_assert_less(shellface, 2);
1153 
1154  // Clear out any previous contents
1155  vec_to_fill.clear();
1156 
1157  // Only level-0 elements store BCs. If this is not a level-0
1158  // element get its level-0 parent and infer the BCs.
1159  const Elem * searched_elem = elem;
1160 #ifdef LIBMESH_ENABLE_AMR
1161  if (elem->level() != 0)
1162  {
1163  while (searched_elem->parent() != libmesh_nullptr)
1164  {
1165  const Elem * parent = searched_elem->parent();
1166  searched_elem = parent;
1167  }
1168  }
1169 #endif
1170 
1171  // Check each element in the range to see if its shellface matches the requested shellface.
1172  for (const auto & pr : as_range(_boundary_shellface_id.equal_range(searched_elem)))
1173  if (pr.second.first == shellface)
1174  vec_to_fill.push_back(pr.second.second);
1175 }
1176 
1177 
1178 
1179 unsigned int BoundaryInfo::n_shellface_boundary_ids (const Elem * const elem,
1180  const unsigned short int shellface) const
1181 {
1182  std::vector<boundary_id_type> ids;
1183  this->shellface_boundary_ids(elem, shellface, ids);
1184  return ids.size();
1185 }
1186 
1187 
1188 
1190  const unsigned short int shellface,
1191  std::vector<boundary_id_type> & vec_to_fill) const
1192 {
1193  libmesh_assert(elem);
1194 
1195  // Shells only have 2 faces
1196  libmesh_assert_less(shellface, 2);
1197 
1198  // Clear out any previous contents
1199  vec_to_fill.clear();
1200 
1201  // Only level-0 elements store BCs.
1202  if (elem->parent())
1203  return;
1204 
1205  // Check each element in the range to see if its shellface matches the requested shellface.
1206  for (const auto & pr : as_range(_boundary_shellface_id.equal_range(elem)))
1207  if (pr.second.first == shellface)
1208  vec_to_fill.push_back(pr.second.second);
1209 }
1210 
1211 
1212 #ifdef LIBMESH_ENABLE_DEPRECATED
1214  const unsigned short int side) const
1215 {
1216  libmesh_deprecated();
1217 
1218  std::vector<boundary_id_type> ids;
1219  this->boundary_ids(elem, side, ids);
1220 
1221  // If the set is empty, return invalid_id
1222  if (ids.empty())
1223  return invalid_id;
1224 
1225  // Otherwise, just return the first id we came across for this
1226  // element on this side.
1227  return *(ids.begin());
1228 }
1229 #endif
1230 
1231 
1232 
1233 bool BoundaryInfo::has_boundary_id(const Elem * const elem,
1234  const unsigned short int side,
1235  const boundary_id_type id) const
1236 {
1237  std::vector<boundary_id_type> ids;
1238  this->boundary_ids(elem, side, ids);
1239  return (std::find(ids.begin(), ids.end(), id) != ids.end());
1240 }
1241 
1242 
1243 
1244 #ifdef LIBMESH_ENABLE_DEPRECATED
1245 std::vector<boundary_id_type> BoundaryInfo::boundary_ids (const Elem * const elem,
1246  const unsigned short int side) const
1247 {
1248  libmesh_deprecated();
1249 
1250  std::vector<boundary_id_type> ids;
1251  this->boundary_ids(elem, side, ids);
1252  return ids;
1253 }
1254 #endif
1255 
1256 
1257 
1258 void BoundaryInfo::boundary_ids (const Elem * const elem,
1259  const unsigned short int side,
1260  std::vector<boundary_id_type> & vec_to_fill) const
1261 {
1262  libmesh_assert(elem);
1263 
1264  // Clear out any previous contents
1265  vec_to_fill.clear();
1266 
1267  // Only level-0 elements store BCs. If this is not a level-0
1268  // element get its level-0 parent and infer the BCs.
1269  const Elem * searched_elem = elem;
1270  if (elem->level() != 0)
1271  {
1272  if (elem->neighbor_ptr(side) == libmesh_nullptr)
1273  searched_elem = elem->top_parent ();
1274 #ifdef LIBMESH_ENABLE_AMR
1275  else
1276  while (searched_elem->parent() != libmesh_nullptr)
1277  {
1278  const Elem * parent = searched_elem->parent();
1279  if (parent->is_child_on_side(parent->which_child_am_i(searched_elem), side) == false)
1280  return;
1281  searched_elem = parent;
1282  }
1283 #endif
1284  }
1285 
1286  // Check each element in the range to see if its side matches the requested side.
1287  for (const auto & pr : as_range(_boundary_side_id.equal_range(searched_elem)))
1288  if (pr.second.first == side)
1289  vec_to_fill.push_back(pr.second.second);
1290 }
1291 
1292 
1293 
1294 
1295 unsigned int BoundaryInfo::n_boundary_ids (const Elem * const elem,
1296  const unsigned short int side) const
1297 {
1298  std::vector<boundary_id_type> ids;
1299  this->boundary_ids(elem, side, ids);
1300  return ids.size();
1301 }
1302 
1303 
1304 
1305 #ifdef LIBMESH_ENABLE_DEPRECATED
1306 std::vector<boundary_id_type> BoundaryInfo::raw_boundary_ids (const Elem * const elem,
1307  const unsigned short int side) const
1308 {
1309  libmesh_deprecated();
1310 
1311  std::vector<boundary_id_type> ids;
1312  this->raw_boundary_ids(elem, side, ids);
1313  return ids;
1314 }
1315 #endif
1316 
1317 
1318 
1319 void BoundaryInfo::raw_boundary_ids (const Elem * const elem,
1320  const unsigned short int side,
1321  std::vector<boundary_id_type> & vec_to_fill) const
1322 {
1323  libmesh_assert(elem);
1324 
1325  // Clear out any previous contents
1326  vec_to_fill.clear();
1327 
1328  // Only level-0 elements store BCs.
1329  if (elem->parent())
1330  return;
1331 
1332  // Check each element in the range to see if its side matches the requested side.
1333  for (const auto & pr : as_range(_boundary_side_id.equal_range(elem)))
1334  if (pr.second.first == side)
1335  vec_to_fill.push_back(pr.second.second);
1336 }
1337 
1338 
1339 
1340 void BoundaryInfo::copy_boundary_ids (const BoundaryInfo & old_boundary_info,
1341  const Elem * const old_elem,
1342  const Elem * const new_elem)
1343 {
1344  libmesh_assert_equal_to (old_elem->n_sides(), new_elem->n_sides());
1345  libmesh_assert_equal_to (old_elem->n_edges(), new_elem->n_edges());
1346 
1347  std::vector<boundary_id_type> bndry_ids;
1348 
1349  for (auto s : old_elem->side_index_range())
1350  {
1351  old_boundary_info.raw_boundary_ids (old_elem, s, bndry_ids);
1352  this->add_side (new_elem, s, bndry_ids);
1353  }
1354 
1355  for (auto e : old_elem->edge_index_range())
1356  {
1357  old_boundary_info.raw_edge_boundary_ids (old_elem, e, bndry_ids);
1358  this->add_edge (new_elem, e, bndry_ids);
1359  }
1360 
1361  for (unsigned short sf=0; sf != 2; sf++)
1362  {
1363  old_boundary_info.raw_shellface_boundary_ids (old_elem, sf, bndry_ids);
1364  this->add_shellface (new_elem, sf, bndry_ids);
1365  }
1366 }
1367 
1368 
1369 
1370 void BoundaryInfo::remove (const Node * node)
1371 {
1372  libmesh_assert(node);
1373 
1374  // Erase everything associated with node
1375  _boundary_node_id.erase (node);
1376 }
1377 
1378 
1379 
1380 void BoundaryInfo::remove (const Elem * elem)
1381 {
1382  libmesh_assert(elem);
1383 
1384  // Erase everything associated with elem
1385  _boundary_edge_id.erase (elem);
1386  _boundary_side_id.erase (elem);
1387  _boundary_shellface_id.erase (elem);
1388 }
1389 
1390 
1391 
1392 void BoundaryInfo::remove_edge (const Elem * elem,
1393  const unsigned short int edge)
1394 {
1395  libmesh_assert(elem);
1396 
1397  // The user shouldn't be trying to remove only one child's boundary
1398  // id
1399  libmesh_assert_equal_to (elem->level(), 0);
1400 
1401  // Some older compilers don't support erasing from a map with
1402  // const_iterators, so we explicitly use non-const iterators here.
1403  std::pair<erase_iter, erase_iter>
1404  e = _boundary_edge_id.equal_range(elem);
1405 
1406  // elem may be there, maybe multiple occurrences
1407  while (e.first != e.second)
1408  {
1409  // if this is true we found the requested edge
1410  // of the element and want to erase the id
1411  if (e.first->second.first == edge)
1412  {
1413  // (postfix++ - increment the iterator before it's invalid)
1414  _boundary_edge_id.erase(e.first++);
1415  }
1416  else
1417  ++e.first;
1418  }
1419 }
1420 
1421 
1422 
1423 void BoundaryInfo::remove_edge (const Elem * elem,
1424  const unsigned short int edge,
1425  const boundary_id_type id)
1426 {
1427  libmesh_assert(elem);
1428 
1429  // The user shouldn't be trying to remove only one child's boundary
1430  // id
1431  libmesh_assert_equal_to (elem->level(), 0);
1432 
1433  // Some older compilers don't support erasing from a map with
1434  // const_iterators, so we explicitly use non-const iterators here.
1435  std::pair<erase_iter, erase_iter>
1436  e = _boundary_edge_id.equal_range(elem);
1437 
1438  // elem may be there, maybe multiple occurrences
1439  while (e.first != e.second)
1440  {
1441  // if this is true we found the requested edge
1442  // of the element and want to erase the requested id
1443  if (e.first->second.first == edge &&
1444  e.first->second.second == id)
1445  {
1446  // (postfix++ - increment the iterator before it's invalid)
1447  _boundary_edge_id.erase(e.first++);
1448  }
1449  else
1450  ++e.first;
1451  }
1452 }
1453 
1454 
1456  const unsigned short int shellface)
1457 {
1458  libmesh_assert(elem);
1459 
1460  // The user shouldn't be trying to remove only one child's boundary
1461  // id
1462  libmesh_assert_equal_to (elem->level(), 0);
1463 
1464  // Shells only have 2 faces
1465  libmesh_assert_less(shellface, 2);
1466 
1467  // Some older compilers don't support erasing from a map with
1468  // const_iterators, so we explicitly use non-const iterators here.
1469  std::pair<erase_iter, erase_iter>
1470  e = _boundary_shellface_id.equal_range(elem);
1471 
1472  // elem may be there, maybe multiple occurrences
1473  while (e.first != e.second)
1474  {
1475  // if this is true we found the requested shellface
1476  // of the element and want to erase the id
1477  if (e.first->second.first == shellface)
1478  {
1479  // (postfix++ - increment the iterator before it's invalid)
1480  _boundary_shellface_id.erase(e.first++);
1481  }
1482  else
1483  ++e.first;
1484  }
1485 }
1486 
1487 
1488 
1490  const unsigned short int shellface,
1491  const boundary_id_type id)
1492 {
1493  libmesh_assert(elem);
1494 
1495  // The user shouldn't be trying to remove only one child's boundary
1496  // id
1497  libmesh_assert_equal_to (elem->level(), 0);
1498 
1499  // Shells only have 2 faces
1500  libmesh_assert_less(shellface, 2);
1501 
1502  // Some older compilers don't support erasing from a map with
1503  // const_iterators, so we explicitly use non-const iterators here.
1504  std::pair<erase_iter, erase_iter>
1505  e = _boundary_shellface_id.equal_range(elem);
1506 
1507  // elem may be there, maybe multiple occurrences
1508  while (e.first != e.second)
1509  {
1510  // if this is true we found the requested shellface
1511  // of the element and want to erase the requested id
1512  if (e.first->second.first == shellface &&
1513  e.first->second.second == id)
1514  {
1515  // (postfix++ - increment the iterator before it's invalid)
1516  _boundary_shellface_id.erase(e.first++);
1517  }
1518  else
1519  ++e.first;
1520  }
1521 }
1522 
1523 void BoundaryInfo::remove_side (const Elem * elem,
1524  const unsigned short int side)
1525 {
1526  libmesh_assert(elem);
1527 
1528  // The user shouldn't be trying to remove only one child's boundary
1529  // id
1530  libmesh_assert_equal_to (elem->level(), 0);
1531 
1532  // Some older compilers don't support erasing from a map with
1533  // const_iterators, so we explicitly use non-const iterators here.
1534  std::pair<erase_iter, erase_iter>
1535  e = _boundary_side_id.equal_range(elem);
1536 
1537  // elem may be there, maybe multiple occurrences
1538  while (e.first != e.second)
1539  {
1540  // if this is true we found the requested side
1541  // of the element and want to erase the id
1542  if (e.first->second.first == side)
1543  {
1544  // (postfix++ - increment the iterator before it's invalid)
1545  _boundary_side_id.erase(e.first++);
1546  }
1547  else
1548  ++e.first;
1549  }
1550 }
1551 
1552 
1553 
1554 void BoundaryInfo::remove_side (const Elem * elem,
1555  const unsigned short int side,
1556  const boundary_id_type id)
1557 {
1558  libmesh_assert(elem);
1559 
1560  // Some older compilers don't support erasing from a map with
1561  // const_iterators, so we explicitly use non-const iterators here.
1562  std::pair<erase_iter, erase_iter>
1563  e = _boundary_side_id.equal_range(elem);
1564 
1565  // elem may be there, maybe multiple occurrences
1566  while (e.first != e.second)
1567  {
1568  // if this is true we found the requested side
1569  // of the element and want to erase the requested id
1570  if (e.first->second.first == side &&
1571  e.first->second.second == id)
1572  {
1573  // (postfix++ - increment the iterator before it's invalid)
1574  _boundary_side_id.erase(e.first++);
1575  }
1576  else
1577  ++e.first;
1578  }
1579 }
1580 
1581 
1582 
1584 {
1585  // Erase id from ids containers
1586  _boundary_ids.erase(id);
1587  _side_boundary_ids.erase(id);
1588  _edge_boundary_ids.erase(id);
1589  _shellface_boundary_ids.erase(id);
1590  _node_boundary_ids.erase(id);
1591  _ss_id_to_name.erase(id);
1592  _ns_id_to_name.erase(id);
1593 
1594  // Erase pointers to geometric entities with this id.
1595  for (boundary_node_erase_iter it = _boundary_node_id.begin(); it != _boundary_node_id.end(); /*below*/)
1596  {
1597  if (it->second == id)
1598  _boundary_node_id.erase(it++);
1599  else
1600  ++it;
1601  }
1602 
1603  for (erase_iter it = _boundary_edge_id.begin(); it != _boundary_edge_id.end(); /*below*/)
1604  {
1605  if (it->second.second == id)
1606  _boundary_edge_id.erase(it++);
1607  else
1608  ++it;
1609  }
1610 
1611  for (erase_iter it = _boundary_shellface_id.begin(); it != _boundary_shellface_id.end(); /*below*/)
1612  {
1613  if (it->second.second == id)
1614  _boundary_shellface_id.erase(it++);
1615  else
1616  ++it;
1617  }
1618 
1619  for (erase_iter it = _boundary_side_id.begin(); it != _boundary_side_id.end(); /*below*/)
1620  {
1621  if (it->second.second == id)
1622  _boundary_side_id.erase(it++);
1623  else
1624  ++it;
1625  }
1626 }
1627 
1628 
1629 
1630 unsigned int BoundaryInfo::side_with_boundary_id(const Elem * const elem,
1631  const boundary_id_type boundary_id_in) const
1632 {
1633  const Elem * searched_elem = elem;
1634  if (elem->level() != 0)
1635  searched_elem = elem->top_parent();
1636 
1637  // elem may have zero or multiple occurrences
1638  for (const auto & pr : as_range(_boundary_side_id.equal_range(searched_elem)))
1639  {
1640  // if this is true we found the requested boundary_id
1641  // of the element and want to return the side
1642  if (pr.second.second == boundary_id_in)
1643  {
1644  unsigned int side = pr.second.first;
1645 
1646  // If we're on this external boundary then we share this
1647  // external boundary id
1648  if (elem->neighbor_ptr(side) == libmesh_nullptr)
1649  return side;
1650 
1651  // If we're on an internal boundary then we need to be sure
1652  // it's the same internal boundary as our top_parent
1653  const Elem * p = elem;
1654 
1655 #ifdef LIBMESH_ENABLE_AMR
1656 
1657  while (p != libmesh_nullptr)
1658  {
1659  const Elem * parent = p->parent();
1660  if (!parent->is_child_on_side(parent->which_child_am_i(p), side))
1661  break;
1662  p = parent;
1663  }
1664 #endif
1665  // We're on that side of our top_parent; return it
1666  if (!p)
1667  return side;
1668  }
1669  }
1670 
1671  // if we get here, we found elem in the data structure but not
1672  // the requested boundary id, so return the default value
1673  return libMesh::invalid_uint;
1674 }
1675 
1676 void
1677 BoundaryInfo::build_node_boundary_ids(std::vector<boundary_id_type> & b_ids) const
1678 {
1679  b_ids.clear();
1680 
1681  boundary_node_iter pos = _boundary_node_id.begin();
1682  for (; pos != _boundary_node_id.end(); ++pos)
1683  {
1684  boundary_id_type id = pos->second;
1685 
1686  if (std::find(b_ids.begin(),b_ids.end(),id) == b_ids.end())
1687  b_ids.push_back(id);
1688  }
1689 }
1690 
1691 void
1692 BoundaryInfo::build_side_boundary_ids(std::vector<boundary_id_type> & b_ids) const
1693 {
1694  b_ids.clear();
1695 
1696  boundary_side_iter pos = _boundary_side_id.begin();
1697  for (; pos != _boundary_side_id.end(); ++pos)
1698  {
1699  boundary_id_type id = pos->second.second;
1700 
1701  if (std::find(b_ids.begin(),b_ids.end(),id) == b_ids.end())
1702  b_ids.push_back(id);
1703  }
1704 }
1705 
1706 void
1707 BoundaryInfo::build_shellface_boundary_ids(std::vector<boundary_id_type> & b_ids) const
1708 {
1709  b_ids.clear();
1710 
1712  for (; pos != _boundary_shellface_id.end(); ++pos)
1713  {
1714  boundary_id_type id = pos->second.second;
1715 
1716  if (std::find(b_ids.begin(),b_ids.end(),id) == b_ids.end())
1717  b_ids.push_back(id);
1718  }
1719 }
1720 
1722 {
1723  // in serial we know the number of bcs from the
1724  // size of the container
1725  if (_mesh.is_serial())
1726  return _boundary_side_id.size();
1727 
1728  // in parallel we need to sum the number of local bcs
1729  parallel_object_only();
1730 
1731  std::size_t nbcs=0;
1732 
1733  boundary_side_iter pos = _boundary_side_id.begin();
1734  for (; pos != _boundary_side_id.end(); ++pos)
1735  if (pos->first->processor_id() == this->processor_id())
1736  nbcs++;
1737 
1738  this->comm().sum (nbcs);
1739 
1740  return nbcs;
1741 }
1742 
1743 std::size_t BoundaryInfo::n_edge_conds () const
1744 {
1745  // in serial we know the number of nodesets from the
1746  // size of the container
1747  if (_mesh.is_serial())
1748  return _boundary_edge_id.size();
1749 
1750  // in parallel we need to sum the number of local nodesets
1751  parallel_object_only();
1752 
1753  std::size_t n_edge_bcs=0;
1754 
1755  boundary_edge_iter pos = _boundary_edge_id.begin();
1756  for (; pos != _boundary_edge_id.end(); ++pos)
1757  if (pos->first->processor_id() == this->processor_id())
1758  n_edge_bcs++;
1759 
1760  this->comm().sum (n_edge_bcs);
1761 
1762  return n_edge_bcs;
1763 }
1764 
1765 
1767 {
1768  // in serial we know the number of nodesets from the
1769  // size of the container
1770  if (_mesh.is_serial())
1771  return _boundary_shellface_id.size();
1772 
1773  // in parallel we need to sum the number of local nodesets
1774  parallel_object_only();
1775 
1776  std::size_t n_shellface_bcs=0;
1777 
1779  for (; pos != _boundary_shellface_id.end(); ++pos)
1780  if (pos->first->processor_id() == this->processor_id())
1781  n_shellface_bcs++;
1782 
1783  this->comm().sum (n_shellface_bcs);
1784 
1785  return n_shellface_bcs;
1786 }
1787 
1788 
1789 std::size_t BoundaryInfo::n_nodeset_conds () const
1790 {
1791  // in serial we know the number of nodesets from the
1792  // size of the container
1793  if (_mesh.is_serial())
1794  return _boundary_node_id.size();
1795 
1796  // in parallel we need to sum the number of local nodesets
1797  parallel_object_only();
1798 
1799  std::size_t n_nodesets=0;
1800 
1801  boundary_node_iter pos = _boundary_node_id.begin();
1802  for (; pos != _boundary_node_id.end(); ++pos)
1803  if (pos->first->processor_id() == this->processor_id())
1804  n_nodesets++;
1805 
1806  this->comm().sum (n_nodesets);
1807 
1808  return n_nodesets;
1809 }
1810 
1811 
1812 
1813 void BoundaryInfo::build_node_list (std::vector<dof_id_type> & nl,
1814  std::vector<boundary_id_type> & il) const
1815 {
1816  // Clear the input vectors, just in case they were used for
1817  // something else recently...
1818  nl.clear();
1819  il.clear();
1820 
1821  // Reserve the size, then use push_back
1822  nl.reserve (_boundary_node_id.size());
1823  il.reserve (_boundary_node_id.size());
1824 
1825  boundary_node_iter pos = _boundary_node_id.begin();
1826  for (; pos != _boundary_node_id.end(); ++pos)
1827  {
1828  nl.push_back (pos->first->id());
1829  il.push_back (pos->second);
1830  }
1831 }
1832 
1833 
1834 void
1836 {
1837  // If we're on a distributed mesh, even the owner of a node is not
1838  // guaranteed to be able to properly assign its new boundary id(s)!
1839  // Nodal neighbors are not always ghosted, and a nodal neighbor
1840  // might have a boundary side.
1841  const bool mesh_is_serial = _mesh.is_serial();
1842 
1843  typedef std::set<std::pair<dof_id_type, boundary_id_type>> set_type;
1844 
1845  const processor_id_type n_proc = this->n_processors();
1846  const processor_id_type my_proc_id = this->processor_id();
1847  std::vector<set_type> nodes_to_push(n_proc);
1848 
1849  // Loop over the side list
1850  boundary_side_iter pos = _boundary_side_id.begin();
1851  for (; pos != _boundary_side_id.end(); ++pos)
1852  {
1853  // Don't add remote sides
1854  if (pos->first->is_remote())
1855  continue;
1856 
1857  // Need to loop over the sides of any possible children
1858  std::vector<const Elem *> family;
1859 #ifdef LIBMESH_ENABLE_AMR
1860  pos->first->active_family_tree_by_side (family, pos->second.first);
1861 #else
1862  family.push_back(pos->first);
1863 #endif
1864 
1865  for (std::size_t elem_it=0; elem_it < family.size(); elem_it++)
1866  {
1867  const Elem * cur_elem = family[elem_it];
1868 
1869  UniquePtr<const Elem> side = cur_elem->build_side_ptr(pos->second.first);
1870 
1871  // Add each node node on the side with the side's boundary id
1872  for (auto i : side->node_index_range())
1873  {
1874  const boundary_id_type bcid = pos->second.second;
1875  this->add_node(side->node_ptr(i), bcid);
1876  if (!mesh_is_serial)
1877  {
1878  const processor_id_type proc_id =
1879  side->node_ptr(i)->processor_id();
1880  if (proc_id != my_proc_id)
1881  nodes_to_push[proc_id].insert
1882  (std::make_pair(side->node_id(i), bcid));
1883  }
1884  }
1885  }
1886  }
1887 
1888  // If we're on a serial mesh then we're done.
1889  if (mesh_is_serial)
1890  return;
1891 
1892  // Otherwise we need to push ghost node bcids to their owners, then
1893  // pull ghost node bcids from their owners.
1895  node_pushes_tag = this->comm().get_unique_tag(31337),
1896  node_pulls_tag = this->comm().get_unique_tag(31338),
1897  node_responses_tag = this->comm().get_unique_tag(31339);
1898 
1899  std::vector<Parallel::Request> node_push_requests(n_proc-1);
1900 
1901  for (processor_id_type p = 0; p != n_proc; ++p)
1902  {
1903  if (p == my_proc_id)
1904  continue;
1905 
1907  node_push_requests[p - (p > my_proc_id)];
1908 
1909  this->comm().send
1910  (p, nodes_to_push[p], request, node_pushes_tag);
1911  }
1912 
1913  for (processor_id_type p = 1; p != n_proc; ++p)
1914  {
1915  set_type received_nodes;
1916 
1917  this->comm().receive
1918  (Parallel::any_source, received_nodes, node_pushes_tag);
1919 
1920  for (set_type::const_iterator it = received_nodes.begin(),
1921  end = received_nodes.end(); it != end; ++it)
1922  this->add_node(_mesh.node_ptr(it->first), it->second);
1923  }
1924 
1925  // At this point we should know all the BCs for our own nodes; now
1926  // we need BCs for ghost nodes.
1927  //
1928  // FIXME - parallel_ghost_sync.h doesn't work here because it
1929  // assumes a fixed size datum on each node.
1930  std::vector<std::vector<dof_id_type>> node_ids_requested(n_proc);
1931 
1932  // Determine what nodes we need to request
1933  for (const auto & node : _mesh.node_ptr_range())
1934  {
1935  const processor_id_type pid = node->processor_id();
1936  if (pid != my_proc_id)
1937  node_ids_requested[pid].push_back(node->id());
1938  }
1939 
1940  typedef std::vector<std::pair<dof_id_type, boundary_id_type>> vec_type;
1941 
1942  std::vector<Parallel::Request>
1943  node_pull_requests(n_proc-1),
1944  node_response_requests(n_proc-1);
1945 
1946  // Make all requests
1947  for (processor_id_type p = 0; p != n_proc; ++p)
1948  {
1949  if (p == my_proc_id)
1950  continue;
1951 
1953  node_pull_requests[p - (p > my_proc_id)];
1954 
1955  this->comm().send
1956  (p, node_ids_requested[p], request, node_pulls_tag);
1957  }
1958 
1959  // Process all incoming requests
1960  std::vector<vec_type> responses(n_proc-1);
1961 
1962  for (processor_id_type p = 1; p != n_proc; ++p)
1963  {
1964  std::vector<dof_id_type> requested_nodes;
1965 
1967  status(this->comm().probe (Parallel::any_source, node_pulls_tag));
1968  const processor_id_type
1969  source_pid = cast_int<processor_id_type>(status.source());
1970 
1971  this->comm().receive
1972  (source_pid, requested_nodes, node_pulls_tag);
1973 
1975  node_response_requests[p-1];
1976 
1977  std::vector<boundary_id_type> bcids;
1978 
1979  for (std::vector<dof_id_type>::const_iterator
1980  it = requested_nodes.begin(),
1981  end = requested_nodes.end(); it != end; ++it)
1982  {
1983  this->boundary_ids(_mesh.node_ptr(*it), bcids);
1984 
1985  for (std::size_t i=0; i != bcids.size(); ++i)
1986  {
1987  const boundary_id_type b = bcids[i];
1988  responses[p-1].push_back(std::make_pair(*it, b));
1989  }
1990  }
1991 
1992  this->comm().send
1993  (source_pid, responses[p-1], request, node_responses_tag);
1994  }
1995 
1996  // Process all incoming responses
1997  for (processor_id_type p = 1; p != n_proc; ++p)
1998  {
2000  status(this->comm().probe (Parallel::any_source, node_responses_tag));
2001  const processor_id_type
2002  source_pid = cast_int<processor_id_type>(status.source());
2003 
2004  vec_type response;
2005 
2006  this->comm().receive
2007  (source_pid, response, node_responses_tag);
2008 
2009  for (vec_type::const_iterator
2010  it = response.begin(),
2011  end = response.end(); it != end; ++it)
2012  {
2013  this->add_node(_mesh.node_ptr(it->first), it->second);
2014  }
2015  }
2016 
2017  Parallel::wait (node_push_requests);
2018  Parallel::wait (node_pull_requests);
2019  Parallel::wait (node_response_requests);
2020 }
2021 
2022 
2023 
2024 
2026 {
2027  // Check for early return
2028  if (_boundary_node_id.empty())
2029  {
2030  libMesh::out << "No boundary node IDs have been added: cannot build side list!" << std::endl;
2031  return;
2032  }
2033 
2034  for (const auto & elem : _mesh.active_element_ptr_range())
2035  for (auto side : elem->side_index_range())
2036  {
2037  UniquePtr<const Elem> side_elem = elem->build_side_ptr(side);
2038 
2039  const unsigned short n_nodes = side_elem->n_nodes();
2040 
2041  // map from nodeset_id to count for that ID
2042  std::map<boundary_id_type, unsigned> nodesets_node_count;
2043  for (unsigned node_num=0; node_num < n_nodes; ++node_num)
2044  {
2045  const Node * node = side_elem->node_ptr(node_num);
2046 
2047  // For each nodeset that this node is a member of, increment the associated
2048  // nodeset ID count
2049  for (const auto & pr : as_range(_boundary_node_id.equal_range(node)))
2050  nodesets_node_count[pr.second]++;
2051  }
2052 
2053  // Now check to see what nodeset_counts have the correct
2054  // number of nodes in them. For any that do, add this side to
2055  // the sideset, making sure the sideset inherits the
2056  // nodeset's name, if there is one.
2057  std::map<boundary_id_type, unsigned>::const_iterator nodesets = nodesets_node_count.begin();
2058  for (; nodesets != nodesets_node_count.end(); ++nodesets)
2059  if (nodesets->second == n_nodes)
2060  {
2061  add_side(elem, side, nodesets->first);
2062 
2063  // Let the sideset inherit any non-empty name from the nodeset
2064  std::string & nset_name = nodeset_name(nodesets->first);
2065 
2066  if (nset_name != "")
2067  sideset_name(nodesets->first) = nset_name;
2068  }
2069  } // end for side
2070 }
2071 
2072 
2073 
2074 
2075 void BoundaryInfo::build_side_list (std::vector<dof_id_type> & el,
2076  std::vector<unsigned short int> & sl,
2077  std::vector<boundary_id_type> & il) const
2078 {
2079  // Clear the input vectors, just in case they were used for
2080  // something else recently...
2081  el.clear();
2082  sl.clear();
2083  il.clear();
2084 
2085  // Reserve the size, then use push_back
2086  el.reserve (_boundary_side_id.size());
2087  sl.reserve (_boundary_side_id.size());
2088  il.reserve (_boundary_side_id.size());
2089 
2090  boundary_side_iter pos = _boundary_side_id.begin();
2091  for (; pos != _boundary_side_id.end(); ++pos)
2092  {
2093  el.push_back (pos->first->id());
2094  sl.push_back (pos->second.first);
2095  il.push_back (pos->second.second);
2096  }
2097 }
2098 
2099 void BoundaryInfo::build_active_side_list (std::vector<dof_id_type> & el,
2100  std::vector<unsigned short int> & sl,
2101  std::vector<boundary_id_type> & il) const
2102 {
2103  // Clear the input vectors, just in case they were used for
2104  // something else recently...
2105  el.clear();
2106  sl.clear();
2107  il.clear();
2108 
2109  boundary_side_iter pos = _boundary_side_id.begin();
2110  for (; pos != _boundary_side_id.end(); ++pos)
2111  {
2112  // Don't add remote sides
2113  if (pos->first->is_remote())
2114  continue;
2115 
2116  // Loop over the sides of possible children
2117  std::vector<const Elem *> family;
2118 #ifdef LIBMESH_ENABLE_AMR
2119  pos->first->active_family_tree_by_side(family, pos->second.first);
2120 #else
2121  family.push_back(pos->first);
2122 #endif
2123 
2124  // Populate the list items
2125  for (std::vector<const Elem *>::iterator elem_it = family.begin(); elem_it != family.end(); elem_it++)
2126  {
2127  el.push_back ((*elem_it)->id());
2128  sl.push_back (pos->second.first);
2129  il.push_back (pos->second.second);
2130  }
2131  }
2132 }
2133 
2134 
2135 void BoundaryInfo::build_edge_list (std::vector<dof_id_type> & el,
2136  std::vector<unsigned short int> & sl,
2137  std::vector<boundary_id_type> & il) const
2138 {
2139  // Clear the input vectors, just in case they were used for
2140  // something else recently...
2141  el.clear();
2142  sl.clear();
2143  il.clear();
2144 
2145  // Reserve the size, then use push_back
2146  el.reserve (_boundary_edge_id.size());
2147  sl.reserve (_boundary_edge_id.size());
2148  il.reserve (_boundary_edge_id.size());
2149 
2150  boundary_edge_iter pos = _boundary_edge_id.begin();
2151  for (; pos != _boundary_edge_id.end(); ++pos)
2152  {
2153  el.push_back (pos->first->id());
2154  sl.push_back (pos->second.first);
2155  il.push_back (pos->second.second);
2156  }
2157 }
2158 
2159 
2160 void BoundaryInfo::build_shellface_list (std::vector<dof_id_type> & el,
2161  std::vector<unsigned short int> & sl,
2162  std::vector<boundary_id_type> & il) const
2163 {
2164  // Clear the input vectors, just in case they were used for
2165  // something else recently...
2166  el.clear();
2167  sl.clear();
2168  il.clear();
2169 
2170  // Reserve the size, then use push_back
2171  el.reserve (_boundary_shellface_id.size());
2172  sl.reserve (_boundary_shellface_id.size());
2173  il.reserve (_boundary_shellface_id.size());
2174 
2176  for (; pos != _boundary_shellface_id.end(); ++pos)
2177  {
2178  el.push_back (pos->first->id());
2179  sl.push_back (pos->second.first);
2180  il.push_back (pos->second.second);
2181  }
2182 }
2183 
2184 
2185 void BoundaryInfo::print_info(std::ostream & out_stream) const
2186 {
2187  // Print out the nodal BCs
2188  if (!_boundary_node_id.empty())
2189  {
2190  out_stream << "Nodal Boundary conditions:" << std::endl
2191  << "--------------------------" << std::endl
2192  << " (Node No., ID) " << std::endl;
2193 
2196 
2197  for (; it != end; ++it)
2198  out_stream << " (" << (*it).first->id()
2199  << ", " << (*it).second
2200  << ")" << std::endl;
2201  }
2202 
2203  // Print out the element edge BCs
2204  if (!_boundary_edge_id.empty())
2205  {
2206  out_stream << std::endl
2207  << "Edge Boundary conditions:" << std::endl
2208  << "-------------------------" << std::endl
2209  << " (Elem No., Edge No., ID) " << std::endl;
2210 
2213 
2214  for (; it != end; ++it)
2215  out_stream << " (" << (*it).first->id()
2216  << ", " << (*it).second.first
2217  << ", " << (*it).second.second
2218  << ")" << std::endl;
2219  }
2220 
2221  // Print out the element shellface BCs
2222  if (!_boundary_shellface_id.empty())
2223  {
2224  out_stream << std::endl
2225  << "Shell-face Boundary conditions:" << std::endl
2226  << "-------------------------" << std::endl
2227  << " (Elem No., Shell-face No., ID) " << std::endl;
2228 
2231 
2232  for (; it != end; ++it)
2233  out_stream << " (" << (*it).first->id()
2234  << ", " << (*it).second.first
2235  << ", " << (*it).second.second
2236  << ")" << std::endl;
2237  }
2238 
2239  // Print out the element side BCs
2240  if (!_boundary_side_id.empty())
2241  {
2242  out_stream << std::endl
2243  << "Side Boundary conditions:" << std::endl
2244  << "-------------------------" << std::endl
2245  << " (Elem No., Side No., ID) " << std::endl;
2246 
2249 
2250  for (; it != end; ++it)
2251  out_stream << " (" << (*it).first->id()
2252  << ", " << (*it).second.first
2253  << ", " << (*it).second.second
2254  << ")" << std::endl;
2255  }
2256 }
2257 
2258 
2259 
2260 void BoundaryInfo::print_summary(std::ostream & out_stream) const
2261 {
2262  // Print out the nodal BCs
2263  if (!_boundary_node_id.empty())
2264  {
2265  out_stream << "Nodal Boundary conditions:" << std::endl
2266  << "--------------------------" << std::endl
2267  << " (ID, number of nodes) " << std::endl;
2268 
2269  std::map<boundary_id_type, std::size_t> ID_counts;
2270 
2273 
2274  for (; it != end; ++it)
2275  ID_counts[(*it).second]++;
2276 
2277  std::map<boundary_id_type, std::size_t>::const_iterator ID_it = ID_counts.begin();
2278  const std::map<boundary_id_type, std::size_t>::const_iterator ID_end = ID_counts.end();
2279 
2280  for (; ID_it != ID_end; ++ID_it)
2281  out_stream << " (" << (*ID_it).first
2282  << ", " << (*ID_it).second
2283  << ")" << std::endl;
2284  }
2285 
2286  // Print out the element edge BCs
2287  if (!_boundary_edge_id.empty())
2288  {
2289  out_stream << std::endl
2290  << "Edge Boundary conditions:" << std::endl
2291  << "-------------------------" << std::endl
2292  << " (ID, number of edges) " << std::endl;
2293 
2294  std::map<boundary_id_type, std::size_t> ID_counts;
2295 
2298 
2299  for (; it != end; ++it)
2300  ID_counts[(*it).second.second]++;
2301 
2302  std::map<boundary_id_type, std::size_t>::const_iterator ID_it = ID_counts.begin();
2303  const std::map<boundary_id_type, std::size_t>::const_iterator ID_end = ID_counts.end();
2304 
2305  for (; ID_it != ID_end; ++ID_it)
2306  out_stream << " (" << (*ID_it).first
2307  << ", " << (*ID_it).second
2308  << ")" << std::endl;
2309  }
2310 
2311 
2312  // Print out the element edge BCs
2313  if (!_boundary_shellface_id.empty())
2314  {
2315  out_stream << std::endl
2316  << "Shell-face Boundary conditions:" << std::endl
2317  << "-------------------------" << std::endl
2318  << " (ID, number of shellfaces) " << std::endl;
2319 
2320  std::map<boundary_id_type, std::size_t> ID_counts;
2321 
2324 
2325  for (; it != end; ++it)
2326  ID_counts[(*it).second.second]++;
2327 
2328  std::map<boundary_id_type, std::size_t>::const_iterator ID_it = ID_counts.begin();
2329  const std::map<boundary_id_type, std::size_t>::const_iterator ID_end = ID_counts.end();
2330 
2331  for (; ID_it != ID_end; ++ID_it)
2332  out_stream << " (" << (*ID_it).first
2333  << ", " << (*ID_it).second
2334  << ")" << std::endl;
2335  }
2336 
2337  // Print out the element side BCs
2338  if (!_boundary_side_id.empty())
2339  {
2340  out_stream << std::endl
2341  << "Side Boundary conditions:" << std::endl
2342  << "-------------------------" << std::endl
2343  << " (ID, number of sides) " << std::endl;
2344 
2345  std::map<boundary_id_type, std::size_t> ID_counts;
2346 
2349 
2350  for (; it != end; ++it)
2351  ID_counts[(*it).second.second]++;
2352 
2353  std::map<boundary_id_type, std::size_t>::const_iterator ID_it = ID_counts.begin();
2354  const std::map<boundary_id_type, std::size_t>::const_iterator ID_end = ID_counts.end();
2355 
2356  for (; ID_it != ID_end; ++ID_it)
2357  out_stream << " (" << (*ID_it).first
2358  << ", " << (*ID_it).second
2359  << ")" << std::endl;
2360  }
2361 }
2362 
2363 
2365 {
2366  static const std::string empty_string;
2367  std::map<boundary_id_type, std::string>::const_iterator it =
2368  _ss_id_to_name.find(id);
2369  if (it == _ss_id_to_name.end())
2370  return empty_string;
2371  else
2372  return it->second;
2373 }
2374 
2375 
2377 {
2378  return _ss_id_to_name[id];
2379 }
2380 
2382 {
2383  static const std::string empty_string;
2384  std::map<boundary_id_type, std::string>::const_iterator it =
2385  _ns_id_to_name.find(id);
2386  if (it == _ns_id_to_name.end())
2387  return empty_string;
2388  else
2389  return it->second;
2390 }
2391 
2393 {
2394  return _ns_id_to_name[id];
2395 }
2396 
2398 {
2399  // Search sidesets
2400  std::map<boundary_id_type, std::string>::const_iterator
2401  iter = _ss_id_to_name.begin(),
2402  end_iter = _ss_id_to_name.end();
2403 
2404  for (; iter != end_iter; ++iter)
2405  if (iter->second == name)
2406  return iter->first;
2407 
2408  // Search nodesets
2409  iter = _ns_id_to_name.begin();
2410  end_iter = _ns_id_to_name.end();
2411  for (; iter != end_iter; ++iter)
2412  if (iter->second == name)
2413  return iter->first;
2414 
2415  // If we made it here without returning, we don't have a sideset or
2416  // nodeset by the requested name, so return invalid_id
2417  return invalid_id;
2418 }
2419 
2420 
2421 
2422 void BoundaryInfo::_find_id_maps(const std::set<boundary_id_type> & requested_boundary_ids,
2423  dof_id_type first_free_node_id,
2424  std::map<dof_id_type, dof_id_type> * node_id_map,
2425  dof_id_type first_free_elem_id,
2426  std::map<std::pair<dof_id_type, unsigned char>, dof_id_type> * side_id_map,
2427  const std::set<subdomain_id_type> & subdomains_relative_to)
2428 {
2429  // We'll do the same modulus trick that DistributedMesh uses to avoid
2430  // id conflicts between different processors
2431  dof_id_type
2432  next_node_id = first_free_node_id + this->processor_id(),
2433  next_elem_id = first_free_elem_id + this->processor_id();
2434 
2435  // We'll pass through the mesh once first to build
2436  // the maps and count boundary nodes and elements.
2437  // To find local boundary nodes, we have to examine all elements
2438  // here rather than just local elements, because it's possible to
2439  // have a local boundary node that's not on a local boundary
2440  // element, e.g. at the tip of a triangle.
2441 
2442  // We'll loop through two different ranges here: first all elements,
2443  // looking for local nodes, and second through unpartitioned
2444  // elements, looking for all remaining nodes.
2446  bool hit_end_el = false;
2447  const MeshBase::const_element_iterator end_unpartitioned_el =
2449 
2451  !hit_end_el || (el != end_unpartitioned_el); ++el)
2452  {
2453  if ((el == end_el) && !hit_end_el)
2454  {
2455  // Note that we're done with local nodes and just looking
2456  // for remaining unpartitioned nodes
2457  hit_end_el = true;
2458 
2459  // Join up the local results from other processors
2460  if (side_id_map)
2461  this->comm().set_union(*side_id_map);
2462  if (node_id_map)
2463  this->comm().set_union(*node_id_map);
2464 
2465  // Finally we'll pass through any unpartitioned elements to add them
2466  // to the maps and counts.
2467  next_node_id = first_free_node_id + this->n_processors();
2468  next_elem_id = first_free_elem_id + this->n_processors();
2469 
2471  if (el == end_unpartitioned_el)
2472  break;
2473  }
2474 
2475  const Elem * elem = *el;
2476 
2477  // If the subdomains_relative_to container has the
2478  // invalid_subdomain_id, we fall back on the "old" behavior of
2479  // adding sides regardless of this Elem's subdomain. Otherwise,
2480  // if the subdomains_relative_to container doesn't contain the
2481  // current Elem's subdomain_id(), we won't add any sides from
2482  // it.
2483  if (!subdomains_relative_to.count(Elem::invalid_subdomain_id) &&
2484  !subdomains_relative_to.count(elem->subdomain_id()))
2485  continue;
2486 
2487  // Get the top-level parent for this element. This is used for
2488  // searching for boundary sides on this element.
2489  const Elem * top_parent = elem->top_parent();
2490 
2491  // Find all the boundary side ids for this Elem.
2492  const std::pair<boundary_side_iter, boundary_side_iter>
2493  bounds = _boundary_side_id.equal_range(top_parent);
2494 
2495  for (auto s : elem->side_index_range())
2496  {
2497  bool add_this_side = false;
2498  boundary_id_type this_bcid = invalid_id;
2499 
2500  for (const auto & pr : as_range(bounds))
2501  {
2502  this_bcid = pr.second.second;
2503 
2504  // if this side is flagged with a boundary condition
2505  // and the user wants this id
2506  if ((pr.second.first == s) &&
2507  (requested_boundary_ids.count(this_bcid)))
2508  {
2509  add_this_side = true;
2510  break;
2511  }
2512  }
2513 
2514  // We may still want to add this side if the user called
2515  // sync() with no requested_boundary_ids. This corresponds
2516  // to the "old" style of calling sync() in which the entire
2517  // boundary was copied to the BoundaryMesh, and handles the
2518  // case where elements on the geometric boundary are not in
2519  // any sidesets.
2520  if (bounds.first == bounds.second &&
2521  requested_boundary_ids.count(invalid_id) &&
2522  elem->neighbor_ptr(s) == libmesh_nullptr)
2523  add_this_side = true;
2524 
2525  if (add_this_side)
2526  {
2527  // We only assign ids for our own and for
2528  // unpartitioned elements
2529  if (side_id_map &&
2530  ((elem->processor_id() == this->processor_id()) ||
2531  (elem->processor_id() ==
2533  {
2534  std::pair<dof_id_type, unsigned char> side_pair(elem->id(), s);
2535  libmesh_assert (!side_id_map->count(side_pair));
2536  (*side_id_map)[side_pair] = next_elem_id;
2537  next_elem_id += this->n_processors() + 1;
2538  }
2539 
2540  // Use a proxy element for the side to query nodes
2542  for (auto n : side->node_index_range())
2543  {
2544  const Node & node = side->node_ref(n);
2545 
2546  // In parallel we don't know enough to number
2547  // others' nodes ourselves.
2548  if (!hit_end_el &&
2549  (node.processor_id() != this->processor_id()))
2550  continue;
2551 
2552  dof_id_type node_id = node.id();
2553  if (node_id_map && !node_id_map->count(node_id))
2554  {
2555  (*node_id_map)[node_id] = next_node_id;
2556  next_node_id += this->n_processors() + 1;
2557  }
2558  }
2559  }
2560  }
2561  }
2562 
2563  // FIXME: ought to renumber side/node_id_map image to be contiguous
2564  // to save memory, also ought to reserve memory
2565 }
2566 
2567 
2568 } // namespace libMesh
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > >::const_iterator boundary_shellface_iter
Typedef for iterators into the _boundary_shellface_id container.
virtual bool is_node_on_side(const unsigned int n, const unsigned int s) const =0
void add_elements(const std::set< boundary_id_type > &requested_boundary_ids, UnstructuredMesh &boundary_mesh)
Generates elements along the boundary of our _mesh, which use pre-existing nodes on the boundary_mesh...
virtual unique_id_type parallel_max_unique_id() const =0
std::map< boundary_id_type, std::string > _ns_id_to_name
This structure maintains the mapping of named node sets for file formats that support named blocks...
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:117
std::set< boundary_id_type > _node_boundary_ids
Set of user-specified boundary IDs for nodes only.
Status wait(Request &r)
Wait for a non-blocking send or receive to finish.
Definition: parallel.h:565
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > >::iterator erase_iter
Some older compilers don&#39;t support erasing from a map with const_iterators, so we need to use a non-c...
virtual UniquePtr< Elem > build_side_ptr(const unsigned int i, bool proxy=true)=0
void make_node_unique_ids_parallel_consistent(MeshBase &)
Assuming all unique_ids on local nodes are globally unique, and assuming all processor ids are parall...
void remove_edge(const Elem *elem, const unsigned short int edge)
Removes all boundary conditions associated with edge edge of element elem, if any exist...
A Node is like a Point, but with more information.
Definition: node.h:52
virtual bool is_serial() const
Definition: mesh_base.h:140
std::string & nodeset_name(boundary_id_type id)
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value.
Definition: libmesh.h:184
void set_parent(Elem *p)
Sets the pointer to the element&#39;s parent.
Definition: elem.h:2362
std::set< boundary_id_type > _edge_boundary_ids
Set of user-specified boundary IDs for edges only.
void sync(UnstructuredMesh &boundary_mesh)
Generates boundary_mesh data structures corresponding to the mesh data structures.
void build_node_list_from_side_list()
Adds nodes with boundary ids based on the side&#39;s boundary ids they are connected to.
~BoundaryInfo()
Destructor.
const unsigned int any_source
Accept from any source.
Definition: parallel.h:204
Encapsulates the MPI_Status struct.
Definition: parallel.h:461
virtual bool is_edge_on_side(const unsigned int e, const unsigned int s) const =0
virtual bool is_child_on_edge(const unsigned int c, const unsigned int e) const
Definition: elem.C:1665
void print_summary(std::ostream &out=libMesh::out) const
Prints a summary of the boundary information.
void remove_id(boundary_id_type id)
Removes all entities (nodes, sides, edges, shellfaces) with boundary id id from their respective cont...
The definition of the const_element_iterator struct.
Definition: mesh_base.h:1494
BoundaryInfo(MeshBase &m)
Constructor.
Definition: boundary_info.C:49
IntRange< unsigned short > side_index_range() const
Definition: elem.h:2083
virtual unsigned int n_edges() const =0
std::multimap< const Node *, boundary_id_type >::iterator boundary_node_erase_iter
Some older compilers don&#39;t support erasing from a map with const_iterators, so we need to use a non-c...
void remove_shellface(const Elem *elem, const unsigned short int shellface)
Removes all boundary conditions associated with shell face shellface of element elem, if any exist.
static void set_node_processor_ids(MeshBase &mesh)
This function is called after partitioning to set the processor IDs for the nodes.
Definition: partitioner.C:416
unsigned int n_partitions() const
Definition: mesh_base.h:833
void raw_shellface_boundary_ids(const Elem *const elem, const unsigned short int shellface, std::vector< boundary_id_type > &vec_to_fill) const
std::set< boundary_id_type > _boundary_ids
A collection of user-specified boundary ids for sides, edges, nodes, and shell faces.
const Elem * parent() const
Definition: elem.h:2346
void remove(const Node *node)
Removes the boundary conditions associated with node node, if any exist.
unsigned short int side
Definition: xdr_io.C:49
std::multimap< const Node *, boundary_id_type >::const_iterator boundary_node_iter
Typedef for iterators into the _boundary_node_id container.
std::size_t n_edge_conds() const
boundary_id_type get_id_by_name(const std::string &name) const
processor_id_type n_processors() const
const Elem * interior_parent() const
Definition: elem.C:951
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
uint8_t processor_id_type
Definition: id_types.h:99
MPI_Request request
Request object for non-blocking I/O.
Definition: parallel.h:171
const class libmesh_nullptr_t libmesh_nullptr
virtual const Node * node_ptr(const dof_id_type i) const =0
std::set< boundary_id_type > _side_boundary_ids
Set of user-specified boundary IDs for sides only.
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > >::const_iterator boundary_edge_iter
Typedef for iterators into the _boundary_edge_id container.
std::size_t n_boundary_conds() const
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
unsigned int n_shellface_boundary_ids(const Elem *const elem, const unsigned short int shellface) const
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > >::const_iterator boundary_side_iter
Typedef for iterators into the _boundary_side_id container.
The libMesh namespace provides an interface to certain functionality in the library.
void build_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of element numbers, sides, and ids for those sides.
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > > _boundary_side_id
Data structure that maps sides of elements to boundary ids.
Real distance(const Point &p)
virtual bool is_child_on_side(const unsigned int c, const unsigned int s) const =0
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id)=0
Add a new Node at Point p to the end of the vertex array, with processor_id procid.
void set_interior_parent(Elem *p)
Sets the pointer to the element&#39;s interior_parent.
Definition: elem.C:1003
IntRange< unsigned short > edge_index_range() const
Definition: elem.h:2074
This is the MeshBase class.
Definition: mesh_base.h:68
const std::string & get_nodeset_name(boundary_id_type id) const
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual unsigned int n_nodes() const =0
std::map< boundary_id_type, std::string > _ss_id_to_name
This structure maintains the mapping of named side sets for file formats that support named blocks...
const Elem * neighbor_ptr(unsigned int i) const
Definition: elem.h:1967
std::vector< boundary_id_type > boundary_ids(const Node *node) const
std::size_t n_shellface_conds() const
virtual element_iterator elements_begin()=0
Iterate over all the elements in the Mesh.
void clear_boundary_node_ids()
Clears all the boundary information from all of the nodes in the mesh.
virtual dof_id_type max_elem_id() const =0
static const subdomain_id_type invalid_subdomain_id
A static integral constant representing an invalid subdomain id.
Definition: elem.h:237
void add_node(const Node *node, const boundary_id_type id)
Add Node node with boundary id id to the boundary information data structures.
const dof_id_type n_nodes
Definition: tecplot_io.C:67
unsigned int side_with_boundary_id(const Elem *const elem, const boundary_id_type boundary_id) const
int8_t boundary_id_type
Definition: id_types.h:51
void shellface_boundary_ids(const Elem *const elem, const unsigned short int shellface, std::vector< boundary_id_type > &vec_to_fill) const
void build_side_boundary_ids(std::vector< boundary_id_type > &b_ids) const
Builds the list of unique side boundary ids.
static const boundary_id_type invalid_id
Number used for internal use.
This is the MeshCommunication class.
virtual SimpleRange< element_iterator > element_ptr_range()=0
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:1874
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:335
std::multimap< const Node *, boundary_id_type > _boundary_node_id
Data structure that maps nodes in the mesh to boundary ids.
The UnstructuredMesh class is derived from the MeshBase class.
void build_side_list_from_node_list()
Adds sides to a sideset if every node on that side are in the same sideset.
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
Encapsulates the MPI tag integers.
Definition: parallel.h:227
unsigned int n_edge_boundary_ids(const Elem *const elem, const unsigned short int edge) const
virtual element_iterator elements_end()=0
virtual SimpleRange< node_iterator > node_ptr_range()=0
The BoundaryInfo class contains information relevant to boundary conditions including storing faces...
Definition: boundary_info.h:56
SimpleRange< I > as_range(const std::pair< I, I > &p)
Helper function that allows us to treat a homogenous pair as a range.
Definition: simple_range.h:57
MPI_Status status
Status object for querying messages.
Definition: parallel.h:176
void prepare_for_use(const bool skip_renumber_nodes_and_elements=false, const bool skip_find_neighbors=false)
Prepare a newly created (or read) mesh for use.
Definition: mesh_base.C:174
boundary_id_type boundary_id(const Elem *const elem, const unsigned short int side) const
void send(const unsigned int dest_processor_id, const T &buf, const MessageTag &tag=no_tag) const
Blocking-send to one processor with data-defined type.
subdomain_id_type subdomain_id() const
Definition: elem.h:1951
virtual unsigned int n_sides() const =0
void set_neighbor(const unsigned int i, Elem *n)
Assigns n as the neighbor.
Definition: elem.h:2000
void build_node_list(std::vector< dof_id_type > &node_id_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of nodes and ids for those nodes.
The DistributedMesh class is derived from the MeshBase class, and is intended to provide identical fu...
This class forms the base class for all other classes that are expected to be implemented in parallel...
void regenerate_id_sets()
Clears and regenerates the cached sets of ids.
virtual void clear()
Deletes all the data that are currently stored.
Definition: mesh_base.C:285
void clear()
Clears the underlying data structures and restores the object to a pristine state with no data stored...
std::string & sideset_name(boundary_id_type id)
BoundaryInfo & operator=(const BoundaryInfo &other_boundary_info)
Actual copying operation.
Definition: boundary_info.C:55
virtual element_iterator pid_elements_begin(processor_id_type proc_id)=0
Iterate over all elements with a specified processor id.
MeshBase & _mesh
The Mesh this boundary info pertains to.
unsigned int which_child_am_i(const Elem *e) const
Definition: elem.h:2487
void remove_side(const Elem *elem, const unsigned short int side)
Removes all boundary conditions associated with side side of element elem, if any exist...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Point & point(const unsigned int i) const
Definition: elem.h:1809
void build_active_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of active element numbers, sides, and ids for those sides.
virtual UniquePtr< Partitioner > & partitioner()
A partitioner to use at each prepare_for_use()
Definition: mesh_base.h:112
void copy_boundary_ids(const BoundaryInfo &old_boundary_info, const Elem *const old_elem, const Elem *const new_elem)
void build_edge_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &edge_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of element numbers, edges, and boundary ids for those edges.
std::size_t n_boundary_ids() const
Temporarily serialize a DistributedMesh for output; a distributed mesh is allgathered by the MeshSeri...
Encapsulates the MPI_Request.
Definition: parallel.h:517
void print_info(std::ostream &out=libMesh::out) const
Prints the boundary information data structure.
const Parallel::Communicator & comm() const
OStreamProxy out
std::vector< boundary_id_type > edge_boundary_ids(const Elem *const elem, const unsigned short int edge) const
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
void add_side(const dof_id_type elem, const unsigned short int side, const boundary_id_type id)
Add side side of element number elem with boundary id id to the boundary information data structure...
void add_shellface(const dof_id_type elem, const unsigned short int shellface, const boundary_id_type id)
Add shell face shellface of element number elem with boundary id id to the boundary information data ...
virtual unsigned int n_vertices() const =0
unsigned int level() const
Definition: elem.h:2388
Status receive(const unsigned int dest_processor_id, T &buf, const MessageTag &tag=any_tag) const
Blocking-receive from one processor with data-defined type.
std::vector< boundary_id_type > raw_edge_boundary_ids(const Elem *const elem, const unsigned short int edge) const
std::size_t n_nodeset_conds() const
std::set< boundary_id_type > _shellface_boundary_ids
Set of user-specified boundary IDs for shellfaces only.
const Elem * top_parent() const
Definition: elem.h:2370
void _find_id_maps(const std::set< boundary_id_type > &requested_boundary_ids, dof_id_type first_free_node_id, std::map< dof_id_type, dof_id_type > *node_id_map, dof_id_type first_free_elem_id, std::map< std::pair< dof_id_type, unsigned char >, dof_id_type > *side_id_map, const std::set< subdomain_id_type > &subdomains_relative_to)
Helper method for finding consistent maps of interior to boundary dof_object ids. ...
void build_shellface_boundary_ids(std::vector< boundary_id_type > &b_ids) const
Builds the list of unique shellface boundary ids.
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > > _boundary_shellface_id
Data structure that maps faces of shell elements to boundary ids.
virtual void delete_remote_elements()
When supported, deletes all nonlocal elements of the mesh except for "ghosts" which touch a local ele...
Definition: mesh_base.h:182
dof_id_type id() const
Definition: dof_object.h:632
virtual const Elem * elem_ptr(const dof_id_type i) const =0
void sum(T &r) const
Take a local variable and replace it with the sum of it&#39;s values on all processors.
const std::string & get_sideset_name(boundary_id_type id) const
void build_shellface_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &shellface_list, std::vector< boundary_id_type > &bc_id_list) const
Creates a list of element numbers, shellfaces, and boundary ids for those shellfaces.
uint8_t unique_id_type
Definition: id_types.h:79
unsigned int & set_n_partitions()
Definition: mesh_base.h:1329
void build_node_boundary_ids(std::vector< boundary_id_type > &b_ids) const
Builds the list of unique node boundary ids.
std::vector< boundary_id_type > raw_boundary_ids(const Elem *const elem, const unsigned short int side) const
MessageTag get_unique_tag(int tagvalue) const
Get a tag that is unique to this Communicator.
std::multimap< const Elem *, std::pair< unsigned short int, boundary_id_type > > _boundary_edge_id
Data structure that maps edges of elements to boundary ids.
status probe(const unsigned int src_processor_id, const MessageTag &tag=any_tag, const Communicator &comm=Communicator_World)
void get_side_and_node_maps(UnstructuredMesh &boundary_mesh, std::map< dof_id_type, dof_id_type > &node_id_map, std::map< dof_id_type, unsigned char > &side_id_map, Real tolerance=1.e-6)
Suppose we have used sync to create boundary_mesh.
processor_id_type processor_id() const
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type processor_id() const
Definition: dof_object.h:694
virtual void libmesh_assert_valid_parallel_ids() const libmesh_override
Verify id and processor_id consistency of our elements and nodes containers.
void add_edge(const dof_id_type elem, const unsigned short int edge, const boundary_id_type id)
Add edge edge of element number elem with boundary id id to the boundary information data structure...
void set_union(T &data, const unsigned int root_id) const
Take a container of local variables on each processor, and collect their union over all processors...
const RemoteElem * remote_elem
Definition: remote_elem.C:57
virtual element_iterator pid_elements_end(processor_id_type proc_id)=0