www.mooseframework.org
MooseMesh.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "MooseMesh.h"
11 #include "Factory.h"
13 #include "Assembly.h"
14 #include "MooseUtils.h"
15 #include "MooseApp.h"
16 #include "RelationshipManager.h"
17 #include "PointListAdaptor.h"
18 
19 #include <utility>
20 
21 // libMesh
22 #include "libmesh/bounding_box.h"
23 #include "libmesh/boundary_info.h"
24 #include "libmesh/mesh_tools.h"
25 #include "libmesh/parallel.h"
26 #include "libmesh/mesh_communication.h"
27 #include "libmesh/parallel_mesh.h"
28 #include "libmesh/periodic_boundary_base.h"
29 #include "libmesh/fe_base.h"
30 #include "libmesh/fe_interface.h"
31 #include "libmesh/serial_mesh.h"
32 #include "libmesh/mesh_inserter_iterator.h"
33 #include "libmesh/mesh_communication.h"
34 #include "libmesh/mesh_inserter_iterator.h"
35 #include "libmesh/mesh_tools.h"
36 #include "libmesh/parallel.h"
37 #include "libmesh/parallel_elem.h"
38 #include "libmesh/parallel_mesh.h"
39 #include "libmesh/parallel_node.h"
40 #include "libmesh/parallel_ghost_sync.h"
41 #include "libmesh/utility.h"
42 #include "libmesh/remote_elem.h"
43 #include "libmesh/linear_partitioner.h"
44 #include "libmesh/centroid_partitioner.h"
45 #include "libmesh/parmetis_partitioner.h"
46 #include "libmesh/hilbert_sfc_partitioner.h"
47 #include "libmesh/morton_sfc_partitioner.h"
48 #include "libmesh/edge_edge2.h"
49 #include "libmesh/mesh_refinement.h"
50 #include "libmesh/quadrature.h"
51 #include "libmesh/boundary_info.h"
52 #include "libmesh/periodic_boundaries.h"
53 #include "libmesh/quadrature_gauss.h"
54 #include "libmesh/point_locator_base.h"
55 #include "libmesh/default_coupling.h"
56 #include "libmesh/ghost_point_neighbors.h"
57 
58 static const int GRAIN_SIZE =
59  1; // the grain_size does not have much influence on our execution speed
60 
61 template <>
64 {
66 
67  MooseEnum parallel_type("DEFAULT REPLICATED DISTRIBUTED", "DEFAULT");
68  params.addParam<MooseEnum>("parallel_type",
69  parallel_type,
70  "DEFAULT: Use libMesh::ReplicatedMesh unless --distributed-mesh is "
71  "specified on the command line "
72  "REPLICATED: Always use libMesh::ReplicatedMesh "
73  "DISTRIBUTED: Always use libMesh::DistributedMesh");
74 
75  params.addParam<bool>(
76  "allow_renumbering",
77  true,
78  "If allow_renumbering=false, node and element numbers are kept fixed until deletion");
79 
80  params.addParam<bool>("nemesis",
81  false,
82  "If nemesis=true and file=foo.e, actually reads "
83  "foo.e.N.0, foo.e.N.1, ... foo.e.N.N-1, "
84  "where N = # CPUs, with NemesisIO.");
85 
86  MooseEnum dims("1=1 2 3", "1");
87  params.addParam<MooseEnum>("dim",
88  dims,
89  "This is only required for certain mesh formats where "
90  "the dimension of the mesh cannot be autodetected. "
91  "In particular you must supply this for GMSH meshes. "
92  "Note: This is completely ignored for ExodusII meshes!");
93 
94  MooseEnum partitioning("default=-3 metis=-2 parmetis=-1 linear=0 centroid hilbert_sfc morton_sfc",
95  "default");
96  params.addParam<MooseEnum>(
97  "partitioner",
98  partitioning,
99  "Specifies a mesh partitioner to use when splitting the mesh for a parallel computation.");
100  MooseEnum direction("x y z radial");
101  params.addParam<MooseEnum>("centroid_partitioner_direction",
102  direction,
103  "Specifies the sort direction if using the centroid partitioner. "
104  "Available options: x, y, z, radial");
105 
106  MooseEnum patch_update_strategy("never always auto iteration", "never");
107  params.addParam<MooseEnum>(
108  "patch_update_strategy",
109  patch_update_strategy,
110  "How often to update the geometric search 'patch'. The default is to "
111  "never update it (which is the most efficient but could be a problem "
112  "with lots of relative motion). 'always' will update the patch for all "
113  "slave nodes at the beginning of every timestep which might be time "
114  "consuming. 'auto' will attempt to determine at the start of which "
115  "timesteps the patch for all slave nodes needs to be updated automatically."
116  "'iteration' updates the patch at every nonlinear iteration for a "
117  "subset of slave nodes for which penetration is not detected. If there "
118  "can be substantial relative motion between the master and slave surfaces "
119  "during the nonlinear iterations within a timestep, it is advisable to use "
120  "'iteration' option to ensure accurate contact detection.");
121 
122  // Note: This parameter is named to match 'construct_side_list_from_node_list' in SetupMeshAction
123  params.addParam<bool>(
124  "construct_node_list_from_side_list",
125  true,
126  "Whether or not to generate nodesets from the sidesets (usually a good idea).");
127  params.addParam<unsigned int>(
128  "patch_size", 40, "The number of nodes to consider in the NearestNode neighborhood.");
129  params.addParam<unsigned int>("ghosting_patch_size",
130  "The number of nearest neighbors considered "
131  "for ghosting purposes when 'iteration' "
132  "patch update strategy is used. Default is "
133  "5 * patch_size.");
134  params.addParam<unsigned int>("max_leaf_size",
135  10,
136  "The maximum number of points in each leaf of the KDTree used in "
137  "the nearest neighbor search. As the leaf size becomes larger,"
138  "KDTree construction becomes faster but the nearest neighbor search"
139  "becomes slower.");
140 
141  params.registerBase("MooseMesh");
142 
143  // groups
144  params.addParamNamesToGroup(
145  "dim nemesis patch_update_strategy construct_node_list_from_side_list patch_size",
146  "Advanced");
147  params.addParamNamesToGroup("partitioner centroid_partitioner_direction", "Partitioning");
148 
149  return params;
150 }
151 
153  : MooseObject(parameters),
154  Restartable(this, "Mesh"),
155  PerfGraphInterface(this),
156  _parallel_type(getParam<MooseEnum>("parallel_type").getEnum<MooseMesh::ParallelType>()),
157  _use_distributed_mesh(false),
158  _distribution_overridden(false),
159  _parallel_type_overridden(false),
160  _partitioner_name(getParam<MooseEnum>("partitioner")),
161  _partitioner_overridden(false),
162  _custom_partitioner_requested(false),
163  _uniform_refine_level(0),
164  _is_nemesis(getParam<bool>("nemesis")),
165  _is_prepared(false),
166  _needs_prepare_for_use(false),
167  _node_to_elem_map_built(false),
168  _node_to_active_semilocal_elem_map_built(false),
169  _patch_size(getParam<unsigned int>("patch_size")),
170  _ghosting_patch_size(isParamValid("ghosting_patch_size")
171  ? getParam<unsigned int>("ghosting_patch_size")
172  : 5 * _patch_size),
173  _max_leaf_size(getParam<unsigned int>("max_leaf_size")),
174  _patch_update_strategy(
175  getParam<MooseEnum>("patch_update_strategy").getEnum<Moose::PatchUpdateType>()),
176  _regular_orthogonal_mesh(false),
177  _allow_recovery(true),
178  _construct_node_list_from_side_list(getParam<bool>("construct_node_list_from_side_list")),
179  _prepare_timer(registerTimedSection("prepare", 2)),
180  _update_timer(registerTimedSection("update", 3)),
181  _mesh_changed_timer(registerTimedSection("meshChanged", 3)),
182  _cache_changed_lists_timer(registerTimedSection("cacheChangedLists", 5)),
183  _update_active_semi_local_node_range_timer(
184  registerTimedSection("updateActiveSemiLocalNodeRange", 5)),
185  _build_node_list_timer(registerTimedSection("buildNodeList", 5)),
186  _build_bnd_elem_list_timer(registerTimedSection("buildBndElemList", 5)),
187  _node_to_elem_map_timer(registerTimedSection("nodeToElemMap", 5)),
188  _node_to_active_semilocal_elem_map_timer(
189  registerTimedSection("nodeToActiveSemilocalElemMap", 5)),
190  _get_active_local_element_range_timer(registerTimedSection("getActiveLocalElementRange", 5)),
191  _get_active_node_range_timer(registerTimedSection("getActiveNodeRange", 5)),
192  _get_local_node_range_timer(registerTimedSection("getLocalNodeRange", 5)),
193  _get_boundary_node_range_timer(registerTimedSection("getBoundaryNodeRange", 5)),
194  _get_boundary_element_range_timer(registerTimedSection("getBoundaryElementRange", 5)),
195  _cache_info_timer(registerTimedSection("cacheInfo", 3)),
196  _build_periodic_node_map_timer(registerTimedSection("buildPeriodicNodeMap", 5)),
197  _build_periodic_node_sets_timer(registerTimedSection("buildPeriodicNodeSets", 5)),
198  _detect_orthogonal_dim_ranges_timer(registerTimedSection("detectOrthogonalDimRanges", 5)),
199  _detect_paired_sidesets_timer(registerTimedSection("detectPairedSidesets", 5)),
200  _build_refinement_map_timer(registerTimedSection("buildRefinementMap", 5)),
201  _build_coarsening_map_timer(registerTimedSection("buildCoarseningMap", 5)),
202  _find_adaptivity_qp_maps_timer(registerTimedSection("findAdaptivityQpMaps", 5)),
203  _build_refinement_and_coarsening_maps_timer(
204  registerTimedSection("buildRefinementAndCoarseningMaps", 5)),
205  _change_boundary_id_timer(registerTimedSection("changeBoundaryId", 6)),
206  _init_timer(registerTimedSection("init", 2)),
207  _read_recovered_mesh_timer(registerTimedSection("readRecoveredMesh", 2)),
208  _ghost_ghosted_boundaries_timer(registerTimedSection("GhostGhostedBoundaries", 3)),
209  _need_delete(false)
210 {
211  if (isParamValid("ghosting_patch_size") && (_patch_update_strategy != Moose::Iteration))
212  mooseError("Ghosting patch size parameter has to be set in the mesh block "
213  "only when 'iteration' patch update strategy is used.");
214 }
215 
216 MooseMesh::MooseMesh(const MooseMesh & other_mesh)
217  : MooseObject(other_mesh._pars),
218  Restartable(this, "Mesh"),
219  PerfGraphInterface(this, "CopiedMesh"),
220  _parallel_type(other_mesh._parallel_type),
221  _use_distributed_mesh(other_mesh._use_distributed_mesh),
222  _distribution_overridden(other_mesh._distribution_overridden),
223  _mesh(other_mesh.getMesh().clone()),
224  _partitioner_name(other_mesh._partitioner_name),
225  _partitioner_overridden(other_mesh._partitioner_overridden),
226  _uniform_refine_level(other_mesh.uniformRefineLevel()),
227  _is_nemesis(false),
228  _is_prepared(false),
229  _needs_prepare_for_use(false),
230  _node_to_elem_map_built(false),
231  _patch_size(other_mesh._patch_size),
232  _ghosting_patch_size(other_mesh._ghosting_patch_size),
233  _max_leaf_size(other_mesh._max_leaf_size),
234  _patch_update_strategy(other_mesh._patch_update_strategy),
235  _regular_orthogonal_mesh(false),
236  _construct_node_list_from_side_list(other_mesh._construct_node_list_from_side_list),
237  _prepare_timer(registerTimedSection("prepare", 2)),
238  _update_timer(registerTimedSection("update", 2)),
239  _mesh_changed_timer(registerTimedSection("meshChanged", 3)),
240  _cache_changed_lists_timer(registerTimedSection("cacheChangedLists", 5)),
241  _update_active_semi_local_node_range_timer(
242  registerTimedSection("updateActiveSemiLocalNodeRange", 5)),
243  _build_node_list_timer(registerTimedSection("buildNodeList", 5)),
244  _build_bnd_elem_list_timer(registerTimedSection("buildBndElemList", 5)),
245  _node_to_elem_map_timer(registerTimedSection("nodeToElemMap", 5)),
246  _node_to_active_semilocal_elem_map_timer(
247  registerTimedSection("nodeToActiveSemilocalElemMap", 5)),
248  _get_active_local_element_range_timer(registerTimedSection("getActiveLocalElementRange", 5)),
249  _get_active_node_range_timer(registerTimedSection("getActiveNodeRange", 5)),
250  _get_local_node_range_timer(registerTimedSection("getLocalNodeRange", 5)),
251  _get_boundary_node_range_timer(registerTimedSection("getBoundaryNodeRange", 5)),
252  _get_boundary_element_range_timer(registerTimedSection("getBoundaryElementRange", 5)),
253  _cache_info_timer(registerTimedSection("cacheInfo", 3)),
254  _build_periodic_node_map_timer(registerTimedSection("buildPeriodicNodeMap", 5)),
255  _build_periodic_node_sets_timer(registerTimedSection("buildPeriodicNodeSets", 5)),
256  _detect_orthogonal_dim_ranges_timer(registerTimedSection("detectOrthogonalDimRanges", 5)),
257  _detect_paired_sidesets_timer(registerTimedSection("detectPairedSidesets", 5)),
258  _build_refinement_map_timer(registerTimedSection("buildRefinementMap", 5)),
259  _build_coarsening_map_timer(registerTimedSection("buildCoarseningMap", 5)),
260  _find_adaptivity_qp_maps_timer(registerTimedSection("findAdaptivityQpMaps", 5)),
261  _build_refinement_and_coarsening_maps_timer(
262  registerTimedSection("buildRefinementAndCoarseningMaps", 5)),
263  _change_boundary_id_timer(registerTimedSection("changeBoundaryId", 5)),
264  _init_timer(registerTimedSection("init", 2)),
265  _read_recovered_mesh_timer(registerTimedSection("readRecoveredMesh", 2)),
266  _ghost_ghosted_boundaries_timer(registerTimedSection("GhostGhostedBoundaries", 3))
267 {
268  // Note: this calls BoundaryInfo::operator= without changing the
269  // ownership semantics of either Mesh's BoundaryInfo object.
270  getMesh().get_boundary_info() = other_mesh.getMesh().get_boundary_info();
271 
272  const std::set<SubdomainID> & subdomains = other_mesh.meshSubdomains();
273  for (const auto & sbd_id : subdomains)
274  setSubdomainName(sbd_id, other_mesh.getMesh().subdomain_name(sbd_id));
275 
276  // Get references to BoundaryInfo objects to make the code below cleaner...
277  const BoundaryInfo & other_boundary_info = other_mesh.getMesh().get_boundary_info();
278  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
279 
280  // Use the other BoundaryInfo object to build the list of side boundary ids
281  std::vector<BoundaryID> side_boundaries;
282  other_boundary_info.build_side_boundary_ids(side_boundaries);
283 
284  // Assign those boundary ids in our BoundaryInfo object
285  for (const auto & side_bnd_id : side_boundaries)
286  boundary_info.sideset_name(side_bnd_id) = other_boundary_info.get_sideset_name(side_bnd_id);
287 
288  // Do the same thing for node boundary ids
289  std::vector<BoundaryID> node_boundaries;
290  other_boundary_info.build_node_boundary_ids(node_boundaries);
291 
292  for (const auto & node_bnd_id : node_boundaries)
293  boundary_info.nodeset_name(node_bnd_id) = other_boundary_info.get_nodeset_name(node_bnd_id);
294 
295  _bounds.resize(other_mesh._bounds.size());
296  for (std::size_t i = 0; i < _bounds.size(); ++i)
297  {
298  _bounds[i].resize(other_mesh._bounds[i].size());
299  for (std::size_t j = 0; j < _bounds[i].size(); ++j)
300  _bounds[i][j] = other_mesh._bounds[i][j];
301  }
302 }
303 
305 {
306  freeBndNodes();
307  freeBndElems();
309 }
310 
311 void
313 {
314  // free memory
315  for (auto & bnode : _bnd_nodes)
316  delete bnode;
317 
318  for (auto & it : _node_set_nodes)
319  it.second.clear();
320 
321  _node_set_nodes.clear();
322 
323  for (auto & it : _bnd_node_ids)
324  it.second.clear();
325 
326  _bnd_node_ids.clear();
327 }
328 
329 void
331 {
332  // free memory
333  for (auto & belem : _bnd_elems)
334  delete belem;
335 
336  for (auto & it : _bnd_elem_ids)
337  it.second.clear();
338 
339  _bnd_elem_ids.clear();
340 }
341 
342 void
344 {
345  TIME_SECTION(_prepare_timer);
346 
347  mooseAssert(_mesh, "The MeshBase has not been constructed");
348 
349  if (dynamic_cast<DistributedMesh *>(&getMesh()) && !_is_nemesis)
350  {
351  // Call prepare_for_use() and don't mess with the renumbering
352  // setting
353  if (force || _needs_prepare_for_use)
354  getMesh().prepare_for_use();
355  }
356  else
357  {
358  // Call prepare_for_use() and DO NOT allow renumbering
359  getMesh().allow_renumbering(false);
360  if (force || _needs_prepare_for_use)
361  getMesh().prepare_for_use();
362  }
363 
364  // Collect (local) subdomain IDs
365  _mesh_subdomains.clear();
366  for (const auto & elem : getMesh().element_ptr_range())
367  _mesh_subdomains.insert(elem->subdomain_id());
368 
369  // Make sure nodesets have been generated
371 
372  // Collect (local) boundary IDs
373  const std::set<BoundaryID> & local_bids = getMesh().get_boundary_info().get_boundary_ids();
374  _mesh_boundary_ids.insert(local_bids.begin(), local_bids.end());
375 
376  const std::set<BoundaryID> & local_node_bids =
377  getMesh().get_boundary_info().get_node_boundary_ids();
378  _mesh_nodeset_ids.insert(local_node_bids.begin(), local_node_bids.end());
379 
380  const std::set<BoundaryID> & local_side_bids =
381  getMesh().get_boundary_info().get_side_boundary_ids();
382  _mesh_sideset_ids.insert(local_side_bids.begin(), local_side_bids.end());
383 
384  // Communicate subdomain and boundary IDs if this is a parallel mesh
385  if (!getMesh().is_serial())
386  {
387  _communicator.set_union(_mesh_subdomains);
388  _communicator.set_union(_mesh_boundary_ids);
389  _communicator.set_union(_mesh_nodeset_ids);
390  _communicator.set_union(_mesh_sideset_ids);
391  }
392 
394 
395  update();
396 
397  // Prepared has been called
398  _is_prepared = true;
399  _needs_prepare_for_use = false;
400 }
401 
402 void
404 {
405  TIME_SECTION(_update_timer);
406 
407  // Rebuild the boundary conditions
409 
410  // Update the node to elem map
411  _node_to_elem_map.clear();
412  _node_to_elem_map_built = false;
415 
416  buildNodeList();
418  cacheInfo();
419 }
420 
421 const Node &
422 MooseMesh::node(const dof_id_type i) const
423 {
424  mooseDeprecated("MooseMesh::node() is deprecated, please use MooseMesh::nodeRef() instead");
425  return nodeRef(i);
426 }
427 
428 Node &
429 MooseMesh::node(const dof_id_type i)
430 {
431  mooseDeprecated("MooseMesh::node() is deprecated, please use MooseMesh::nodeRef() instead");
432  return nodeRef(i);
433 }
434 
435 const Node &
436 MooseMesh::nodeRef(const dof_id_type i) const
437 {
438  if (i > getMesh().max_node_id())
439  return *(*_quadrature_nodes.find(i)).second;
440 
441  return getMesh().node_ref(i);
442 }
443 
444 Node &
445 MooseMesh::nodeRef(const dof_id_type i)
446 {
447  if (i > getMesh().max_node_id())
448  return *_quadrature_nodes[i];
449 
450  return getMesh().node_ref(i);
451 }
452 
453 const Node *
454 MooseMesh::nodePtr(const dof_id_type i) const
455 {
456  if (i > getMesh().max_node_id())
457  return (*_quadrature_nodes.find(i)).second;
458 
459  return getMesh().node_ptr(i);
460 }
461 
462 Node *
463 MooseMesh::nodePtr(const dof_id_type i)
464 {
465  if (i > getMesh().max_node_id())
466  return _quadrature_nodes[i];
467 
468  return getMesh().node_ptr(i);
469 }
470 
471 const Node *
472 MooseMesh::queryNodePtr(const dof_id_type i) const
473 {
474  if (i > getMesh().max_node_id())
475  return (*_quadrature_nodes.find(i)).second;
476 
477  return getMesh().query_node_ptr(i);
478 }
479 
480 Node *
481 MooseMesh::queryNodePtr(const dof_id_type i)
482 {
483  if (i > getMesh().max_node_id())
484  return _quadrature_nodes[i];
485 
486  return getMesh().query_node_ptr(i);
487 }
488 
489 void
491 {
492  TIME_SECTION(_mesh_changed_timer);
493 
494  update();
495 
496  // Delete all of the cached ranges
497  _active_local_elem_range.reset();
498  _active_node_range.reset();
500  _local_node_range.reset();
501  _bnd_node_range.reset();
502  _bnd_elem_range.reset();
503 
504  // Rebuild the ranges
510 
511  // Call the callback function onMeshChanged
512  onMeshChanged();
513 }
514 
515 void
517 {
518 }
519 
520 void
522 {
523  TIME_SECTION(_cache_changed_lists_timer);
524 
525  ConstElemRange elem_range(getMesh().local_elements_begin(), getMesh().local_elements_end(), 1);
526  CacheChangedListsThread cclt(*this);
527  Threads::parallel_reduce(elem_range, cclt);
528 
530 
531  _refined_elements = libmesh_make_unique<ConstElemPointerRange>(cclt._refined_elements.begin(),
532  cclt._refined_elements.end());
533  _coarsened_elements = libmesh_make_unique<ConstElemPointerRange>(cclt._coarsened_elements.begin(),
534  cclt._coarsened_elements.end());
536 }
537 
540 {
541  return _refined_elements.get();
542 }
543 
546 {
547  return _coarsened_elements.get();
548 }
549 
550 const std::vector<const Elem *> &
551 MooseMesh::coarsenedElementChildren(const Elem * elem) const
552 {
553  auto elem_to_child_pair = _coarsened_element_children.find(elem);
554  mooseAssert(elem_to_child_pair != _coarsened_element_children.end(), "Missing element in map");
555  return elem_to_child_pair->second;
556 }
557 
558 void
559 MooseMesh::updateActiveSemiLocalNodeRange(std::set<dof_id_type> & ghosted_elems)
560 {
562 
563  _semilocal_node_list.clear();
564 
565  // First add the nodes connected to local elems
566  ConstElemRange * active_local_elems = getActiveLocalElementRange();
567  for (const auto & elem : *active_local_elems)
568  {
569  for (unsigned int n = 0; n < elem->n_nodes(); ++n)
570  {
571  // Since elem is const here but we require a non-const Node * to
572  // store in the _semilocal_node_list (otherwise things like
573  // UpdateDisplacedMeshThread don't work), we are using a
574  // const_cast. A more long-term fix would be to have
575  // getActiveLocalElementRange return a non-const ElemRange.
576  Node * node = const_cast<Node *>(elem->node_ptr(n));
577 
578  _semilocal_node_list.insert(node);
579  }
580  }
581 
582  // Now add the nodes connected to ghosted_elems
583  for (const auto & ghost_elem_id : ghosted_elems)
584  {
585  Elem * elem = getMesh().elem_ptr(ghost_elem_id);
586  for (unsigned int n = 0; n < elem->n_nodes(); n++)
587  {
588  Node * node = elem->node_ptr(n);
589 
590  _semilocal_node_list.insert(node);
591  }
592  }
593 
594  // Now create the actual range
595  _active_semilocal_node_range = libmesh_make_unique<SemiLocalNodeRange>(
597 }
598 
599 bool
600 MooseMesh::isSemiLocal(Node * const node) const
601 {
602  return _semilocal_node_list.find(node) != _semilocal_node_list.end();
603 }
604 
610 {
611 public:
613 
614  bool operator()(const BndNode * const & lhs, const BndNode * const & rhs)
615  {
616  if (lhs->_bnd_id < rhs->_bnd_id)
617  return true;
618 
619  if (lhs->_bnd_id > rhs->_bnd_id)
620  return false;
621 
622  if (lhs->_node->id() < rhs->_node->id())
623  return true;
624 
625  if (lhs->_node->id() > rhs->_node->id())
626  return false;
627 
628  return false;
629  }
630 };
631 
632 void
634 {
635  TIME_SECTION(_build_node_list_timer);
636 
637  freeBndNodes();
638 
639  auto bc_tuples = getMesh().get_boundary_info().build_node_list();
640 
641  int n = bc_tuples.size();
642  _bnd_nodes.clear();
643  _bnd_nodes.reserve(n);
644  for (const auto & t : bc_tuples)
645  {
646  auto node_id = std::get<0>(t);
647  auto bc_id = std::get<1>(t);
648 
649  _bnd_nodes.push_back(new BndNode(getMesh().node_ptr(node_id), bc_id));
650  _node_set_nodes[bc_id].push_back(node_id);
651  _bnd_node_ids[bc_id].insert(node_id);
652  }
653 
654  _bnd_nodes.reserve(_bnd_nodes.size() + _extra_bnd_nodes.size());
655  for (unsigned int i = 0; i < _extra_bnd_nodes.size(); i++)
656  {
657  BndNode * bnode = new BndNode(_extra_bnd_nodes[i]._node, _extra_bnd_nodes[i]._bnd_id);
658  _bnd_nodes.push_back(bnode);
659  _bnd_node_ids[std::get<1>(bc_tuples[i])].insert(_extra_bnd_nodes[i]._node->id());
660  }
661 
662  // This sort is here so that boundary conditions are always applied in the same order
663  std::sort(_bnd_nodes.begin(), _bnd_nodes.end(), BndNodeCompare());
664 }
665 
666 void
668 {
669  TIME_SECTION(_build_bnd_elem_list_timer);
670 
671  freeBndElems();
672 
673  auto bc_tuples = getMesh().get_boundary_info().build_active_side_list();
674 
675  int n = bc_tuples.size();
676  _bnd_elems.clear();
677  _bnd_elems.reserve(n);
678  for (const auto & t : bc_tuples)
679  {
680  auto elem_id = std::get<0>(t);
681  auto side_id = std::get<1>(t);
682  auto bc_id = std::get<2>(t);
683 
684  _bnd_elems.push_back(new BndElement(getMesh().elem_ptr(elem_id), side_id, bc_id));
685  _bnd_elem_ids[bc_id].insert(elem_id);
686  }
687 }
688 
689 const std::map<dof_id_type, std::vector<dof_id_type>> &
691 {
692  if (!_node_to_elem_map_built) // Guard the creation with a double checked lock
693  {
694  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
695 
697  {
698  TIME_SECTION(_node_to_elem_map_timer);
699 
700  for (const auto & elem : getMesh().active_element_ptr_range())
701  for (unsigned int n = 0; n < elem->n_nodes(); n++)
702  _node_to_elem_map[elem->node_id(n)].push_back(elem->id());
703 
704  _node_to_elem_map_built = true; // MUST be set at the end for double-checked locking to work!
705  }
706  }
707 
708  return _node_to_elem_map;
709 }
710 
711 const std::map<dof_id_type, std::vector<dof_id_type>> &
713 {
714  if (!_node_to_active_semilocal_elem_map_built) // Guard the creation with a double checked lock
715  {
716  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
717 
719 
721  {
722  for (const auto & elem :
723  as_range(getMesh().semilocal_elements_begin(), getMesh().semilocal_elements_end()))
724  if (elem->active())
725  for (unsigned int n = 0; n < elem->n_nodes(); n++)
726  _node_to_active_semilocal_elem_map[elem->node_id(n)].push_back(elem->id());
727 
729  true; // MUST be set at the end for double-checked locking to work!
730  }
731  }
732 
734 }
735 
736 ConstElemRange *
738 {
740  {
742 
743  _active_local_elem_range = libmesh_make_unique<ConstElemRange>(
744  getMesh().active_local_elements_begin(), getMesh().active_local_elements_end(), GRAIN_SIZE);
745  }
746 
747  return _active_local_elem_range.get();
748 }
749 
750 NodeRange *
752 {
753  if (!_active_node_range)
754  {
755  TIME_SECTION(_get_active_node_range_timer);
756 
757  _active_node_range = libmesh_make_unique<NodeRange>(
758  getMesh().active_nodes_begin(), getMesh().active_nodes_end(), GRAIN_SIZE);
759  }
760 
761  return _active_node_range.get();
762 }
763 
766 {
767  mooseAssert(_active_semilocal_node_range,
768  "_active_semilocal_node_range has not been created yet!");
769 
770  return _active_semilocal_node_range.get();
771 }
772 
773 ConstNodeRange *
775 {
776  if (!_local_node_range)
777  {
778  TIME_SECTION(_get_local_node_range_timer);
779 
780  _local_node_range = libmesh_make_unique<ConstNodeRange>(
781  getMesh().local_nodes_begin(), getMesh().local_nodes_end(), GRAIN_SIZE);
782  }
783 
784  return _local_node_range.get();
785 }
786 
789 {
790  if (!_bnd_node_range)
791  {
792  TIME_SECTION(_get_boundary_node_range_timer);
794  libmesh_make_unique<ConstBndNodeRange>(bndNodesBegin(), bndNodesEnd(), GRAIN_SIZE);
795  }
796 
797  return _bnd_node_range.get();
798 }
799 
802 {
803  if (!_bnd_elem_range)
804  {
805  TIME_SECTION(_get_boundary_element_range_timer);
807  libmesh_make_unique<ConstBndElemRange>(bndElemsBegin(), bndElemsEnd(), GRAIN_SIZE);
808  }
809 
810  return _bnd_elem_range.get();
811 }
812 
813 const std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>> &
815 {
816  return _bnd_elem_ids;
817 }
818 
819 void
821 {
822  TIME_SECTION(_cache_info_timer);
823 
824  // TODO: Thread this!
825  for (const auto & elem : getMesh().element_ptr_range())
826  {
827  SubdomainID subdomain_id = elem->subdomain_id();
828 
829  for (unsigned int side = 0; side < elem->n_sides(); side++)
830  {
831  std::vector<BoundaryID> boundaryids = getBoundaryIDs(elem, side);
832 
833  std::set<BoundaryID> & subdomain_set = _subdomain_boundary_ids[subdomain_id];
834 
835  subdomain_set.insert(boundaryids.begin(), boundaryids.end());
836  }
837 
838  for (unsigned int nd = 0; nd < elem->n_nodes(); ++nd)
839  {
840  Node & node = *elem->node_ptr(nd);
841  _block_node_list[node.id()].insert(elem->subdomain_id());
842  }
843  }
844 }
845 
846 const std::set<SubdomainID> &
847 MooseMesh::getNodeBlockIds(const Node & node) const
848 {
849  std::map<dof_id_type, std::set<SubdomainID>>::const_iterator it =
850  _block_node_list.find(node.id());
851 
852  if (it == _block_node_list.end())
853  mooseError("Unable to find node: ", node.id(), " in any block list.");
854 
855  return it->second;
856 }
857 
858 // default begin() accessor
861 {
862  Predicates::NotNull<bnd_node_iterator_imp> p;
863  return bnd_node_iterator(_bnd_nodes.begin(), _bnd_nodes.end(), p);
864 }
865 
866 // default end() accessor
869 {
870  Predicates::NotNull<bnd_node_iterator_imp> p;
871  return bnd_node_iterator(_bnd_nodes.end(), _bnd_nodes.end(), p);
872 }
873 
874 // default begin() accessor
877 {
878  Predicates::NotNull<bnd_elem_iterator_imp> p;
879  return bnd_elem_iterator(_bnd_elems.begin(), _bnd_elems.end(), p);
880 }
881 
882 // default end() accessor
885 {
886  Predicates::NotNull<bnd_elem_iterator_imp> p;
887  return bnd_elem_iterator(_bnd_elems.end(), _bnd_elems.end(), p);
888 }
889 
890 const Node *
891 MooseMesh::addUniqueNode(const Point & p, Real tol)
892 {
897  if (getMesh().n_nodes() != _node_map.size())
898  {
899  _node_map.clear();
900  _node_map.reserve(getMesh().n_nodes());
901  for (const auto & node : getMesh().node_ptr_range())
902  _node_map.push_back(node);
903  }
904 
905  Node * node = nullptr;
906  for (unsigned int i = 0; i < _node_map.size(); ++i)
907  {
908  if (p.relative_fuzzy_equals(*_node_map[i], tol))
909  {
910  node = _node_map[i];
911  break;
912  }
913  }
914  if (node == nullptr)
915  {
916  node = getMesh().add_node(new Node(p));
917  _node_map.push_back(node);
918  }
919 
920  mooseAssert(node != nullptr, "Node is NULL");
921  return node;
922 }
923 
924 Node *
926  const unsigned short int side,
927  const unsigned int qp,
928  BoundaryID bid,
929  const Point & point)
930 {
931  Node * qnode;
932 
933  if (_elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side].find(qp) ==
935  {
936  // Create a new node id starting from the max node id and counting down. This will be the least
937  // likely to collide with an existing node id.
938  // Note that we are using numeric_limits<unsigned>::max even
939  // though max_id is stored as a dof_id_type. I tried this with
940  // numeric_limits<dof_id_type>::max and it broke several tests in
941  // MOOSE. So, this is some kind of a magic number that we will
942  // just continue to use...
943  dof_id_type max_id = std::numeric_limits<unsigned int>::max() - 100;
944  dof_id_type new_id = max_id - _quadrature_nodes.size();
945 
946  if (new_id <= getMesh().max_node_id())
947  mooseError("Quadrature node id collides with existing node id!");
948 
949  qnode = new Node(point, new_id);
950 
951  // Keep track of this new node in two different ways for easy lookup
952  _quadrature_nodes[new_id] = qnode;
953  _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp] = qnode;
954 
955  if (elem->active())
956  {
957  _node_to_elem_map[new_id].push_back(elem->id());
958  _node_to_active_semilocal_elem_map[new_id].push_back(elem->id());
959  }
960  }
961  else
962  qnode = _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp];
963 
964  BndNode * bnode = new BndNode(qnode, bid);
965  _bnd_nodes.push_back(bnode);
966  _bnd_node_ids[bid].insert(qnode->id());
967 
968  _extra_bnd_nodes.push_back(*bnode);
969 
970  // Do this so the range will be regenerated next time it is accessed
971  _bnd_node_range.reset();
972 
973  return qnode;
974 }
975 
976 Node *
978  const unsigned short int side,
979  const unsigned int qp)
980 {
981  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes.find(elem->id()) !=
983  "Elem has no quadrature nodes!");
984  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes[elem->id()].find(side) !=
986  "Side has no quadrature nodes!");
987  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side].find(qp) !=
989  "qp not found on side!");
990 
991  return _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp];
992 }
993 
994 void
996 {
997  // Delete all the quadrature nodes
998  for (auto & it : _quadrature_nodes)
999  delete it.second;
1000 
1001  _quadrature_nodes.clear();
1003  _extra_bnd_nodes.clear();
1004 }
1005 
1006 BoundaryID
1007 MooseMesh::getBoundaryID(const BoundaryName & boundary_name) const
1008 {
1009  if (boundary_name == "ANY_BOUNDARY_ID")
1010  mooseError("Please use getBoundaryIDs() when passing \"ANY_BOUNDARY_ID\"");
1011 
1013  std::istringstream ss(boundary_name);
1014 
1015  if (!(ss >> id))
1016  id = getMesh().get_boundary_info().get_id_by_name(boundary_name);
1017 
1018  return id;
1019 }
1020 
1021 std::vector<BoundaryID>
1022 MooseMesh::getBoundaryIDs(const std::vector<BoundaryName> & boundary_name,
1023  bool generate_unknown) const
1024 {
1025  const BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1026  const std::map<BoundaryID, std::string> & sideset_map = boundary_info.get_sideset_name_map();
1027  const std::map<BoundaryID, std::string> & nodeset_map = boundary_info.get_nodeset_name_map();
1028 
1029  std::set<BoundaryID> boundary_ids = boundary_info.get_boundary_ids();
1030 
1031  // On a distributed mesh we may have boundary ids that only exist on
1032  // other processors.
1033  if (!this->getMesh().is_replicated())
1034  _communicator.set_union(boundary_ids);
1035 
1036  BoundaryID max_boundary_id = boundary_ids.empty() ? 0 : *(boundary_ids.rbegin());
1037 
1038  std::vector<BoundaryID> ids(boundary_name.size());
1039  for (unsigned int i = 0; i < boundary_name.size(); i++)
1040  {
1041  if (boundary_name[i] == "ANY_BOUNDARY_ID")
1042  {
1043  ids.assign(_mesh_boundary_ids.begin(), _mesh_boundary_ids.end());
1044  if (i)
1045  mooseWarning("You passed \"ANY_BOUNDARY_ID\" in addition to other boundary_names. This "
1046  "may be a logic error.");
1047  break;
1048  }
1049 
1050  BoundaryID id;
1051  std::istringstream ss(boundary_name[i]);
1052 
1053  if (!(ss >> id) || !ss.eof())
1054  {
1060  if (generate_unknown &&
1061  !MooseUtils::doesMapContainValue(sideset_map, std::string(boundary_name[i])) &&
1062  !MooseUtils::doesMapContainValue(nodeset_map, std::string(boundary_name[i])))
1063  id = ++max_boundary_id;
1064  else
1065  id = boundary_info.get_id_by_name(boundary_name[i]);
1066  }
1067 
1068  ids[i] = id;
1069  }
1070 
1071  return ids;
1072 }
1073 
1075 MooseMesh::getSubdomainID(const SubdomainName & subdomain_name) const
1076 {
1077  if (subdomain_name == "ANY_BLOCK_ID")
1078  mooseError("Please use getSubdomainIDs() when passing \"ANY_BLOCK_ID\"");
1079 
1081  std::istringstream ss(subdomain_name);
1082 
1083  if (!(ss >> id) || !ss.eof())
1084  id = getMesh().get_id_by_name(subdomain_name);
1085 
1086  return id;
1087 }
1088 
1089 std::vector<SubdomainID>
1090 MooseMesh::getSubdomainIDs(const std::vector<SubdomainName> & subdomain_name) const
1091 {
1092  std::vector<SubdomainID> ids(subdomain_name.size());
1093 
1094  for (unsigned int i = 0; i < subdomain_name.size(); i++)
1095  {
1096  if (subdomain_name[i] == "ANY_BLOCK_ID")
1097  {
1098  ids.assign(_mesh_subdomains.begin(), _mesh_subdomains.end());
1099  if (i)
1100  mooseWarning("You passed \"ANY_BLOCK_ID\" in addition to other block names. This may be a "
1101  "logic error.");
1102  break;
1103  }
1104 
1106  std::istringstream ss(subdomain_name[i]);
1107 
1108  if (!(ss >> id) || !ss.eof())
1109  id = getMesh().get_id_by_name(subdomain_name[i]);
1110 
1111  ids[i] = id;
1112  }
1113 
1114  return ids;
1115 }
1116 
1117 void
1118 MooseMesh::setSubdomainName(SubdomainID subdomain_id, SubdomainName name)
1119 {
1120  getMesh().subdomain_name(subdomain_id) = name;
1121 }
1122 
1123 const std::string &
1125 {
1126  return getMesh().subdomain_name(subdomain_id);
1127 }
1128 
1129 void
1130 MooseMesh::setBoundaryName(BoundaryID boundary_id, BoundaryName name)
1131 {
1132  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1133 
1134  std::vector<BoundaryID> side_boundaries;
1135  boundary_info.build_side_boundary_ids(side_boundaries);
1136 
1137  // We need to figure out if this boundary is a sideset or nodeset
1138  if (std::find(side_boundaries.begin(), side_boundaries.end(), boundary_id) !=
1139  side_boundaries.end())
1140  boundary_info.sideset_name(boundary_id) = name;
1141  else
1142  boundary_info.nodeset_name(boundary_id) = name;
1143 }
1144 
1145 const std::string &
1147 {
1148  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1149 
1150  std::vector<BoundaryID> side_boundaries;
1151  boundary_info.build_side_boundary_ids(side_boundaries);
1152 
1153  // We need to figure out if this boundary is a sideset or nodeset
1154  if (std::find(side_boundaries.begin(), side_boundaries.end(), boundary_id) !=
1155  side_boundaries.end())
1156  return boundary_info.get_sideset_name(boundary_id);
1157  else
1158  return boundary_info.get_nodeset_name(boundary_id);
1159 }
1160 
1161 // specialization for PointListAdaptor<MooseMesh::PeriodicNodeInfo>
1162 template <>
1163 inline const Point &
1165  const MooseMesh::PeriodicNodeInfo & item) const
1166 {
1167  return *(item.first);
1168 }
1169 
1170 void
1171 MooseMesh::buildPeriodicNodeMap(std::multimap<dof_id_type, dof_id_type> & periodic_node_map,
1172  unsigned int var_number,
1173  PeriodicBoundaries * pbs) const
1174 {
1175  TIME_SECTION(_build_periodic_node_map_timer);
1176 
1177  // clear existing map
1178  periodic_node_map.clear();
1179 
1180  // get periodic nodes
1181  std::vector<PeriodicNodeInfo> periodic_nodes;
1182  for (const auto & t : getMesh().get_boundary_info().build_node_list())
1183  {
1184  // unfortunately libMesh does not give us a pointer, so we have to look it up ourselves
1185  auto node = _mesh->node_ptr(std::get<0>(t));
1186  mooseAssert(node != nullptr,
1187  "libMesh::BoundaryInfo::build_node_list() returned an ID for a non-existing node");
1188  auto bc_id = std::get<1>(t);
1189  periodic_nodes.emplace_back(node, bc_id);
1190  }
1191 
1192  // sort by boundary id
1193  std::sort(periodic_nodes.begin(),
1194  periodic_nodes.end(),
1195  [](const PeriodicNodeInfo & a, const PeriodicNodeInfo & b) -> bool {
1196  return a.second > b.second;
1197  });
1198 
1199  // build kd-tree
1200  using KDTreeType = nanoflann::KDTreeSingleIndexAdaptor<
1201  nanoflann::L2_Simple_Adaptor<Real, PointListAdaptor<PeriodicNodeInfo>>,
1203  LIBMESH_DIM>;
1204  const unsigned int max_leaf_size = 20; // slightly affects runtime
1205  auto point_list =
1206  PointListAdaptor<PeriodicNodeInfo>(periodic_nodes.begin(), periodic_nodes.end());
1207  auto kd_tree = libmesh_make_unique<KDTreeType>(
1208  LIBMESH_DIM, point_list, nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size));
1209  mooseAssert(kd_tree != nullptr, "KDTree was not properly initialized.");
1210  kd_tree->buildIndex();
1211 
1212  // data structures for kd-tree search
1213  nanoflann::SearchParams search_params;
1214  std::vector<std::pair<std::size_t, Real>> ret_matches;
1215 
1216  // iterate over periodic nodes (boundary ids are in contiguous blocks)
1217  PeriodicBoundaryBase * periodic = nullptr;
1218  BoundaryID current_bc_id = BoundaryInfo::invalid_id;
1219  for (auto & pair : periodic_nodes)
1220  {
1221  // entering a new block of boundary IDs
1222  if (pair.second != current_bc_id)
1223  {
1224  current_bc_id = pair.second;
1225  periodic = pbs->boundary(current_bc_id);
1226  if (periodic && !periodic->is_my_variable(var_number))
1227  periodic = nullptr;
1228  }
1229 
1230  // variable is not periodic at this node, skip
1231  if (!periodic)
1232  continue;
1233 
1234  // clear result buffer
1235  ret_matches.clear();
1236 
1237  // id of the current node
1238  const auto id = pair.first->id();
1239 
1240  // position where we expect a periodic partner for the current node and boundary
1241  Point search_point = periodic->get_corresponding_pos(*pair.first);
1242 
1243  // search at the expected point
1244  kd_tree->radiusSearch(&(search_point)(0), libMesh::TOLERANCE, ret_matches, search_params);
1245  for (auto & match_pair : ret_matches)
1246  {
1247  const auto & match = periodic_nodes[match_pair.first];
1248  // add matched node if the boundary id is the corresponding id in the periodic pair
1249  if (match.second == periodic->pairedboundary)
1250  periodic_node_map.emplace(id, match.first->id());
1251  }
1252  }
1253 }
1254 
1255 void
1256 MooseMesh::buildPeriodicNodeSets(std::map<BoundaryID, std::set<dof_id_type>> & periodic_node_sets,
1257  unsigned int var_number,
1258  PeriodicBoundaries * pbs) const
1259 {
1260  TIME_SECTION(_build_periodic_node_sets_timer);
1261 
1262  periodic_node_sets.clear();
1263 
1264  // Loop over all the boundary nodes adding the periodic nodes to the appropriate set
1265  for (const auto & t : getMesh().get_boundary_info().build_node_list())
1266  {
1267  auto node_id = std::get<0>(t);
1268  auto bc_id = std::get<1>(t);
1269 
1270  // Is this current node on a known periodic boundary?
1271  if (periodic_node_sets.find(bc_id) != periodic_node_sets.end())
1272  periodic_node_sets[bc_id].insert(node_id);
1273  else // This still might be a periodic node but we just haven't seen this boundary_id yet
1274  {
1275  const PeriodicBoundaryBase * periodic = pbs->boundary(bc_id);
1276  if (periodic && periodic->is_my_variable(var_number))
1277  periodic_node_sets[bc_id].insert(node_id);
1278  }
1279  }
1280 }
1281 
1282 bool
1284 {
1286 
1288  return true;
1289 
1290  std::vector<Real> min(3, std::numeric_limits<Real>::max());
1291  std::vector<Real> max(3, std::numeric_limits<Real>::min());
1292  unsigned int dim = getMesh().mesh_dimension();
1293 
1294  // Find the bounding box of our mesh
1295  for (const auto & node : getMesh().node_ptr_range())
1296  // Check all coordinates, we don't know if this mesh might be lying in a higher dim even if the
1297  // mesh dimension is lower.
1298  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1299  {
1300  if ((*node)(i) < min[i])
1301  min[i] = (*node)(i);
1302  if ((*node)(i) > max[i])
1303  max[i] = (*node)(i);
1304  }
1305 
1306  this->comm().max(max);
1307  this->comm().min(min);
1308 
1309  _extreme_nodes.resize(8); // 2^LIBMESH_DIM
1310  // Now make sure that there are actual nodes at all of the extremes
1311  std::vector<bool> extreme_matches(8, false);
1312  std::vector<unsigned int> comp_map(3);
1313  for (const auto & node : getMesh().node_ptr_range())
1314  {
1315  // See if the current node is located at one of the extremes
1316  unsigned int coord_match = 0;
1317 
1318  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1319  {
1320  if (std::abs((*node)(i)-min[i]) < tol)
1321  {
1322  comp_map[i] = MIN;
1323  ++coord_match;
1324  }
1325  else if (std::abs((*node)(i)-max[i]) < tol)
1326  {
1327  comp_map[i] = MAX;
1328  ++coord_match;
1329  }
1330  }
1331 
1332  if (coord_match == LIBMESH_DIM) // Found a coordinate at one of the extremes
1333  {
1334  _extreme_nodes[comp_map[X] * 4 + comp_map[Y] * 2 + comp_map[Z]] = node;
1335  extreme_matches[comp_map[X] * 4 + comp_map[Y] * 2 + comp_map[Z]] = true;
1336  }
1337  }
1338 
1339  // See if we matched all of the extremes for the mesh dimension
1340  this->comm().max(extreme_matches);
1341  if (std::count(extreme_matches.begin(), extreme_matches.end(), true) == (1 << dim))
1342  _regular_orthogonal_mesh = true;
1343 
1344  // Set the bounds
1345  _bounds.resize(LIBMESH_DIM);
1346  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1347  {
1348  _bounds[i].resize(2);
1349  _bounds[i][MIN] = min[i];
1350  _bounds[i][MAX] = max[i];
1351  }
1352 
1353  return _regular_orthogonal_mesh;
1354 }
1355 
1356 void
1358 {
1359  TIME_SECTION(_detect_paired_sidesets_timer);
1360 
1361  // Loop over level-0 elements (since boundary condition information
1362  // is only directly stored for them) and find sidesets with normals
1363  // that point in the -x, +x, -y, +y, and -z, +z direction. If there
1364  // is a unique sideset id for each direction, then the paired
1365  // sidesets consist of (-x,+x), (-y,+y), (-z,+z). If there are
1366  // multiple sideset ids for a given direction, then we can't pick a
1367  // single pair for that direction. In that case, we'll just return
1368  // as was done in the original algorithm.
1369 
1370  // Points used for direction comparison
1371  const Point minus_x(-1, 0, 0), plus_x(1, 0, 0), minus_y(0, -1, 0), plus_y(0, 1, 0),
1372  minus_z(0, 0, -1), plus_z(0, 0, 1);
1373 
1374  // we need to test all element dimensions from dim down to 1
1375  const unsigned int dim = getMesh().mesh_dimension();
1376 
1377  // boundary id sets for elements of different dimensions
1378  std::vector<std::set<BoundaryID>> minus_x_ids(dim), plus_x_ids(dim), minus_y_ids(dim),
1379  plus_y_ids(dim), minus_z_ids(dim), plus_z_ids(dim);
1380 
1381  std::vector<std::unique_ptr<FEBase>> fe_faces(dim);
1382  std::vector<std::unique_ptr<QGauss>> qfaces(dim);
1383  for (unsigned side_dim = 0; side_dim < dim; ++side_dim)
1384  {
1385  // Face is assumed to be flat, therefore normal is assumed to be
1386  // constant over the face, therefore only compute it at 1 qp.
1387  qfaces[side_dim] = std::unique_ptr<QGauss>(new QGauss(side_dim, CONSTANT));
1388 
1389  // A first-order Lagrange FE for the face.
1390  fe_faces[side_dim] = FEBase::build(side_dim + 1, FEType(FIRST, LAGRANGE));
1391  fe_faces[side_dim]->attach_quadrature_rule(qfaces[side_dim].get());
1392  }
1393 
1394  // We need this to get boundary ids for each boundary face we encounter.
1395  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1396  std::vector<boundary_id_type> face_ids;
1397 
1398  for (auto & elem : as_range(getMesh().level_elements_begin(0), getMesh().level_elements_end(0)))
1399  {
1400  // dimension of the current element and its normals
1401  unsigned int side_dim = elem->dim() - 1;
1402  const std::vector<Point> & normals = fe_faces[side_dim]->get_normals();
1403 
1404  // loop over element sides
1405  for (unsigned int s = 0; s < elem->n_sides(); s++)
1406  {
1407  // If side is on the boundary
1408  if (elem->neighbor_ptr(s) == nullptr)
1409  {
1410  std::unique_ptr<Elem> side = elem->build_side_ptr(s);
1411 
1412  fe_faces[side_dim]->reinit(elem, s);
1413 
1414  // Get the boundary ID(s) for this side. If there is more
1415  // than 1 boundary id, then we already can't determine a
1416  // unique pairing of sides in this direction, but we'll just
1417  // keep going to keep the logic simple.
1418  boundary_info.boundary_ids(elem, s, face_ids);
1419 
1420  // x-direction faces
1421  if (normals[0].absolute_fuzzy_equals(minus_x))
1422  minus_x_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1423  else if (normals[0].absolute_fuzzy_equals(plus_x))
1424  plus_x_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1425 
1426  // y-direction faces
1427  else if (normals[0].absolute_fuzzy_equals(minus_y))
1428  minus_y_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1429  else if (normals[0].absolute_fuzzy_equals(plus_y))
1430  plus_y_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1431 
1432  // z-direction faces
1433  else if (normals[0].absolute_fuzzy_equals(minus_z))
1434  minus_z_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1435  else if (normals[0].absolute_fuzzy_equals(plus_z))
1436  plus_z_ids[side_dim].insert(face_ids.begin(), face_ids.end());
1437  }
1438  }
1439  }
1440 
1441  for (unsigned side_dim = 0; side_dim < dim; ++side_dim)
1442  {
1443  // If unique pairings were found, fill up the _paired_boundary data
1444  // structure with that information.
1445  if (minus_x_ids[side_dim].size() == 1 && plus_x_ids[side_dim].size() == 1)
1446  _paired_boundary.emplace_back(
1447  std::make_pair(*(minus_x_ids[side_dim].begin()), *(plus_x_ids[side_dim].begin())));
1448 
1449  if (minus_y_ids[side_dim].size() == 1 && plus_y_ids[side_dim].size() == 1)
1450  _paired_boundary.emplace_back(
1451  std::make_pair(*(minus_y_ids[side_dim].begin()), *(plus_y_ids[side_dim].begin())));
1452 
1453  if (minus_z_ids[side_dim].size() == 1 && plus_z_ids[side_dim].size() == 1)
1454  _paired_boundary.emplace_back(
1455  std::make_pair(*(minus_z_ids[side_dim].begin()), *(plus_z_ids[side_dim].begin())));
1456  }
1457 }
1458 
1459 Real
1460 MooseMesh::dimensionWidth(unsigned int component) const
1461 {
1462  return getMaxInDimension(component) - getMinInDimension(component);
1463 }
1464 
1465 Real
1466 MooseMesh::getMinInDimension(unsigned int component) const
1467 {
1468  mooseAssert(_mesh, "The MeshBase has not been constructed");
1469  mooseAssert(component < _bounds.size(), "Requested dimension out of bounds");
1470 
1471  return _bounds[component][MIN];
1472 }
1473 
1474 Real
1475 MooseMesh::getMaxInDimension(unsigned int component) const
1476 {
1477  mooseAssert(_mesh, "The MeshBase has not been constructed");
1478  mooseAssert(component < _bounds.size(), "Requested dimension out of bounds");
1479 
1480  return _bounds[component][MAX];
1481 }
1482 
1483 void
1484 MooseMesh::addPeriodicVariable(unsigned int var_num, BoundaryID primary, BoundaryID secondary)
1485 {
1487  return;
1488 
1489  _periodic_dim[var_num].resize(dimension());
1490 
1491  _half_range = Point(dimensionWidth(0) / 2.0, dimensionWidth(1) / 2.0, dimensionWidth(2) / 2.0);
1492 
1493  for (unsigned int component = 0; component < dimension(); ++component)
1494  {
1495  const std::pair<BoundaryID, BoundaryID> * boundary_ids = getPairedBoundaryMapping(component);
1496 
1497  if (boundary_ids != nullptr &&
1498  ((boundary_ids->first == primary && boundary_ids->second == secondary) ||
1499  (boundary_ids->first == secondary && boundary_ids->second == primary)))
1500  _periodic_dim[var_num][component] = true;
1501  }
1502 }
1503 
1504 bool
1505 MooseMesh::isTranslatedPeriodic(unsigned int nonlinear_var_num, unsigned int component) const
1506 {
1507  mooseAssert(component < dimension(), "Requested dimension out of bounds");
1508 
1509  if (_periodic_dim.find(nonlinear_var_num) != _periodic_dim.end())
1510  return _periodic_dim.at(nonlinear_var_num)[component];
1511  else
1512  return false;
1513 }
1514 
1516 MooseMesh::minPeriodicVector(unsigned int nonlinear_var_num, Point p, Point q) const
1517 {
1518  for (unsigned int i = 0; i < dimension(); ++i)
1519  {
1520  // check to see if we're closer in real or periodic space in x, y, and z
1521  if (isTranslatedPeriodic(nonlinear_var_num, i))
1522  {
1523  // Need to test order before differencing
1524  if (p(i) > q(i))
1525  {
1526  if (p(i) - q(i) > _half_range(i))
1527  p(i) -= _half_range(i) * 2;
1528  }
1529  else
1530  {
1531  if (q(i) - p(i) > _half_range(i))
1532  p(i) += _half_range(i) * 2;
1533  }
1534  }
1535  }
1536 
1537  return q - p;
1538 }
1539 
1540 Real
1541 MooseMesh::minPeriodicDistance(unsigned int nonlinear_var_num, Point p, Point q) const
1542 {
1543  return minPeriodicVector(nonlinear_var_num, p, q).norm();
1544 }
1545 
1546 const std::pair<BoundaryID, BoundaryID> *
1548 {
1550  mooseError("Trying to retrieve automatic paired mapping for a mesh that is not regular and "
1551  "orthogonal");
1552 
1553  mooseAssert(component < dimension(), "Requested dimension out of bounds");
1554 
1555  if (_paired_boundary.empty())
1557 
1558  if (component < _paired_boundary.size())
1559  return &_paired_boundary[component];
1560  else
1561  return nullptr;
1562 }
1563 
1564 void
1566 {
1568 
1569  std::map<ElemType, Elem *> canonical_elems;
1570 
1571  // First, loop over all elements and find a canonical element for each type
1572  // Doing it this way guarantees that this is going to work in parallel
1573  for (const auto & elem : getMesh().element_ptr_range()) // TODO: Thread this
1574  {
1575  ElemType type = elem->type();
1576 
1577  if (canonical_elems.find(type) ==
1578  canonical_elems.end()) // If we haven't seen this type of elem before save it
1579  canonical_elems[type] = elem;
1580  else
1581  {
1582  Elem * stored = canonical_elems[type];
1583  if (elem->id() < stored->id()) // Arbitrarily keep the one with a lower id
1584  canonical_elems[type] = elem;
1585  }
1586  }
1587  // Now build the maps using these templates
1588  // Note: This MUST be done NOT threaded!
1589  for (const auto & can_it : canonical_elems)
1590  {
1591  Elem * elem = can_it.second;
1592 
1593  // Need to do this just once to get the right qrules put in place
1594  assembly->setCurrentSubdomainID(elem->subdomain_id());
1595  assembly->reinit(elem);
1596  assembly->reinit(elem, 0);
1597  auto && qrule = assembly->writeableQRule();
1598  auto && qrule_face = assembly->writeableQRuleFace();
1599 
1600  // Volume to volume projection for refinement
1601  buildRefinementMap(*elem, *qrule, *qrule_face, -1, -1, -1);
1602 
1603  // Volume to volume projection for coarsening
1604  buildCoarseningMap(*elem, *qrule, *qrule_face, -1);
1605 
1606  // Map the sides of children
1607  for (unsigned int side = 0; side < elem->n_sides(); side++)
1608  {
1609  // Side to side for sides that match parent's sides
1610  buildRefinementMap(*elem, *qrule, *qrule_face, side, -1, side);
1611  buildCoarseningMap(*elem, *qrule, *qrule_face, side);
1612  }
1613 
1614  // Child side to parent volume mapping for "internal" child sides
1615  for (unsigned int child = 0; child < elem->n_children(); ++child)
1616  for (unsigned int side = 0; side < elem->n_sides();
1617  ++side) // Assume children have the same number of sides!
1618  if (!elem->is_child_on_side(child, side)) // Otherwise we already computed that map
1619  buildRefinementMap(*elem, *qrule, *qrule_face, -1, child, side);
1620  }
1621 }
1622 
1623 void
1625  QBase & qrule,
1626  QBase & qrule_face,
1627  int parent_side,
1628  int child,
1629  int child_side)
1630 {
1631  TIME_SECTION(_build_refinement_map_timer);
1632 
1633  if (child == -1) // Doing volume mapping or parent side mapping
1634  {
1635  mooseAssert(parent_side == child_side,
1636  "Parent side must match child_side if not passing a specific child!");
1637 
1638  std::pair<int, ElemType> the_pair(parent_side, elem.type());
1639 
1640  if (_elem_type_to_refinement_map.find(the_pair) != _elem_type_to_refinement_map.end())
1641  mooseError("Already built a qp refinement map!");
1642 
1643  std::vector<std::pair<unsigned int, QpMap>> coarsen_map;
1644  std::vector<std::vector<QpMap>> & refinement_map = _elem_type_to_refinement_map[the_pair];
1646  &elem, qrule, qrule_face, refinement_map, coarsen_map, parent_side, child, child_side);
1647  }
1648  else // Need to map a child side to parent volume qps
1649  {
1650  std::pair<int, int> child_pair(child, child_side);
1651 
1652  if (_elem_type_to_child_side_refinement_map.find(elem.type()) !=
1654  _elem_type_to_child_side_refinement_map[elem.type()].find(child_pair) !=
1656  mooseError("Already built a qp refinement map!");
1657 
1658  std::vector<std::pair<unsigned int, QpMap>> coarsen_map;
1659  std::vector<std::vector<QpMap>> & refinement_map =
1660  _elem_type_to_child_side_refinement_map[elem.type()][child_pair];
1662  &elem, qrule, qrule_face, refinement_map, coarsen_map, parent_side, child, child_side);
1663  }
1664 }
1665 
1666 const std::vector<std::vector<QpMap>> &
1667 MooseMesh::getRefinementMap(const Elem & elem, int parent_side, int child, int child_side)
1668 {
1669  if (child == -1) // Doing volume mapping or parent side mapping
1670  {
1671  mooseAssert(parent_side == child_side,
1672  "Parent side must match child_side if not passing a specific child!");
1673 
1674  std::pair<int, ElemType> the_pair(parent_side, elem.type());
1675 
1676  if (_elem_type_to_refinement_map.find(the_pair) == _elem_type_to_refinement_map.end())
1677  mooseError("Could not find a suitable qp refinement map!");
1678 
1679  return _elem_type_to_refinement_map[the_pair];
1680  }
1681  else // Need to map a child side to parent volume qps
1682  {
1683  std::pair<int, int> child_pair(child, child_side);
1684 
1685  if (_elem_type_to_child_side_refinement_map.find(elem.type()) ==
1687  _elem_type_to_child_side_refinement_map[elem.type()].find(child_pair) ==
1689  mooseError("Could not find a suitable qp refinement map!");
1690 
1691  return _elem_type_to_child_side_refinement_map[elem.type()][child_pair];
1692  }
1693 
1700 }
1701 
1702 void
1703 MooseMesh::buildCoarseningMap(const Elem & elem, QBase & qrule, QBase & qrule_face, int input_side)
1704 {
1705  TIME_SECTION(_build_coarsening_map_timer);
1706 
1707  std::pair<int, ElemType> the_pair(input_side, elem.type());
1708 
1709  if (_elem_type_to_coarsening_map.find(the_pair) != _elem_type_to_coarsening_map.end())
1710  mooseError("Already built a qp coarsening map!");
1711 
1712  std::vector<std::vector<QpMap>> refinement_map;
1713  std::vector<std::pair<unsigned int, QpMap>> & coarsen_map =
1714  _elem_type_to_coarsening_map[the_pair];
1715 
1716  // The -1 here is for a specific child. We don't do that for coarsening maps
1717  // Also note that we're always mapping the same side to the same side (which is guaranteed by
1718  // libMesh).
1720  &elem, qrule, qrule_face, refinement_map, coarsen_map, input_side, -1, input_side);
1721 
1728 }
1729 
1730 const std::vector<std::pair<unsigned int, QpMap>> &
1731 MooseMesh::getCoarseningMap(const Elem & elem, int input_side)
1732 {
1733  std::pair<int, ElemType> the_pair(input_side, elem.type());
1734 
1735  if (_elem_type_to_coarsening_map.find(the_pair) == _elem_type_to_coarsening_map.end())
1736  mooseError("Could not find a suitable qp refinement map!");
1737 
1738  return _elem_type_to_coarsening_map[the_pair];
1739 }
1740 
1741 void
1742 MooseMesh::mapPoints(const std::vector<Point> & from,
1743  const std::vector<Point> & to,
1744  std::vector<QpMap> & qp_map)
1745 {
1746  unsigned int n_from = from.size();
1747  unsigned int n_to = to.size();
1748 
1749  qp_map.resize(n_from);
1750 
1751  for (unsigned int i = 0; i < n_from; ++i)
1752  {
1753  const Point & from_point = from[i];
1754 
1755  QpMap & current_map = qp_map[i];
1756 
1757  for (unsigned int j = 0; j < n_to; ++j)
1758  {
1759  const Point & to_point = to[j];
1760  Real distance = (from_point - to_point).norm();
1761 
1762  if (distance < current_map._distance)
1763  {
1764  current_map._distance = distance;
1765  current_map._from = i;
1766  current_map._to = j;
1767  }
1768  }
1769  }
1770 }
1771 
1772 void
1773 MooseMesh::findAdaptivityQpMaps(const Elem * template_elem,
1774  QBase & qrule,
1775  QBase & qrule_face,
1776  std::vector<std::vector<QpMap>> & refinement_map,
1777  std::vector<std::pair<unsigned int, QpMap>> & coarsen_map,
1778  int parent_side,
1779  int child,
1780  int child_side)
1781 {
1782  TIME_SECTION(_find_adaptivity_qp_maps_timer);
1783 
1784  ReplicatedMesh mesh(_communicator);
1785  mesh.skip_partitioning(true);
1786 
1787  unsigned int dim = template_elem->dim();
1788  mesh.set_mesh_dimension(dim);
1789 
1790  for (unsigned int i = 0; i < template_elem->n_nodes(); ++i)
1791  mesh.add_point(template_elem->point(i));
1792 
1793  Elem * elem = mesh.add_elem(Elem::build(template_elem->type()).release());
1794 
1795  for (unsigned int i = 0; i < template_elem->n_nodes(); ++i)
1796  elem->set_node(i) = mesh.node_ptr(i);
1797 
1798  std::unique_ptr<FEBase> fe(FEBase::build(dim, FEType()));
1799  fe->get_phi();
1800  const std::vector<Point> & q_points_volume = fe->get_xyz();
1801 
1802  std::unique_ptr<FEBase> fe_face(FEBase::build(dim, FEType()));
1803  fe_face->get_phi();
1804  const std::vector<Point> & q_points_face = fe_face->get_xyz();
1805 
1806  fe->attach_quadrature_rule(&qrule);
1807  fe_face->attach_quadrature_rule(&qrule_face);
1808 
1809  // The current q_points
1810  const std::vector<Point> * q_points;
1811 
1812  if (parent_side != -1)
1813  {
1814  fe_face->reinit(elem, parent_side);
1815  q_points = &q_points_face;
1816  }
1817  else
1818  {
1819  fe->reinit(elem);
1820  q_points = &q_points_volume;
1821  }
1822 
1823  std::vector<Point> parent_ref_points;
1824 
1825  FEInterface::inverse_map(elem->dim(), FEType(), elem, *q_points, parent_ref_points);
1826  MeshRefinement mesh_refinement(mesh);
1827  mesh_refinement.uniformly_refine(1);
1828 
1829  std::map<unsigned int, std::vector<Point>> child_to_ref_points;
1830 
1831  unsigned int n_children = elem->n_children();
1832 
1833  refinement_map.resize(n_children);
1834 
1835  std::vector<unsigned int> children;
1836 
1837  if (child != -1) // Passed in a child explicitly
1838  children.push_back(child);
1839  else
1840  {
1841  children.resize(n_children);
1842  for (unsigned int child = 0; child < n_children; ++child)
1843  children[child] = child;
1844  }
1845 
1846  for (unsigned int i = 0; i < children.size(); ++i)
1847  {
1848  unsigned int child = children[i];
1849 
1850  if ((parent_side != -1 && !elem->is_child_on_side(child, parent_side)))
1851  continue;
1852 
1853  const Elem * child_elem = elem->child_ptr(child);
1854 
1855  if (child_side != -1)
1856  {
1857  fe_face->reinit(child_elem, child_side);
1858  q_points = &q_points_face;
1859  }
1860  else
1861  {
1862  fe->reinit(child_elem);
1863  q_points = &q_points_volume;
1864  }
1865 
1866  std::vector<Point> child_ref_points;
1867 
1868  FEInterface::inverse_map(elem->dim(), FEType(), elem, *q_points, child_ref_points);
1869  child_to_ref_points[child] = child_ref_points;
1870 
1871  std::vector<QpMap> & qp_map = refinement_map[child];
1872 
1873  // Find the closest parent_qp to each child_qp
1874  mapPoints(child_ref_points, parent_ref_points, qp_map);
1875  }
1876 
1877  coarsen_map.resize(parent_ref_points.size());
1878 
1879  // For each parent qp find the closest child qp
1880  for (unsigned int child = 0; child < n_children; child++)
1881  {
1882  if (parent_side != -1 && !elem->is_child_on_side(child, child_side))
1883  continue;
1884 
1885  std::vector<Point> & child_ref_points = child_to_ref_points[child];
1886 
1887  std::vector<QpMap> qp_map;
1888 
1889  // Find all of the closest points from parent_qp to _THIS_ child's qp
1890  mapPoints(parent_ref_points, child_ref_points, qp_map);
1891 
1892  // Check those to see if they are closer than what we currently have for each point
1893  for (unsigned int parent_qp = 0; parent_qp < parent_ref_points.size(); ++parent_qp)
1894  {
1895  std::pair<unsigned int, QpMap> & child_and_map = coarsen_map[parent_qp];
1896  unsigned int & closest_child = child_and_map.first;
1897  QpMap & closest_map = child_and_map.second;
1898 
1899  QpMap & current_map = qp_map[parent_qp];
1900 
1901  if (current_map._distance < closest_map._distance)
1902  {
1903  closest_child = child;
1904  closest_map = current_map;
1905  }
1906  }
1907  }
1908 }
1909 
1910 void
1911 MooseMesh::changeBoundaryId(const boundary_id_type old_id,
1912  const boundary_id_type new_id,
1913  bool delete_prev)
1914 {
1915  TIME_SECTION(_change_boundary_id_timer);
1916 
1917  // Get a reference to our BoundaryInfo object, we will use it several times below...
1918  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1919 
1920  // Container to catch ids passed back from BoundaryInfo
1921  std::vector<boundary_id_type> old_ids;
1922 
1923  // Only level-0 elements store BCs. Loop over them.
1924  for (auto & elem : as_range(getMesh().level_elements_begin(0), getMesh().level_elements_end(0)))
1925  {
1926  unsigned int n_sides = elem->n_sides();
1927  for (unsigned int s = 0; s != n_sides; ++s)
1928  {
1929  boundary_info.boundary_ids(elem, s, old_ids);
1930  if (std::find(old_ids.begin(), old_ids.end(), old_id) != old_ids.end())
1931  {
1932  std::vector<boundary_id_type> new_ids(old_ids);
1933  std::replace(new_ids.begin(), new_ids.end(), old_id, new_id);
1934  if (delete_prev)
1935  {
1936  boundary_info.remove_side(elem, s);
1937  boundary_info.add_side(elem, s, new_ids);
1938  }
1939  else
1940  boundary_info.add_side(elem, s, new_ids);
1941  }
1942  }
1943  }
1944 
1945  // Remove any remaining references to the old ID from the
1946  // BoundaryInfo object. This prevents things like empty sidesets
1947  // from showing up when printing information, etc.
1948  if (delete_prev)
1949  boundary_info.remove_id(old_id);
1950 }
1951 
1952 const RealVectorValue &
1954 {
1955  mooseAssert(_boundary_to_normal_map.get() != nullptr, "Boundary To Normal Map not built!");
1956 
1957  // Note: Boundaries that are not in the map (existing boundaries) will default
1958  // construct a new RealVectorValue - (x,y,z)=(0, 0, 0)
1959  return (*_boundary_to_normal_map)[id];
1960 }
1961 
1962 MooseMesh &
1964 {
1965  mooseError("MooseMesh::clone() is no longer supported, use MooseMesh::safeClone() instead.");
1966 }
1967 
1968 std::unique_ptr<MeshBase>
1970 {
1971  switch (_parallel_type)
1972  {
1973  case ParallelType::DEFAULT:
1974  // The user did not specify 'parallel_type = XYZ' in the input file,
1975  // so we allow the --distributed-mesh command line arg to possibly turn
1976  // on DistributedMesh. If the command line arg is not present, we pick ReplicatedMesh.
1978  _use_distributed_mesh = true;
1979  break;
1983  break;
1985  _use_distributed_mesh = true;
1986  break;
1987  }
1988 
1989  // If the user specifies 'nemesis = true' in the Mesh block, or they are using --use-split,
1990  // we must use DistributedMesh.
1991  if (_is_nemesis || _app.isUseSplit())
1992  _use_distributed_mesh = true;
1993 
1994  unsigned dim = getParam<MooseEnum>("dim");
1995 
1996  std::unique_ptr<MeshBase> mesh;
1998  {
1999  if (override_type == ParallelType::REPLICATED)
2000  mooseError("The requested override_type of \"Replicated\" may not be used when MOOSE is "
2001  "running with a DistributedMesh");
2002 
2003  mesh = libmesh_make_unique<DistributedMesh>(_communicator, dim);
2004  if (_partitioner_name != "default" && _partitioner_name != "parmetis")
2005  {
2006  _partitioner_name = "parmetis";
2007  _partitioner_overridden = true;
2008  }
2009  }
2010  else
2011  {
2012  if (override_type == ParallelType::DISTRIBUTED)
2013  mooseError("The requested override_type of \"Distributed\" may not be used when MOOSE is "
2014  "running with a ReplicatedMesh");
2015 
2016  mesh = libmesh_make_unique<ReplicatedMesh>(_communicator, dim);
2017  }
2018 
2019  if (!getParam<bool>("allow_renumbering"))
2020  mesh->allow_renumbering(false);
2021 
2022  return mesh;
2023 }
2024 
2025 void
2026 MooseMesh::setMeshBase(std::unique_ptr<MeshBase> mesh_base)
2027 {
2028  _mesh = std::move(mesh_base);
2029 }
2030 
2031 void
2033 {
2040  if (!_mesh)
2042 
2044  mooseError("You cannot use the mesh splitter capability with DistributedMesh!");
2045 
2046  TIME_SECTION(_init_timer);
2047 
2049  {
2050  // Check of partitioner is supplied (not allowed if custom partitioner is used)
2051  if (!parameters().isParamSetByAddParam("partitioner"))
2052  mooseError("If partitioner block is provided, partitioner keyword cannot be used!");
2053  // Set custom partitioner
2054  if (!_custom_partitioner.get())
2055  mooseError("Custom partitioner requested but not set!");
2056  getMesh().partitioner().reset(_custom_partitioner.release());
2057  }
2058  else
2059  {
2060  // Set standard partitioner
2061  // Set the partitioner based on partitioner name
2062  switch (_partitioner_name)
2063  {
2064  case -3: // default
2065  // We'll use the default partitioner, but notify the user of which one is being used...
2067  _partitioner_name = "parmetis";
2068  else
2069  _partitioner_name = "metis";
2070  break;
2071 
2072  // No need to explicitily create the metis or parmetis partitioners,
2073  // They are the default for serial and parallel mesh respectively
2074  case -2: // metis
2075  case -1: // parmetis
2076  break;
2077 
2078  case 0: // linear
2079  getMesh().partitioner().reset(new LinearPartitioner);
2080  break;
2081  case 1: // centroid
2082  {
2083  if (!isParamValid("centroid_partitioner_direction"))
2084  mooseError("If using the centroid partitioner you _must_ specify "
2085  "centroid_partitioner_direction!");
2086 
2087  MooseEnum direction = getParam<MooseEnum>("centroid_partitioner_direction");
2088 
2089  if (direction == "x")
2090  getMesh().partitioner().reset(new CentroidPartitioner(CentroidPartitioner::X));
2091  else if (direction == "y")
2092  getMesh().partitioner().reset(new CentroidPartitioner(CentroidPartitioner::Y));
2093  else if (direction == "z")
2094  getMesh().partitioner().reset(new CentroidPartitioner(CentroidPartitioner::Z));
2095  else if (direction == "radial")
2096  getMesh().partitioner().reset(new CentroidPartitioner(CentroidPartitioner::RADIAL));
2097  break;
2098  }
2099  case 2: // hilbert_sfc
2100  getMesh().partitioner().reset(new HilbertSFCPartitioner);
2101  break;
2102  case 3: // morton_sfc
2103  getMesh().partitioner().reset(new MortonSFCPartitioner);
2104  break;
2105  }
2106  }
2107 
2109  {
2110  // Some partitioners are not idempotent. Some recovery data
2111  // files require partitioning to match mesh partitioning. This
2112  // means that, when recovering, we can't safely repartition.
2113  const bool skip_partitioning_later = getMesh().skip_partitioning();
2114  getMesh().skip_partitioning(true);
2115  const bool allow_renumbering_later = getMesh().allow_renumbering();
2116  getMesh().allow_renumbering(false);
2117 
2118  // For now, only read the recovery mesh on the Ultimate Master..
2119  // sub-apps need to just build their mesh like normal
2120  {
2121  TIME_SECTION(_read_recovered_mesh_timer);
2122  getMesh().read(_app.getRecoverFileBase() + "_mesh." + _app.getRecoverFileSuffix());
2123  }
2124 
2125  getMesh().allow_renumbering(allow_renumbering_later);
2126  getMesh().skip_partitioning(skip_partitioning_later);
2127  }
2128  else // Normally just build the mesh
2129  buildMesh();
2130 }
2131 
2132 unsigned int
2134 {
2135  return getMesh().mesh_dimension();
2136 }
2137 
2138 unsigned int
2140 {
2141  const Real abs_zero = 1e-12;
2142 
2143  // See if the mesh is completely containd in the z and y planes to calculate effective spatial dim
2144  for (unsigned int dim = LIBMESH_DIM; dim >= 1; --dim)
2145  if (dimensionWidth(dim - 1) >= abs_zero)
2146  return dim;
2147 
2148  // If we get here, we have a 1D mesh on the x-axis.
2149  return 1;
2150 }
2151 
2152 std::vector<BoundaryID>
2153 MooseMesh::getBoundaryIDs(const Elem * const elem, const unsigned short int side) const
2154 {
2155  std::vector<BoundaryID> ids;
2156  getMesh().get_boundary_info().boundary_ids(elem, side, ids);
2157  return ids;
2158 }
2159 
2160 const std::set<BoundaryID> &
2162 {
2163  return getMesh().get_boundary_info().get_boundary_ids();
2164 }
2165 
2166 void
2168 {
2170  getMesh().get_boundary_info().build_node_list_from_side_list();
2171 }
2172 
2173 void
2174 MooseMesh::buildSideList(std::vector<dof_id_type> & el,
2175  std::vector<unsigned short int> & sl,
2176  std::vector<boundary_id_type> & il)
2177 {
2178 #ifdef LIBMESH_ENABLE_DEPRECATED
2179  mooseDeprecated("The version of MooseMesh::buildSideList() taking three arguments is "
2180  "deprecated, call the version that returns a vector of tuples instead.");
2181  getMesh().get_boundary_info().build_side_list(el, sl, il);
2182 #else
2183  libmesh_ignore(el);
2184  libmesh_ignore(sl);
2185  libmesh_ignore(il);
2186  mooseError("The version of MooseMesh::buildSideList() taking three "
2187  "arguments is not available in your version of libmesh, call the "
2188  "version that returns a vector of tuples instead.");
2189 #endif
2190 }
2191 
2192 std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>>
2194 {
2195  return getMesh().get_boundary_info().build_side_list();
2196 }
2197 
2198 unsigned int
2199 MooseMesh::sideWithBoundaryID(const Elem * const elem, const BoundaryID boundary_id) const
2200 {
2201  return getMesh().get_boundary_info().side_with_boundary_id(elem, boundary_id);
2202 }
2203 
2204 MeshBase::const_node_iterator
2206 {
2207  return getMesh().local_nodes_begin();
2208 }
2209 
2210 MeshBase::const_node_iterator
2212 {
2213  return getMesh().local_nodes_end();
2214 }
2215 
2216 MeshBase::const_element_iterator
2218 {
2219  return getMesh().active_local_elements_begin();
2220 }
2221 
2222 const MeshBase::const_element_iterator
2224 {
2225  return getMesh().active_local_elements_end();
2226 }
2227 
2228 dof_id_type
2230 {
2231  return getMesh().n_nodes();
2232 }
2233 
2234 dof_id_type
2236 {
2237  return getMesh().n_elem();
2238 }
2239 
2240 dof_id_type
2242 {
2243  return getMesh().max_node_id();
2244 }
2245 
2246 dof_id_type
2248 {
2249  return getMesh().max_elem_id();
2250 }
2251 
2252 Elem *
2253 MooseMesh::elem(const dof_id_type i)
2254 {
2255  mooseDeprecated("MooseMesh::elem() is deprecated, please use MooseMesh::elemPtr() instead");
2256  return elemPtr(i);
2257 }
2258 
2259 const Elem *
2260 MooseMesh::elem(const dof_id_type i) const
2261 {
2262  mooseDeprecated("MooseMesh::elem() is deprecated, please use MooseMesh::elemPtr() instead");
2263  return elemPtr(i);
2264 }
2265 
2266 Elem *
2267 MooseMesh::elemPtr(const dof_id_type i)
2268 {
2269  return getMesh().elem_ptr(i);
2270 }
2271 
2272 const Elem *
2273 MooseMesh::elemPtr(const dof_id_type i) const
2274 {
2275  return getMesh().elem_ptr(i);
2276 }
2277 
2278 Elem *
2279 MooseMesh::queryElemPtr(const dof_id_type i)
2280 {
2281  return getMesh().query_elem_ptr(i);
2282 }
2283 
2284 const Elem *
2285 MooseMesh::queryElemPtr(const dof_id_type i) const
2286 {
2287  return getMesh().query_elem_ptr(i);
2288 }
2289 
2290 bool
2292 {
2293  return _is_prepared;
2294 }
2295 
2296 void
2298 {
2299  _is_prepared = state;
2300 
2305  if (!state)
2306  _regular_orthogonal_mesh = false;
2307 }
2308 
2309 void
2311 {
2312  _needs_prepare_for_use = true;
2313 }
2314 
2315 const std::set<SubdomainID> &
2317 {
2318  return _mesh_subdomains;
2319 }
2320 
2321 const std::set<BoundaryID> &
2323 {
2324  return _mesh_boundary_ids;
2325 }
2326 
2327 const std::set<BoundaryID> &
2329 {
2330  return _mesh_sideset_ids;
2331 }
2332 
2333 const std::set<BoundaryID> &
2335 {
2336  return _mesh_nodeset_ids;
2337 }
2338 
2339 void
2340 MooseMesh::setMeshBoundaryIDs(std::set<BoundaryID> boundary_IDs)
2341 {
2342  _mesh_boundary_ids = boundary_IDs;
2343 }
2344 
2345 void
2347  std::unique_ptr<std::map<BoundaryID, RealVectorValue>> boundary_map)
2348 {
2349  _boundary_to_normal_map = std::move(boundary_map);
2350 }
2351 
2352 void
2353 MooseMesh::setBoundaryToNormalMap(std::map<BoundaryID, RealVectorValue> * boundary_map)
2354 {
2355  mooseDeprecated("setBoundaryToNormalMap(std::map<BoundaryID, RealVectorValue> * boundary_map) is "
2356  "deprecated, use the unique_ptr version instead");
2357  _boundary_to_normal_map.reset(boundary_map);
2358 }
2359 
2360 unsigned int
2362 {
2363  return _uniform_refine_level;
2364 }
2365 
2366 void
2368 {
2369  _uniform_refine_level = level;
2370 }
2371 
2372 void
2374 {
2375  _ghosted_boundaries.insert(boundary_id);
2376 }
2377 
2378 void
2379 MooseMesh::setGhostedBoundaryInflation(const std::vector<Real> & inflation)
2380 {
2381  _ghosted_boundaries_inflation = inflation;
2382 }
2383 
2384 const std::set<unsigned int> &
2386 {
2387  return _ghosted_boundaries;
2388 }
2389 
2390 const std::vector<Real> &
2392 {
2394 }
2395 
2396 namespace // Anonymous namespace for helpers
2397 {
2398 // A class for templated methods that expect output iterator
2399 // arguments, which adds objects to the Mesh.
2400 // Although any mesh_inserter_iterator can add any object, we
2401 // template it around object type so that type inference and
2402 // iterator_traits will work.
2403 // This object specifically is used to insert extra ghost elems into the mesh
2404 template <typename T>
2405 struct extra_ghost_elem_inserter : std::iterator<std::output_iterator_tag, T>
2406 {
2407  extra_ghost_elem_inserter(DistributedMesh & m) : mesh(m) {}
2408 
2409  void operator=(const Elem * e) { mesh.add_extra_ghost_elem(const_cast<Elem *>(e)); }
2410 
2411  void operator=(Node * n) { mesh.add_node(n); }
2412 
2413  void operator=(Point * p) { mesh.add_point(*p); }
2414 
2415  extra_ghost_elem_inserter & operator++() { return *this; }
2416 
2417  extra_ghost_elem_inserter operator++(int) { return extra_ghost_elem_inserter(*this); }
2418 
2419  // We don't return a reference-to-T here because we don't want to
2420  // construct one or have any of its methods called. We just want
2421  // to allow the returned object to be able to do mesh insertions
2422  // with operator=().
2423  extra_ghost_elem_inserter & operator*() { return *this; }
2424 
2425 private:
2426  DistributedMesh & mesh;
2427 };
2428 
2439 struct CompareElemsByLevel
2440 {
2441  bool operator()(const Elem * a, const Elem * b) const
2442  {
2443  libmesh_assert(a);
2444  libmesh_assert(b);
2445  const unsigned int al = a->level(), bl = b->level();
2446  const dof_id_type aid = a->id(), bid = b->id();
2447 
2448  return (al == bl) ? aid < bid : al < bl;
2449  }
2450 };
2451 
2452 } // anonymous namespace
2453 
2454 void
2456 {
2457  // No need to do this if using a serial mesh
2458  if (!_use_distributed_mesh)
2459  return;
2460 
2461  TIME_SECTION(_ghost_ghosted_boundaries_timer);
2462 
2463  DistributedMesh & mesh = dynamic_cast<DistributedMesh &>(getMesh());
2464 
2465  // We would like to clear ghosted elements that were added by
2466  // previous invocations of this method; however we can't do so
2467  // simply without also clearing ghosted elements that were added by
2468  // other code; e.g. OversampleOutput. So for now we'll just
2469  // swallow the inefficiency that can come from leaving unnecessary
2470  // elements ghosted after AMR.
2471  // mesh.clear_extra_ghost_elems();
2472 
2473  std::set<const Elem *, CompareElemsByLevel> boundary_elems_to_ghost;
2474  std::set<Node *> connected_nodes_to_ghost;
2475 
2476  std::vector<const Elem *> family_tree;
2477 
2478  for (const auto & t : mesh.get_boundary_info().build_side_list())
2479  {
2480  auto elem_id = std::get<0>(t);
2481  auto bc_id = std::get<2>(t);
2482 
2483  if (_ghosted_boundaries.find(bc_id) != _ghosted_boundaries.end())
2484  {
2485  Elem * elem = mesh.elem_ptr(elem_id);
2486 
2487 #ifdef LIBMESH_ENABLE_AMR
2488  elem->family_tree(family_tree);
2489  Elem * parent = elem->parent();
2490  while (parent)
2491  {
2492  family_tree.push_back(parent);
2493  parent = parent->parent();
2494  }
2495 #else
2496  family_tree.clear();
2497  family_tree.push_back(elem);
2498 #endif
2499  for (const auto & felem : family_tree)
2500  {
2501  boundary_elems_to_ghost.insert(felem);
2502 
2503  // The entries of connected_nodes_to_ghost need to be
2504  // non-constant, so that they will work in things like
2505  // UpdateDisplacedMeshThread. The container returned by
2506  // family_tree contains const Elems even when the Elem
2507  // it is called on is non-const, so once that interface
2508  // gets fixed we can remove this const_cast.
2509  for (unsigned int n = 0; n < felem->n_nodes(); ++n)
2510  connected_nodes_to_ghost.insert(const_cast<Node *>(felem->node_ptr(n)));
2511  }
2512  }
2513  }
2514 
2515  mesh.comm().allgather_packed_range(&mesh,
2516  connected_nodes_to_ghost.begin(),
2517  connected_nodes_to_ghost.end(),
2518  extra_ghost_elem_inserter<Node>(mesh));
2519  mesh.comm().allgather_packed_range(&mesh,
2520  boundary_elems_to_ghost.begin(),
2521  boundary_elems_to_ghost.end(),
2522  extra_ghost_elem_inserter<Elem>(mesh));
2523 }
2524 
2525 unsigned int
2527 {
2528  return _patch_size;
2529 }
2530 
2531 void
2533 {
2534  _patch_update_strategy = patch_update_strategy;
2535 }
2536 
2537 const Moose::PatchUpdateType &
2539 {
2540  return _patch_update_strategy;
2541 }
2542 
2543 BoundingBox
2544 MooseMesh::getInflatedProcessorBoundingBox(Real inflation_multiplier) const
2545 {
2546  // Grab a bounding box to speed things up. Note that
2547  // local_bounding_box is *not* equivalent to processor_bounding_box
2548  // with processor_id() except in serial.
2549  BoundingBox bbox = MeshTools::create_local_bounding_box(getMesh());
2550 
2551  // Inflate the bbox just a bit to deal with roundoff
2552  // Adding 1% of the diagonal size in each direction on each end
2553  Real inflation_amount = inflation_multiplier * (bbox.max() - bbox.min()).norm();
2554  Point inflation(inflation_amount, inflation_amount, inflation_amount);
2555 
2556  bbox.first -= inflation; // min
2557  bbox.second += inflation; // max
2558 
2559  return bbox;
2560 }
2561 
2562 MooseMesh::operator libMesh::MeshBase &() { return getMesh(); }
2563 
2564 MooseMesh::operator const libMesh::MeshBase &() const { return getMesh(); }
2565 
2566 MeshBase &
2568 {
2569  mooseAssert(_mesh, "Mesh hasn't been created");
2570  return *_mesh;
2571 }
2572 
2573 const MeshBase &
2575 {
2576  mooseAssert(_mesh, "Mesh hasn't been created");
2577  return *_mesh;
2578 }
2579 
2580 ExodusII_IO *
2582 {
2583  // TODO: Implement or remove
2584  return nullptr;
2585 }
2586 
2587 void
2588 MooseMesh::printInfo(std::ostream & os) const
2589 {
2590  getMesh().print_info(os);
2591 }
2592 
2593 const std::vector<dof_id_type> &
2594 MooseMesh::getNodeList(boundary_id_type nodeset_id) const
2595 {
2596  std::map<boundary_id_type, std::vector<dof_id_type>>::const_iterator it =
2597  _node_set_nodes.find(nodeset_id);
2598 
2599  if (it == _node_set_nodes.end())
2600  {
2601  // On a distributed mesh we might not know about a remote nodeset,
2602  // so we'll return an empty vector and hope the nodeset exists
2603  // elsewhere.
2604  if (!getMesh().is_serial())
2605  {
2606  static const std::vector<dof_id_type> empty_vec;
2607  return empty_vec;
2608  }
2609  // On a replicated mesh we should know about every nodeset and if
2610  // we're asked for one that doesn't exist then it must be a bug.
2611  else
2612  {
2613  mooseError("Unable to nodeset ID: ", nodeset_id, '.');
2614  }
2615  }
2616 
2617  return it->second;
2618 }
2619 
2620 const std::set<BoundaryID> &
2622 {
2623  std::map<SubdomainID, std::set<BoundaryID>>::const_iterator it =
2624  _subdomain_boundary_ids.find(subdomain_id);
2625 
2626  if (it == _subdomain_boundary_ids.end())
2627  mooseError("Unable to find subdomain ID: ", subdomain_id, '.');
2628 
2629  return it->second;
2630 }
2631 
2632 bool
2633 MooseMesh::isBoundaryNode(dof_id_type node_id) const
2634 {
2635  bool found_node = false;
2636  for (const auto & it : _bnd_node_ids)
2637  {
2638  if (it.second.find(node_id) != it.second.end())
2639  {
2640  found_node = true;
2641  break;
2642  }
2643  }
2644  return found_node;
2645 }
2646 
2647 bool
2648 MooseMesh::isBoundaryNode(dof_id_type node_id, BoundaryID bnd_id) const
2649 {
2650  bool found_node = false;
2651  std::map<boundary_id_type, std::set<dof_id_type>>::const_iterator it = _bnd_node_ids.find(bnd_id);
2652  if (it != _bnd_node_ids.end())
2653  if (it->second.find(node_id) != it->second.end())
2654  found_node = true;
2655  return found_node;
2656 }
2657 
2658 bool
2659 MooseMesh::isBoundaryElem(dof_id_type elem_id) const
2660 {
2661  bool found_elem = false;
2662  for (const auto & it : _bnd_elem_ids)
2663  {
2664  if (it.second.find(elem_id) != it.second.end())
2665  {
2666  found_elem = true;
2667  break;
2668  }
2669  }
2670  return found_elem;
2671 }
2672 
2673 bool
2674 MooseMesh::isBoundaryElem(dof_id_type elem_id, BoundaryID bnd_id) const
2675 {
2676  bool found_elem = false;
2677  auto it = _bnd_elem_ids.find(bnd_id);
2678  if (it != _bnd_elem_ids.end())
2679  if (it->second.find(elem_id) != it->second.end())
2680  found_elem = true;
2681  return found_elem;
2682 }
2683 
2684 void
2685 MooseMesh::errorIfDistributedMesh(std::string name) const
2686 {
2688  mooseError("Cannot use ",
2689  name,
2690  " with DistributedMesh!\n",
2691  "Consider specifying parallel_type = 'replicated' in your input file\n",
2692  "to prevent it from being run with DistributedMesh.");
2693 }
2694 
2695 void
2696 MooseMesh::setCustomPartitioner(Partitioner * partitioner)
2697 {
2698  _custom_partitioner = partitioner->clone();
2699 }
2700 
2701 bool
2703 {
2705 }
2706 
2707 bool
2709 {
2710  bool mesh_has_second_order_elements = false;
2711  for (auto it = activeLocalElementsBegin(), end = activeLocalElementsEnd(); it != end; ++it)
2712  if ((*it)->default_order() == SECOND)
2713  {
2714  mesh_has_second_order_elements = true;
2715  break;
2716  }
2717 
2718  // We checked our local elements, so take the max over all processors.
2719  comm().max(mesh_has_second_order_elements);
2720  return mesh_has_second_order_elements;
2721 }
2722 
2723 void
2725 {
2727 }
2728 
2729 std::unique_ptr<PointLocatorBase>
2731 {
2732  return getMesh().sub_point_locator();
2733 }
ParallelType _parallel_type
Can be set to DISTRIBUTED, REPLICATED, or DEFAULT.
Definition: MooseMesh.h:859
virtual bnd_node_iterator bndNodesEnd()
Definition: MooseMesh.C:868
virtual bnd_elem_iterator bndElemsEnd()
Definition: MooseMesh.C:884
PerfID _build_refinement_map_timer
Definition: MooseMesh.h:1155
std::vector< std::vector< Real > > _bounds
The bounds in each dimension of the mesh for regular orthogonal meshes.
Definition: MooseMesh.h:1012
std::set< Node * > _semilocal_node_list
Used for generating the semilocal node range.
Definition: MooseMesh.h:921
std::map< dof_id_type, Node * > _quadrature_nodes
Definition: MooseMesh.h:979
virtual Real getMaxInDimension(unsigned int component) const
Definition: MooseMesh.C:1475
bool _node_to_elem_map_built
Definition: MooseMesh.h:938
std::vector< Node * > _extreme_nodes
A vector containing the nodes at the corners of a regular orthogonal mesh.
Definition: MooseMesh.h:1034
Node * addQuadratureNode(const Elem *elem, const unsigned short int side, const unsigned int qp, BoundaryID bid, const Point &point)
Adds a fictitious "QuadratureNode".
Definition: MooseMesh.C:925
std::vector< std::tuple< dof_id_type, unsigned short int, boundary_id_type > > buildSideList()
As above, but uses the non-deprecated std::tuple interface.
Definition: MooseMesh.C:2193
MeshBase::const_element_iterator activeLocalElementsBegin()
Calls active_local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2217
PerfID _build_coarsening_map_timer
Definition: MooseMesh.h:1156
const std::set< BoundaryID > & meshNodesetIds() const
Returns a read-only reference to the set of nodesets currently present in the Mesh.
Definition: MooseMesh.C:2334
virtual Real getMinInDimension(unsigned int component) const
Returns the min or max of the requested dimension respectively.
Definition: MooseMesh.C:1466
ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
Definition: MooseMesh.C:737
virtual void onMeshChanged()
Declares a callback function that is executed at the conclusion of meshChanged(). ...
Definition: MooseMesh.C:516
bool prepared() const
Setter/getter for the _is_prepared flag.
Definition: MooseMesh.C:2291
void needsPrepareForUse()
If this method is called, we will call libMesh&#39;s prepare_for_use method when we call Moose&#39;s prepare ...
Definition: MooseMesh.C:2310
const std::set< BoundaryID > & getBoundaryIDs() const
Returns a const reference to a set of all user-specified boundary IDs.
Definition: MooseMesh.C:2161
bool _is_nemesis
True if a Nemesis Mesh was read in.
Definition: MooseMesh.h:899
virtual MooseMesh & clone() const
Clone method.
Definition: MooseMesh.C:1963
MetaPhysicL::DualNumber< T, D > abs(const MetaPhysicL::DualNumber< T, D > &in)
bool isCustomPartitionerRequested() const
Setter and getter for _custom_partitioner_requested.
Definition: MooseMesh.C:2702
A class for creating restricted objects.
Definition: Restartable.h:29
unsigned int _uniform_refine_level
The level of uniform refinement requested (set to zero if AMR is disabled)
Definition: MooseMesh.h:893
const std::set< BoundaryID > & getSubdomainBoundaryIds(SubdomainID subdomain_id) const
Get the list of boundary ids associated with the given subdomain id.
Definition: MooseMesh.C:2621
Helper class for sorting Boundary Nodes so that we always get the same order of application for bound...
Definition: MooseMesh.C:609
RealVectorValue _half_range
A convenience vector used to hold values in each dimension representing half of the range...
Definition: MooseMesh.h:1031
Keeps track of stuff related to assembling.
Definition: Assembly.h:62
virtual ~MooseMesh()
Definition: MooseMesh.C:304
void freeBndElems()
Definition: MooseMesh.C:330
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:2267
PerfID _build_node_list_timer
Definition: MooseMesh.h:1141
The definition of the bnd_elem_iterator struct.
Definition: MooseMesh.h:1214
void mooseWarning(Args &&... args) const
Definition: MooseObject.h:155
bool isBoundaryNode(dof_id_type node_id) const
Returns true if the requested node is in the list of boundary nodes, false otherwise.
Definition: MooseMesh.C:2633
const std::map< dof_id_type, std::vector< dof_id_type > > & nodeToActiveSemilocalElemMap()
If not already created, creates a map from every node to all active semilocal elements to which they ...
Definition: MooseMesh.C:712
bool _is_prepared
True if prepare has been called on the mesh.
Definition: MooseMesh.h:902
const std::set< BoundaryID > & meshSidesetIds() const
Returns a read-only reference to the set of sidesets currently present in the Mesh.
Definition: MooseMesh.C:2328
VectorValue< Real > RealVectorValue
Definition: Assembly.h:31
const BoundaryID INVALID_BOUNDARY_ID
Definition: MooseTypes.C:18
bool _custom_partitioner_requested
Definition: MooseMesh.h:877
const std::unordered_map< boundary_id_type, std::unordered_set< dof_id_type > > & getBoundariesToElems() const
Returns a map of boundaries to elements.
Definition: MooseMesh.C:814
std::map< dof_id_type, std::vector< dof_id_type > > _node_to_elem_map
A map of all of the current nodes to the elements that they are connected to.
Definition: MooseMesh.h:937
bool isSemiLocal(Node *const node) const
Returns true if the node is semi-local.
Definition: MooseMesh.C:600
std::string getRecoverFileBase()
The file_base for the recovery file.
Definition: MooseApp.h:363
RealVectorValue minPeriodicVector(unsigned int nonlinear_var_num, Point p, Point q) const
This function returns the minimum vector between two points on the mesh taking into account periodici...
Definition: MooseMesh.C:1516
StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > ConstBndElemRange
Definition: MooseMesh.h:1259
const std::vector< std::pair< unsigned int, QpMap > > & getCoarseningMap(const Elem &elem, int input_side)
Get the coarsening map for a given element type.
Definition: MooseMesh.C:1731
bool doesMapContainValue(const std::map< T1, T2 > &the_map, const T2 &value)
This routine is a simple helper function for searching a map by values instead of keys...
Definition: MooseUtils.h:208
MeshBase::const_node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2205
std::unordered_map< boundary_id_type, std::unordered_set< dof_id_type > > _bnd_elem_ids
Map of set of elem IDs connected to each boundary.
Definition: MooseMesh.h:977
const std::vector< std::vector< QpMap > > & getRefinementMap(const Elem &elem, int parent_side, int child, int child_side)
Get the refinement map for a given element type.
Definition: MooseMesh.C:1667
const std::string & getBoundaryName(BoundaryID boundary_id)
Return the name of the boundary given the id.
Definition: MooseMesh.C:1146
void cacheChangedLists()
Cache information about what elements were refined and coarsened in the previous step.
Definition: MooseMesh.C:521
bool isTranslatedPeriodic(unsigned int nonlinear_var_num, unsigned int component) const
Returns whether this generated mesh is periodic in the given dimension for the given variable...
Definition: MooseMesh.C:1505
const Point & getPoint(const PointObject &item) const
get a Point reference from the PointData object at index idx in the list
The definition of the bnd_node_iterator struct.
Definition: MooseMesh.h:1171
virtual const Node * queryNodePtr(const dof_id_type i) const
Definition: MooseMesh.C:472
Helper object for holding qp mapping info.
Definition: MooseMesh.h:55
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const std::vector< Real > & getGhostedBoundaryInflation() const
Return a writable reference to the _ghosted_boundaries_inflation vector.
Definition: MooseMesh.C:2391
std::unique_ptr< ConstElemPointerRange > _refined_elements
The elements that were just refined.
Definition: MooseMesh.h:908
virtual dof_id_type maxElemId() const
Definition: MooseMesh.C:2247
virtual bnd_elem_iterator bndElemsBegin()
Return iterators to the beginning/end of the boundary elements list.
Definition: MooseMesh.C:876
MooseEnum _partitioner_name
The partitioner used on this mesh.
Definition: MooseMesh.h:872
std::map< dof_id_type, std::set< SubdomainID > > _block_node_list
list of nodes that belongs to a specified block (domain)
Definition: MooseMesh.h:985
bool _node_to_active_semilocal_elem_map_built
Definition: MooseMesh.h:942
std::map< boundary_id_type, std::set< dof_id_type > > _bnd_node_ids
Map of sets of node IDs in each boundary.
Definition: MooseMesh.h:969
virtual void init()
Initialize the Mesh object.
Definition: MooseMesh.C:2032
ConstElemPointerRange * refinedElementRange() const
Return a range that is suitable for threaded execution over elements that were just refined...
Definition: MooseMesh.C:539
PerfID _get_boundary_node_range_timer
Definition: MooseMesh.h:1148
void detectPairedSidesets()
This routine detects paired sidesets of a regular orthogonal mesh (.i.e.
Definition: MooseMesh.C:1357
const std::set< SubdomainID > & getNodeBlockIds(const Node &node) const
Return list of blocks to which the given node belongs.
Definition: MooseMesh.C:847
MooseMesh(const InputParameters &parameters)
Typical "Moose-style" constructor and copy constructor.
Definition: MooseMesh.C:152
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
void setMeshBoundaryIDs(std::set< BoundaryID > boundary_IDs)
Sets the set of BoundaryIDs Is called by AddAllSideSetsByNormals.
Definition: MooseMesh.C:2340
std::map< boundary_id_type, std::vector< dof_id_type > > _node_set_nodes
list of nodes that belongs to a specified nodeset: indexing [nodeset_id] -> [array of node ids] ...
Definition: MooseMesh.h:988
void cacheInfo()
Definition: MooseMesh.C:820
void changeBoundaryId(const boundary_id_type old_id, const boundary_id_type new_id, bool delete_prev)
Change all the boundary IDs for a given side from old_id to new_id.
Definition: MooseMesh.C:1911
std::map< std::pair< int, ElemType >, std::vector< std::pair< unsigned int, QpMap > > > _elem_type_to_coarsening_map
Holds mappings for volume to volume and parent side to child side.
Definition: MooseMesh.h:1124
virtual void buildMesh()=0
Must be overridden by child classes.
bool isSplitMesh() const
Whether or not this is a split mesh operation.
Definition: MooseApp.C:871
unsigned int _to
The qp to map to.
Definition: MooseMesh.h:64
BoundaryID _bnd_id
boundary id for the node
Definition: BndNode.h:26
ConstNodeRange * getLocalNodeRange()
Definition: MooseMesh.C:774
Node * _node
pointer to the node
Definition: BndNode.h:24
virtual const Node & nodeRef(const dof_id_type i) const
Definition: MooseMesh.C:436
bool _allow_recovery
Whether or not this Mesh is allowed to read a recovery file.
Definition: MooseMesh.h:1130
bool operator()(const BndNode *const &lhs, const BndNode *const &rhs)
Definition: MooseMesh.C:614
void buildNodeListFromSideList()
Calls BoundaryInfo::build_node_list_from_side_list().
Definition: MooseMesh.C:2167
void setPatchUpdateStrategy(Moose::PatchUpdateType patch_update_strategy)
Set the patch size update strategy.
Definition: MooseMesh.C:2532
void reinit(const Elem *elem)
Reinitialize objects (JxW, q_points, ...) for an elements.
Definition: Assembly.C:1567
std::unique_ptr< NodeRange > _active_node_range
Definition: MooseMesh.h:930
const std::set< BoundaryID > & meshBoundaryIds() const
Returns a read-only reference to the set of boundary IDs currently present in the Mesh...
Definition: MooseMesh.C:2322
auto operator*(const T1 &a, const RankFourTensorTempl< T2 > &b) -> typename std::enable_if< ScalarTraits< T1 >::value, RankFourTensorTempl< decltype(T1() *T2())>>::type
virtual Elem * queryElemPtr(const dof_id_type i)
Definition: MooseMesh.C:2279
PerfID _detect_paired_sidesets_timer
Definition: MooseMesh.h:1154
void setIsCustomPartitionerRequested(bool cpr)
Definition: MooseMesh.C:2724
MeshBase::const_node_iterator localNodesEnd()
Definition: MooseMesh.C:2211
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
std::unique_ptr< StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > > _bnd_elem_range
Definition: MooseMesh.h:934
virtual bnd_node_iterator bndNodesBegin()
Return iterators to the beginning/end of the boundary nodes list.
Definition: MooseMesh.C:860
const SubdomainID INVALID_BLOCK_ID
Definition: MooseTypes.C:16
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
void mapPoints(const std::vector< Point > &from, const std::vector< Point > &to, std::vector< QpMap > &qp_map)
Find the closest points that map "from" to "to" and fill up "qp_map".
Definition: MooseMesh.C:1742
StoredRange< std::vector< const Elem * >::iterator, const Elem * > ConstElemPointerRange
Definition: MooseTypes.h:166
PerfID _cache_changed_lists_timer
Definition: MooseMesh.h:1139
QBase *const & writeableQRuleFace()
Returns the reference to the current quadrature being used on a current face.
Definition: Assembly.h:241
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2685
bool getDistributedMeshOnCommandLine() const
Returns true if the user specified –distributed-mesh (or –parallel-mesh, for backwards compatibilit...
Definition: MooseApp.h:323
std::map< ElemType, std::map< std::pair< int, int >, std::vector< std::vector< QpMap > > > > _elem_type_to_child_side_refinement_map
Holds mappings for "internal" child sides to parent volume. The second key is (child, child_side).
Definition: MooseMesh.h:1120
static const int GRAIN_SIZE
Definition: MooseMesh.C:58
bool _use_distributed_mesh
False by default.
Definition: MooseMesh.h:864
std::unique_ptr< std::map< BoundaryID, RealVectorValue > > _boundary_to_normal_map
The boundary to normal map - valid only when AddAllSideSetsByNormals is active.
Definition: MooseMesh.h:962
SemiLocalNodeRange * getActiveSemiLocalNodeRange() const
Definition: MooseMesh.C:765
void clearQuadratureNodes()
Clear out any existing quadrature nodes.
Definition: MooseMesh.C:995
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2567
void setSubdomainName(SubdomainID subdomain_id, SubdomainName name)
This method returns a writable reference to a subdomain name based on the id parameter.
Definition: MooseMesh.C:1118
PerfID _find_adaptivity_qp_maps_timer
Definition: MooseMesh.h:1157
std::vector< BndNode * > _bnd_nodes
array of boundary nodes
Definition: MooseMesh.h:965
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:42
PetscInt m
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimsension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh m...
Definition: MooseMesh.C:2133
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
PerfID _build_periodic_node_sets_timer
Definition: MooseMesh.h:1152
unsigned int _from
The qp to map from.
Definition: MooseMesh.h:61
std::pair< const Node *, BoundaryID > PeriodicNodeInfo
Helper type for building periodic node maps.
Definition: MooseMesh.h:838
QBase *const & writeableQRule()
Returns the reference to the current quadrature being used.
Definition: Assembly.h:174
PerfID _update_timer
Definition: MooseMesh.h:1137
const std::set< unsigned int > & getGhostedBoundaries() const
Return a writable reference to the set of ghosted boundary IDs.
Definition: MooseMesh.C:2385
bool _construct_node_list_from_side_list
Whether or not to allow generation of nodesets from sidesets.
Definition: MooseMesh.h:1133
boundary_id_type BoundaryID
PerfID _init_timer
Definition: MooseMesh.h:1160
void buildCoarseningMap(const Elem &elem, QBase &qrule, QBase &qrule_face, int input_side)
Build the coarsening map for a given element type.
Definition: MooseMesh.C:1703
unsigned int sideWithBoundaryID(const Elem *const elem, const BoundaryID boundary_id) const
Calls BoundaryInfo::side_with_boundary_id().
Definition: MooseMesh.C:2199
const std::vector< dof_id_type > & getNodeList(boundary_id_type nodeset_id) const
Return a writable reference to a vector of node IDs that belong to nodeset_id.
Definition: MooseMesh.C:2594
virtual const Node * nodePtr(const dof_id_type i) const
Definition: MooseMesh.C:454
bool isUseSplit() const
Whether or not we are running with pre-split (distributed mesh)
Definition: MooseApp.C:877
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
PerfID _read_recovered_mesh_timer
Definition: MooseMesh.h:1161
void update()
Calls buildNodeListFromSideList(), buildNodeList(), and buildBndElemList().
Definition: MooseMesh.C:403
const std::pair< BoundaryID, BoundaryID > * getPairedBoundaryMapping(unsigned int component)
This function attempts to return the paired boundary ids for the given component. ...
Definition: MooseMesh.C:1547
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
std::set< BoundaryID > _mesh_nodeset_ids
Definition: MooseMesh.h:958
PerfID _mesh_changed_timer
Definition: MooseMesh.h:1138
InputParameters validParams< MooseMesh >()
Definition: MooseMesh.C:63
PerfID _prepare_timer
Timers.
Definition: MooseMesh.h:1136
PerfID _cache_info_timer
Definition: MooseMesh.h:1150
std::unique_ptr< ConstElemPointerRange > _coarsened_elements
The elements that were just coarsened.
Definition: MooseMesh.h:911
subdomain_id_type SubdomainID
std::set< BoundaryID > _mesh_boundary_ids
A set of boundary IDs currently present in the mesh.
Definition: MooseMesh.h:956
const Node * addUniqueNode(const Point &p, Real tol=1e-6)
Add a new node to the mesh.
Definition: MooseMesh.C:891
std::vector< BndNode > _extra_bnd_nodes
Definition: MooseMesh.h:982
void printInfo(std::ostream &os=libMesh::out) const
Calls print_info() on the underlying Mesh.
Definition: MooseMesh.C:2588
virtual ExodusII_IO * exReader() const
Not implemented – always returns NULL.
Definition: MooseMesh.C:2581
void prepare(bool force=false)
Calls prepare_for_use() if force=true on the underlying Mesh object, then communicates various bounda...
Definition: MooseMesh.C:343
unsigned int uniformRefineLevel() const
Returns the level of uniform refinement requested (zero if AMR is disabled).
Definition: MooseMesh.C:2361
std::vector< BndElement * > _bnd_elems
array of boundary elems
Definition: MooseMesh.h:972
Real minPeriodicDistance(unsigned int nonlinear_var_num, Point p, Point q) const
This function returns the distance between two points on the mesh taking into account periodicity for...
Definition: MooseMesh.C:1541
bool hasSecondOrderElements()
check if the mesh has SECOND order elements
Definition: MooseMesh.C:2708
unsigned int getPatchSize() const
Getter for the patch_size parameter.
Definition: MooseMesh.C:2526
void buildRefinementAndCoarseningMaps(Assembly *assembly)
Create the refinement and coarsening maps necessary for projection of stateful material properties wh...
Definition: MooseMesh.C:1565
void findAdaptivityQpMaps(const Elem *template_elem, QBase &qrule, QBase &qrule_face, std::vector< std::vector< QpMap >> &refinement_map, std::vector< std::pair< unsigned int, QpMap >> &coarsen_map, int parent_side, int child, int child_side)
Given an elem type, get maps that tell us what qp&#39;s are closest to each other between a parent and it...
Definition: MooseMesh.C:1773
bool _parallel_type_overridden
Definition: MooseMesh.h:866
PerfID _build_periodic_node_map_timer
Definition: MooseMesh.h:1151
Interface for objects that needs transient capabilities.
std::vector< Node * > _node_map
Vector of all the Nodes in the mesh for determining when to add a new point.
Definition: MooseMesh.h:1006
std::map< dof_id_type, std::map< unsigned int, std::map< dof_id_type, Node * > > > _elem_to_side_to_qp_to_quadrature_nodes
Definition: MooseMesh.h:981
std::unique_ptr< StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > > _bnd_node_range
Definition: MooseMesh.h:932
PerfID _get_local_node_range_timer
Definition: MooseMesh.h:1147
Node * getQuadratureNode(const Elem *elem, const unsigned short int side, const unsigned int qp)
Get a specified quadrature node.
Definition: MooseMesh.C:977
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:25
virtual dof_id_type nNodes() const
Calls n_nodes/elem() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2229
const std::vector< const Elem * > & coarsenedElementChildren(const Elem *elem) const
Get the newly removed children element ids for an element that was just coarsened.
Definition: MooseMesh.C:551
void setBoundaryName(BoundaryID boundary_id, BoundaryName name)
This method returns a writable reference to a boundary name based on the id parameter.
Definition: MooseMesh.C:1130
PerfID _change_boundary_id_timer
Definition: MooseMesh.h:1159
std::unique_ptr< MeshBase > buildMeshBaseObject(ParallelType override_type=ParallelType::DEFAULT)
Method to construct a libMesh::MeshBase object that is normally set and used by the MooseMesh object ...
Definition: MooseMesh.C:1969
Moose::PatchUpdateType _patch_update_strategy
The patch update strategy.
Definition: MooseMesh.h:1003
void addPeriodicVariable(unsigned int var_num, BoundaryID primary, BoundaryID secondary)
For "regular orthogonal" meshes, determine if variable var_num is periodic with respect to the primar...
Definition: MooseMesh.C:1484
StoredRange< std::set< Node * >::iterator, Node * > SemiLocalNodeRange
Definition: MooseMesh.h:47
void setCurrentSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:293
PerfID _build_refinement_and_coarsening_maps_timer
Definition: MooseMesh.h:1158
std::set< BoundaryID > _mesh_sideset_ids
Definition: MooseMesh.h:957
void setBoundaryToNormalMap(std::unique_ptr< std::map< BoundaryID, RealVectorValue >> boundary_map)
Sets the mapping between BoundaryID and normal vector Is called by AddAllSideSetsByNormals.
Definition: MooseMesh.C:2346
bool _partitioner_overridden
Definition: MooseMesh.h:873
bool detectOrthogonalDimRanges(Real tol=1e-6)
This routine determines whether the Mesh is a regular orthogonal mesh (i.e.
Definition: MooseMesh.C:1283
void updateActiveSemiLocalNodeRange(std::set< dof_id_type > &ghosted_elems)
Clears the "semi-local" node list and rebuilds it.
Definition: MooseMesh.C:559
std::map< const Elem *, std::vector< const Elem * > > _coarsened_element_children
Map of Parent elements to children elements for elements that were just coarsened.
std::unique_ptr< libMesh::MeshBase > _mesh
Pointer to underlying libMesh mesh object.
Definition: MooseMesh.h:869
NodeRange * getActiveNodeRange()
Definition: MooseMesh.C:751
void buildPeriodicNodeSets(std::map< BoundaryID, std::set< dof_id_type >> &periodic_node_sets, unsigned int var_number, PeriodicBoundaries *pbs) const
This routine builds a datastructure of node ids organized by periodic boundary ids.
Definition: MooseMesh.C:1256
void buildRefinementMap(const Elem &elem, QBase &qrule, QBase &qrule_face, int parent_side, int child, int child_side)
Build the refinement map for a given element type.
Definition: MooseMesh.C:1624
PetscInt n
void setGhostedBoundaryInflation(const std::vector< Real > &inflation)
This sets the inflation amount for the bounding box for each partition for use in ghosting boundaries...
Definition: MooseMesh.C:2379
void freeBndNodes()
Definition: MooseMesh.C:312
Real dimensionWidth(unsigned int component) const
Returns the width of the requested dimension.
Definition: MooseMesh.C:1460
PerfID _get_boundary_element_range_timer
Definition: MooseMesh.h:1149
PatchUpdateType
Type of patch update strategy for modeling node-face constraints or contact.
Definition: MooseTypes.h:705
void mooseDeprecated(Args &&... args) const
Definition: MooseObject.h:161
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:177
PerfID _build_bnd_elem_list_timer
Definition: MooseMesh.h:1142
std::map< SubdomainID, std::set< BoundaryID > > _subdomain_boundary_ids
Holds a map from subomdain ids to the boundary ids that are attached to it.
Definition: MooseMesh.h:1127
std::vector< const Elem * > _refined_elements
The elements that were just refined.
std::unique_ptr< ConstElemRange > _active_local_elem_range
A range for use with threading.
Definition: MooseMesh.h:927
MPI_Comm comm
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
Definition: MooseMesh.C:422
BoundingBox getInflatedProcessorBoundingBox(Real inflation_multiplier=0.01) const
Get a (slightly inflated) processor bounding box.
Definition: MooseMesh.C:2544
std::map< unsigned int, std::vector< bool > > _periodic_dim
A map of vectors indicating which dimensions are periodic in a regular orthogonal mesh for the specif...
Definition: MooseMesh.h:1026
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:59
void setMeshBase(std::unique_ptr< MeshBase > mesh_base)
Method to set the mesh_base object.
Definition: MooseMesh.C:2026
void buildBndElemList()
Definition: MooseMesh.C:667
std::vector< Real > _ghosted_boundaries_inflation
Definition: MooseMesh.h:991
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
unsigned int _patch_size
The number of nodes to consider in the NearestNode neighborhood.
Definition: MooseMesh.h:994
std::string getRecoverFileSuffix()
The suffix for the recovery file.
Definition: MooseApp.h:379
PerfID _ghost_ghosted_boundaries_timer
Definition: MooseMesh.h:1162
std::vector< SubdomainID > getSubdomainIDs(const std::vector< SubdomainName > &subdomain_name) const
Get the associated subdomainIDs for the subdomain names that are passed in.
Definition: MooseMesh.C:1090
void addGhostedBoundary(BoundaryID boundary_id)
This will add the boundary ids to be ghosted to this processor.
Definition: MooseMesh.C:2373
PerfID _update_active_semi_local_node_range_timer
Definition: MooseMesh.h:1140
virtual dof_id_type maxNodeId() const
Calls max_node/elem_id() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2241
virtual Elem * elem(const dof_id_type i)
Various accessors (pointers/references) for Elem "i".
Definition: MooseMesh.C:2253
StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > * getBoundaryElementRange()
Definition: MooseMesh.C:801
Definition: Moose.h:112
const RealVectorValue & getNormalByBoundaryID(BoundaryID id) const
Returns the normal vector associated with a given BoundaryID.
Definition: MooseMesh.C:1953
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > ConstBndNodeRange
Some useful StoredRange typedefs.
Definition: MooseMesh.h:1258
std::set< SubdomainID > _mesh_subdomains
A set of subdomain IDs currently present in the mesh.
Definition: MooseMesh.h:948
ConstElemPointerRange * coarsenedElementRange() const
Return a range that is suitable for threaded execution over elements that were just coarsened...
Definition: MooseMesh.C:545
const Moose::PatchUpdateType & getPatchUpdateStrategy() const
Get the current patch update strategy.
Definition: MooseMesh.C:2538
std::vector< std::pair< BoundaryID, BoundaryID > > _paired_boundary
A vector holding the paired boundaries for a regular orthogonal mesh.
Definition: MooseMesh.h:1015
void ghostGhostedBoundaries()
Actually do the ghosting of boundaries that need to be ghosted to this processor. ...
Definition: MooseMesh.C:2455
std::set< unsigned int > _ghosted_boundaries
Definition: MooseMesh.h:990
PerfID _node_to_active_semilocal_elem_map_timer
Definition: MooseMesh.h:1144
void buildNodeList()
Calls BoundaryInfo::build_node_list()/build_side_list() and makes separate copies of Nodes/Elems in t...
Definition: MooseMesh.C:633
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
bool isBoundaryElem(dof_id_type elem_id) const
Returns true if the requested element is in the list of boundary elements, false otherwise.
Definition: MooseMesh.C:2659
std::map< const Elem *, std::vector< const Elem * > > _coarsened_element_children
Map of Parent elements to child elements for elements that were just coarsened.
Definition: MooseMesh.h:918
bool isUltimateMaster()
Whether or not this app is the ultimate master app.
Definition: MooseApp.h:528
bool isRecovering() const
Whether or not this is a "recover" calculation.
Definition: MooseApp.C:859
const MeshBase::const_element_iterator activeLocalElementsEnd()
Definition: MooseMesh.C:2223
std::vector< const Elem * > _coarsened_elements
The elements that were just coarsened.
virtual std::unique_ptr< PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libMesh mesh (default)...
Definition: MooseMesh.C:2730
bool _needs_prepare_for_use
True if prepare_for_use should be called when Mesh is prepared.
Definition: MooseMesh.h:905
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
Definition: MooseMesh.C:788
std::map< std::pair< int, ElemType >, std::vector< std::vector< QpMap > > > _elem_type_to_refinement_map
Holds mappings for volume to volume and parent side to child side.
Definition: MooseMesh.h:1116
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
Definition: MooseMesh.h:1009
void buildPeriodicNodeMap(std::multimap< dof_id_type, dof_id_type > &periodic_node_map, unsigned int var_number, PeriodicBoundaries *pbs) const
This routine builds a multimap of boundary ids to matching boundary ids across all periodic boundarie...
Definition: MooseMesh.C:1171
std::map< dof_id_type, std::vector< dof_id_type > > _node_to_active_semilocal_elem_map
A map of all of the current nodes to the active elements that they are connected to.
Definition: MooseMesh.h:941
Real _distance
The distance between them.
Definition: MooseMesh.h:67
void setCustomPartitioner(Partitioner *partitioner)
Setter for custom partitioner.
Definition: MooseMesh.C:2696
PerfID _detect_orthogonal_dim_ranges_timer
Definition: MooseMesh.h:1153
std::unique_ptr< Partitioner > _custom_partitioner
The custom partitioner.
Definition: MooseMesh.h:876
void meshChanged()
Declares that the MooseMesh has changed, invalidates cached data and rebuilds caches.
Definition: MooseMesh.C:490
PerfID _get_active_local_element_range_timer
Definition: MooseMesh.h:1145
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
Definition: MooseMesh.C:1007
virtual dof_id_type nElem() const
Definition: MooseMesh.C:2235
const std::string & getSubdomainName(SubdomainID subdomain_id)
Return the name of a block given an id.
Definition: MooseMesh.C:1124
const std::map< dof_id_type, std::vector< dof_id_type > > & nodeToElemMap()
If not already created, creates a map from every node to all elements to which they are connected...
Definition: MooseMesh.C:690
std::unique_ptr< SemiLocalNodeRange > _active_semilocal_node_range
Definition: MooseMesh.h:929
PerfID _get_active_node_range_timer
Definition: MooseMesh.h:1146
void setUniformRefineLevel(unsigned int)
Set uniform refinement level.
Definition: MooseMesh.C:2367
const std::set< SubdomainID > & meshSubdomains() const
Returns a read-only reference to the set of subdomains currently present in the Mesh.
Definition: MooseMesh.C:2316
SubdomainID getSubdomainID(const SubdomainName &subdomain_name) const
Get the associated subdomain ID for the subdomain name.
Definition: MooseMesh.C:1075
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
std::unique_ptr< ConstNodeRange > _local_node_range
Definition: MooseMesh.h:931
PerfID _node_to_elem_map_timer
Definition: MooseMesh.h:1143
virtual unsigned int effectiveSpatialDimension() const
Returns the effective spatial dimension determined by the coordinates actually used by the mesh...
Definition: MooseMesh.C:2139