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