https://mooseframework.inl.gov
MooseMesh.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "MooseError.h"
11 #include "MooseMesh.h"
12 #include "Factory.h"
14 #include "MooseUtils.h"
15 #include "MooseApp.h"
16 #include "RelationshipManager.h"
17 #include "PointListAdaptor.h"
18 #include "Executioner.h"
19 #include "NonlinearSystemBase.h"
20 #include "LinearSystem.h"
21 #include "AuxiliarySystem.h"
22 #include "Assembly.h"
23 #include "SubProblem.h"
24 #include "MooseVariableBase.h"
25 #include "MooseMeshUtils.h"
26 #include "MooseAppCoordTransform.h"
27 #include "FEProblemBase.h"
28 
29 #include <utility>
30 
31 // libMesh
32 #include "libmesh/bounding_box.h"
33 #include "libmesh/boundary_info.h"
34 #include "libmesh/mesh_tools.h"
35 #include "libmesh/parallel.h"
36 #include "libmesh/mesh_communication.h"
37 #include "libmesh/periodic_boundary_base.h"
38 #include "libmesh/fe_base.h"
39 #include "libmesh/fe_interface.h"
40 #include "libmesh/mesh_inserter_iterator.h"
41 #include "libmesh/mesh_communication.h"
42 #include "libmesh/mesh_inserter_iterator.h"
43 #include "libmesh/mesh_tools.h"
44 #include "libmesh/parallel.h"
45 #include "libmesh/parallel_elem.h"
46 #include "libmesh/parallel_node.h"
47 #include "libmesh/parallel_ghost_sync.h"
48 #include "libmesh/utility.h"
49 #include "libmesh/remote_elem.h"
50 #include "libmesh/linear_partitioner.h"
51 #include "libmesh/centroid_partitioner.h"
52 #include "libmesh/parmetis_partitioner.h"
53 #include "libmesh/hilbert_sfc_partitioner.h"
54 #include "libmesh/morton_sfc_partitioner.h"
55 #include "libmesh/edge_edge2.h"
56 #include "libmesh/mesh_refinement.h"
57 #include "libmesh/quadrature.h"
58 #include "libmesh/boundary_info.h"
59 #include "libmesh/periodic_boundaries.h"
60 #include "libmesh/quadrature_gauss.h"
61 #include "libmesh/point_locator_base.h"
62 #include "libmesh/default_coupling.h"
63 #include "libmesh/ghost_point_neighbors.h"
64 #include "libmesh/fe_type.h"
65 #include "libmesh/enum_to_string.h"
66 
67 static const int GRAIN_SIZE =
68  1; // the grain_size does not have much influence on our execution speed
69 
70 using namespace libMesh;
71 
72 // Make newer nanoflann API compatible with older nanoflann versions
73 #if NANOFLANN_VERSION < 0x150
74 namespace nanoflann
75 {
76 typedef SearchParams SearchParameters;
77 
78 template <typename T, typename U>
79 using ResultItem = std::pair<T, U>;
80 }
81 #endif
82 
85 {
87 
88  MooseEnum parallel_type("DEFAULT REPLICATED DISTRIBUTED", "DEFAULT");
89  params.addParam<MooseEnum>("parallel_type",
90  parallel_type,
91  "DEFAULT: Use libMesh::ReplicatedMesh unless --distributed-mesh is "
92  "specified on the command line "
93  "REPLICATED: Always use libMesh::ReplicatedMesh "
94  "DISTRIBUTED: Always use libMesh::DistributedMesh");
95 
96  params.addParam<bool>(
97  "allow_renumbering",
98  true,
99  "If allow_renumbering=false, node and element numbers are kept fixed until deletion");
100 
101  // TODO: this parameter does not belong here, it's only for FileMesh
102  params.addParam<bool>("nemesis",
103  false,
104  "If nemesis=true and file=foo.e, actually reads "
105  "foo.e.N.0, foo.e.N.1, ... foo.e.N.N-1, "
106  "where N = # CPUs, with NemesisIO.");
107 
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 iteration", "never");
119  params.addParam<MooseEnum>(
120  "patch_update_strategy",
121  patch_update_strategy,
122  "How often to update the geometric search 'patch'. The default is to "
123  "never update it (which is the most efficient but could be a problem "
124  "with lots of relative motion). 'always' will update the patch for all "
125  "secondary nodes at the beginning of every timestep which might be time "
126  "consuming. 'auto' will attempt to determine at the start of which "
127  "timesteps the patch for all secondary nodes needs to be updated automatically."
128  "'iteration' updates the patch at every nonlinear iteration for a "
129  "subset of secondary nodes for which penetration is not detected. If there "
130  "can be substantial relative motion between the primary and secondary surfaces "
131  "during the nonlinear iterations within a timestep, it is advisable to use "
132  "'iteration' option to ensure accurate contact detection.");
133 
134  // Note: This parameter is named to match 'construct_side_list_from_node_list' in SetupMeshAction
135  params.addParam<bool>(
136  "construct_node_list_from_side_list",
137  true,
138  "Whether or not to generate nodesets from the sidesets (usually a good idea).");
139  params.addParam<unsigned int>(
140  "patch_size", 40, "The number of nodes to consider in the NearestNode neighborhood.");
141  params.addParam<unsigned int>("ghosting_patch_size",
142  "The number of nearest neighbors considered "
143  "for ghosting purposes when 'iteration' "
144  "patch update strategy is used. Default is "
145  "5 * patch_size.");
146  params.addParam<unsigned int>("max_leaf_size",
147  10,
148  "The maximum number of points in each leaf of the KDTree used in "
149  "the nearest neighbor search. As the leaf size becomes larger,"
150  "KDTree construction becomes faster but the nearest neighbor search"
151  "becomes slower.");
152 
153  params.addParam<bool>("build_all_side_lowerd_mesh",
154  false,
155  "True to build the lower-dimensional mesh for all sides.");
156 
157  params.addParam<bool>("skip_refine_when_use_split",
158  true,
159  "True to skip uniform refinements when using a pre-split mesh.");
160 
161  params.addParam<std::vector<SubdomainID>>(
162  "add_subdomain_ids",
163  "The listed subdomain ids will be assumed valid for the mesh. This permits setting up "
164  "subdomain restrictions for subdomains initially containing no elements, which can occur, "
165  "for example, in additive manufacturing simulations which dynamically add and remove "
166  "elements. Names for this subdomains may be provided using add_subdomain_names. In this case "
167  "this list and add_subdomain_names must contain the same number of items.");
168  params.addParam<std::vector<SubdomainName>>(
169  "add_subdomain_names",
170  "The listed subdomain names will be assumed valid for the mesh. This permits setting up "
171  "subdomain restrictions for subdomains initially containing no elements, which can occur, "
172  "for example, in additive manufacturing simulations which dynamically add and remove "
173  "elements. IDs for this subdomains may be provided using add_subdomain_ids. Otherwise IDs "
174  "are automatically assigned. In case add_subdomain_ids is set too, both lists must contain "
175  "the same number of items.");
176 
177  params.addParam<std::vector<BoundaryID>>(
178  "add_sideset_ids",
179  "The listed sideset ids will be assumed valid for the mesh. This permits setting up boundary "
180  "restrictions for sidesets initially containing no sides. Names for this sidesets may be "
181  "provided using add_sideset_names. In this case this list and add_sideset_names must contain "
182  "the same number of items.");
183  params.addParam<std::vector<BoundaryName>>(
184  "add_sideset_names",
185  "The listed sideset names will be assumed valid for the mesh. This permits setting up "
186  "boundary restrictions for sidesets initially containing no sides. Ids for this sidesets may "
187  "be provided using add_sideset_ids. In this case this list and add_sideset_ids must contain "
188  "the same number of items.");
189 
190  params.addParam<std::vector<BoundaryID>>(
191  "add_nodeset_ids",
192  "The listed nodeset ids will be assumed valid for the mesh. This permits setting up boundary "
193  "restrictions for node initially containing no sides. Names for this nodesets may be "
194  "provided using add_nodeset_names. In this case this list and add_nodeset_names must contain "
195  "the same number of items.");
196  params.addParam<std::vector<BoundaryName>>(
197  "add_nodeset_names",
198  "The listed nodeset names will be assumed valid for the mesh. This permits setting up "
199  "boundary restrictions for nodesets initially containing no sides. Ids for this nodesets may "
200  "be provided using add_nodesets_ids. In this case this list and add_nodesets_ids must "
201  "contain the same number of items.");
202 
204 
205  // This indicates that the derived mesh type accepts a MeshGenerator, and should be set to true in
206  // derived types that do so.
207  params.addPrivateParam<bool>("_mesh_generator_mesh", false);
208 
209  // Whether or not the mesh is pre split
210  params.addPrivateParam<bool>("_is_split", false);
211 
212  params.registerBase("MooseMesh");
213 
214  // groups
215  params.addParamNamesToGroup("patch_update_strategy patch_size max_leaf_size", "Geometric search");
216  params.addParamNamesToGroup("nemesis", "Advanced");
217  params.addParamNamesToGroup("add_subdomain_ids add_subdomain_names add_sideset_ids "
218  "add_sideset_names add_nodeset_ids add_nodeset_names",
219  "Pre-declaration of future mesh sub-entities");
220  params.addParamNamesToGroup("construct_node_list_from_side_list build_all_side_lowerd_mesh",
221  "Automatic definition of mesh element sides entities");
222  params.addParamNamesToGroup("partitioner centroid_partitioner_direction", "Partitioning");
223 
224  return params;
225 }
226 
228  : MooseObject(parameters),
229  Restartable(this, "Mesh"),
230  PerfGraphInterface(this),
231  _parallel_type(getParam<MooseEnum>("parallel_type").getEnum<MooseMesh::ParallelType>()),
232  _use_distributed_mesh(false),
233  _distribution_overridden(false),
234  _parallel_type_overridden(false),
235  _mesh(nullptr),
236  _partitioner_name(getParam<MooseEnum>("partitioner")),
237  _partitioner_overridden(false),
238  _custom_partitioner_requested(false),
239  _uniform_refine_level(0),
240  _skip_refine_when_use_split(getParam<bool>("skip_refine_when_use_split")),
241  _skip_deletion_repartition_after_refine(false),
242  _is_nemesis(getParam<bool>("nemesis")),
243  _node_to_elem_map_built(false),
244  _node_to_active_semilocal_elem_map_built(false),
245  _patch_size(getParam<unsigned int>("patch_size")),
246  _ghosting_patch_size(isParamValid("ghosting_patch_size")
247  ? getParam<unsigned int>("ghosting_patch_size")
248  : 5 * _patch_size),
249  _max_leaf_size(getParam<unsigned int>("max_leaf_size")),
250  _patch_update_strategy(
251  getParam<MooseEnum>("patch_update_strategy").getEnum<Moose::PatchUpdateType>()),
252  _regular_orthogonal_mesh(false),
253  _is_split(getParam<bool>("_is_split")),
254  _has_lower_d(false),
255  _allow_recovery(true),
256  _construct_node_list_from_side_list(getParam<bool>("construct_node_list_from_side_list")),
257  _need_delete(false),
258  _allow_remote_element_removal(true),
259  _need_ghost_ghosted_boundaries(true),
260  _is_displaced(false),
261  _coord_sys(
262  declareRestartableData<std::map<SubdomainID, Moose::CoordinateSystemType>>("coord_sys")),
263  _rz_coord_axis(getParam<MooseEnum>("rz_coord_axis")),
264  _coord_system_set(false),
265  _doing_p_refinement(false)
266 {
267  if (isParamValid("ghosting_patch_size") && (_patch_update_strategy != Moose::Iteration))
268  mooseError("Ghosting patch size parameter has to be set in the mesh block "
269  "only when 'iteration' patch update strategy is used.");
270 
271  if (isParamValid("coord_block"))
272  {
273  if (isParamValid("block"))
274  paramWarning("block",
275  "You set both 'Mesh/block' and 'Mesh/coord_block'. The value of "
276  "'Mesh/coord_block' will be used.");
277 
278  _provided_coord_blocks = getParam<std::vector<SubdomainName>>("coord_block");
279  }
280  else if (isParamValid("block"))
281  _provided_coord_blocks = getParam<std::vector<SubdomainName>>("block");
282 
283  if (getParam<bool>("build_all_side_lowerd_mesh"))
284  // Do not initially allow removal of remote elements
286 
288 }
289 
290 MooseMesh::MooseMesh(const MooseMesh & other_mesh)
291  : MooseObject(other_mesh._pars),
292  Restartable(this, "Mesh"),
293  PerfGraphInterface(this, "CopiedMesh"),
294  _built_from_other_mesh(true),
295  _parallel_type(other_mesh._parallel_type),
296  _use_distributed_mesh(other_mesh._use_distributed_mesh),
297  _distribution_overridden(other_mesh._distribution_overridden),
298  _parallel_type_overridden(other_mesh._parallel_type_overridden),
299  _mesh(other_mesh.getMesh().clone()),
300  _partitioner_name(other_mesh._partitioner_name),
301  _partitioner_overridden(other_mesh._partitioner_overridden),
302  _custom_partitioner_requested(other_mesh._custom_partitioner_requested),
303  _uniform_refine_level(other_mesh.uniformRefineLevel()),
304  _skip_refine_when_use_split(other_mesh._skip_refine_when_use_split),
305  _skip_deletion_repartition_after_refine(other_mesh._skip_deletion_repartition_after_refine),
306  _is_nemesis(false),
307  _node_to_elem_map_built(false),
308  _node_to_active_semilocal_elem_map_built(false),
309  _patch_size(other_mesh._patch_size),
310  _ghosting_patch_size(other_mesh._ghosting_patch_size),
311  _max_leaf_size(other_mesh._max_leaf_size),
312  _patch_update_strategy(other_mesh._patch_update_strategy),
313  _regular_orthogonal_mesh(false),
314  _is_split(other_mesh._is_split),
315  _lower_d_interior_blocks(other_mesh._lower_d_interior_blocks),
316  _lower_d_boundary_blocks(other_mesh._lower_d_boundary_blocks),
317  _has_lower_d(other_mesh._has_lower_d),
318  _allow_recovery(other_mesh._allow_recovery),
319  _construct_node_list_from_side_list(other_mesh._construct_node_list_from_side_list),
320  _need_delete(other_mesh._need_delete),
321  _allow_remote_element_removal(other_mesh._allow_remote_element_removal),
322  _need_ghost_ghosted_boundaries(other_mesh._need_ghost_ghosted_boundaries),
323  _coord_sys(other_mesh._coord_sys),
324  _rz_coord_axis(other_mesh._rz_coord_axis),
325  _subdomain_id_to_rz_coord_axis(other_mesh._subdomain_id_to_rz_coord_axis),
326  _coord_system_set(other_mesh._coord_system_set),
327  _provided_coord_blocks(other_mesh._provided_coord_blocks),
328  _doing_p_refinement(other_mesh._doing_p_refinement)
329 {
330  // Note: this calls BoundaryInfo::operator= without changing the
331  // ownership semantics of either Mesh's BoundaryInfo object.
332  getMesh().get_boundary_info() = other_mesh.getMesh().get_boundary_info();
333 
334  const std::set<SubdomainID> & subdomains = other_mesh.meshSubdomains();
335  for (const auto & sbd_id : subdomains)
336  setSubdomainName(sbd_id, other_mesh.getMesh().subdomain_name(sbd_id));
337 
338  // Get references to BoundaryInfo objects to make the code below cleaner...
339  const BoundaryInfo & other_boundary_info = other_mesh.getMesh().get_boundary_info();
340  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
341 
342  // Use the other BoundaryInfo object to build the list of side boundary ids
343  std::vector<BoundaryID> side_boundaries;
344  other_boundary_info.build_side_boundary_ids(side_boundaries);
345 
346  // Assign those boundary ids in our BoundaryInfo object
347  for (const auto & side_bnd_id : side_boundaries)
348  boundary_info.sideset_name(side_bnd_id) = other_boundary_info.get_sideset_name(side_bnd_id);
349 
350  // Do the same thing for node boundary ids
351  std::vector<BoundaryID> node_boundaries;
352  other_boundary_info.build_node_boundary_ids(node_boundaries);
353 
354  for (const auto & node_bnd_id : node_boundaries)
355  boundary_info.nodeset_name(node_bnd_id) = other_boundary_info.get_nodeset_name(node_bnd_id);
356 
357  _bounds.resize(other_mesh._bounds.size());
358  for (std::size_t i = 0; i < _bounds.size(); ++i)
359  {
360  _bounds[i].resize(other_mesh._bounds[i].size());
361  for (std::size_t j = 0; j < _bounds[i].size(); ++j)
362  _bounds[i][j] = other_mesh._bounds[i][j];
363  }
364 
366 }
367 
369 {
370  freeBndNodes();
371  freeBndElems();
373 }
374 
375 void
377 {
378  // free memory
379  for (auto & bnode : _bnd_nodes)
380  delete bnode;
381 
382  for (auto & it : _node_set_nodes)
383  it.second.clear();
384 
385  _node_set_nodes.clear();
386 
387  for (auto & it : _bnd_node_ids)
388  it.second.clear();
389 
390  _bnd_node_ids.clear();
391  _bnd_node_range.reset();
392 }
393 
394 void
396 {
397  // free memory
398  for (auto & belem : _bnd_elems)
399  delete belem;
400 
401  for (auto & it : _bnd_elem_ids)
402  it.second.clear();
403 
404  _bnd_elem_ids.clear();
405  _bnd_elem_range.reset();
406 }
407 
408 bool
409 MooseMesh::prepare(const MeshBase * const mesh_to_clone)
410 {
411  TIME_SECTION("prepare", 2, "Preparing Mesh", true);
412 
413  bool called_prepare_for_use = false;
414 
415  mooseAssert(_mesh, "The MeshBase has not been constructed");
416 
417  if (!dynamic_cast<DistributedMesh *>(&getMesh()) || _is_nemesis)
418  // For whatever reason we do not want to allow renumbering here nor ever in the future?
419  getMesh().allow_renumbering(false);
420 
421  if (mesh_to_clone)
422  {
423  mooseAssert(mesh_to_clone->is_prepared(),
424  "The mesh we wish to clone from must already be prepared");
425  _mesh = mesh_to_clone->clone();
426  _moose_mesh_prepared = false;
427  }
428  else if (!_mesh->is_prepared())
429  {
430  _mesh->prepare_for_use();
431  _moose_mesh_prepared = false;
432  called_prepare_for_use = true;
433  }
434 
436  return called_prepare_for_use;
437 
438  // Collect (local) subdomain IDs
439  _mesh_subdomains.clear();
440  for (const auto & elem : getMesh().element_ptr_range())
442 
443  // add explicitly requested subdomains
444  if (isParamValid("add_subdomain_ids") && !isParamValid("add_subdomain_names"))
445  {
446  // only subdomain ids are explicitly given
447  const auto & add_subdomain_id = getParam<std::vector<SubdomainID>>("add_subdomain_ids");
448  _mesh_subdomains.insert(add_subdomain_id.begin(), add_subdomain_id.end());
449  }
450  else if (isParamValid("add_subdomain_ids") && isParamValid("add_subdomain_names"))
451  {
452  const auto add_subdomain =
453  getParam<SubdomainID, SubdomainName>("add_subdomain_ids", "add_subdomain_names");
454  for (const auto & [sub_id, sub_name] : add_subdomain)
455  {
456  // add subdomain id
457  _mesh_subdomains.insert(sub_id);
458  // set name of the subdomain just added
459  setSubdomainName(sub_id, sub_name);
460  }
461  }
462  else if (isParamValid("add_subdomain_names"))
463  {
464  // the user has defined add_subdomain_names, but not add_subdomain_ids
465  const auto & add_subdomain_names = getParam<std::vector<SubdomainName>>("add_subdomain_names");
466 
467  // to define subdomain ids, we need the largest subdomain id defined yet.
468  subdomain_id_type offset = 0;
469  if (!_mesh_subdomains.empty())
470  offset = *_mesh_subdomains.rbegin();
471 
472  // add all subdomains (and auto-assign ids)
473  for (const SubdomainName & sub_name : add_subdomain_names)
474  {
475  // to avoid two subdomains with the same ID (notably on recover)
477  continue;
478  const auto sub_id = ++offset;
479  // add subdomain id
480  _mesh_subdomains.insert(sub_id);
481  // set name of the subdomain just added
482  setSubdomainName(sub_id, sub_name);
483  }
484  }
485 
486  // Make sure nodesets have been generated
488 
489  // Collect (local) boundary IDs
490  const std::set<BoundaryID> & local_bids = getMesh().get_boundary_info().get_boundary_ids();
491  _mesh_boundary_ids.insert(local_bids.begin(), local_bids.end());
492 
493  const std::set<BoundaryID> & local_node_bids =
495  _mesh_nodeset_ids.insert(local_node_bids.begin(), local_node_bids.end());
496 
497  const std::set<BoundaryID> & local_side_bids =
499  _mesh_sideset_ids.insert(local_side_bids.begin(), local_side_bids.end());
500 
501  // Add explicitly requested sidesets/nodesets
502  // This is done *after* the side boundaries (e.g. "right", ...) have been generated.
503  auto add_sets = [this](const bool sidesets, auto & set_ids)
504  {
505  const std::string type = sidesets ? "sideset" : "nodeset";
506  const std::string id_param = "add_" + type + "_ids";
507  const std::string name_param = "add_" + type + "_names";
508 
509  if (isParamValid(id_param))
510  {
511  const auto & add_ids = getParam<std::vector<BoundaryID>>(id_param);
512  _mesh_boundary_ids.insert(add_ids.begin(), add_ids.end());
513  set_ids.insert(add_ids.begin(), add_ids.end());
515  {
516  const auto & add_names = getParam<std::vector<BoundaryName>>(name_param);
517  mooseAssert(add_names.size() == add_ids.size(),
518  "Id and name sets must be the same size when adding.");
519  for (const auto i : index_range(add_ids))
520  setBoundaryName(add_ids[i], add_names[i]);
521  }
522  }
523  else if (isParamValid(name_param))
524  {
525  // the user has defined names, but not ids
526  const auto & add_names = getParam<std::vector<BoundaryName>>(name_param);
527 
528  auto & mesh_ids = sidesets ? _mesh_sideset_ids : _mesh_nodeset_ids;
529 
530  // to define ids, we need the largest id defined yet.
531  boundary_id_type offset = 0;
532  if (!mesh_ids.empty())
533  offset = *mesh_ids.rbegin();
534  if (!_mesh_boundary_ids.empty())
535  offset = std::max(offset, *_mesh_boundary_ids.rbegin());
536 
537  // add all sidesets/nodesets (and auto-assign ids)
538  for (const auto & name : add_names)
539  {
540  // to avoid two sets with the same ID (notably on recover)
542  continue;
543  const auto id = ++offset;
544  // add sideset id
545  _mesh_boundary_ids.insert(id);
546  set_ids.insert(id);
547  // set name of the sideset just added
548  setBoundaryName(id, name);
549  }
550  }
551  };
552 
553  add_sets(true, _mesh_sideset_ids);
554  add_sets(false, _mesh_nodeset_ids);
555 
556  // Communicate subdomain and boundary IDs if this is a parallel mesh
557  if (!getMesh().is_serial())
558  {
563  }
564 
566  {
567  if (!_coord_system_set)
568  setCoordSystem(_provided_coord_blocks, getParam<MultiMooseEnum>("coord_type"));
569  else if (_pars.isParamSetByUser("coord_type"))
570  mooseError(
571  "Trying to set coordinate system type information based on the user input file, but "
572  "the coordinate system type information has already been set programmatically! "
573  "Either remove your coordinate system type information from the input file, or contact "
574  "your application developer");
575  }
576 
577  // Set general axisymmetric axes if provided
578  if (isParamValid("rz_coord_blocks") && isParamValid("rz_coord_origins") &&
579  isParamValid("rz_coord_directions"))
580  {
581  const auto rz_coord_blocks = getParam<std::vector<SubdomainName>>("rz_coord_blocks");
582  const auto rz_coord_origins = getParam<std::vector<Point>>("rz_coord_origins");
583  const auto rz_coord_directions = getParam<std::vector<RealVectorValue>>("rz_coord_directions");
584  if (rz_coord_origins.size() == rz_coord_blocks.size() &&
585  rz_coord_directions.size() == rz_coord_blocks.size())
586  {
587  std::vector<std::pair<Point, RealVectorValue>> rz_coord_axes;
588  for (unsigned int i = 0; i < rz_coord_origins.size(); ++i)
589  rz_coord_axes.push_back(std::make_pair(rz_coord_origins[i], rz_coord_directions[i]));
590 
591  setGeneralAxisymmetricCoordAxes(rz_coord_blocks, rz_coord_axes);
592 
593  if (isParamSetByUser("rz_coord_axis"))
594  mooseError("The parameter 'rz_coord_axis' may not be provided if 'rz_coord_blocks', "
595  "'rz_coord_origins', and 'rz_coord_directions' are provided.");
596  }
597  else
598  mooseError("The parameters 'rz_coord_blocks', 'rz_coord_origins', and "
599  "'rz_coord_directions' must all have the same size.");
600  }
601  else if (isParamValid("rz_coord_blocks") || isParamValid("rz_coord_origins") ||
602  isParamValid("rz_coord_directions"))
603  mooseError("If any of the parameters 'rz_coord_blocks', 'rz_coord_origins', and "
604  "'rz_coord_directions' are provided, then all must be provided.");
605 
607 
608  update();
609 
610  // Check if there is subdomain name duplication for the same subdomain ID
612 
613  _moose_mesh_prepared = true;
614 
615  return called_prepare_for_use;
616 }
617 
618 void
620 {
621  TIME_SECTION("update", 3, "Updating Mesh", true);
622 
623  // Rebuild the boundary conditions
625 
626  // Update the node to elem map
627  _node_to_elem_map.clear();
628  _node_to_elem_map_built = false;
631 
632  buildNodeList();
634  cacheInfo();
635  buildElemIDInfo();
636 
637  // this will make moose mesh aware of p-refinement added by mesh generators including
638  // a file mesh generator loading a restart checkpoint file
639  _max_p_level = 0;
640  _max_h_level = 0;
641  for (const auto & elem : getMesh().active_local_element_ptr_range())
642  {
643  if (elem->p_level() > _max_p_level)
645  if (elem->level() > _max_h_level)
646  _max_h_level = elem->level();
647  }
648  comm().max(_max_p_level);
649  comm().max(_max_h_level);
650 
651  // the flag might have been set by calling doingPRefinement(true)
653 
655 }
656 
657 void
659 {
660  auto & mesh = getMesh();
661 
662  if (!mesh.is_serial())
663  mooseError(
664  "Hybrid finite element method must use replicated mesh.\nCurrently lower-dimensional mesh "
665  "does not support mesh re-partitioning and a debug assertion being hit related with "
666  "neighbors of lower-dimensional element, with distributed mesh.");
667 
668  // Lower-D element build requires neighboring element information
669  if (!mesh.is_prepared())
671 
672  // maximum number of sides of all elements
673  unsigned int max_n_sides = 0;
674 
675  // remove existing lower-d element first
676  std::set<Elem *> deleteable_elems;
677  for (auto & elem : mesh.element_ptr_range())
680  deleteable_elems.insert(elem);
681  else if (elem->n_sides() > max_n_sides)
682  max_n_sides = elem->n_sides();
683 
684  for (auto & elem : deleteable_elems)
686  for (const auto & id : _lower_d_interior_blocks)
687  _mesh_subdomains.erase(id);
688  for (const auto & id : _lower_d_boundary_blocks)
689  _mesh_subdomains.erase(id);
690  _lower_d_interior_blocks.clear();
691  _lower_d_boundary_blocks.clear();
692 
693  mesh.comm().max(max_n_sides);
694 
695  deleteable_elems.clear();
696 
697  // get all side types
698  std::set<int> interior_side_types;
699  std::set<int> boundary_side_types;
700  for (const auto & elem : mesh.active_element_ptr_range())
701  for (const auto side : elem->side_index_range())
702  {
703  Elem * neig = elem->neighbor_ptr(side);
704  std::unique_ptr<Elem> side_elem(elem->build_side_ptr(side));
705  if (neig)
706  interior_side_types.insert(side_elem->type());
707  else
708  boundary_side_types.insert(side_elem->type());
709  }
710  mesh.comm().set_union(interior_side_types);
711  mesh.comm().set_union(boundary_side_types);
712 
713  // assign block ids for different side types
714  std::map<ElemType, SubdomainID> interior_block_ids;
715  std::map<ElemType, SubdomainID> boundary_block_ids;
716  // we assume this id is not used by the mesh
718  for (const auto & tpid : interior_side_types)
719  {
720  const auto type = ElemType(tpid);
721  mesh.subdomain_name(id) = "INTERNAL_SIDE_LOWERD_SUBDOMAIN_" + Utility::enum_to_string(type);
722  interior_block_ids[type] = id;
723  _lower_d_interior_blocks.insert(id);
724  if (_mesh_subdomains.count(id) > 0)
725  mooseError("Trying to add a mesh block with id ", id, " that has existed in the mesh");
726  _mesh_subdomains.insert(id);
727  --id;
728  }
729  for (const auto & tpid : boundary_side_types)
730  {
731  const auto type = ElemType(tpid);
732  mesh.subdomain_name(id) = "BOUNDARY_SIDE_LOWERD_SUBDOMAIN_" + Utility::enum_to_string(type);
733  boundary_block_ids[type] = id;
734  _lower_d_boundary_blocks.insert(id);
735  if (_mesh_subdomains.count(id) > 0)
736  mooseError("Trying to add a mesh block with id ", id, " that has existed in the mesh");
737  _mesh_subdomains.insert(id);
738  --id;
739  }
740 
741  dof_id_type max_elem_id = mesh.max_elem_id();
742  unique_id_type max_unique_id = mesh.parallel_max_unique_id();
743 
744  std::vector<Elem *> side_elems;
746  for (const auto & elem : mesh.active_element_ptr_range())
747  {
748  // skip existing lower-d elements
749  if (elem->interior_parent())
750  continue;
751 
752  for (const auto side : elem->side_index_range())
753  {
754  Elem * neig = elem->neighbor_ptr(side);
755 
756  bool build_side = false;
757  if (!neig)
758  build_side = true;
759  else
760  {
761  mooseAssert(!neig->is_remote(), "We error if the mesh is not serial");
762  if (!neig->active())
763  build_side = true;
764  else if (neig->level() == elem->level() && elem->id() < neig->id())
765  build_side = true;
766  }
767 
768  if (build_side)
769  {
770  std::unique_ptr<Elem> side_elem(elem->build_side_ptr(side));
771 
772  // The side will be added with the same processor id as the parent.
773  side_elem->processor_id() = elem->processor_id();
774 
775  // Add subdomain ID
776  if (neig)
777  side_elem->subdomain_id() = interior_block_ids.at(side_elem->type());
778  else
779  side_elem->subdomain_id() = boundary_block_ids.at(side_elem->type());
780 
781  // set ids consistently across processors (these ids will be temporary)
782  side_elem->set_id(max_elem_id + elem->id() * max_n_sides + side);
783  side_elem->set_unique_id(max_unique_id + elem->id() * max_n_sides + side);
784 
785  // Also assign the side's interior parent, so it is always
786  // easy to figure out the Elem we came from.
787  // Note: the interior parent could be a ghost element.
788  side_elem->set_interior_parent(elem);
789 
790  side_elems.push_back(side_elem.release());
791 
792  // add link between higher d element to lower d element
793  auto pair = std::make_pair(elem, side);
794  auto link = std::make_pair(pair, side_elems.back());
795  auto ilink = std::make_pair(side_elems.back(), side);
798  }
799  }
800  }
801 
802  // finally, add the lower-dimensional element to the mesh
803  // Note: lower-d interior element will exist on a processor if its associated interior
804  // parent exists on a processor whether or not being a ghost. Lower-d elements will
805  // get its interior parent's processor id.
806  for (auto & elem : side_elems)
807  mesh.add_elem(elem);
808 
809  // we do all the stuff in prepare_for_use such as renumber_nodes_and_elements(),
810  // update_parallel_id_counts(), cache_elem_dims(), etc. except partitioning here.
811  const bool skip_partitioning_old = mesh.skip_partitioning();
812  mesh.skip_partitioning(true);
813  // Finding neighbors is ambiguous for lower-dimensional elements on interior faces
814  mesh.allow_find_neighbors(false);
816  mesh.skip_partitioning(skip_partitioning_old);
817 }
818 
819 const Node &
821 {
822  mooseDeprecated("MooseMesh::node() is deprecated, please use MooseMesh::nodeRef() instead");
823  return nodeRef(i);
824 }
825 
826 Node &
828 {
829  mooseDeprecated("MooseMesh::node() is deprecated, please use MooseMesh::nodeRef() instead");
830  return nodeRef(i);
831 }
832 
833 const Node &
835 {
836  const auto node_ptr = queryNodePtr(i);
837  mooseAssert(node_ptr, "Missing node");
838  return *node_ptr;
839 }
840 
841 Node &
843 {
844  return const_cast<Node &>(const_cast<const MooseMesh *>(this)->nodeRef(i));
845 }
846 
847 const Node *
849 {
850  return &nodeRef(i);
851 }
852 
853 Node *
855 {
856  return &nodeRef(i);
857 }
858 
859 const Node *
861 {
862  if (i > getMesh().max_node_id())
863  {
864  auto it = _quadrature_nodes.find(i);
865  if (it == _quadrature_nodes.end())
866  return nullptr;
867  auto & node_ptr = it->second;
868  mooseAssert(node_ptr, "Uninitialized quadrature node");
869  return node_ptr;
870  }
871 
872  return getMesh().query_node_ptr(i);
873 }
874 
875 Node *
877 {
878  return const_cast<Node *>(const_cast<const MooseMesh *>(this)->queryNodePtr(i));
879 }
880 
881 void
883 {
884  TIME_SECTION("meshChanged", 3, "Updating Because Mesh Changed");
885 
886  update();
887 
888  // Delete all of the cached ranges
889  _active_local_elem_range.reset();
890  _active_node_range.reset();
892  _local_node_range.reset();
893  _bnd_node_range.reset();
894  _bnd_elem_range.reset();
895 
896  // Rebuild the ranges
902 
903  // Call the callback function onMeshChanged
904  onMeshChanged();
905 }
906 
907 void
909 {
910 }
911 
912 void
914 {
915  TIME_SECTION("cacheChangedLists", 5, "Caching Changed Lists");
916 
917  ConstElemRange elem_range(getMesh().local_elements_begin(), getMesh().local_elements_end(), 1);
918  CacheChangedListsThread cclt(*this);
919  Threads::parallel_reduce(elem_range, cclt);
920 
922 
923  _refined_elements = std::make_unique<ConstElemPointerRange>(cclt._refined_elements.begin(),
924  cclt._refined_elements.end());
925  _coarsened_elements = std::make_unique<ConstElemPointerRange>(cclt._coarsened_elements.begin(),
926  cclt._coarsened_elements.end());
928 }
929 
932 {
933  return _refined_elements.get();
934 }
935 
938 {
939  return _coarsened_elements.get();
940 }
941 
942 const std::vector<const Elem *> &
944 {
945  auto elem_to_child_pair = _coarsened_element_children.find(elem);
946  mooseAssert(elem_to_child_pair != _coarsened_element_children.end(), "Missing element in map");
947  return elem_to_child_pair->second;
948 }
949 
950 void
951 MooseMesh::updateActiveSemiLocalNodeRange(std::set<dof_id_type> & ghosted_elems)
952 {
953  TIME_SECTION("updateActiveSemiLocalNodeRange", 5, "Updating ActiveSemiLocalNode Range");
954 
955  _semilocal_node_list.clear();
956 
957  // First add the nodes connected to local elems
958  ConstElemRange * active_local_elems = getActiveLocalElementRange();
959  for (const auto & elem : *active_local_elems)
960  {
961  for (unsigned int n = 0; n < elem->n_nodes(); ++n)
962  {
963  // Since elem is const here but we require a non-const Node * to
964  // store in the _semilocal_node_list (otherwise things like
965  // UpdateDisplacedMeshThread don't work), we are using a
966  // const_cast. A more long-term fix would be to have
967  // getActiveLocalElementRange return a non-const ElemRange.
968  Node * node = const_cast<Node *>(elem->node_ptr(n));
969 
970  _semilocal_node_list.insert(node);
971  }
972  }
973 
974  // Now add the nodes connected to ghosted_elems
975  for (const auto & ghost_elem_id : ghosted_elems)
976  {
977  Elem * elem = getMesh().elem_ptr(ghost_elem_id);
978  for (unsigned int n = 0; n < elem->n_nodes(); n++)
979  {
980  Node * node = elem->node_ptr(n);
981 
982  _semilocal_node_list.insert(node);
983  }
984  }
985 
986  // Now create the actual range
987  _active_semilocal_node_range = std::make_unique<SemiLocalNodeRange>(_semilocal_node_list.begin(),
988  _semilocal_node_list.end());
989 }
990 
991 bool
992 MooseMesh::isSemiLocal(Node * const node) const
993 {
994  return _semilocal_node_list.find(node) != _semilocal_node_list.end();
995 }
996 
1002 {
1003 public:
1005 
1006  bool operator()(const BndNode * const & lhs, const BndNode * const & rhs)
1007  {
1008  if (lhs->_bnd_id < rhs->_bnd_id)
1009  return true;
1010 
1011  if (lhs->_bnd_id > rhs->_bnd_id)
1012  return false;
1013 
1014  if (lhs->_node->id() < rhs->_node->id())
1015  return true;
1016 
1017  if (lhs->_node->id() > rhs->_node->id())
1018  return false;
1019 
1020  return false;
1021  }
1022 };
1023 
1024 void
1026 {
1027  TIME_SECTION("buildNodeList", 5, "Building Node List");
1028 
1029  freeBndNodes();
1030 
1031  auto bc_tuples = getMesh().get_boundary_info().build_node_list();
1032 
1033  int n = bc_tuples.size();
1034  _bnd_nodes.clear();
1035  _bnd_nodes.reserve(n);
1036  for (const auto & t : bc_tuples)
1037  {
1038  auto node_id = std::get<0>(t);
1039  auto bc_id = std::get<1>(t);
1040 
1041  _bnd_nodes.push_back(new BndNode(getMesh().node_ptr(node_id), bc_id));
1042  _node_set_nodes[bc_id].push_back(node_id);
1043  _bnd_node_ids[bc_id].insert(node_id);
1044  }
1045 
1046  _bnd_nodes.reserve(_bnd_nodes.size() + _extra_bnd_nodes.size());
1047  for (unsigned int i = 0; i < _extra_bnd_nodes.size(); i++)
1048  {
1049  BndNode * bnode = new BndNode(_extra_bnd_nodes[i]._node, _extra_bnd_nodes[i]._bnd_id);
1050  _bnd_nodes.push_back(bnode);
1051  _bnd_node_ids[std::get<1>(bc_tuples[i])].insert(_extra_bnd_nodes[i]._node->id());
1052  }
1053 
1054  // This sort is here so that boundary conditions are always applied in the same order
1055  std::sort(_bnd_nodes.begin(), _bnd_nodes.end(), BndNodeCompare());
1056 }
1057 
1058 void
1060 {
1061  unsigned int n = getMesh().n_elem_integers() + 1;
1062 
1063  _block_id_mapping.clear();
1064  _max_ids.clear();
1065  _min_ids.clear();
1066  _id_identical_flag.clear();
1067 
1068  _block_id_mapping.resize(n);
1071  _id_identical_flag.resize(n, std::vector<bool>(n, true));
1072  for (const auto & elem : getMesh().active_local_element_ptr_range())
1073  for (unsigned int i = 0; i < n; ++i)
1074  {
1075  auto id = (i == n - 1 ? elem->subdomain_id() : elem->get_extra_integer(i));
1076  _block_id_mapping[i][elem->subdomain_id()].insert(id);
1077  if (id > _max_ids[i])
1078  _max_ids[i] = id;
1079  if (id < _min_ids[i])
1080  _min_ids[i] = id;
1081  for (unsigned int j = 0; j < n; ++j)
1082  {
1083  auto idj = (j == n - 1 ? elem->subdomain_id() : elem->get_extra_integer(j));
1084  if (i != j && _id_identical_flag[i][j] && id != idj)
1085  _id_identical_flag[i][j] = false;
1086  }
1087  }
1088 
1089  for (unsigned int i = 0; i < n; ++i)
1090  {
1091  for (auto & blk : meshSubdomains())
1092  comm().set_union(_block_id_mapping[i][blk]);
1093  comm().min(_id_identical_flag[i]);
1094  }
1095  comm().max(_max_ids);
1096  comm().min(_min_ids);
1097 }
1098 
1099 std::unordered_map<dof_id_type, std::set<dof_id_type>>
1100 MooseMesh::getElemIDMapping(const std::string & from_id_name, const std::string & to_id_name) const
1101 {
1102  auto & mesh_base = getMesh();
1103 
1104  if (!mesh_base.has_elem_integer(from_id_name))
1105  mooseError("Mesh does not have the element integer name '", from_id_name, "'");
1106  if (!mesh_base.has_elem_integer(to_id_name))
1107  mooseError("Mesh does not have the element integer name '", to_id_name, "'");
1108 
1109  const auto id1 = mesh_base.get_elem_integer_index(from_id_name);
1110  const auto id2 = mesh_base.get_elem_integer_index(to_id_name);
1111 
1112  std::unordered_map<dof_id_type, std::set<dof_id_type>> id_map;
1113  for (const auto id : getAllElemIDs(id1))
1114  id_map[id] = std::set<dof_id_type>();
1115 
1116  for (const auto & elem : mesh_base.active_local_element_ptr_range())
1117  id_map[elem->get_extra_integer(id1)].insert(elem->get_extra_integer(id2));
1118 
1119  for (auto & [id, ids] : id_map)
1120  {
1121  libmesh_ignore(id); // avoid overzealous gcc 9.4 unused var warning
1122  comm().set_union(ids);
1123  }
1124 
1125  return id_map;
1126 }
1127 
1128 std::set<dof_id_type>
1129 MooseMesh::getAllElemIDs(unsigned int elem_id_index) const
1130 {
1131  std::set<dof_id_type> unique_ids;
1132  for (auto & pair : _block_id_mapping[elem_id_index])
1133  for (auto & id : pair.second)
1134  unique_ids.insert(id);
1135  return unique_ids;
1136 }
1137 
1138 std::set<dof_id_type>
1139 MooseMesh::getElemIDsOnBlocks(unsigned int elem_id_index, const std::set<SubdomainID> & blks) const
1140 {
1141  std::set<dof_id_type> unique_ids;
1142  for (auto & blk : blks)
1143  {
1144  auto it = _block_id_mapping[elem_id_index].find(blk);
1145  if (it == _block_id_mapping[elem_id_index].end())
1146  mooseError("Block ", blk, " is not available on the mesh");
1147 
1148  for (auto & mid : it->second)
1149  unique_ids.insert(mid);
1150  }
1151  return unique_ids;
1152 }
1153 
1154 void
1156 {
1157  TIME_SECTION("buildBndElemList", 5, "Building Boundary Elements List");
1158 
1159  freeBndElems();
1160 
1161  auto bc_tuples = getMesh().get_boundary_info().build_active_side_list();
1162 
1163  int n = bc_tuples.size();
1164  _bnd_elems.clear();
1165  _bnd_elems.reserve(n);
1166  for (const auto & t : bc_tuples)
1167  {
1168  auto elem_id = std::get<0>(t);
1169  auto side_id = std::get<1>(t);
1170  auto bc_id = std::get<2>(t);
1171 
1172  _bnd_elems.push_back(new BndElement(getMesh().elem_ptr(elem_id), side_id, bc_id));
1173  _bnd_elem_ids[bc_id].insert(elem_id);
1174  }
1175 }
1176 
1177 const std::map<dof_id_type, std::vector<dof_id_type>> &
1179 {
1180  if (!_node_to_elem_map_built) // Guard the creation with a double checked lock
1181  {
1182  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1183 
1185  {
1186  // This is allowing the timing to be run even with threads
1187  // This is safe because all threads will be waiting on this section when it runs
1188  // NOTE: Do not copy this construction to other places without thinking REALLY hard about it
1189  // The PerfGraph is NOT threadsafe and will cause all kinds of havok if care isn't taken
1190  auto in_threads = Threads::in_threads;
1191  Threads::in_threads = false;
1192  TIME_SECTION("nodeToElemMap", 5, "Building Node To Elem Map");
1193  Threads::in_threads = in_threads;
1194 
1195  for (const auto & elem : getMesh().active_element_ptr_range())
1196  for (unsigned int n = 0; n < elem->n_nodes(); n++)
1197  _node_to_elem_map[elem->node_id(n)].push_back(elem->id());
1198 
1199  _node_to_elem_map_built = true; // MUST be set at the end for double-checked locking to work!
1200  }
1201  }
1202  return _node_to_elem_map;
1203 }
1204 
1205 const std::map<dof_id_type, std::vector<dof_id_type>> &
1207 {
1208  if (!_node_to_active_semilocal_elem_map_built) // Guard the creation with a double checked lock
1209  {
1210  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1211 
1212  // This is allowing the timing to be run even with threads
1213  // This is safe because all threads will be waiting on this section when it runs
1214  // NOTE: Do not copy this construction to other places without thinking REALLY hard about it
1215  // The PerfGraph is NOT threadsafe and will cause all kinds of havok if care isn't taken
1216  auto in_threads = Threads::in_threads;
1217  Threads::in_threads = false;
1218  TIME_SECTION("nodeToActiveSemilocalElemMap", 5, "Building SemiLocalElemMap");
1219  Threads::in_threads = in_threads;
1220 
1222  {
1223  for (const auto & elem :
1224  as_range(getMesh().semilocal_elements_begin(), getMesh().semilocal_elements_end()))
1225  if (elem->active())
1226  for (unsigned int n = 0; n < elem->n_nodes(); n++)
1228 
1230  true; // MUST be set at the end for double-checked locking to work!
1231  }
1232  }
1233 
1235 }
1236 
1239 {
1241  {
1242  TIME_SECTION("getActiveLocalElementRange", 5);
1243 
1244  _active_local_elem_range = std::make_unique<ConstElemRange>(
1245  getMesh().active_local_elements_begin(), getMesh().active_local_elements_end(), GRAIN_SIZE);
1246  }
1247 
1248  return _active_local_elem_range.get();
1249 }
1250 
1251 NodeRange *
1253 {
1254  if (!_active_node_range)
1255  {
1256  TIME_SECTION("getActiveNodeRange", 5);
1257 
1258  _active_node_range = std::make_unique<NodeRange>(
1259  getMesh().active_nodes_begin(), getMesh().active_nodes_end(), GRAIN_SIZE);
1260  }
1261 
1262  return _active_node_range.get();
1263 }
1264 
1267 {
1268  mooseAssert(_active_semilocal_node_range,
1269  "_active_semilocal_node_range has not been created yet!");
1270 
1271  return _active_semilocal_node_range.get();
1272 }
1273 
1276 {
1277  if (!_local_node_range)
1278  {
1279  TIME_SECTION("getLocalNodeRange", 5);
1280 
1281  _local_node_range = std::make_unique<ConstNodeRange>(
1282  getMesh().local_nodes_begin(), getMesh().local_nodes_end(), GRAIN_SIZE);
1283  }
1284 
1285  return _local_node_range.get();
1286 }
1287 
1290 {
1291  if (!_bnd_node_range)
1292  {
1293  TIME_SECTION("getBoundaryNodeRange", 5);
1294 
1295  _bnd_node_range =
1296  std::make_unique<ConstBndNodeRange>(bndNodesBegin(), bndNodesEnd(), GRAIN_SIZE);
1297  }
1298 
1299  return _bnd_node_range.get();
1300 }
1301 
1304 {
1305  if (!_bnd_elem_range)
1306  {
1307  TIME_SECTION("getBoundaryElementRange", 5);
1308 
1309  _bnd_elem_range =
1310  std::make_unique<ConstBndElemRange>(bndElemsBegin(), bndElemsEnd(), GRAIN_SIZE);
1311  }
1312 
1313  return _bnd_elem_range.get();
1314 }
1315 
1316 const std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>> &
1318 {
1319  mooseDeprecated("MooseMesh::getBoundariesToElems is deprecated, "
1320  "use MooseMesh::getBoundariesToActiveSemiLocalElemIds");
1322 }
1323 
1324 const std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>> &
1326 {
1327  return _bnd_elem_ids;
1328 }
1329 
1330 std::unordered_set<dof_id_type>
1332 {
1333  // The boundary to element map is computed on every mesh update
1334  const auto it = _bnd_elem_ids.find(bid);
1335  if (it == _bnd_elem_ids.end())
1336  // Boundary is not local to this domain, return an empty set
1337  return std::unordered_set<dof_id_type>{};
1338  return it->second;
1339 }
1340 
1341 std::unordered_set<dof_id_type>
1343 {
1344  // Vector of boundary elems is updated every mesh update
1345  std::unordered_set<dof_id_type> neighbor_elems;
1346  for (const auto & bnd_elem : _bnd_elems)
1347  {
1348  const auto & [elem_ptr, elem_side, elem_bid] = *bnd_elem;
1349  if (elem_bid == bid)
1350  {
1351  const auto * neighbor = elem_ptr->neighbor_ptr(elem_side);
1352  // Dont add fully remote elements, ghosted is fine
1353  if (neighbor && neighbor != libMesh::remote_elem)
1354  {
1355  // handle mesh refinement, only return active elements near the boundary
1356  if (neighbor->active())
1357  neighbor_elems.insert(neighbor->id());
1358  else
1359  {
1360  std::vector<const Elem *> family;
1361  neighbor->active_family_tree_by_neighbor(family, elem_ptr);
1362  for (const auto & child_neighbor : family)
1363  neighbor_elems.insert(child_neighbor->id());
1364  }
1365  }
1366  }
1367  }
1368 
1369  return neighbor_elems;
1370 }
1371 
1372 bool
1374  const std::set<SubdomainID> & blk_group) const
1375 {
1376  mooseAssert(_bnd_elem_range, "Boundary element range is not initialized");
1377 
1378  // Loop over all side elements of the mesh, select those on the boundary
1379  for (const auto & bnd_elem : *_bnd_elem_range)
1380  {
1381  const auto & [elem_ptr, elem_side, elem_bid] = *bnd_elem;
1382  if (elem_bid == bid)
1383  {
1384  // If an element is internal to the group of subdomain, check the neighbor
1385  if (blk_group.find(elem_ptr->subdomain_id()) != blk_group.end())
1386  {
1387  const auto * const neighbor = elem_ptr->neighbor_ptr(elem_side);
1388 
1389  // If we did not ghost the neighbor, we cannot decide
1390  if (neighbor == libMesh::remote_elem)
1391  mooseError("Insufficient level of geometrical ghosting to determine "
1392  "if a boundary is internal to the mesh");
1393  // If the neighbor does not exist, then we are on the edge of the mesh
1394  if (!neighbor)
1395  continue;
1396  // If the neighbor is also in the group of subdomain,
1397  // then the boundary cuts the subdomains
1398  if (blk_group.find(neighbor->subdomain_id()) != blk_group.end())
1399  return false;
1400  }
1401  }
1402  }
1403  return true;
1404 }
1405 
1406 void
1408 {
1409  TIME_SECTION("cacheInfo", 3);
1410 
1411  _has_lower_d = false;
1412  _sub_to_data.clear();
1414  _block_node_list.clear();
1417  _lower_d_interior_blocks.clear();
1418  _lower_d_boundary_blocks.clear();
1419 
1420  auto & mesh = getMesh();
1421 
1422  // TODO: Thread this!
1423  for (const auto & elem : mesh.element_ptr_range())
1424  {
1425  const Elem * ip_elem = elem->interior_parent();
1426 
1427  if (ip_elem)
1428  {
1429  if (elem->active())
1430  _sub_to_data[elem->subdomain_id()].is_lower_d = true;
1431  unsigned int ip_side = ip_elem->which_side_am_i(elem);
1432 
1433  // For some grid sequencing tests: ip_side == libMesh::invalid_uint
1434  if (ip_side != libMesh::invalid_uint)
1435  {
1436  auto pair = std::make_pair(ip_elem, ip_side);
1438  std::pair<std::pair<const Elem *, unsigned short int>, const Elem *>(pair, elem));
1440  std::pair<const Elem *, unsigned short int>(elem, ip_side));
1441 
1442  auto id = elem->subdomain_id();
1443  if (ip_elem->neighbor_ptr(ip_side))
1444  {
1445  if (mesh.subdomain_name(id).find("INTERNAL_SIDE_LOWERD_SUBDOMAIN_") != std::string::npos)
1446  _lower_d_interior_blocks.insert(id);
1447  }
1448  else
1449  {
1450  if (mesh.subdomain_name(id).find("BOUNDARY_SIDE_LOWERD_SUBDOMAIN_") != std::string::npos)
1451  _lower_d_boundary_blocks.insert(id);
1452  }
1453  }
1454  }
1455 
1456  for (unsigned int nd = 0; nd < elem->n_nodes(); ++nd)
1457  {
1458  Node & node = *elem->node_ptr(nd);
1459  _block_node_list[node.id()].insert(elem->subdomain_id());
1460  }
1461  }
1464 
1465  for (const auto & elem : mesh.active_local_element_ptr_range())
1466  {
1467  SubdomainID subdomain_id = elem->subdomain_id();
1468  auto & sub_data = _sub_to_data[subdomain_id];
1469  for (unsigned int side = 0; side < elem->n_sides(); side++)
1470  {
1471  std::vector<BoundaryID> boundary_ids = getBoundaryIDs(elem, side);
1472  sub_data.boundary_ids.insert(boundary_ids.begin(), boundary_ids.end());
1473 
1474  Elem * neig = elem->neighbor_ptr(side);
1475  if (neig)
1476  {
1477  _neighbor_subdomain_boundary_ids[neig->subdomain_id()].insert(boundary_ids.begin(),
1478  boundary_ids.end());
1479  SubdomainID neighbor_subdomain_id = neig->subdomain_id();
1480  if (neighbor_subdomain_id != subdomain_id)
1481  sub_data.neighbor_subs.insert(neighbor_subdomain_id);
1482  }
1483  }
1484  }
1485 
1486  for (const auto blk_id : _mesh_subdomains)
1487  {
1488  auto & sub_data = _sub_to_data[blk_id];
1489  _communicator.set_union(sub_data.neighbor_subs);
1490  _communicator.set_union(sub_data.boundary_ids);
1491  _communicator.max(sub_data.is_lower_d);
1492  if (sub_data.is_lower_d)
1493  _has_lower_d = true;
1495  }
1496 }
1497 
1498 const std::set<SubdomainID> &
1500 {
1501  auto it = _block_node_list.find(node.id());
1502 
1503  if (it == _block_node_list.end())
1504  mooseError("Unable to find node: ", node.id(), " in any block list.");
1505 
1506  return it->second;
1507 }
1508 
1511 {
1512  return face_info_iterator(
1513  _face_info.begin(),
1514  _face_info.end(),
1516 }
1517 
1520 {
1521  return face_info_iterator(
1522  _face_info.end(),
1523  _face_info.end(),
1525 }
1526 
1529 {
1530  return elem_info_iterator(_elem_info.begin(),
1531  _elem_info.end(),
1533 }
1534 
1537 {
1538  return elem_info_iterator(_elem_info.end(),
1539  _elem_info.end(),
1541 }
1542 
1543 // default begin() accessor
1546 {
1548  return bnd_node_iterator(_bnd_nodes.begin(), _bnd_nodes.end(), p);
1549 }
1550 
1551 // default end() accessor
1554 {
1556  return bnd_node_iterator(_bnd_nodes.end(), _bnd_nodes.end(), p);
1557 }
1558 
1559 // default begin() accessor
1562 {
1564  return bnd_elem_iterator(_bnd_elems.begin(), _bnd_elems.end(), p);
1565 }
1566 
1567 // default end() accessor
1570 {
1572  return bnd_elem_iterator(_bnd_elems.end(), _bnd_elems.end(), p);
1573 }
1574 
1575 const Node *
1576 MooseMesh::addUniqueNode(const Point & p, Real tol)
1577 {
1582  if (getMesh().n_nodes() != _node_map.size())
1583  {
1584  _node_map.clear();
1585  _node_map.reserve(getMesh().n_nodes());
1586  for (const auto & node : getMesh().node_ptr_range())
1587  _node_map.push_back(node);
1588  }
1589 
1590  Node * node = nullptr;
1591  for (unsigned int i = 0; i < _node_map.size(); ++i)
1592  {
1593  if (p.relative_fuzzy_equals(*_node_map[i], tol))
1594  {
1595  node = _node_map[i];
1596  break;
1597  }
1598  }
1599  if (node == nullptr)
1600  {
1601  node = getMesh().add_node(new Node(p));
1602  _node_map.push_back(node);
1603  }
1604 
1605  mooseAssert(node != nullptr, "Node is NULL");
1606  return node;
1607 }
1608 
1609 Node *
1611  const unsigned short int side,
1612  const unsigned int qp,
1613  BoundaryID bid,
1614  const Point & point)
1615 {
1616  Node * qnode;
1617 
1618  if (_elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side].find(qp) ==
1620  {
1621  // Create a new node id starting from the max node id and counting down. This will be the least
1622  // likely to collide with an existing node id.
1623  // Note that we are using numeric_limits<unsigned>::max even
1624  // though max_id is stored as a dof_id_type. I tried this with
1625  // numeric_limits<dof_id_type>::max and it broke several tests in
1626  // MOOSE. So, this is some kind of a magic number that we will
1627  // just continue to use...
1629  dof_id_type new_id = max_id - _quadrature_nodes.size();
1630 
1631  if (new_id <= getMesh().max_node_id())
1632  mooseError("Quadrature node id collides with existing node id!");
1633 
1634  qnode = new Node(point, new_id);
1635 
1636  // Keep track of this new node in two different ways for easy lookup
1637  _quadrature_nodes[new_id] = qnode;
1638  _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp] = qnode;
1639 
1640  if (elem->active())
1641  {
1642  _node_to_elem_map[new_id].push_back(elem->id());
1643  _node_to_active_semilocal_elem_map[new_id].push_back(elem->id());
1644  }
1645  }
1646  else
1647  qnode = _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp];
1648 
1649  BndNode * bnode = new BndNode(qnode, bid);
1650  _bnd_nodes.push_back(bnode);
1651  _bnd_node_ids[bid].insert(qnode->id());
1652 
1653  _extra_bnd_nodes.push_back(*bnode);
1654 
1655  // Do this so the range will be regenerated next time it is accessed
1656  _bnd_node_range.reset();
1657 
1658  return qnode;
1659 }
1660 
1661 Node *
1663  const unsigned short int side,
1664  const unsigned int qp)
1665 {
1666  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes.find(elem->id()) !=
1668  "Elem has no quadrature nodes!");
1669  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes[elem->id()].find(side) !=
1671  "Side has no quadrature nodes!");
1672  mooseAssert(_elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side].find(qp) !=
1674  "qp not found on side!");
1675 
1676  return _elem_to_side_to_qp_to_quadrature_nodes[elem->id()][side][qp];
1677 }
1678 
1679 void
1681 {
1682  // Delete all the quadrature nodes
1683  for (auto & it : _quadrature_nodes)
1684  delete it.second;
1685 
1686  _quadrature_nodes.clear();
1688  _extra_bnd_nodes.clear();
1689 }
1690 
1691 BoundaryID
1692 MooseMesh::getBoundaryID(const BoundaryName & boundary_name) const
1693 {
1694  if (boundary_name == "ANY_BOUNDARY_ID")
1695  mooseError("Please use getBoundaryIDs() when passing \"ANY_BOUNDARY_ID\"");
1696 
1697  return MooseMeshUtils::getBoundaryID(boundary_name, getMesh());
1698 }
1699 
1700 const Elem *
1701 MooseMesh::getLowerDElem(const Elem * elem, unsigned short int side) const
1702 {
1703  auto it = _higher_d_elem_side_to_lower_d_elem.find(std::make_pair(elem, side));
1704 
1705  if (it != _higher_d_elem_side_to_lower_d_elem.end())
1706  return it->second;
1707  else
1708  return nullptr;
1709 }
1710 
1711 unsigned int
1712 MooseMesh::getHigherDSide(const Elem * elem) const
1713 {
1714  auto it = _lower_d_elem_to_higher_d_elem_side.find(elem);
1715 
1716  if (it != _lower_d_elem_to_higher_d_elem_side.end())
1717  return it->second;
1718  else
1719  return libMesh::invalid_uint;
1720 }
1721 
1722 std::vector<BoundaryID>
1723 MooseMesh::getBoundaryIDs(const std::vector<BoundaryName> & boundary_name,
1724  bool generate_unknown) const
1725 {
1727  getMesh(), boundary_name, generate_unknown, _mesh_boundary_ids);
1728 }
1729 
1731 MooseMesh::getSubdomainID(const SubdomainName & subdomain_name) const
1732 {
1733  return MooseMeshUtils::getSubdomainID(subdomain_name, getMesh());
1734 }
1735 
1736 std::vector<SubdomainID>
1737 MooseMesh::getSubdomainIDs(const std::vector<SubdomainName> & subdomain_name) const
1738 {
1739  return MooseMeshUtils::getSubdomainIDs(getMesh(), subdomain_name);
1740 }
1741 
1742 std::set<SubdomainID>
1743 MooseMesh::getSubdomainIDs(const std::set<SubdomainName> & subdomain_name) const
1744 {
1745  return MooseMeshUtils::getSubdomainIDs(getMesh(), subdomain_name);
1746 }
1747 
1748 void
1749 MooseMesh::setSubdomainName(SubdomainID subdomain_id, const SubdomainName & name)
1750 {
1751  mooseAssert(name != "ANY_BLOCK_ID", "Cannot set subdomain name to 'ANY_BLOCK_ID'");
1752  getMesh().subdomain_name(subdomain_id) = name;
1753 }
1754 
1755 void
1756 MooseMesh::setSubdomainName(MeshBase & mesh, SubdomainID subdomain_id, const SubdomainName & name)
1757 {
1758  mooseAssert(name != "ANY_BLOCK_ID", "Cannot set subdomain name to 'ANY_BLOCK_ID'");
1759  mesh.subdomain_name(subdomain_id) = name;
1760 }
1761 
1762 const std::string &
1764 {
1765  return getMesh().subdomain_name(subdomain_id);
1766 }
1767 
1768 std::vector<SubdomainName>
1769 MooseMesh::getSubdomainNames(const std::vector<SubdomainID> & subdomain_ids) const
1770 {
1771  std::vector<SubdomainName> names(subdomain_ids.size());
1772 
1773  for (unsigned int i = 0; i < subdomain_ids.size(); i++)
1774  names[i] = getSubdomainName(subdomain_ids[i]);
1775 
1776  return names;
1777 }
1778 
1779 void
1780 MooseMesh::setBoundaryName(BoundaryID boundary_id, BoundaryName name)
1781 {
1782  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1783 
1784  // We need to figure out if this boundary is a sideset or nodeset
1785  if (boundary_info.get_side_boundary_ids().count(boundary_id))
1786  boundary_info.sideset_name(boundary_id) = name;
1787  else
1788  boundary_info.nodeset_name(boundary_id) = name;
1789 }
1790 
1791 const std::string &
1793 {
1794  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
1795 
1796  // We need to figure out if this boundary is a sideset or nodeset
1797  if (boundary_info.get_side_boundary_ids().count(boundary_id))
1798  return boundary_info.get_sideset_name(boundary_id);
1799  else
1800  return boundary_info.get_nodeset_name(boundary_id);
1801 }
1802 
1803 // specialization for PointListAdaptor<MooseMesh::PeriodicNodeInfo>
1804 template <>
1805 inline const Point &
1807  const MooseMesh::PeriodicNodeInfo & item) const
1808 {
1809  return *(item.first);
1810 }
1811 
1812 void
1813 MooseMesh::buildPeriodicNodeMap(std::multimap<dof_id_type, dof_id_type> & periodic_node_map,
1814  unsigned int var_number,
1815  libMesh::PeriodicBoundaries * pbs) const
1816 {
1817  TIME_SECTION("buildPeriodicNodeMap", 5);
1818 
1819  // clear existing map
1820  periodic_node_map.clear();
1821 
1822  // get periodic nodes
1823  std::vector<PeriodicNodeInfo> periodic_nodes;
1824  for (const auto & t : getMesh().get_boundary_info().build_node_list())
1825  {
1826  // unfortunately libMesh does not give us a pointer, so we have to look it up ourselves
1827  auto node = _mesh->node_ptr(std::get<0>(t));
1828  mooseAssert(node != nullptr,
1829  "libMesh::BoundaryInfo::build_node_list() returned an ID for a non-existing node");
1830  auto bc_id = std::get<1>(t);
1831  periodic_nodes.emplace_back(node, bc_id);
1832  }
1833 
1834  // sort by boundary id
1835  std::sort(periodic_nodes.begin(),
1836  periodic_nodes.end(),
1837  [](const PeriodicNodeInfo & a, const PeriodicNodeInfo & b) -> bool
1838  { return a.second > b.second; });
1839 
1840  // build kd-tree
1841  using KDTreeType = nanoflann::KDTreeSingleIndexAdaptor<
1842  nanoflann::L2_Simple_Adaptor<Real, PointListAdaptor<PeriodicNodeInfo>, Real, std::size_t>,
1844  LIBMESH_DIM,
1845  std::size_t>;
1846  const unsigned int max_leaf_size = 20; // slightly affects runtime
1847  auto point_list =
1848  PointListAdaptor<PeriodicNodeInfo>(periodic_nodes.begin(), periodic_nodes.end());
1849  auto kd_tree = std::make_unique<KDTreeType>(
1850  LIBMESH_DIM, point_list, nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size));
1851  mooseAssert(kd_tree != nullptr, "KDTree was not properly initialized.");
1852  kd_tree->buildIndex();
1853 
1854  // data structures for kd-tree search
1855  nanoflann::SearchParameters search_params;
1856  std::vector<nanoflann::ResultItem<std::size_t, Real>> ret_matches;
1857 
1858  // iterate over periodic nodes (boundary ids are in contiguous blocks)
1859  libMesh::PeriodicBoundaryBase * periodic = nullptr;
1860  BoundaryID current_bc_id = BoundaryInfo::invalid_id;
1861  for (auto & pair : periodic_nodes)
1862  {
1863  // entering a new block of boundary IDs
1864  if (pair.second != current_bc_id)
1865  {
1866  current_bc_id = pair.second;
1867  periodic = pbs->boundary(current_bc_id);
1868  if (periodic && !periodic->is_my_variable(var_number))
1869  periodic = nullptr;
1870  }
1871 
1872  // variable is not periodic at this node, skip
1873  if (!periodic)
1874  continue;
1875 
1876  // clear result buffer
1877  ret_matches.clear();
1878 
1879  // id of the current node
1880  const auto id = pair.first->id();
1881 
1882  // position where we expect a periodic partner for the current node and boundary
1883  Point search_point = periodic->get_corresponding_pos(*pair.first);
1884 
1885  // search at the expected point
1886  kd_tree->radiusSearch(&(search_point)(0), libMesh::TOLERANCE, ret_matches, search_params);
1887  for (auto & match_pair : ret_matches)
1888  {
1889  const auto & match = periodic_nodes[match_pair.first];
1890  // add matched node if the boundary id is the corresponding id in the periodic pair
1891  if (match.second == periodic->pairedboundary)
1892  periodic_node_map.emplace(id, match.first->id());
1893  }
1894  }
1895 }
1896 
1897 void
1898 MooseMesh::buildPeriodicNodeSets(std::map<BoundaryID, std::set<dof_id_type>> & periodic_node_sets,
1899  unsigned int var_number,
1900  libMesh::PeriodicBoundaries * pbs) const
1901 {
1902  TIME_SECTION("buildPeriodicNodeSets", 5);
1903 
1904  periodic_node_sets.clear();
1905 
1906  // Loop over all the boundary nodes adding the periodic nodes to the appropriate set
1907  for (const auto & t : getMesh().get_boundary_info().build_node_list())
1908  {
1909  auto node_id = std::get<0>(t);
1910  auto bc_id = std::get<1>(t);
1911 
1912  // Is this current node on a known periodic boundary?
1913  if (periodic_node_sets.find(bc_id) != periodic_node_sets.end())
1914  periodic_node_sets[bc_id].insert(node_id);
1915  else // This still might be a periodic node but we just haven't seen this boundary_id yet
1916  {
1917  const libMesh::PeriodicBoundaryBase * periodic = pbs->boundary(bc_id);
1918  if (periodic && periodic->is_my_variable(var_number))
1919  periodic_node_sets[bc_id].insert(node_id);
1920  }
1921  }
1922 }
1923 
1924 bool
1926 {
1927  TIME_SECTION("detectOrthogonalDimRanges", 5);
1928 
1930  return true;
1931 
1932  std::vector<Real> min(3, std::numeric_limits<Real>::max());
1933  std::vector<Real> max(3, std::numeric_limits<Real>::min());
1934  unsigned int dim = getMesh().mesh_dimension();
1935 
1936  // Find the bounding box of our mesh
1937  for (const auto & node : getMesh().node_ptr_range())
1938  // Check all coordinates, we don't know if this mesh might be lying in a higher dim even if the
1939  // mesh dimension is lower.
1940  for (const auto i : make_range(Moose::dim))
1941  {
1942  if ((*node)(i) < min[i])
1943  min[i] = (*node)(i);
1944  if ((*node)(i) > max[i])
1945  max[i] = (*node)(i);
1946  }
1947 
1948  this->comm().max(max);
1949  this->comm().min(min);
1950 
1951  _extreme_nodes.resize(8); // 2^LIBMESH_DIM
1952  // Now make sure that there are actual nodes at all of the extremes
1953  std::vector<bool> extreme_matches(8, false);
1954  std::vector<unsigned int> comp_map(3);
1955  for (const auto & node : getMesh().node_ptr_range())
1956  {
1957  // See if the current node is located at one of the extremes
1958  unsigned int coord_match = 0;
1959 
1960  for (const auto i : make_range(Moose::dim))
1961  {
1962  if (std::abs((*node)(i)-min[i]) < tol)
1963  {
1964  comp_map[i] = MIN;
1965  ++coord_match;
1966  }
1967  else if (std::abs((*node)(i)-max[i]) < tol)
1968  {
1969  comp_map[i] = MAX;
1970  ++coord_match;
1971  }
1972  }
1973 
1974  if (coord_match == LIBMESH_DIM) // Found a coordinate at one of the extremes
1975  {
1976  _extreme_nodes[comp_map[X] * 4 + comp_map[Y] * 2 + comp_map[Z]] = node;
1977  extreme_matches[comp_map[X] * 4 + comp_map[Y] * 2 + comp_map[Z]] = true;
1978  }
1979  }
1980 
1981  // See if we matched all of the extremes for the mesh dimension
1982  this->comm().max(extreme_matches);
1983  if (std::count(extreme_matches.begin(), extreme_matches.end(), true) == (1 << dim))
1984  _regular_orthogonal_mesh = true;
1985 
1986  // Set the bounds
1987  _bounds.resize(LIBMESH_DIM);
1988  for (const auto i : make_range(Moose::dim))
1989  {
1990  _bounds[i].resize(2);
1991  _bounds[i][MIN] = min[i];
1992  _bounds[i][MAX] = max[i];
1993  }
1994 
1995  return _regular_orthogonal_mesh;
1996 }
1997 
1998 void
2000 {
2001  TIME_SECTION("detectPairedSidesets", 5);
2002 
2003  // Loop over level-0 elements (since boundary condition information
2004  // is only directly stored for them) and find sidesets with normals
2005  // that point in the -x, +x, -y, +y, and -z, +z direction. If there
2006  // is a unique sideset id for each direction, then the paired
2007  // sidesets consist of (-x,+x), (-y,+y), (-z,+z). If there are
2008  // multiple sideset ids for a given direction, then we can't pick a
2009  // single pair for that direction. In that case, we'll just return
2010  // as was done in the original algorithm.
2011 
2012  // Points used for direction comparison
2013  const Point minus_x(-1, 0, 0), plus_x(1, 0, 0), minus_y(0, -1, 0), plus_y(0, 1, 0),
2014  minus_z(0, 0, -1), plus_z(0, 0, 1);
2015 
2016  // we need to test all element dimensions from dim down to 1
2017  const unsigned int dim = getMesh().mesh_dimension();
2018 
2019  // boundary id sets for elements of different dimensions
2020  std::vector<std::set<BoundaryID>> minus_x_ids(dim), plus_x_ids(dim), minus_y_ids(dim),
2021  plus_y_ids(dim), minus_z_ids(dim), plus_z_ids(dim);
2022 
2023  std::vector<std::unique_ptr<FEBase>> fe_faces(dim);
2024  std::vector<std::unique_ptr<libMesh::QGauss>> qfaces(dim);
2025  for (unsigned side_dim = 0; side_dim < dim; ++side_dim)
2026  {
2027  // Face is assumed to be flat, therefore normal is assumed to be
2028  // constant over the face, therefore only compute it at 1 qp.
2029  qfaces[side_dim] = std::unique_ptr<libMesh::QGauss>(new libMesh::QGauss(side_dim, CONSTANT));
2030 
2031  // A first-order Lagrange FE for the face.
2032  fe_faces[side_dim] = FEBase::build(side_dim + 1, FEType(FIRST, libMesh::LAGRANGE));
2033  fe_faces[side_dim]->attach_quadrature_rule(qfaces[side_dim].get());
2034  }
2035 
2036  // We need this to get boundary ids for each boundary face we encounter.
2037  BoundaryInfo & boundary_info = getMesh().get_boundary_info();
2038  std::vector<boundary_id_type> face_ids;
2039 
2040  for (auto & elem : as_range(getMesh().level_elements_begin(0), getMesh().level_elements_end(0)))
2041  {
2042  // dimension of the current element and its normals
2043  unsigned int side_dim = elem->dim() - 1;
2044  const std::vector<Point> & normals = fe_faces[side_dim]->get_normals();
2045 
2046  // loop over element sides
2047  for (unsigned int s = 0; s < elem->n_sides(); s++)
2048  {
2049  // If side is on the boundary
2050  if (elem->neighbor_ptr(s) == nullptr)
2051  {
2052  std::unique_ptr<Elem> side = elem->build_side_ptr(s);
2053 
2054  fe_faces[side_dim]->reinit(elem, s);
2055 
2056  // Get the boundary ID(s) for this side. If there is more
2057  // than 1 boundary id, then we already can't determine a
2058  // unique pairing of sides in this direction, but we'll just
2059  // keep going to keep the logic simple.
2060  boundary_info.boundary_ids(elem, s, face_ids);
2061 
2062  // x-direction faces
2063  if (normals[0].absolute_fuzzy_equals(minus_x))
2064  minus_x_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2065  else if (normals[0].absolute_fuzzy_equals(plus_x))
2066  plus_x_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2067 
2068  // y-direction faces
2069  else if (normals[0].absolute_fuzzy_equals(minus_y))
2070  minus_y_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2071  else if (normals[0].absolute_fuzzy_equals(plus_y))
2072  plus_y_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2073 
2074  // z-direction faces
2075  else if (normals[0].absolute_fuzzy_equals(minus_z))
2076  minus_z_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2077  else if (normals[0].absolute_fuzzy_equals(plus_z))
2078  plus_z_ids[side_dim].insert(face_ids.begin(), face_ids.end());
2079  }
2080  }
2081  }
2082 
2083  // For a distributed mesh, boundaries may be distributed as well. We therefore collect information
2084  // from everyone. If the mesh is already serial, then there is no need to do an allgather. Note
2085  // that this is just going to gather information about what the periodic bc ids are. We are not
2086  // gathering any remote elements or anything like that. It's just that the GhostPointNeighbors
2087  // ghosting functor currently relies on the fact that every process agrees on whether we have
2088  // periodic boundaries; every process that thinks there are periodic boundaries will call
2089  // MeshBase::sub_point_locator which makes a parallel_object_only() assertion (right or wrong). So
2090  // we all need to go there (or not go there)
2091  if (_use_distributed_mesh && !_mesh->is_serial())
2092  {
2093  // Pack all data together so that we send them via one communication
2094  // pair: boundary side --> boundary ids.
2095  std::vector<std::pair<boundary_id_type, boundary_id_type>> vecdata;
2096  // We check boundaries on all dimensions
2097  for (unsigned side_dim = 0; side_dim < dim; ++side_dim)
2098  {
2099  // "6" means: we have at most 6 boundaries. It is true for generated simple mesh
2100  // "detectPairedSidesets" is designed for only simple meshes
2101  for (auto bd = minus_x_ids[side_dim].begin(); bd != minus_x_ids[side_dim].end(); bd++)
2102  vecdata.emplace_back(side_dim * 6 + 0, *bd);
2103 
2104  for (auto bd = plus_x_ids[side_dim].begin(); bd != plus_x_ids[side_dim].end(); bd++)
2105  vecdata.emplace_back(side_dim * 6 + 1, *bd);
2106 
2107  for (auto bd = minus_y_ids[side_dim].begin(); bd != minus_y_ids[side_dim].end(); bd++)
2108  vecdata.emplace_back(side_dim * 6 + 2, *bd);
2109 
2110  for (auto bd = plus_y_ids[side_dim].begin(); bd != plus_y_ids[side_dim].end(); bd++)
2111  vecdata.emplace_back(side_dim * 6 + 3, *bd);
2112 
2113  for (auto bd = minus_z_ids[side_dim].begin(); bd != minus_z_ids[side_dim].end(); bd++)
2114  vecdata.emplace_back(side_dim * 6 + 4, *bd);
2115 
2116  for (auto bd = plus_z_ids[side_dim].begin(); bd != plus_z_ids[side_dim].end(); bd++)
2117  vecdata.emplace_back(side_dim * 6 + 5, *bd);
2118  }
2119 
2120  _communicator.allgather(vecdata, false);
2121 
2122  // Unpack data, and add them into minus/plus_x/y_ids
2123  for (auto pair = vecdata.begin(); pair != vecdata.end(); pair++)
2124  {
2125  // Convert data from the long vector, and add data to separated sets
2126  auto side_dim = pair->first / 6;
2127  auto side = pair->first % 6;
2128 
2129  switch (side)
2130  {
2131  case 0:
2132  minus_x_ids[side_dim].insert(pair->second);
2133  break;
2134  case 1:
2135  plus_x_ids[side_dim].insert(pair->second);
2136  break;
2137  case 2:
2138  minus_y_ids[side_dim].insert(pair->second);
2139  break;
2140  case 3:
2141  plus_y_ids[side_dim].insert(pair->second);
2142  break;
2143  case 4:
2144  minus_z_ids[side_dim].insert(pair->second);
2145  break;
2146  case 5:
2147  plus_z_ids[side_dim].insert(pair->second);
2148  break;
2149  default:
2150  mooseError("Unknown boundary side ", side);
2151  }
2152  }
2153 
2154  } // end if (_use_distributed_mesh && !_need_ghost_ghosted_boundaries)
2155 
2156  for (unsigned side_dim = 0; side_dim < dim; ++side_dim)
2157  {
2158  // If unique pairings were found, fill up the _paired_boundary data
2159  // structure with that information.
2160  if (minus_x_ids[side_dim].size() == 1 && plus_x_ids[side_dim].size() == 1)
2161  _paired_boundary.emplace_back(
2162  std::make_pair(*(minus_x_ids[side_dim].begin()), *(plus_x_ids[side_dim].begin())));
2163  else
2165  "For side dimension " + std::to_string(side_dim) +
2166  " we did not find paired boundaries (sidesets) in X due to the presence of " +
2167  std::to_string(minus_x_ids[side_dim].size()) + " -X normal and " +
2168  std::to_string(plus_x_ids[side_dim].size()) + " +X normal boundaries.");
2169 
2170  if (minus_y_ids[side_dim].size() == 1 && plus_y_ids[side_dim].size() == 1)
2171  _paired_boundary.emplace_back(
2172  std::make_pair(*(minus_y_ids[side_dim].begin()), *(plus_y_ids[side_dim].begin())));
2173  else
2175  "For side dimension " + std::to_string(side_dim) +
2176  " we did not find paired boundaries (sidesets) in Y due to the presence of " +
2177  std::to_string(minus_y_ids[side_dim].size()) + " -Y normal and " +
2178  std::to_string(plus_y_ids[side_dim].size()) + " +Y normal boundaries.");
2179 
2180  if (minus_z_ids[side_dim].size() == 1 && plus_z_ids[side_dim].size() == 1)
2181  _paired_boundary.emplace_back(
2182  std::make_pair(*(minus_z_ids[side_dim].begin()), *(plus_z_ids[side_dim].begin())));
2183  else
2185  "For side dimension " + std::to_string(side_dim) +
2186  " we did not find paired boundaries (sidesets) in Z due to the presence of " +
2187  std::to_string(minus_z_ids[side_dim].size()) + " -Z normal and " +
2188  std::to_string(plus_z_ids[side_dim].size()) + " +Z normal boundaries.");
2189  }
2190 }
2191 
2192 Real
2193 MooseMesh::dimensionWidth(unsigned int component) const
2194 {
2195  return getMaxInDimension(component) - getMinInDimension(component);
2196 }
2197 
2198 Real
2199 MooseMesh::getMinInDimension(unsigned int component) const
2200 {
2201  mooseAssert(_mesh, "The MeshBase has not been constructed");
2202  mooseAssert(component < _bounds.size(), "Requested dimension out of bounds");
2203 
2204  return _bounds[component][MIN];
2205 }
2206 
2207 Real
2208 MooseMesh::getMaxInDimension(unsigned int component) const
2209 {
2210  mooseAssert(_mesh, "The MeshBase has not been constructed");
2211  mooseAssert(component < _bounds.size(), "Requested dimension out of bounds");
2212 
2213  return _bounds[component][MAX];
2214 }
2215 
2216 void
2217 MooseMesh::addPeriodicVariable(unsigned int var_num, BoundaryID primary, BoundaryID secondary)
2218 {
2220  return;
2221 
2222  _periodic_dim[var_num].resize(dimension());
2223 
2224  _half_range = Point(dimensionWidth(0) / 2.0, dimensionWidth(1) / 2.0, dimensionWidth(2) / 2.0);
2225 
2226  bool component_found = false;
2227  for (unsigned int component = 0; component < dimension(); ++component)
2228  {
2229  const std::pair<BoundaryID, BoundaryID> * boundary_ids = getPairedBoundaryMapping(component);
2230 
2231  if (boundary_ids != nullptr &&
2232  ((boundary_ids->first == primary && boundary_ids->second == secondary) ||
2233  (boundary_ids->first == secondary && boundary_ids->second == primary)))
2234  {
2235  _periodic_dim[var_num][component] = true;
2236  component_found = true;
2237  }
2238  }
2239  if (!component_found)
2240  mooseWarning("Could not find a match between boundary '",
2241  getBoundaryName(primary),
2242  "' and '",
2243  getBoundaryName(secondary),
2244  "' to set periodic boundary conditions for variable (index:",
2245  var_num,
2246  ") in either the X, Y or Z direction. The periodic dimension of the mesh for this "
2247  "variable will not be stored.");
2248 }
2249 
2250 bool
2251 MooseMesh::isTranslatedPeriodic(unsigned int nonlinear_var_num, unsigned int component) const
2252 {
2253  mooseAssert(component < dimension(), "Requested dimension out of bounds");
2254 
2255  if (_periodic_dim.find(nonlinear_var_num) != _periodic_dim.end())
2256  return _periodic_dim.at(nonlinear_var_num)[component];
2257  else
2258  return false;
2259 }
2260 
2262 MooseMesh::minPeriodicVector(unsigned int nonlinear_var_num, Point p, Point q) const
2263 {
2264  for (unsigned int i = 0; i < dimension(); ++i)
2265  {
2266  // check to see if we're closer in real or periodic space in x, y, and z
2267  if (isTranslatedPeriodic(nonlinear_var_num, i))
2268  {
2269  // Need to test order before differencing
2270  if (p(i) > q(i))
2271  {
2272  if (p(i) - q(i) > _half_range(i))
2273  p(i) -= _half_range(i) * 2;
2274  }
2275  else
2276  {
2277  if (q(i) - p(i) > _half_range(i))
2278  p(i) += _half_range(i) * 2;
2279  }
2280  }
2281  }
2282 
2283  return q - p;
2284 }
2285 
2286 Real
2287 MooseMesh::minPeriodicDistance(unsigned int nonlinear_var_num, Point p, Point q) const
2288 {
2289  return minPeriodicVector(nonlinear_var_num, p, q).norm();
2290 }
2291 
2292 const std::pair<BoundaryID, BoundaryID> *
2294 {
2296  mooseError("Trying to retrieve automatic paired mapping for a mesh that is not regular and "
2297  "orthogonal");
2298 
2299  mooseAssert(component < dimension(), "Requested dimension out of bounds");
2300 
2301  if (_paired_boundary.empty())
2303 
2304  if (component < _paired_boundary.size())
2305  return &_paired_boundary[component];
2306  else
2307  return nullptr;
2308 }
2309 
2310 void
2312 {
2313  std::map<ElemType, Elem *> canonical_elems;
2314 
2315  // First, loop over all elements and find a canonical element for each type
2316  // Doing it this way guarantees that this is going to work in parallel
2317  for (const auto & elem : getMesh().element_ptr_range()) // TODO: Thread this
2318  {
2319  ElemType type = elem->type();
2320 
2321  if (canonical_elems.find(type) ==
2322  canonical_elems.end()) // If we haven't seen this type of elem before save it
2323  canonical_elems[type] = elem;
2324  else
2325  {
2326  Elem * stored = canonical_elems[type];
2327  if (elem->id() < stored->id()) // Arbitrarily keep the one with a lower id
2328  canonical_elems[type] = elem;
2329  }
2330  }
2331  // Now build the maps using these templates
2332  // Note: This MUST be done NOT threaded!
2333  for (const auto & can_it : canonical_elems)
2334  {
2335  Elem * elem = can_it.second;
2336 
2337  // Need to do this just once to get the right qrules put in place
2338  assembly->setCurrentSubdomainID(elem->subdomain_id());
2339  assembly->reinit(elem);
2340  assembly->reinit(elem, 0);
2341  auto && qrule = assembly->writeableQRule();
2342  auto && qrule_face = assembly->writeableQRuleFace();
2343 
2344  // Volume to volume projection for refinement
2345  buildRefinementMap(*elem, *qrule, *qrule_face, -1, -1, -1);
2346 
2347  // Volume to volume projection for coarsening
2348  buildCoarseningMap(*elem, *qrule, *qrule_face, -1);
2349 
2350  // Map the sides of children
2351  for (unsigned int side = 0; side < elem->n_sides(); side++)
2352  {
2353  // Side to side for sides that match parent's sides
2354  buildRefinementMap(*elem, *qrule, *qrule_face, side, -1, side);
2355  buildCoarseningMap(*elem, *qrule, *qrule_face, side);
2356  }
2357 
2358  // Child side to parent volume mapping for "internal" child sides
2359  for (unsigned int child = 0; child < elem->n_children(); ++child)
2360  for (unsigned int side = 0; side < elem->n_sides();
2361  ++side) // Assume children have the same number of sides!
2362  if (!elem->is_child_on_side(child, side)) // Otherwise we already computed that map
2363  buildRefinementMap(*elem, *qrule, *qrule_face, -1, child, side);
2364  }
2365 }
2366 
2367 void
2369 {
2374 
2375  std::map<ElemType, std::pair<Elem *, unsigned int>> elems_and_max_p_level;
2376 
2377  for (const auto & elem : getMesh().active_element_ptr_range())
2378  {
2379  const auto type = elem->type();
2380  auto & [picked_elem, max_p_level] = elems_and_max_p_level[type];
2381  if (!picked_elem)
2382  picked_elem = elem;
2383  max_p_level = std::max(max_p_level, elem->p_level());
2384  }
2385 
2386  // The only requirement on the FEType is that it can be arbitrarily p-refined
2387  const FEType p_refinable_fe_type(CONSTANT, libMesh::MONOMIAL);
2388  std::vector<Point> volume_ref_points_coarse, volume_ref_points_fine, face_ref_points_coarse,
2389  face_ref_points_fine;
2390  std::vector<unsigned int> p_levels;
2391 
2392  for (auto & [elem_type, elem_p_level_pair] : elems_and_max_p_level)
2393  {
2394  auto & [moose_elem, max_p_level] = elem_p_level_pair;
2395  const auto dim = moose_elem->dim();
2396  // Need to do this just once to get the right qrules put in place
2397  assembly->setCurrentSubdomainID(moose_elem->subdomain_id());
2398  assembly->reinit(moose_elem);
2399  assembly->reinit(moose_elem, 0);
2400  auto & qrule = assembly->writeableQRule();
2401  auto & qrule_face = assembly->writeableQRuleFace();
2402 
2403  libMesh::Parallel::Communicator self_comm{};
2404  ReplicatedMesh mesh(self_comm);
2406  for (const auto & nd : moose_elem->node_ref_range())
2407  mesh.add_point(nd);
2408 
2409  Elem * const elem = mesh.add_elem(Elem::build(elem_type).release());
2410  for (const auto i : elem->node_index_range())
2411  elem->set_node(i, mesh.node_ptr(i));
2412 
2413  std::unique_ptr<FEBase> fe_face(FEBase::build(dim, p_refinable_fe_type));
2414  fe_face->get_phi();
2415  const auto & face_phys_points = fe_face->get_xyz();
2416  fe_face->attach_quadrature_rule(qrule_face);
2417 
2418  qrule->init(*elem);
2419  volume_ref_points_coarse = qrule->get_points();
2420  fe_face->reinit(elem, (unsigned int)0);
2421  libMesh::FEMap::inverse_map(dim, elem, face_phys_points, face_ref_points_coarse);
2422 
2423  p_levels.resize(max_p_level + 1);
2424  std::iota(p_levels.begin(), p_levels.end(), 0);
2425  libMesh::MeshRefinement mesh_refinement(mesh);
2426 
2427  for (const auto p_level : p_levels)
2428  {
2429  mesh_refinement.uniformly_p_refine(1);
2430  qrule->init(*elem);
2431  volume_ref_points_fine = qrule->get_points();
2432  fe_face->reinit(elem, (unsigned int)0);
2433  libMesh::FEMap::inverse_map(dim, elem, face_phys_points, face_ref_points_fine);
2434 
2435  const auto map_key = std::make_pair(elem_type, p_level);
2436  auto & volume_refine_map = _elem_type_to_p_refinement_map[map_key];
2437  auto & face_refine_map = _elem_type_to_p_refinement_side_map[map_key];
2438  auto & volume_coarsen_map = _elem_type_to_p_coarsening_map[map_key];
2439  auto & face_coarsen_map = _elem_type_to_p_coarsening_side_map[map_key];
2440 
2441  auto fill_maps = [this](const auto & coarse_ref_points,
2442  const auto & fine_ref_points,
2443  auto & coarsen_map,
2444  auto & refine_map)
2445  {
2446  mapPoints(fine_ref_points, coarse_ref_points, refine_map);
2447  mapPoints(coarse_ref_points, fine_ref_points, coarsen_map);
2448  };
2449 
2450  fill_maps(
2451  volume_ref_points_coarse, volume_ref_points_fine, volume_coarsen_map, volume_refine_map);
2452  fill_maps(face_ref_points_coarse, face_ref_points_fine, face_coarsen_map, face_refine_map);
2453 
2454  // With this level's maps filled our fine points now become our coarse points
2455  volume_ref_points_fine.swap(volume_ref_points_coarse);
2456  face_ref_points_fine.swap(face_ref_points_coarse);
2457  }
2458  }
2459 }
2460 
2461 void
2463 {
2464  TIME_SECTION("buildRefinementAndCoarseningMaps", 5, "Building Refinement And Coarsening Maps");
2465  if (doingPRefinement())
2467  else
2469 }
2470 
2471 void
2473  QBase & qrule,
2474  QBase & qrule_face,
2475  int parent_side,
2476  int child,
2477  int child_side)
2478 {
2479  TIME_SECTION("buildRefinementMap", 5, "Building Refinement Map");
2480 
2481  if (child == -1) // Doing volume mapping or parent side mapping
2482  {
2483  mooseAssert(parent_side == child_side,
2484  "Parent side must match child_side if not passing a specific child!");
2485 
2486  std::pair<int, ElemType> the_pair(parent_side, elem.type());
2487 
2488  if (_elem_type_to_refinement_map.find(the_pair) != _elem_type_to_refinement_map.end())
2489  mooseError("Already built a qp refinement map!");
2490 
2491  std::vector<std::pair<unsigned int, QpMap>> coarsen_map;
2492  std::vector<std::vector<QpMap>> & refinement_map = _elem_type_to_refinement_map[the_pair];
2494  &elem, qrule, qrule_face, refinement_map, coarsen_map, parent_side, child, child_side);
2495  }
2496  else // Need to map a child side to parent volume qps
2497  {
2498  std::pair<int, int> child_pair(child, child_side);
2499 
2502  _elem_type_to_child_side_refinement_map[elem.type()].find(child_pair) !=
2504  mooseError("Already built a qp refinement map!");
2505 
2506  std::vector<std::pair<unsigned int, QpMap>> coarsen_map;
2507  std::vector<std::vector<QpMap>> & refinement_map =
2510  &elem, qrule, qrule_face, refinement_map, coarsen_map, parent_side, child, child_side);
2511  }
2512 }
2513 
2514 const std::vector<std::vector<QpMap>> &
2515 MooseMesh::getRefinementMap(const Elem & elem, int parent_side, int child, int child_side)
2516 {
2517  if (child == -1) // Doing volume mapping or parent side mapping
2518  {
2519  mooseAssert(parent_side == child_side,
2520  "Parent side must match child_side if not passing a specific child!");
2521 
2522  std::pair<int, ElemType> the_pair(parent_side, elem.type());
2523 
2524  if (_elem_type_to_refinement_map.find(the_pair) == _elem_type_to_refinement_map.end())
2525  mooseError("Could not find a suitable qp refinement map!");
2526 
2527  return _elem_type_to_refinement_map[the_pair];
2528  }
2529  else // Need to map a child side to parent volume qps
2530  {
2531  std::pair<int, int> child_pair(child, child_side);
2532 
2535  _elem_type_to_child_side_refinement_map[elem.type()].find(child_pair) ==
2537  mooseError("Could not find a suitable qp refinement map!");
2538 
2539  return _elem_type_to_child_side_refinement_map[elem.type()][child_pair];
2540  }
2541 
2548 }
2549 
2550 void
2551 MooseMesh::buildCoarseningMap(const Elem & elem, QBase & qrule, QBase & qrule_face, int input_side)
2552 {
2553  TIME_SECTION("buildCoarseningMap", 5, "Building Coarsening Map");
2554 
2555  std::pair<int, ElemType> the_pair(input_side, elem.type());
2556 
2557  if (_elem_type_to_coarsening_map.find(the_pair) != _elem_type_to_coarsening_map.end())
2558  mooseError("Already built a qp coarsening map!");
2559 
2560  std::vector<std::vector<QpMap>> refinement_map;
2561  std::vector<std::pair<unsigned int, QpMap>> & coarsen_map =
2562  _elem_type_to_coarsening_map[the_pair];
2563 
2564  // The -1 here is for a specific child. We don't do that for coarsening maps
2565  // Also note that we're always mapping the same side to the same side (which is guaranteed by
2566  // libMesh).
2568  &elem, qrule, qrule_face, refinement_map, coarsen_map, input_side, -1, input_side);
2569 
2576 }
2577 
2578 const std::vector<std::pair<unsigned int, QpMap>> &
2579 MooseMesh::getCoarseningMap(const Elem & elem, int input_side)
2580 {
2581  std::pair<int, ElemType> the_pair(input_side, elem.type());
2582 
2583  if (_elem_type_to_coarsening_map.find(the_pair) == _elem_type_to_coarsening_map.end())
2584  mooseError("Could not find a suitable qp refinement map!");
2585 
2586  return _elem_type_to_coarsening_map[the_pair];
2587 }
2588 
2589 void
2590 MooseMesh::mapPoints(const std::vector<Point> & from,
2591  const std::vector<Point> & to,
2592  std::vector<QpMap> & qp_map)
2593 {
2594  unsigned int n_from = from.size();
2595  unsigned int n_to = to.size();
2596 
2597  qp_map.resize(n_from);
2598 
2599  for (unsigned int i = 0; i < n_from; ++i)
2600  {
2601  const Point & from_point = from[i];
2602 
2603  QpMap & current_map = qp_map[i];
2604 
2605  for (unsigned int j = 0; j < n_to; ++j)
2606  {
2607  const Point & to_point = to[j];
2608  Real distance = (from_point - to_point).norm();
2609 
2610  if (distance < current_map._distance)
2611  {
2612  current_map._distance = distance;
2613  current_map._from = i;
2614  current_map._to = j;
2615  }
2616  }
2617  }
2618 }
2619 
2620 void
2622  QBase & qrule,
2623  QBase & qrule_face,
2624  std::vector<std::vector<QpMap>> & refinement_map,
2625  std::vector<std::pair<unsigned int, QpMap>> & coarsen_map,
2626  int parent_side,
2627  int child,
2628  int child_side)
2629 {
2630  TIME_SECTION("findAdaptivityQpMaps", 5);
2631 
2633  mesh.skip_partitioning(true);
2634 
2635  unsigned int dim = template_elem->dim();
2637 
2638  for (unsigned int i = 0; i < template_elem->n_nodes(); ++i)
2639  mesh.add_point(template_elem->point(i));
2640 
2641  Elem * elem = mesh.add_elem(Elem::build(template_elem->type()).release());
2642 
2643  for (unsigned int i = 0; i < template_elem->n_nodes(); ++i)
2644  elem->set_node(i, mesh.node_ptr(i));
2645 
2646  std::unique_ptr<FEBase> fe(FEBase::build(dim, FEType()));
2647  fe->get_phi();
2648  const std::vector<Point> & q_points_volume = fe->get_xyz();
2649 
2650  std::unique_ptr<FEBase> fe_face(FEBase::build(dim, FEType()));
2651  fe_face->get_phi();
2652  const std::vector<Point> & q_points_face = fe_face->get_xyz();
2653 
2654  fe->attach_quadrature_rule(&qrule);
2655  fe_face->attach_quadrature_rule(&qrule_face);
2656 
2657  // The current q_points (locations in *physical* space)
2658  const std::vector<Point> * q_points;
2659 
2660  if (parent_side != -1)
2661  {
2662  fe_face->reinit(elem, parent_side);
2663  q_points = &q_points_face;
2664  }
2665  else
2666  {
2667  fe->reinit(elem);
2668  q_points = &q_points_volume;
2669  }
2670 
2671  std::vector<Point> parent_ref_points;
2672 
2673  libMesh::FEMap::inverse_map(elem->dim(), elem, *q_points, parent_ref_points);
2674  libMesh::MeshRefinement mesh_refinement(mesh);
2675  mesh_refinement.uniformly_refine(1);
2676 
2677  // A map from the child element index to the locations of all the child's quadrature points in
2678  // *reference* space. Note that we use a map here instead of a vector because the caller can
2679  // pass an explicit child index. We are not guaranteed to have a sequence from [0, n_children)
2680  std::map<unsigned int, std::vector<Point>> child_to_ref_points;
2681 
2682  unsigned int n_children = elem->n_children();
2683 
2684  refinement_map.resize(n_children);
2685 
2686  std::vector<unsigned int> children;
2687 
2688  if (child != -1) // Passed in a child explicitly
2689  children.push_back(child);
2690  else
2691  {
2692  children.resize(n_children);
2693  for (unsigned int child = 0; child < n_children; ++child)
2694  children[child] = child;
2695  }
2696 
2697  for (unsigned int i = 0; i < children.size(); ++i)
2698  {
2699  unsigned int child = children[i];
2700 
2701  if ((parent_side != -1 && !elem->is_child_on_side(child, parent_side)))
2702  continue;
2703 
2704  const Elem * child_elem = elem->child_ptr(child);
2705 
2706  if (child_side != -1)
2707  {
2708  fe_face->reinit(child_elem, child_side);
2709  q_points = &q_points_face;
2710  }
2711  else
2712  {
2713  fe->reinit(child_elem);
2714  q_points = &q_points_volume;
2715  }
2716 
2717  std::vector<Point> child_ref_points;
2718 
2719  libMesh::FEMap::inverse_map(elem->dim(), elem, *q_points, child_ref_points);
2720  child_to_ref_points[child] = child_ref_points;
2721 
2722  std::vector<QpMap> & qp_map = refinement_map[child];
2723 
2724  // Find the closest parent_qp to each child_qp
2725  mapPoints(child_ref_points, parent_ref_points, qp_map);
2726  }
2727 
2728  coarsen_map.resize(parent_ref_points.size());
2729 
2730  // For each parent qp find the closest child qp
2731  for (unsigned int child = 0; child < n_children; child++)
2732  {
2733  if (parent_side != -1 && !elem->is_child_on_side(child, child_side))
2734  continue;
2735 
2736  std::vector<Point> & child_ref_points = child_to_ref_points[child];
2737 
2738  std::vector<QpMap> qp_map;
2739 
2740  // Find all of the closest points from parent_qp to _THIS_ child's qp
2741  mapPoints(parent_ref_points, child_ref_points, qp_map);
2742 
2743  // Check those to see if they are closer than what we currently have for each point
2744  for (unsigned int parent_qp = 0; parent_qp < parent_ref_points.size(); ++parent_qp)
2745  {
2746  std::pair<unsigned int, QpMap> & child_and_map = coarsen_map[parent_qp];
2747  unsigned int & closest_child = child_and_map.first;
2748  QpMap & closest_map = child_and_map.second;
2749 
2750  QpMap & current_map = qp_map[parent_qp];
2751 
2752  if (current_map._distance < closest_map._distance)
2753  {
2754  closest_child = child;
2755  closest_map = current_map;
2756  }
2757  }
2758  }
2759 }
2760 
2761 void
2763  const boundary_id_type new_id,
2764  bool delete_prev)
2765 {
2766  TIME_SECTION("changeBoundaryId", 6);
2767  changeBoundaryId(getMesh(), old_id, new_id, delete_prev);
2768 }
2769 
2770 void
2772  const boundary_id_type old_id,
2773  const boundary_id_type new_id,
2774  bool delete_prev)
2775 {
2776  // Get a reference to our BoundaryInfo object, we will use it several times below...
2777  BoundaryInfo & boundary_info = mesh.get_boundary_info();
2778 
2779  // Container to catch ids passed back from BoundaryInfo
2780  std::vector<boundary_id_type> old_ids;
2781 
2782  // Only level-0 elements store BCs. Loop over them.
2783  for (auto & elem : as_range(mesh.level_elements_begin(0), mesh.level_elements_end(0)))
2784  {
2785  unsigned int n_sides = elem->n_sides();
2786  for (unsigned int s = 0; s != n_sides; ++s)
2787  {
2788  boundary_info.boundary_ids(elem, s, old_ids);
2789  if (std::find(old_ids.begin(), old_ids.end(), old_id) != old_ids.end())
2790  {
2791  std::vector<boundary_id_type> new_ids(old_ids);
2792  std::replace(new_ids.begin(), new_ids.end(), old_id, new_id);
2793  if (delete_prev)
2794  {
2795  boundary_info.remove_side(elem, s);
2796  boundary_info.add_side(elem, s, new_ids);
2797  }
2798  else
2799  boundary_info.add_side(elem, s, new_ids);
2800  }
2801  }
2802  }
2803 
2804  // Remove any remaining references to the old ID from the
2805  // BoundaryInfo object. This prevents things like empty sidesets
2806  // from showing up when printing information, etc.
2807  if (delete_prev)
2808  boundary_info.remove_id(old_id);
2809 }
2810 
2811 const RealVectorValue &
2813 {
2814  mooseAssert(_boundary_to_normal_map.get() != nullptr, "Boundary To Normal Map not built!");
2815 
2816  // Note: Boundaries that are not in the map (existing boundaries) will default
2817  // construct a new RealVectorValue - (x,y,z)=(0, 0, 0)
2818  return (*_boundary_to_normal_map)[id];
2819 }
2820 
2821 MooseMesh &
2823 {
2824  mooseError("MooseMesh::clone() is no longer supported, use MooseMesh::safeClone() instead.");
2825 }
2826 
2827 void
2829 {
2830  switch (_parallel_type)
2831  {
2832  case ParallelType::DEFAULT:
2833  // The user did not specify 'parallel_type = XYZ' in the input file,
2834  // so we allow the --distributed-mesh command line arg to possibly turn
2835  // on DistributedMesh. If the command line arg is not present, we pick ReplicatedMesh.
2837  _use_distributed_mesh = true;
2838  break;
2842  _use_distributed_mesh = false;
2843  break;
2845  _use_distributed_mesh = true;
2846  break;
2847  }
2848 
2849  // If the user specifies 'nemesis = true' in the Mesh block, or they are using --use-split,
2850  // we must use DistributedMesh.
2851  if (_is_nemesis || _is_split)
2852  _use_distributed_mesh = true;
2853 }
2854 
2855 std::unique_ptr<MeshBase>
2857 {
2858  std::unique_ptr<MeshBase> mesh;
2860  mesh = buildTypedMesh<DistributedMesh>(dim);
2861  else
2862  mesh = buildTypedMesh<ReplicatedMesh>(dim);
2863 
2864  return mesh;
2865 }
2866 
2867 void
2868 MooseMesh::setMeshBase(std::unique_ptr<MeshBase> mesh_base)
2869 {
2870  _mesh = std::move(mesh_base);
2871  _mesh->allow_remote_element_removal(_allow_remote_element_removal);
2872 }
2873 
2874 void
2876 {
2883  if (!_mesh)
2885 
2887  mooseError("You cannot use the mesh splitter capability with DistributedMesh!");
2888 
2889  TIME_SECTION("init", 2);
2890 
2892  {
2893  // Some partitioners are not idempotent. Some recovery data
2894  // files require partitioning to match mesh partitioning. This
2895  // means that, when recovering, we can't safely repartition.
2896  const bool skip_partitioning_later = getMesh().skip_partitioning();
2897  getMesh().skip_partitioning(true);
2898  const bool allow_renumbering_later = getMesh().allow_renumbering();
2899  getMesh().allow_renumbering(false);
2900 
2901  // For now, only read the recovery mesh on the Ultimate Master..
2902  // sub-apps need to just build their mesh like normal
2903  {
2904  TIME_SECTION("readRecoveredMesh", 2);
2906  }
2907 
2908  getMesh().allow_renumbering(allow_renumbering_later);
2909  getMesh().skip_partitioning(skip_partitioning_later);
2910  }
2911  else // Normally just build the mesh
2912  {
2913  // Don't allow partitioning during building
2914  if (_app.isSplitMesh())
2915  getMesh().skip_partitioning(true);
2916  buildMesh();
2917 
2918  // Re-enable partitioning so the splitter can partition!
2919  if (_app.isSplitMesh())
2920  getMesh().skip_partitioning(false);
2921 
2922  if (getParam<bool>("build_all_side_lowerd_mesh"))
2923  buildLowerDMesh();
2924  }
2925 }
2926 
2927 unsigned int
2929 {
2930  return getMesh().mesh_dimension();
2931 }
2932 
2933 unsigned int
2935 {
2936  const Real abs_zero = 1e-12;
2937 
2938  // See if the mesh is completely containd in the z and y planes to calculate effective spatial
2939  // dim
2940  for (unsigned int dim = LIBMESH_DIM; dim >= 1; --dim)
2941  if (dimensionWidth(dim - 1) >= abs_zero)
2942  return dim;
2943 
2944  // If we get here, we have a 1D mesh on the x-axis.
2945  return 1;
2946 }
2947 
2948 unsigned int
2949 MooseMesh::getBlocksMaxDimension(const std::vector<SubdomainName> & blocks) const
2950 {
2951  const auto & mesh = getMesh();
2952 
2953  // Take a shortcut if possible
2954  if (const auto & elem_dims = mesh.elem_dimensions(); mesh.is_prepared() && elem_dims.size() == 1)
2955  return *elem_dims.begin();
2956 
2957  unsigned short dim = 0;
2958  const auto subdomain_ids = getSubdomainIDs(blocks);
2959  const std::set<SubdomainID> subdomain_ids_set(subdomain_ids.begin(), subdomain_ids.end());
2960  for (const auto & elem : mesh.active_subdomain_set_elements_ptr_range(subdomain_ids_set))
2961  dim = std::max(dim, elem->dim());
2962 
2963  // Get the maximumal globally
2965  return dim;
2966 }
2967 
2968 std::vector<BoundaryID>
2969 MooseMesh::getBoundaryIDs(const Elem * const elem, const unsigned short int side) const
2970 {
2971  std::vector<BoundaryID> ids;
2972  getMesh().get_boundary_info().boundary_ids(elem, side, ids);
2973  return ids;
2974 }
2975 
2976 const std::set<BoundaryID> &
2978 {
2980 }
2981 
2982 void
2984 {
2987 }
2988 
2989 void
2990 MooseMesh::buildSideList(std::vector<dof_id_type> & el,
2991  std::vector<unsigned short int> & sl,
2992  std::vector<boundary_id_type> & il)
2993 {
2994 #ifdef LIBMESH_ENABLE_DEPRECATED
2995  mooseDeprecated("The version of MooseMesh::buildSideList() taking three arguments is "
2996  "deprecated, call the version that returns a vector of tuples instead.");
2997  getMesh().get_boundary_info().build_side_list(el, sl, il);
2998 #else
2999  libmesh_ignore(el);
3000  libmesh_ignore(sl);
3001  libmesh_ignore(il);
3002  mooseError("The version of MooseMesh::buildSideList() taking three "
3003  "arguments is not available in your version of libmesh, call the "
3004  "version that returns a vector of tuples instead.");
3005 #endif
3006 }
3007 
3008 std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>>
3010 {
3012 }
3013 
3014 std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>>
3016 {
3018 }
3019 
3020 unsigned int
3021 MooseMesh::sideWithBoundaryID(const Elem * const elem, const BoundaryID boundary_id) const
3022 {
3023  return getMesh().get_boundary_info().side_with_boundary_id(elem, boundary_id);
3024 }
3025 
3028 {
3029  return getMesh().local_nodes_begin();
3030 }
3031 
3034 {
3035  return getMesh().local_nodes_end();
3036 }
3037 
3038 MeshBase::const_node_iterator
3040 {
3041  return getMesh().local_nodes_begin();
3042 }
3043 
3044 MeshBase::const_node_iterator
3046 {
3047  return getMesh().local_nodes_end();
3048 }
3049 
3050 MeshBase::element_iterator
3052 {
3053  return getMesh().active_local_elements_begin();
3054 }
3055 
3056 const MeshBase::element_iterator
3058 {
3059  return getMesh().active_local_elements_end();
3060 }
3061 
3062 MeshBase::const_element_iterator
3064 {
3065  return getMesh().active_local_elements_begin();
3066 }
3067 
3068 const MeshBase::const_element_iterator
3070 {
3071  return getMesh().active_local_elements_end();
3072 }
3073 
3076 {
3077  return getMesh().n_nodes();
3078 }
3079 
3082 {
3083  return getMesh().n_elem();
3084 }
3085 
3088 {
3089  return getMesh().max_node_id();
3090 }
3091 
3094 {
3095  return getMesh().max_elem_id();
3096 }
3097 
3098 Elem *
3100 {
3101  mooseDeprecated("MooseMesh::elem() is deprecated, please use MooseMesh::elemPtr() instead");
3102  return elemPtr(i);
3103 }
3104 
3105 const Elem *
3107 {
3108  mooseDeprecated("MooseMesh::elem() is deprecated, please use MooseMesh::elemPtr() instead");
3109  return elemPtr(i);
3110 }
3111 
3112 Elem *
3114 {
3115  return getMesh().elem_ptr(i);
3116 }
3117 
3118 const Elem *
3120 {
3121  return getMesh().elem_ptr(i);
3122 }
3123 
3124 Elem *
3126 {
3127  return getMesh().query_elem_ptr(i);
3128 }
3129 
3130 const Elem *
3132 {
3133  return getMesh().query_elem_ptr(i);
3134 }
3135 
3136 bool
3138 {
3139  return _mesh->is_prepared() && _moose_mesh_prepared;
3140 }
3141 
3142 void
3144 {
3145  if (state)
3146  mooseError("We don't have any right to tell the libmesh mesh that it *is* prepared. Only a "
3147  "call to prepare_for_use should tell us that");
3148 
3149  // Some people may call this even before we have a MeshBase object. This isn't dangerous really
3150  // because when the MeshBase object is born, it knows it's in an unprepared state
3151  if (_mesh)
3152  _mesh->set_isnt_prepared();
3153 
3154  // If the libMesh mesh isn't preparead, then our MooseMesh wrapper is also no longer prepared
3155  _moose_mesh_prepared = false;
3156 
3161  _regular_orthogonal_mesh = false;
3162 }
3163 
3164 void
3166 {
3167  prepared(false);
3168 }
3169 
3170 const std::set<SubdomainID> &
3172 {
3173  return _mesh_subdomains;
3174 }
3175 
3176 const std::set<BoundaryID> &
3178 {
3179  return _mesh_boundary_ids;
3180 }
3181 
3182 const std::set<BoundaryID> &
3184 {
3185  return _mesh_sideset_ids;
3186 }
3187 
3188 const std::set<BoundaryID> &
3190 {
3191  return _mesh_nodeset_ids;
3192 }
3193 
3194 void
3195 MooseMesh::setMeshBoundaryIDs(std::set<BoundaryID> boundary_IDs)
3196 {
3197  _mesh_boundary_ids = boundary_IDs;
3198 }
3199 
3200 void
3202  std::unique_ptr<std::map<BoundaryID, RealVectorValue>> boundary_map)
3203 {
3204  _boundary_to_normal_map = std::move(boundary_map);
3205 }
3206 
3207 void
3208 MooseMesh::setBoundaryToNormalMap(std::map<BoundaryID, RealVectorValue> * boundary_map)
3209 {
3210  mooseDeprecated("setBoundaryToNormalMap(std::map<BoundaryID, RealVectorValue> * boundary_map) is "
3211  "deprecated, use the unique_ptr version instead");
3212  _boundary_to_normal_map.reset(boundary_map);
3213 }
3214 
3215 unsigned int
3217 {
3218  return _uniform_refine_level;
3219 }
3220 
3221 void
3222 MooseMesh::setUniformRefineLevel(unsigned int level, bool deletion)
3223 {
3224  _uniform_refine_level = level;
3226 }
3227 
3228 void
3230 {
3231  _ghosted_boundaries.insert(boundary_id);
3232 }
3233 
3234 void
3235 MooseMesh::setGhostedBoundaryInflation(const std::vector<Real> & inflation)
3236 {
3237  _ghosted_boundaries_inflation = inflation;
3238 }
3239 
3240 const std::set<unsigned int> &
3242 {
3243  return _ghosted_boundaries;
3244 }
3245 
3246 const std::vector<Real> &
3248 {
3250 }
3251 
3252 namespace // Anonymous namespace for helpers
3253 {
3254 // A class for templated methods that expect output iterator
3255 // arguments, which adds objects to the Mesh.
3256 // Although any mesh_inserter_iterator can add any object, we
3257 // template it around object type so that type inference and
3258 // iterator_traits will work.
3259 // This object specifically is used to insert extra ghost elems into the mesh
3260 template <typename T>
3261 struct extra_ghost_elem_inserter
3262 {
3263  using iterator_category = std::output_iterator_tag;
3264  using value_type = T;
3265 
3266  extra_ghost_elem_inserter(DistributedMesh & m) : mesh(m) {}
3267 
3268  void operator=(const Elem * e) { mesh.add_extra_ghost_elem(const_cast<Elem *>(e)); }
3269 
3270  void operator=(Node * n) { mesh.add_node(n); }
3271 
3272  void operator=(Point * p) { mesh.add_point(*p); }
3273 
3274  extra_ghost_elem_inserter & operator++() { return *this; }
3275 
3276  extra_ghost_elem_inserter operator++(int) { return extra_ghost_elem_inserter(*this); }
3277 
3278  // We don't return a reference-to-T here because we don't want to
3279  // construct one or have any of its methods called. We just want
3280  // to allow the returned object to be able to do mesh insertions
3281  // with operator=().
3282  extra_ghost_elem_inserter & operator*() { return *this; }
3283 
3284 private:
3286 };
3287 
3298 struct CompareElemsByLevel
3299 {
3300  bool operator()(const Elem * a, const Elem * b) const
3301  {
3302  libmesh_assert(a);
3303  libmesh_assert(b);
3304  const unsigned int al = a->level(), bl = b->level();
3305  const dof_id_type aid = a->id(), bid = b->id();
3306 
3307  return (al == bl) ? aid < bid : al < bl;
3308  }
3309 };
3310 
3311 } // anonymous namespace
3312 
3313 void
3315 {
3316  // No need to do this if using a serial mesh
3317  // We do not need to ghost boundary elements when _need_ghost_ghosted_boundaries
3318  // is not true. _need_ghost_ghosted_boundaries can be set by a mesh generator
3319  // where boundaries are already ghosted accordingly
3321  return;
3322 
3323  TIME_SECTION("GhostGhostedBoundaries", 3);
3324 
3325  parallel_object_only();
3326 
3327  DistributedMesh & mesh = dynamic_cast<DistributedMesh &>(getMesh());
3328 
3329  // We clear ghosted elements that were added by previous invocations of this
3330  // method but leave ghosted elements that were added by other code, e.g.
3331  // OversampleOutput, untouched
3332  mesh.clear_extra_ghost_elems(_ghost_elems_from_ghost_boundaries);
3334 
3335  std::set<const Elem *, CompareElemsByLevel> boundary_elems_to_ghost;
3336  std::set<Node *> connected_nodes_to_ghost;
3337 
3338  std::vector<const Elem *> family_tree;
3339 
3340  for (const auto & t : mesh.get_boundary_info().build_side_list())
3341  {
3342  auto elem_id = std::get<0>(t);
3343  auto bc_id = std::get<2>(t);
3344 
3345  if (_ghosted_boundaries.find(bc_id) != _ghosted_boundaries.end())
3346  {
3347  Elem * elem = mesh.elem_ptr(elem_id);
3348 
3349 #ifdef LIBMESH_ENABLE_AMR
3351  Elem * parent = elem->parent();
3352  while (parent)
3353  {
3354  family_tree.push_back(parent);
3355  parent = parent->parent();
3356  }
3357 #else
3358  family_tree.clear();
3359  family_tree.push_back(elem);
3360 #endif
3361  for (const auto & felem : family_tree)
3362  {
3363  boundary_elems_to_ghost.insert(felem);
3364 
3365  // The entries of connected_nodes_to_ghost need to be
3366  // non-constant, so that they will work in things like
3367  // UpdateDisplacedMeshThread. The container returned by
3368  // family_tree contains const Elems even when the Elem
3369  // it is called on is non-const, so once that interface
3370  // gets fixed we can remove this const_cast.
3371  for (unsigned int n = 0; n < felem->n_nodes(); ++n)
3372  connected_nodes_to_ghost.insert(const_cast<Node *>(felem->node_ptr(n)));
3373  }
3374  }
3375  }
3376 
3377  // We really do want to store this by value instead of by reference
3378  const auto prior_ghost_elems = mesh.extra_ghost_elems();
3379 
3381  connected_nodes_to_ghost.begin(),
3382  connected_nodes_to_ghost.end(),
3383  extra_ghost_elem_inserter<Node>(mesh));
3384 
3386  boundary_elems_to_ghost.begin(),
3387  boundary_elems_to_ghost.end(),
3388  extra_ghost_elem_inserter<Elem>(mesh));
3389 
3390  const auto & current_ghost_elems = mesh.extra_ghost_elems();
3391 
3392  std::set_difference(current_ghost_elems.begin(),
3393  current_ghost_elems.end(),
3394  prior_ghost_elems.begin(),
3395  prior_ghost_elems.end(),
3396  std::inserter(_ghost_elems_from_ghost_boundaries,
3398 }
3399 
3400 unsigned int
3402 {
3403  return _patch_size;
3404 }
3405 
3406 void
3408 {
3409  _patch_update_strategy = patch_update_strategy;
3410 }
3411 
3412 const Moose::PatchUpdateType &
3414 {
3415  return _patch_update_strategy;
3416 }
3417 
3419 MooseMesh::getInflatedProcessorBoundingBox(Real inflation_multiplier) const
3420 {
3421  // Grab a bounding box to speed things up. Note that
3422  // local_bounding_box is *not* equivalent to processor_bounding_box
3423  // with processor_id() except in serial.
3424  BoundingBox bbox = MeshTools::create_local_bounding_box(getMesh());
3425 
3426  // Inflate the bbox just a bit to deal with roundoff
3427  // Adding 1% of the diagonal size in each direction on each end
3428  Real inflation_amount = inflation_multiplier * (bbox.max() - bbox.min()).norm();
3429  Point inflation(inflation_amount, inflation_amount, inflation_amount);
3430 
3431  bbox.first -= inflation; // min
3432  bbox.second += inflation; // max
3433 
3434  return bbox;
3435 }
3436 
3437 MooseMesh::operator libMesh::MeshBase &() { return getMesh(); }
3438 
3439 MooseMesh::operator const libMesh::MeshBase &() const { return getMesh(); }
3440 
3441 const MeshBase *
3443 {
3444  return _mesh.get();
3445 }
3446 
3447 MeshBase &
3449 {
3450  mooseAssert(_mesh, "Mesh hasn't been created");
3451  return *_mesh;
3452 }
3453 
3454 const MeshBase &
3456 {
3457  mooseAssert(_mesh, "Mesh hasn't been created");
3458  return *_mesh;
3459 }
3460 
3461 void
3462 MooseMesh::printInfo(std::ostream & os, const unsigned int verbosity /* = 0 */) const
3463 {
3464  os << '\n';
3465  getMesh().print_info(os, verbosity);
3466  os << std::flush;
3467 }
3468 
3469 const std::vector<dof_id_type> &
3471 {
3472  std::map<boundary_id_type, std::vector<dof_id_type>>::const_iterator it =
3473  _node_set_nodes.find(nodeset_id);
3474 
3475  if (it == _node_set_nodes.end())
3476  {
3477  // On a distributed mesh we might not know about a remote nodeset,
3478  // so we'll return an empty vector and hope the nodeset exists
3479  // elsewhere.
3480  if (!getMesh().is_serial())
3481  {
3482  static const std::vector<dof_id_type> empty_vec;
3483  return empty_vec;
3484  }
3485  // On a replicated mesh we should know about every nodeset and if
3486  // we're asked for one that doesn't exist then it must be a bug.
3487  else
3488  {
3489  mooseError("Unable to nodeset ID: ", nodeset_id, '.');
3490  }
3491  }
3492 
3493  return it->second;
3494 }
3495 
3496 const std::set<BoundaryID> &
3498 {
3499  const auto it = _sub_to_data.find(subdomain_id);
3500 
3501  if (it == _sub_to_data.end())
3502  mooseError("Unable to find subdomain ID: ", subdomain_id, '.');
3503 
3504  return it->second.boundary_ids;
3505 }
3506 
3507 std::set<BoundaryID>
3509 {
3510  const auto & bnd_ids = getSubdomainBoundaryIds(subdomain_id);
3511  std::set<BoundaryID> boundary_ids(bnd_ids.begin(), bnd_ids.end());
3512  std::unordered_map<SubdomainID, std::set<BoundaryID>>::const_iterator it =
3513  _neighbor_subdomain_boundary_ids.find(subdomain_id);
3514 
3515  boundary_ids.insert(it->second.begin(), it->second.end());
3516 
3517  return boundary_ids;
3518 }
3519 
3520 std::set<SubdomainID>
3522 {
3523  std::set<SubdomainID> subdomain_ids;
3524  for (const auto & [sub_id, data] : _sub_to_data)
3525  if (data.boundary_ids.find(bid) != data.boundary_ids.end())
3526  subdomain_ids.insert(sub_id);
3527 
3528  return subdomain_ids;
3529 }
3530 
3531 std::set<SubdomainID>
3533 {
3534  std::set<SubdomainID> subdomain_ids;
3535  for (const auto & it : _neighbor_subdomain_boundary_ids)
3536  if (it.second.find(bid) != it.second.end())
3537  subdomain_ids.insert(it.first);
3538 
3539  return subdomain_ids;
3540 }
3541 
3542 std::set<SubdomainID>
3544 {
3545  std::set<SubdomainID> subdomain_ids = getBoundaryConnectedBlocks(bid);
3546  for (const auto & it : _neighbor_subdomain_boundary_ids)
3547  if (it.second.find(bid) != it.second.end())
3548  subdomain_ids.insert(it.first);
3549 
3550  return subdomain_ids;
3551 }
3552 
3553 const std::set<SubdomainID> &
3555 {
3556  const auto it = _sub_to_data.find(subdomain_id);
3557 
3558  if (it == _sub_to_data.end())
3559  mooseError("Unable to find subdomain ID: ", subdomain_id, '.');
3560 
3561  return it->second.neighbor_subs;
3562 }
3563 
3564 bool
3566 {
3567  bool found_node = false;
3568  for (const auto & it : _bnd_node_ids)
3569  {
3570  if (it.second.find(node_id) != it.second.end())
3571  {
3572  found_node = true;
3573  break;
3574  }
3575  }
3576  return found_node;
3577 }
3578 
3579 bool
3581 {
3582  bool found_node = false;
3583  std::map<boundary_id_type, std::set<dof_id_type>>::const_iterator it = _bnd_node_ids.find(bnd_id);
3584  if (it != _bnd_node_ids.end())
3585  if (it->second.find(node_id) != it->second.end())
3586  found_node = true;
3587  return found_node;
3588 }
3589 
3590 bool
3592 {
3593  bool found_elem = false;
3594  for (const auto & it : _bnd_elem_ids)
3595  {
3596  if (it.second.find(elem_id) != it.second.end())
3597  {
3598  found_elem = true;
3599  break;
3600  }
3601  }
3602  return found_elem;
3603 }
3604 
3605 bool
3607 {
3608  bool found_elem = false;
3609  auto it = _bnd_elem_ids.find(bnd_id);
3610  if (it != _bnd_elem_ids.end())
3611  if (it->second.find(elem_id) != it->second.end())
3612  found_elem = true;
3613  return found_elem;
3614 }
3615 
3616 void
3617 MooseMesh::errorIfDistributedMesh(std::string name) const
3618 {
3620  mooseError("Cannot use ",
3621  name,
3622  " with DistributedMesh!\n",
3623  "Consider specifying parallel_type = 'replicated' in your input file\n",
3624  "to prevent it from being run with DistributedMesh.");
3625 }
3626 
3627 void
3629 {
3630  if (_use_distributed_mesh && (_partitioner_name != "default" && _partitioner_name != "parmetis"))
3631  {
3632  _partitioner_name = "parmetis";
3633  _partitioner_overridden = true;
3634  }
3635 
3637 }
3638 
3639 void
3641  MooseEnum & partitioner,
3642  bool use_distributed_mesh,
3643  const InputParameters & params,
3644  MooseObject & context_obj)
3645 {
3646  // Set the partitioner based on partitioner name
3647  switch (partitioner)
3648  {
3649  case -3: // default
3650  // We'll use the default partitioner, but notify the user of which one is being used...
3651  if (use_distributed_mesh)
3652  partitioner = "parmetis";
3653  else
3654  partitioner = "metis";
3655  break;
3656 
3657  // No need to explicitily create the metis or parmetis partitioners,
3658  // They are the default for serial and parallel mesh respectively
3659  case -2: // metis
3660  case -1: // parmetis
3661  break;
3662 
3663  case 0: // linear
3664  mesh_base.partitioner().reset(new libMesh::LinearPartitioner);
3665  break;
3666  case 1: // centroid
3667  {
3668  if (!params.isParamValid("centroid_partitioner_direction"))
3669  context_obj.paramError(
3670  "centroid_partitioner_direction",
3671  "If using the centroid partitioner you _must_ specify centroid_partitioner_direction!");
3672 
3673  MooseEnum direction = params.get<MooseEnum>("centroid_partitioner_direction");
3674 
3675  if (direction == "x")
3676  mesh_base.partitioner().reset(
3678  else if (direction == "y")
3679  mesh_base.partitioner().reset(
3681  else if (direction == "z")
3682  mesh_base.partitioner().reset(
3684  else if (direction == "radial")
3685  mesh_base.partitioner().reset(
3687  break;
3688  }
3689  case 2: // hilbert_sfc
3690  mesh_base.partitioner().reset(new libMesh::HilbertSFCPartitioner);
3691  break;
3692  case 3: // morton_sfc
3693  mesh_base.partitioner().reset(new libMesh::MortonSFCPartitioner);
3694  break;
3695  }
3696 }
3697 
3698 void
3700 {
3701  _custom_partitioner = partitioner->clone();
3702 }
3703 
3704 bool
3706 {
3708 }
3709 
3710 bool
3712 {
3713  bool mesh_has_second_order_elements = false;
3714  for (auto it = activeLocalElementsBegin(), end = activeLocalElementsEnd(); it != end; ++it)
3715  if ((*it)->default_order() == SECOND)
3716  {
3717  mesh_has_second_order_elements = true;
3718  break;
3719  }
3720 
3721  // We checked our local elements, so take the max over all processors.
3722  comm().max(mesh_has_second_order_elements);
3723  return mesh_has_second_order_elements;
3724 }
3725 
3726 void
3728 {
3730 }
3731 
3732 std::unique_ptr<libMesh::PointLocatorBase>
3734 {
3735  return getMesh().sub_point_locator();
3736 }
3737 
3738 void
3740 {
3741  mooseAssert(!Threads::in_threads,
3742  "This routine has not been implemented for threads. Please query this routine before "
3743  "a threaded region or contact a MOOSE developer to discuss.");
3744  _finite_volume_info_dirty = false;
3745 
3746  using Keytype = std::pair<const Elem *, unsigned short int>;
3747 
3748  // create a map from elem/side --> boundary ids
3749  std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>> side_list =
3751  std::map<Keytype, std::set<boundary_id_type>> side_map;
3752  for (auto & [elem_id, side, bc_id] : side_list)
3753  {
3754  const Elem * elem = _mesh->elem_ptr(elem_id);
3755  Keytype key(elem, side);
3756  auto & bc_set = side_map[key];
3757  bc_set.insert(bc_id);
3758  }
3759 
3760  _face_info.clear();
3761  _all_face_info.clear();
3762  _elem_side_to_face_info.clear();
3763 
3764  _elem_to_elem_info.clear();
3765  _elem_info.clear();
3766 
3767  // by performing the element ID comparison check in the below loop, we are ensuring that we never
3768  // double count face contributions. If a face lies along a process boundary, the only process that
3769  // will contribute to both sides of the face residuals/Jacobians will be the process that owns the
3770  // element with the lower ID.
3771  auto begin = getMesh().active_elements_begin();
3772  auto end = getMesh().active_elements_end();
3773 
3774  // We prepare a map connecting the Elem* and the corresponding ElemInfo
3775  // for the active elements.
3776  for (const Elem * elem : as_range(begin, end))
3777  _elem_to_elem_info.emplace(elem->id(), elem);
3778 
3779  dof_id_type face_index = 0;
3780  for (const Elem * elem : as_range(begin, end))
3781  {
3782  for (unsigned int side = 0; side < elem->n_sides(); ++side)
3783  {
3784  // get the neighbor element
3785  const Elem * neighbor = elem->neighbor_ptr(side);
3786 
3787  // Check if the FaceInfo shall belong to the element. If yes,
3788  // create and initialize the FaceInfo. We need this to ensure that
3789  // we do not duplicate FaceInfo-s.
3790  if (Moose::FV::elemHasFaceInfo(*elem, neighbor))
3791  {
3792  mooseAssert(!neighbor || (neighbor->level() < elem->level() ? neighbor->active() : true),
3793  "If the neighbor is coarser than the element, we expect that the neighbor must "
3794  "be active.");
3795 
3796  // We construct the faceInfo using the elementinfo and side index
3797  _all_face_info.emplace_back(&_elem_to_elem_info[elem->id()], side, face_index++);
3798 
3799  auto & fi = _all_face_info.back();
3800 
3801  // get all the sidesets that this face is contained in and cache them
3802  // in the face info.
3803  std::set<boundary_id_type> & boundary_ids = fi.boundaryIDs();
3804  boundary_ids.clear();
3805 
3806  // We initialize the weights/other information in faceInfo. If the neighbor does not exist
3807  // or is remote (so when we are on some sort of mesh boundary), we initialize the ghost
3808  // cell and use it to compute the weights corresponding to the faceInfo.
3809  if (!neighbor || neighbor == libMesh::remote_elem)
3810  fi.computeBoundaryCoefficients();
3811  else
3812  fi.computeInternalCoefficients(&_elem_to_elem_info[neighbor->id()]);
3813 
3814  auto lit = side_map.find(Keytype(&fi.elem(), fi.elemSideID()));
3815  if (lit != side_map.end())
3816  boundary_ids.insert(lit->second.begin(), lit->second.end());
3817 
3818  if (fi.neighborPtr())
3819  {
3820  auto rit = side_map.find(Keytype(fi.neighborPtr(), fi.neighborSideID()));
3821  if (rit != side_map.end())
3822  boundary_ids.insert(rit->second.begin(), rit->second.end());
3823  }
3824  }
3825  }
3826  }
3827 
3828  // Build the local face info and elem_side to face info maps. We need to do this after
3829  // _all_face_info is finished being constructed because emplace_back invalidates all iterators and
3830  // references if ever the new size exceeds capacity
3831  for (auto & fi : _all_face_info)
3832  {
3833  const Elem * const elem = &fi.elem();
3834  const auto side = fi.elemSideID();
3835 
3836 #ifndef NDEBUG
3837  auto pair_it =
3838 #endif
3839  _elem_side_to_face_info.emplace(std::make_pair(elem, side), &fi);
3840  mooseAssert(pair_it.second, "We should be adding unique FaceInfo objects.");
3841 
3842  // We will add the faces on processor boundaries to the list of face infos on each
3843  // associated processor.
3844  if (fi.elem().processor_id() == this->processor_id() ||
3845  (fi.neighborPtr() && (fi.neighborPtr()->processor_id() == this->processor_id())))
3846  _face_info.push_back(&fi);
3847  }
3848 
3849  for (auto & ei : _elem_to_elem_info)
3850  if (ei.second.elem()->processor_id() == this->processor_id())
3851  _elem_info.push_back(&ei.second);
3852 }
3853 
3854 const FaceInfo *
3855 MooseMesh::faceInfo(const Elem * elem, unsigned int side) const
3856 {
3857  auto it = _elem_side_to_face_info.find(std::make_pair(elem, side));
3858 
3859  if (it == _elem_side_to_face_info.end())
3860  return nullptr;
3861  else
3862  {
3863  mooseAssert(it->second,
3864  "For some reason, the FaceInfo object is NULL! Try calling "
3865  "`buildFiniteVolumeInfo()` before using this accessor!");
3866  return it->second;
3867  }
3868 }
3869 
3870 const ElemInfo &
3872 {
3873  return libmesh_map_find(_elem_to_elem_info, id);
3874 }
3875 
3876 void
3878 {
3880  mooseError("Trying to compute face- and elem-info coords when the information is dirty");
3881 
3882  for (auto & fi : _all_face_info)
3883  {
3884  // get elem & neighbor elements, and set subdomain ids
3885  const SubdomainID elem_subdomain_id = fi.elemSubdomainID();
3886  const SubdomainID neighbor_subdomain_id = fi.neighborSubdomainID();
3887 
3889  *this, elem_subdomain_id, fi.faceCentroid(), fi.faceCoord(), neighbor_subdomain_id);
3890  }
3891 
3892  for (auto & ei : _elem_to_elem_info)
3894  *this, ei.second.subdomain_id(), ei.second.centroid(), ei.second.coordFactor());
3895 }
3896 
3897 MooseEnum
3899 {
3900  MooseEnum partitioning("default=-3 metis=-2 parmetis=-1 linear=0 centroid hilbert_sfc morton_sfc",
3901  "default");
3902  return partitioning;
3903 }
3904 
3905 MooseEnum
3907 {
3909  "EDGE EDGE2 EDGE3 EDGE4 QUAD QUAD4 QUAD8 QUAD9 TRI3 TRI6 HEX HEX8 HEX20 HEX27 TET4 TET10 "
3910  "PRISM6 PRISM15 PRISM18 PYRAMID5 PYRAMID13 PYRAMID14");
3911  return elemTypes;
3912 }
3913 
3914 void
3915 MooseMesh::allowRemoteElementRemoval(const bool allow_remote_element_removal)
3916 {
3917  _allow_remote_element_removal = allow_remote_element_removal;
3918  if (_mesh)
3919  _mesh->allow_remote_element_removal(allow_remote_element_removal);
3920 
3921  if (!allow_remote_element_removal)
3922  // If we're not allowing remote element removal now, then we will need deletion later after
3923  // late geoemetric ghosting functors have been added (late geometric ghosting functor addition
3924  // happens when algebraic ghosting functors are added)
3925  _need_delete = true;
3926 }
3927 
3928 void
3930 {
3932  if (!_mesh)
3933  mooseError("Cannot delete remote elements because we have not yet attached a MeshBase");
3934 
3935  _mesh->allow_remote_element_removal(true);
3936 
3937  _mesh->delete_remote_elements();
3938 }
3939 
3940 void
3942 {
3943  mooseAssert(
3944  !Threads::in_threads,
3945  "Performing writes to faceInfo variable association maps. This must be done unthreaded!");
3946 
3947  const unsigned int num_eqs = _app.feProblem().es().n_systems();
3948 
3949  auto face_lambda = [this](const SubdomainID elem_subdomain_id,
3950  const SubdomainID neighbor_subdomain_id,
3951  SystemBase & sys,
3952  std::vector<std::vector<FaceInfo::VarFaceNeighbors>> & face_type_vector)
3953  {
3954  face_type_vector[sys.number()].resize(sys.nVariables(), FaceInfo::VarFaceNeighbors::NEITHER);
3955  const auto & variables = sys.getVariables(0);
3956 
3957  for (const auto & var : variables)
3958  {
3959  const unsigned int var_num = var->number();
3960  const unsigned int sys_num = var->sys().number();
3961  std::set<SubdomainID> var_subdomains = var->blockIDs();
3971  bool var_defined_elem = var_subdomains.find(elem_subdomain_id) != var_subdomains.end();
3972  bool var_defined_neighbor =
3973  var_subdomains.find(neighbor_subdomain_id) != var_subdomains.end();
3974  if (var_defined_elem && var_defined_neighbor)
3975  face_type_vector[sys_num][var_num] = FaceInfo::VarFaceNeighbors::BOTH;
3976  else if (!var_defined_elem && !var_defined_neighbor)
3977  face_type_vector[sys_num][var_num] = FaceInfo::VarFaceNeighbors::NEITHER;
3978  else
3979  {
3980  // this is a boundary face for this variable, set elem or neighbor
3981  if (var_defined_elem)
3982  face_type_vector[sys_num][var_num] = FaceInfo::VarFaceNeighbors::ELEM;
3983  else if (var_defined_neighbor)
3984  face_type_vector[sys_num][var_num] = FaceInfo::VarFaceNeighbors::NEIGHBOR;
3985  else
3986  mooseError("Should never get here");
3987  }
3988  }
3989  };
3990 
3991  // We loop through the faces and check if they are internal, boundary or external to
3992  // the variables in the problem
3993  for (FaceInfo & face : _all_face_info)
3994  {
3995  const SubdomainID elem_subdomain_id = face.elemSubdomainID();
3996  const SubdomainID neighbor_subdomain_id = face.neighborSubdomainID();
3997 
3998  auto & face_type_vector = face.faceType();
3999 
4000  face_type_vector.clear();
4001  face_type_vector.resize(num_eqs);
4002 
4003  // First, we check the variables in the solver systems (linear/nonlinear)
4004  for (const auto i : make_range(_app.feProblem().numSolverSystems()))
4005  face_lambda(elem_subdomain_id,
4006  neighbor_subdomain_id,
4008  face_type_vector);
4009 
4010  // Then we check the variables in the auxiliary system
4011  face_lambda(elem_subdomain_id,
4012  neighbor_subdomain_id,
4014  face_type_vector);
4015  }
4016 }
4017 
4018 void
4020 {
4021  mooseAssert(!Threads::in_threads,
4022  "Performing writes to elemInfo dof indices. This must be done unthreaded!");
4023 
4024  auto elem_lambda = [](const ElemInfo & elem_info,
4025  SystemBase & sys,
4026  std::vector<std::vector<dof_id_type>> & dof_vector)
4027  {
4028  if (sys.nFVVariables())
4029  {
4030  dof_vector[sys.number()].resize(sys.nVariables(), libMesh::DofObject::invalid_id);
4031  const auto & variables = sys.getVariables(0);
4032 
4033  for (const auto & var : variables)
4034  if (var->isFV())
4035  {
4036  const auto & var_subdomains = var->blockIDs();
4037 
4038  // We will only cache for FV variables and if they live on the current subdomain
4039  if (var_subdomains.find(elem_info.subdomain_id()) != var_subdomains.end())
4040  {
4041  std::vector<dof_id_type> indices;
4042  var->dofMap().dof_indices(elem_info.elem(), indices, var->number());
4043  mooseAssert(indices.size() == 1, "We expect to have only one dof per element!");
4044  dof_vector[sys.number()][var->number()] = indices[0];
4045  }
4046  }
4047  }
4048  };
4049 
4050  const unsigned int num_eqs = _app.feProblem().es().n_systems();
4051 
4052  // We loop through the elements in the mesh and cache the dof indices
4053  // for the corresponding variables.
4054  for (auto & ei_pair : _elem_to_elem_info)
4055  {
4056  auto & elem_info = ei_pair.second;
4057  auto & dof_vector = elem_info.dofIndices();
4058 
4059  dof_vector.clear();
4060  dof_vector.resize(num_eqs);
4061 
4062  // First, we cache the dof indices for the variables in the solver systems (linear, nonlinear)
4063  for (const auto i : make_range(_app.feProblem().numSolverSystems()))
4064  elem_lambda(elem_info, _app.feProblem().getSolverSystem(i), dof_vector);
4065 
4066  // Then we cache the dof indices for the auxvariables
4067  elem_lambda(elem_info, _app.feProblem().getAuxiliarySystem(), dof_vector);
4068  }
4069 }
4070 
4071 void
4073 {
4078 }
4079 
4080 void
4081 MooseMesh::setCoordSystem(const std::vector<SubdomainName> & blocks,
4082  const MultiMooseEnum & coord_sys)
4083 {
4084  TIME_SECTION("setCoordSystem", 5, "Setting Coordinate System");
4086  {
4087  const std::string param_name = isParamValid("coord_block") ? "coord_block" : "block";
4088  mooseWarning("Supplied blocks in the 'setCoordSystem' method do not match the value of the "
4089  "'Mesh/",
4090  param_name,
4091  "' parameter. Did you provide different parameter values for 'Mesh/",
4092  param_name,
4093  "' and 'Problem/block'?. We will honor the parameter value from 'Mesh/",
4094  param_name,
4095  "'");
4096  mooseAssert(_coord_system_set,
4097  "If we are arriving here due to a bad specification in the Problem block, then we "
4098  "should have already set our coordinate system subdomains from the Mesh block");
4099  return;
4100  }
4101  if (_pars.isParamSetByUser("coord_type") && getParam<MultiMooseEnum>("coord_type") != coord_sys)
4102  mooseError("Supplied coordinate systems in the 'setCoordSystem' method do not match the value "
4103  "of the 'Mesh/coord_type' parameter. Did you provide different parameter values for "
4104  "'coord_type' to 'Mesh' and 'Problem'?");
4105 
4106  // If blocks contain ANY_BLOCK_ID, it should be the only block specified, and coord_sys should
4107  // have one and only one entry. In that case, the same coordinate system will be set for all
4108  // subdomains.
4109  if (blocks.size() == 1 && blocks[0] == "ANY_BLOCK_ID")
4110  {
4111  if (coord_sys.size() > 1)
4112  mooseError("If you specify ANY_BLOCK_ID as the only block, you must also specify a single "
4113  "coordinate system for it.");
4114  if (!_mesh->is_prepared())
4115  mooseError(
4116  "You cannot set the coordinate system for ANY_BLOCK_ID before the mesh is prepared. "
4117  "Please call this method after the mesh is prepared.");
4118  const auto coord_type = coord_sys.size() == 0
4120  : Moose::stringToEnum<Moose::CoordinateSystemType>(coord_sys[0]);
4121  for (const auto sid : meshSubdomains())
4122  _coord_sys[sid] = coord_type;
4123  return;
4124  }
4125 
4126  // If multiple blocks are specified, but one of them is ANY_BLOCK_ID, let's emit a helpful error
4127  if (std::find(blocks.begin(), blocks.end(), "ANY_BLOCK_ID") != blocks.end())
4128  mooseError("You cannot specify ANY_BLOCK_ID together with other blocks in the "
4129  "setCoordSystem() method. If you want to set the same coordinate system for all "
4130  "blocks, use ANY_BLOCK_ID as the only block.");
4131 
4132  auto subdomains = meshSubdomains();
4133  // It's possible that a user has called this API before the mesh is prepared and consequently we
4134  // don't yet have the subdomains in meshSubdomains()
4135  for (const auto & sub_name : blocks)
4136  {
4137  const auto sub_id = getSubdomainID(sub_name);
4138  subdomains.insert(sub_id);
4139  }
4140 
4141  if (coord_sys.size() <= 1)
4142  {
4143  // We will specify the same coordinate system for all blocks
4144  const auto coord_type = coord_sys.size() == 0
4146  : Moose::stringToEnum<Moose::CoordinateSystemType>(coord_sys[0]);
4147  for (const auto sid : subdomains)
4148  _coord_sys[sid] = coord_type;
4149  }
4150  else
4151  {
4152  if (blocks.size() != coord_sys.size())
4153  mooseError("Number of blocks and coordinate systems does not match.");
4154 
4155  for (const auto i : index_range(blocks))
4156  {
4157  SubdomainID sid = getSubdomainID(blocks[i]);
4158  Moose::CoordinateSystemType coord_type =
4159  Moose::stringToEnum<Moose::CoordinateSystemType>(coord_sys[i]);
4160  _coord_sys[sid] = coord_type;
4161  }
4162 
4163  for (const auto & sid : subdomains)
4164  if (_coord_sys.find(sid) == _coord_sys.end())
4165  mooseError("Subdomain '" + Moose::stringify(sid) +
4166  "' does not have a coordinate system specified.");
4167  }
4168 
4169  _coord_system_set = true;
4170 
4172 }
4173 
4176 {
4177  auto it = _coord_sys.find(sid);
4178  if (it != _coord_sys.end())
4179  return (*it).second;
4180  else
4181  mooseError("Requested subdomain ", sid, " does not exist.");
4182 }
4183 
4186 {
4187  const auto unique_system = _coord_sys.find(*meshSubdomains().begin())->second;
4188  // Check that it is actually unique
4189  bool result = std::all_of(
4190  std::next(_coord_sys.begin()),
4191  _coord_sys.end(),
4192  [unique_system](
4193  typename std::unordered_map<SubdomainID, Moose::CoordinateSystemType>::const_reference
4194  item) { return (item.second == unique_system); });
4195  if (!result)
4196  mooseError("The unique coordinate system of the mesh was requested by the mesh contains "
4197  "multiple blocks with different coordinate systems");
4198 
4200  mooseError("General axisymmetric coordinate axes are being used, and it is currently "
4201  "conservatively assumed that in this case there is no unique coordinate system.");
4202 
4203  return unique_system;
4204 }
4205 
4206 const std::map<SubdomainID, Moose::CoordinateSystemType> &
4208 {
4209  return _coord_sys;
4210 }
4211 
4212 void
4214 {
4215  _rz_coord_axis = rz_coord_axis;
4216 
4218 }
4219 
4220 void
4222  const std::vector<SubdomainName> & blocks,
4223  const std::vector<std::pair<Point, RealVectorValue>> & axes)
4224 {
4225  // Set the axes for the given blocks
4226  mooseAssert(blocks.size() == axes.size(), "Blocks and axes vectors must be the same length.");
4227  for (const auto i : index_range(blocks))
4228  {
4229  const auto subdomain_id = getSubdomainID(blocks[i]);
4230  const auto it = _coord_sys.find(subdomain_id);
4231  if (it == _coord_sys.end())
4232  mooseError("The block '",
4233  blocks[i],
4234  "' has not set a coordinate system. Make sure to call setCoordSystem() before "
4235  "setGeneralAxisymmetricCoordAxes().");
4236  else
4237  {
4238  if (it->second == Moose::COORD_RZ)
4239  {
4240  const auto direction = axes[i].second;
4241  if (direction.is_zero())
4242  mooseError("Only nonzero vectors may be supplied for RZ directions.");
4243 
4244  _subdomain_id_to_rz_coord_axis[subdomain_id] =
4245  std::make_pair(axes[i].first, direction.unit());
4246  }
4247  else
4248  mooseError("The block '",
4249  blocks[i],
4250  "' was provided in setGeneralAxisymmetricCoordAxes(), but the coordinate system "
4251  "for this block is not 'RZ'.");
4252  }
4253  }
4254 
4255  // Make sure there are no RZ blocks that still do not have axes
4256  const auto all_subdomain_ids = meshSubdomains();
4257  for (const auto subdomain_id : all_subdomain_ids)
4258  if (getCoordSystem(subdomain_id) == Moose::COORD_RZ &&
4259  !_subdomain_id_to_rz_coord_axis.count(subdomain_id))
4260  mooseError("The block '",
4261  getSubdomainName(subdomain_id),
4262  "' was specified to use the 'RZ' coordinate system but was not given in "
4263  "setGeneralAxisymmetricCoordAxes().");
4264 
4266 }
4267 
4268 const std::pair<Point, RealVectorValue> &
4270 {
4271  auto it = _subdomain_id_to_rz_coord_axis.find(subdomain_id);
4272  if (it != _subdomain_id_to_rz_coord_axis.end())
4273  return (*it).second;
4274  else
4275  mooseError("Requested subdomain ", subdomain_id, " does not exist.");
4276 }
4277 
4278 bool
4280 {
4281  return _subdomain_id_to_rz_coord_axis.size() > 0;
4282 }
4283 
4284 void
4286 {
4287  if (!_coord_transform)
4288  _coord_transform = std::make_unique<MooseAppCoordTransform>(*this);
4289  else
4290  _coord_transform->setCoordinateSystem(*this);
4291 }
4292 
4293 unsigned int
4295 {
4297  mooseError("getAxisymmetricRadialCoord() should not be called if "
4298  "setGeneralAxisymmetricCoordAxes() has been called.");
4299 
4300  if (_rz_coord_axis == 0)
4301  return 1; // if the rotation axis is x (0), then the radial direction is y (1)
4302  else
4303  return 0; // otherwise the radial direction is assumed to be x, i.e., the rotation axis is y
4304 }
4305 
4306 void
4308 {
4309  for (const auto & elem : getMesh().element_ptr_range())
4310  {
4311  SubdomainID sid = elem->subdomain_id();
4312  if (_coord_sys[sid] == Moose::COORD_RZ && elem->dim() == 3)
4313  mooseError("An RZ coordinate system was requested for subdomain " + Moose::stringify(sid) +
4314  " which contains 3D elements.");
4315  if (_coord_sys[sid] == Moose::COORD_RSPHERICAL && elem->dim() > 1)
4316  mooseError("An RSPHERICAL coordinate system was requested for subdomain " +
4317  Moose::stringify(sid) + " which contains 2D or 3D elements.");
4318  }
4319 }
4320 
4321 void
4323 {
4324  _coord_sys = other_mesh._coord_sys;
4325  _rz_coord_axis = other_mesh._rz_coord_axis;
4327 }
4328 
4329 const MooseUnits &
4331 {
4332  mooseAssert(_coord_transform, "This must be non-null");
4333  return _coord_transform->lengthUnit();
4334 }
4335 
4336 void
4338 {
4339  std::map<SubdomainName, SubdomainID> subdomain;
4340  for (const auto & sbd_id : _mesh_subdomains)
4341  {
4342  std::string sub_name = getSubdomainName(sbd_id);
4343  if (!sub_name.empty() && subdomain.count(sub_name) > 0)
4344  mooseError("The subdomain name ",
4345  sub_name,
4346  " is used for both subdomain with ID=",
4347  subdomain[sub_name],
4348  " and ID=",
4349  sbd_id,
4350  ", Please rename one of them!");
4351  else
4352  subdomain[sub_name] = sbd_id;
4353  }
4354 }
4355 
4356 const std::vector<QpMap> &
4358  const Elem & elem,
4359  const std::map<std::pair<ElemType, unsigned int>, std::vector<QpMap>> & map) const
4360 {
4361  // We are actually seeking the map stored with the p_level - 1 key, e.g. the refinement map that
4362  // maps from the previous p_level to this element's p_level
4363  return libmesh_map_find(map,
4364  std::make_pair(elem.type(), cast_int<unsigned int>(elem.p_level() - 1)));
4365 }
4366 
4367 const std::vector<QpMap> &
4369  const Elem & elem,
4370  const std::map<std::pair<ElemType, unsigned int>, std::vector<QpMap>> & map) const
4371 {
4372  mooseAssert(elem.active() && elem.p_refinement_flag() == Elem::JUST_COARSENED,
4373  "These are the conditions that should be met for requesting a coarsening map");
4374  return libmesh_map_find(map, std::make_pair(elem.type(), elem.p_level()));
4375 }
4376 
4377 const std::vector<QpMap> &
4379 {
4381 }
4382 
4383 const std::vector<QpMap> &
4385 {
4387 }
4388 
4389 const std::vector<QpMap> &
4391 {
4393 }
4394 
4395 const std::vector<QpMap> &
4397 {
4399 }
4400 
4401 bool
4403 {
4404  return _mesh->skip_noncritical_partitioning();
4405 }
ParallelType _parallel_type
Can be set to DISTRIBUTED, REPLICATED, or DEFAULT.
Definition: MooseMesh.h:1423
static InputParameters validParams()
Typical "Moose-style" constructor and copy constructor.
Definition: MooseMesh.C:84
virtual bnd_node_iterator bndNodesEnd()
Definition: MooseMesh.C:1553
virtual bnd_elem_iterator bndElemsEnd()
Definition: MooseMesh.C:1569
std::vector< std::vector< Real > > _bounds
The bounds in each dimension of the mesh for regular orthogonal meshes.
Definition: MooseMesh.h:1580
std::set< Node * > _semilocal_node_list
Used for generating the semilocal node range.
Definition: MooseMesh.h:1488
void remove_id(boundary_id_type id, bool global=false)
std::map< dof_id_type, Node * > _quadrature_nodes
Definition: MooseMesh.h:1547
const std::vector< QpMap > & getPCoarseningSideMap(const Elem &elem) const
Get the map describing for each side quadrature point (qp) on the coarse level which qp on the previo...
Definition: MooseMesh.C:4396
virtual Real getMaxInDimension(unsigned int component) const
Definition: MooseMesh.C:2208
static const std::string & checkpointSuffix()
The file suffix for the checkpoint mesh.
Definition: MooseApp.C:3072
bool _node_to_elem_map_built
Definition: MooseMesh.h:1506
std::unique_ptr< libMesh::NodeRange > _active_node_range
Definition: MooseMesh.h:1497
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
std::vector< Node * > _extreme_nodes
A vector containing the nodes at the corners of a regular orthogonal mesh.
Definition: MooseMesh.h:1634
ElemType
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:1610
std::vector< std::tuple< dof_id_type, unsigned short int, boundary_id_type > > buildSideList()
As above, but uses the non-deprecated std::tuple interface.
Definition: MooseMesh.C:3009
const std::set< BoundaryID > & meshNodesetIds() const
Returns a read-only reference to the set of nodesets currently present in the Mesh.
Definition: MooseMesh.C:3189
void buildElemIDInfo()
Build extra data for faster access to the information of extra element integers.
Definition: MooseMesh.C:1059
std::vector< FaceInfo > _all_face_info
FaceInfo object storing information for face based loops.
Definition: MooseMesh.h:1604
std::unique_ptr< FEGenericBase< Real > > build(const unsigned int dim, const FEType &fet)
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
std::vector< const FaceInfo * > _face_info
Holds only those FaceInfo objects that have processor_id equal to this process&#39;s id, e.g.
Definition: MooseMesh.h:1608
bool allowRemoteElementRemoval() const
Whether we are allow remote element removal.
Definition: MooseMesh.h:1094
virtual Real getMinInDimension(unsigned int component) const
Returns the min or max of the requested dimension respectively.
Definition: MooseMesh.C:2199
bool elemHasFaceInfo(const Elem &elem, const Elem *const neighbor)
This function infers based on elements if the faceinfo between them belongs to the element or not...
Definition: FVUtils.C:21
libMesh::ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
Definition: MooseMesh.C:1238
virtual void read(const std::string &name, void *mesh_data=nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false)=0
virtual void onMeshChanged()
Declares a callback function that is executed at the conclusion of meshChanged(). ...
Definition: MooseMesh.C:908
bool prepared() const
Setter/getter for whether the mesh is prepared.
Definition: MooseMesh.C:3137
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:3165
const std::set< boundary_id_type > & get_side_boundary_ids() const
const std::set< BoundaryID > & getBoundaryIDs() const
Returns a const reference to a set of all user-specified boundary IDs.
Definition: MooseMesh.C:2977
bool _is_nemesis
True if a Nemesis Mesh was read in.
Definition: MooseMesh.h:1469
static const std::string name_param
The name of the parameter that contains the object name.
Definition: MooseBase.h:55
std::vector< SubdomainName > _provided_coord_blocks
Set for holding user-provided coordinate system type block names.
Definition: MooseMesh.h:1874
virtual MooseMesh & clone() const
Clone method.
Definition: MooseMesh.C:2822
const Elem * parent() const
bool is_prepared() const
const std::vector< QpMap > & getPCoarseningMapHelper(const Elem &elem, const std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap >> &) const
Definition: MooseMesh.C:4368
bool isCustomPartitionerRequested() const
Setter and getter for _custom_partitioner_requested.
Definition: MooseMesh.C:3705
bool _need_ghost_ghosted_boundaries
A parallel mesh generator such as DistributedRectilinearMeshGenerator already make everything ready...
Definition: MooseMesh.h:1837
A class for creating restricted objects.
Definition: Restartable.h:28
bool isUltimateMaster() const
Whether or not this app is the ultimate master app.
Definition: MooseApp.h:813
unsigned int _uniform_refine_level
The level of uniform refinement requested (set to zero if AMR is disabled)
Definition: MooseMesh.h:1457
virtual Node *& set_node(const unsigned int i)
const std::vector< QpMap > & getPRefinementMapHelper(const Elem &elem, const std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap >> &) const
Definition: MooseMesh.C:4357
unsigned int n_systems() const
std::vector< dof_id_type > _min_ids
Minimum integer ID for each extra element integer.
Definition: MooseMesh.h:1844
Helper class for sorting Boundary Nodes so that we always get the same order of application for bound...
Definition: MooseMesh.C:1001
const InputParameters & _pars
The object&#39;s parameters.
Definition: MooseBase.h:362
libMesh::QBase *const & writeableQRule()
Returns the reference to the current quadrature being used.
Definition: Assembly.h:224
virtual unique_id_type parallel_max_unique_id() const=0
std::unordered_set< dof_id_type > getBoundaryActiveSemiLocalElemIds(BoundaryID bid) const
Return all ids of elements which have a side which is part of a sideset.
Definition: MooseMesh.C:1331
auto norm() const -> decltype(std::norm(Real()))
std::string & nodeset_name(boundary_id_type id)
const MooseUnits & lengthUnit() const
Definition: MooseMesh.C:4330
void checkDuplicateSubdomainNames()
Loop through all subdomain IDs and check if there is name duplication used for the subdomains with sa...
Definition: MooseMesh.C:4337
const unsigned int invalid_uint
const std::set< BoundaryID > & getSubdomainBoundaryIds(const SubdomainID subdomain_id) const
Get the list of boundary ids associated with the given subdomain id.
Definition: MooseMesh.C:3497
std::unique_ptr< PointLocatorBase > sub_point_locator() const
void build_node_list_from_side_list()
RealVectorValue _half_range
A convenience vector used to hold values in each dimension representing half of the range...
Definition: MooseMesh.h:1631
Keeps track of stuff related to assembling.
Definition: Assembly.h:101
std::map< libMesh::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:1762
void setCoordData(const MooseMesh &other_mesh)
Set the coordinate system data to that of other_mesh.
Definition: MooseMesh.C:4322
static InputParameters validParams()
Describes the parameters this object can take to setup transformations.
const Elem * interior_parent() const
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
Definition: MooseBase.h:435
virtual ~MooseMesh()
Definition: MooseMesh.C:368
void freeBndElems()
Definition: MooseMesh.C:395
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
Definition: MooseBase.h:384
bool _finite_volume_info_dirty
Definition: MooseMesh.h:1615
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:3113
void allow_renumbering(bool allow)
char ** blocks
The definition of the bnd_elem_iterator struct.
Definition: MooseMesh.h:2040
std::map< SubdomainID, Moose::CoordinateSystemType > & _coord_sys
Type of coordinate system per subdomain.
Definition: MooseMesh.h:1858
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:3565
face_info_iterator ownedFaceInfoBegin()
Iterators to owned faceInfo objects.
Definition: MooseMesh.C:1510
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:1206
std::vector< std::tuple< dof_id_type, unsigned short int, boundary_id_type > > buildActiveSideList() const
Calls BoundaryInfo::build_active_side_list.
Definition: MooseMesh.C:3015
IntRange< unsigned short > side_index_range() const
static void setPartitioner(MeshBase &mesh_base, MooseEnum &partitioner, bool use_distributed_mesh, const InputParameters &params, MooseObject &context_obj)
Method for setting the partitioner on the passed in mesh_base object.
Definition: MooseMesh.C:3640
void skip_partitioning(bool skip)
void buildLowerDMesh()
Build lower-d mesh for all sides.
Definition: MooseMesh.C:658
const std::set< BoundaryID > & meshSidesetIds() const
Returns a read-only reference to the set of sidesets currently present in the Mesh.
Definition: MooseMesh.C:3183
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i)=0
std::unordered_map< const Elem *, unsigned short int > _lower_d_elem_to_higher_d_elem_side
Definition: MooseMesh.h:1812
const BoundaryID INVALID_BOUNDARY_ID
Definition: MooseTypes.C:22
std::vector< std::pair< R1, R2 > > get(const std::string &param1, const std::string &param2) const
Combine two vector parameters into a single vector of pairs.
void cacheFVElementalDoFs() const
Cache the DoF indices for FV variables on each element.
Definition: MooseMesh.C:4019
static constexpr Real TOLERANCE
void cacheFaceInfoVariableOwnership() const
Cache if variables live on the elements connected by the FaceInfo objects.
Definition: MooseMesh.C:3941
bool _custom_partitioner_requested
Definition: MooseMesh.h:1441
const std::unordered_map< boundary_id_type, std::unordered_set< dof_id_type > > & getBoundariesToElems() const
Returns a map of boundaries to ids of elements on the boundary.
Definition: MooseMesh.C:1317
Moose::CoordinateSystemType getUniqueCoordSystem() const
Get the coordinate system from the mesh, it must be the same in all subdomains otherwise this will er...
Definition: MooseMesh.C:4185
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
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:1505
MooseMesh()=delete
bool isSemiLocal(Node *const node) const
Returns true if the node is semi-local.
Definition: MooseMesh.C:992
unsigned int which_side_am_i(const Elem *e) const
virtual bool is_child_on_side(const unsigned int c, const unsigned int s) const=0
const Elem * getLowerDElem(const Elem *, unsigned short int) const
Returns a const pointer to a lower dimensional element that corresponds to a side of a higher dimensi...
Definition: MooseMesh.C:1701
std::unique_ptr< libMesh::StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > > _bnd_elem_range
Definition: MooseMesh.h:1502
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:2262
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:2579
FIRST
std::unordered_map< boundary_id_type, std::unordered_set< dof_id_type > > _bnd_elem_ids
Map of set of elem IDs connected to each boundary.
Definition: MooseMesh.h:1545
void prepare_for_use(const bool skip_renumber_nodes_and_elements, const bool skip_find_neighbors)
RefinementState p_refinement_flag() const
bool _doing_p_refinement
Whether we have p-refinement (as opposed to h-refinement)
Definition: MooseMesh.h:1877
const Elem * elem() const
Definition: ElemInfo.h:34
void determineUseDistributedMesh()
Determine whether to use a distributed mesh.
Definition: MooseMesh.C:2828
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:2515
const boundary_id_type side_id
const std::string & getBoundaryName(BoundaryID boundary_id)
Return the name of the boundary given the id.
Definition: MooseMesh.C:1792
void cacheChangedLists()
Cache information about what elements were refined and coarsened in the previous step.
Definition: MooseMesh.C:913
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:2251
void coordTransformFactor(const SubProblem &s, SubdomainID sub_id, const P &point, C &factor, SubdomainID neighbor_sub_id=libMesh::Elem::invalid_subdomain_id)
Computes a conversion multiplier for use when computing integraals for the current coordinate system ...
Definition: Assembly.C:41
const Point & getPoint(const PointObject &item) const
get a Point reference from the PointData object at index idx in the list
void family_tree(std::vector< const Elem * > &family, bool reset=true) const
unsigned int side_with_boundary_id(const Elem *const elem, const boundary_id_type boundary_id) const
The definition of the bnd_node_iterator struct.
Definition: MooseMesh.h:1997
unsigned int n_elem_integers() const
std::vector< const ElemInfo * > _elem_info
Holds only those ElemInfo objects that have processor_id equal to this process&#39;s id, e.g.
Definition: MooseMesh.h:1600
void buildHRefinementAndCoarseningMaps(Assembly *assembly)
Definition: MooseMesh.C:2311
MeshBase & mesh
boundary_id_type pairedboundary
PeriodicBoundaryBase * boundary(boundary_id_type id)
bool usingGeneralAxisymmetricCoordAxes() const
Returns true if general axisymmetric coordinate axes are being used.
Definition: MooseMesh.C:4279
std::map< std::pair< int, libMesh::ElemType >, std::vector< std::vector< QpMap > > > _elem_type_to_refinement_map
Holds mappings for volume to volume and parent side to child side Map key:
Definition: MooseMesh.h:1753
const std::set< SubdomainID > & getBlockConnectedBlocks(const SubdomainID subdomain_id) const
Get the list of subdomains neighboring a given subdomain.
Definition: MooseMesh.C:3554
virtual const Node * queryNodePtr(const dof_id_type i) const
Definition: MooseMesh.C:860
virtual std::unique_ptr< Partitioner > & partitioner()
std::unordered_map< std::pair< const Elem *, unsigned short int >, const Elem * > _higher_d_elem_side_to_lower_d_elem
Holds a map from a high-order element side to its corresponding lower-d element.
Definition: MooseMesh.h:1811
Helper object for holding qp mapping info.
Definition: MooseMesh.h:69
unsigned int size() const
Return the number of active items in the MultiMooseEnum.
void mooseInfoRepeated(Args &&... args)
Emit an informational message with the given stringified, concatenated args.
Definition: MooseError.h:408
bool _has_lower_d
Whether there are any lower-dimensional blocks that are manifolds of higher-dimensional block faces...
Definition: MooseMesh.h:1816
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const std::vector< Real > & getGhostedBoundaryInflation() const
Return a writable reference to the _ghosted_boundaries_inflation vector.
Definition: MooseMesh.C:3247
std::unique_ptr< ConstElemPointerRange > _refined_elements
The elements that were just refined.
Definition: MooseMesh.h:1475
std::vector< std::vector< bool > > _id_identical_flag
Flags to indicate whether or not any two extra element integers are the same.
Definition: MooseMesh.h:1846
virtual dof_id_type maxElemId() const
Definition: MooseMesh.C:3093
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:159
const Parallel::Communicator & comm() const
virtual bnd_elem_iterator bndElemsBegin()
Return iterators to the beginning/end of the boundary elements list.
Definition: MooseMesh.C:1561
void setUniformRefineLevel(unsigned int, bool deletion=true)
Set uniform refinement level.
Definition: MooseMesh.C:3222
virtual unsigned int n_children() const=0
unsigned int p_level() const
void allgather_packed_range(Context *context, Iter range_begin, const Iter range_end, OutputIter out, std::size_t approx_buffer_size=1000000) const
MooseEnum _partitioner_name
The partitioner used on this mesh.
Definition: MooseMesh.h:1436
std::unique_ptr< libMesh::ConstElemRange > _active_local_elem_range
A range for use with threading.
Definition: MooseMesh.h:1494
std::map< dof_id_type, std::set< SubdomainID > > _block_node_list
list of nodes that belongs to a specified block (domain)
Definition: MooseMesh.h:1553
bool _node_to_active_semilocal_elem_map_built
Definition: MooseMesh.h:1510
void build_side_boundary_ids(std::vector< boundary_id_type > &b_ids) const
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:1537
virtual void init()
Initialize the Mesh object.
Definition: MooseMesh.C:2875
ConstElemPointerRange * refinedElementRange() const
Return a range that is suitable for threaded execution over elements that were just refined...
Definition: MooseMesh.C:931
void boundary_ids(const Node *node, std::vector< boundary_id_type > &vec_to_fill) const
unsigned int getHigherDSide(const Elem *elem) const
Returns the local side ID of the interior parent aligned with the lower dimensional element...
Definition: MooseMesh.C:1712
std::unordered_map< std::pair< const Elem *, unsigned int >, FaceInfo * > _elem_side_to_face_info
Map from elem-side pair to FaceInfo.
Definition: MooseMesh.h:1612
void detectPairedSidesets()
This routine detects paired sidesets of a regular orthogonal mesh (.i.e.
Definition: MooseMesh.C:1999
const std::set< SubdomainID > & getNodeBlockIds(const Node &node) const
Return list of blocks to which the given node belongs.
Definition: MooseMesh.C:1499
virtual std::unique_ptr< Partitioner > clone() const=0
const Parallel::Communicator & _communicator
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_coarsening_side_map
Definition: MooseMesh.h:1784
void setMeshBoundaryIDs(std::set< BoundaryID > boundary_IDs)
Sets the set of BoundaryIDs Is called by AddAllSideSetsByNormals.
Definition: MooseMesh.C:3195
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:1556
std::vector< subdomain_id_type > getSubdomainIDs(const libMesh::MeshBase &mesh, const std::vector< SubdomainName > &subdomain_name)
Get the associated subdomainIDs for the subdomain names that are passed in.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
std::set< SubdomainID > _lower_d_boundary_blocks
Mesh blocks for boundary lower-d elements in different types.
Definition: MooseMesh.h:1808
void cacheInfo()
Definition: MooseMesh.C:1407
std::basic_ostream< charT, traits > * os
Definition: InfixIterator.h:33
void build_active_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
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:2762
Base class for a system (of equations)
Definition: SystemBase.h:84
const BoundaryInfo & get_boundary_info() const
std::set< Elem * > _ghost_elems_from_ghost_boundaries
Set of elements ghosted by ghostGhostedBoundaries.
Definition: MooseMesh.h:1831
SECOND
virtual void buildMesh()=0
Must be overridden by child classes.
void setPartitionerHelper(MeshBase *mesh=nullptr)
Definition: MooseMesh.C:3628
void deleteRemoteElements()
Delete remote elements.
Definition: MooseMesh.C:3929
bool isSplitMesh() const
Whether or not this is a split mesh operation.
Definition: MooseApp.C:1813
unsigned int _to
The qp to map to.
Definition: MooseMesh.h:78
bool in_threads
BoundaryID _bnd_id
boundary id for the node
Definition: BndNode.h:26
libMesh::ConstNodeRange * getLocalNodeRange()
Definition: MooseMesh.C:1275
virtual std::unique_ptr< MeshBase > clone() const=0
Real distance(const Point &p)
virtual const Node & nodeRef(const dof_id_type i) const
Definition: MooseMesh.C:834
bool _allow_recovery
Whether or not this Mesh is allowed to read a recovery file.
Definition: MooseMesh.h:1819
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id)=0
std::vector< SubdomainID > getSubdomainIDs(const std::vector< SubdomainName > &subdomain_names) const
Get the associated subdomainIDs for the subdomain names that are passed in.
Definition: MooseMesh.C:1737
bool operator()(const BndNode *const &lhs, const BndNode *const &rhs)
Definition: MooseMesh.C:1006
void buildNodeListFromSideList()
Calls BoundaryInfo::build_node_list_from_side_list().
Definition: MooseMesh.C:2983
FEProblemBase & feProblem() const
Definition: MooseApp.C:1960
const std::string & getSubdomainName(SubdomainID subdomain_id) const
Return the name of a block given an id.
Definition: MooseMesh.C:1763
void setPatchUpdateStrategy(Moose::PatchUpdateType patch_update_strategy)
Set the patch size update strategy.
Definition: MooseMesh.C:3407
void buildFiniteVolumeInfo() const
Builds the face and elem info vectors that store meta-data needed for looping over and doing calculat...
Definition: MooseMesh.C:3739
std::unordered_map< SubdomainID, std::set< BoundaryID > > _neighbor_subdomain_boundary_ids
Holds a map from neighbor subomdain ids to the boundary ids that are attached to it.
Definition: MooseMesh.h:1803
const std::pair< Point, RealVectorValue > & getGeneralAxisymmetricCoordAxis(SubdomainID subdomain_id) const
Gets the general axisymmetric coordinate axis for a block.
Definition: MooseMesh.C:4269
void reinit(const Elem *elem)
Reinitialize objects (JxW, q_points, ...) for an elements.
SubdomainID getSubdomainID(const SubdomainName &subdomain_name, const MeshBase &mesh)
Gets the subdomain ID associated with the given SubdomainName.
auto max(const L &left, const R &right)
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:3177
std::set< SubdomainID > _lower_d_interior_blocks
Mesh blocks for interior lower-d elements in different types.
Definition: MooseMesh.h:1806
virtual Elem * queryElemPtr(const dof_id_type i)
Definition: MooseMesh.C:3125
elem_info_iterator ownedElemInfoBegin()
Iterators to owned faceInfo objects.
Definition: MooseMesh.C:1528
void build_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
void setIsCustomPartitionerRequested(bool cpr)
Definition: MooseMesh.C:3727
const std::set< boundary_id_type > & get_node_boundary_ids() const
unsigned int _max_h_level
Maximum h-refinement level of all elements.
Definition: MooseMesh.h:1881
virtual bnd_node_iterator bndNodesBegin()
Return iterators to the beginning/end of the boundary nodes list.
Definition: MooseMesh.C:1545
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:2590
BoundaryID getBoundaryID(const BoundaryName &boundary_name, const MeshBase &mesh)
Gets the boundary ID associated with the given BoundaryName.
virtual Point get_corresponding_pos(const Point &pt) const=0
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:3617
const MeshBase * getMeshPtr() const
Definition: MooseMesh.C:3442
bool getDistributedMeshOnCommandLine() const
Returns true if the user specified –distributed-mesh (or –parallel-mesh, for backwards compatibilit...
Definition: MooseApp.h:453
This data structure is used to store geometric and variable related metadata about each cell face in ...
Definition: FaceInfo.h:36
const std::vector< const FaceInfo * > & faceInfo() const
Accessor for local FaceInfo objects.
Definition: MooseMesh.h:2173
void updateCoordTransform()
Update the coordinate transformation object based on our coordinate system data.
Definition: MooseMesh.C:4285
static const int GRAIN_SIZE
Definition: MooseMesh.C:67
bool _use_distributed_mesh
False by default.
Definition: MooseMesh.h:1428
virtual bool is_serial() const
void libmesh_ignore(const Args &...)
void build_node_list(std::vector< dof_id_type > &node_id_list, std::vector< boundary_id_type > &bc_id_list) const
static const subdomain_id_type invalid_subdomain_id
CONSTANT
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:1530
const dof_id_type n_nodes
const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:99
virtual void find_neighbors(const bool reset_remote_elements=false, const bool reset_current_list=true)=0
bool _built_from_other_mesh
Whether or not this mesh was built from another mesh.
Definition: MooseMesh.h:1419
bool _allow_remote_element_removal
Whether to allow removal of remote elements.
Definition: MooseMesh.h:1828
virtual bool skipNoncriticalPartitioning() const
Definition: MooseMesh.C:4402
SemiLocalNodeRange * getActiveSemiLocalNodeRange() const
Definition: MooseMesh.C:1266
std::unordered_set< dof_id_type > getBoundaryActiveNeighborElemIds(BoundaryID bid) const
Return all ids of neighbors of elements which have a side which is part of a sideset.
Definition: MooseMesh.C:1342
int8_t boundary_id_type
void setSubdomainName(SubdomainID subdomain_id, const SubdomainName &name)
This method sets the name for subdomain_id to name.
Definition: MooseMesh.C:1749
virtual void delete_elem(Elem *e)=0
void clearQuadratureNodes()
Clear out any existing quadrature nodes.
Definition: MooseMesh.C:1680
dof_id_type id() const
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:3448
void min(const T &r, T &o, Request &req) const
const std::map< SubdomainID, Moose::CoordinateSystemType > & getCoordSystem() const
Get the map from subdomain ID to coordinate system type, e.g.
Definition: MooseMesh.C:4207
std::vector< BndNode * > _bnd_nodes
array of boundary nodes
Definition: MooseMesh.h:1533
virtual unsigned int n_nodes() const=0
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:27
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh me...
Definition: MooseMesh.C:2928
void mooseDeprecated(Args &&... args) const
Definition: MooseBase.h:310
unsigned int _from
The qp to map from.
Definition: MooseMesh.h:75
std::pair< const Node *, BoundaryID > PeriodicNodeInfo
Helper type for building periodic node maps.
Definition: MooseMesh.h:1074
std::unique_ptr< MeshBase > buildMeshBaseObject(unsigned int dim=libMesh::invalid_uint)
Method to construct a libMesh::MeshBase object that is normally set and used by the MooseMesh object ...
Definition: MooseMesh.C:2856
unsigned int getAxisymmetricRadialCoord() const
Returns the desired radial direction for RZ coordinate transformation.
Definition: MooseMesh.C:4294
const Point & min() const
virtual Elem * add_elem(Elem *e)=0
const std::set< unsigned int > & getGhostedBoundaries() const
Return a writable reference to the set of ghosted boundary IDs.
Definition: MooseMesh.C:3241
bool _construct_node_list_from_side_list
Whether or not to allow generation of nodesets from sidesets.
Definition: MooseMesh.h:1822
boundary_id_type BoundaryID
void print_info(std::ostream &os=libMesh::out, const unsigned int verbosity=0, const bool global=true) const
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
unsigned int sideWithBoundaryID(const Elem *const elem, const BoundaryID boundary_id) const
Calls BoundaryInfo::side_with_boundary_id().
Definition: MooseMesh.C:3021
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:3470
virtual const Node * nodePtr(const dof_id_type i) const
Definition: MooseMesh.C:848
virtual const Node * query_node_ptr(const dof_id_type i) const=0
virtual libMesh::EquationSystems & es() override
std::vector< dof_id_type > _max_ids
Maximum integer ID for each extra element integer.
Definition: MooseMesh.h:1842
virtual dof_id_type max_elem_id() const=0
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
void family_tree(T elem, std::vector< T > &family, bool reset=true)
libmesh_assert(ctx)
void update()
Calls buildNodeListFromSideList(), buildNodeList(), and buildBndElemList().
Definition: MooseMesh.C:619
static const dof_id_type invalid_id
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:89
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:2293
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:33
std::set< BoundaryID > _mesh_nodeset_ids
Definition: MooseMesh.h:1526
virtual bool is_remote() const
std::unique_ptr< MooseAppCoordTransform > _coord_transform
A coordinate transformation object that describes how to transform this problem&#39;s coordinate system i...
Definition: MooseMesh.h:1868
std::set< SubdomainID > getBoundaryConnectedBlocks(const BoundaryID bid) const
Get the list of subdomains associated with the given boundary.
Definition: MooseMesh.C:3521
void checkCoordinateSystems()
Performs a sanity check for every element in the mesh.
Definition: MooseMesh.C:4307
std::vector< BoundaryID > getBoundaryIDs(const libMesh::MeshBase &mesh, const std::vector< BoundaryName > &boundary_name, bool generate_unknown, const std::set< BoundaryID > &mesh_boundary_ids)
Gets the boundary IDs with their names.
std::string & subdomain_name(subdomain_id_type id)
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_refinement_side_map
Definition: MooseMesh.h:1758
const std::set< unsigned char > & elem_dimensions() const
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_coarsening_map
Definition: MooseMesh.h:1782
const bool _is_split
Whether or not we are using a (pre-)split mesh (automatically DistributedMesh)
Definition: MooseMesh.h:1586
void setCoordSystem(const std::vector< SubdomainName > &blocks, const MultiMooseEnum &coord_sys)
Set the coordinate system for the provided blocks to coord_sys.
Definition: MooseMesh.C:4081
std::unique_ptr< ConstElemPointerRange > _coarsened_elements
The elements that were just coarsened.
Definition: MooseMesh.h:1478
void allow_find_neighbors(bool allow)
void set_mesh_dimension(unsigned char d)
std::set< dof_id_type > getAllElemIDs(unsigned int elem_id_index) const
Return all the unique element IDs for an extra element integer with its index.
Definition: MooseMesh.C:1129
std::set< BoundaryID > _mesh_boundary_ids
A set of boundary IDs currently present in the mesh.
Definition: MooseMesh.h:1524
const Node * addUniqueNode(const Point &p, Real tol=1e-6)
Add a new node to the mesh.
Definition: MooseMesh.C:1576
std::vector< BndNode > _extra_bnd_nodes
Definition: MooseMesh.h:1550
MooseApp & _app
The MOOSE application this is associated with.
Definition: MooseBase.h:353
bool _moose_mesh_prepared
True if prepare has been called on the mesh.
Definition: MooseMesh.h:1472
void buildRefinementMap(const Elem &elem, libMesh::QBase &qrule, libMesh::QBase &qrule_face, int parent_side, int child, int child_side)
Build the refinement map for a given element type.
Definition: MooseMesh.C:2472
auto norm(const T &a) -> decltype(std::abs(a))
The definition of the face_info_iterator struct.
Definition: MooseMesh.h:1909
unsigned int uniformRefineLevel() const
Returns the level of uniform refinement requested (zero if AMR is disabled).
Definition: MooseMesh.C:3216
std::vector< BndElement * > _bnd_elems
array of boundary elems
Definition: MooseMesh.h:1540
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
bool _coord_system_set
Whether the coordinate system has been set.
Definition: MooseMesh.h:1871
std::vector< SubdomainName > getSubdomainNames(const std::vector< SubdomainID > &subdomain_ids) const
Get the associated subdomainNames for the subdomain ids that are passed in.
Definition: MooseMesh.C:1769
const std::vector< QpMap > & getPRefinementMap(const Elem &elem) const
Get the map describing for each volumetric quadrature point (qp) on the refined level which qp on the...
Definition: MooseMesh.C:4378
AuxiliarySystem & getAuxiliarySystem()
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:2287
bool hasSecondOrderElements()
check if the mesh has SECOND order elements
Definition: MooseMesh.C:3711
unsigned int getPatchSize() const
Getter for the patch_size parameter.
Definition: MooseMesh.C:3401
void buildRefinementAndCoarseningMaps(Assembly *assembly)
Create the refinement and coarsening maps necessary for projection of stateful material properties wh...
Definition: MooseMesh.C:2462
bool _parallel_type_overridden
Definition: MooseMesh.h:1430
const std::string & get_nodeset_name(boundary_id_type id) const
Interface for objects interacting with the PerfGraph.
MeshBase::element_iterator activeLocalElementsBegin()
Calls active_local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:3051
std::vector< Node * > _node_map
Vector of all the Nodes in the mesh for determining when to add a new point.
Definition: MooseMesh.h:1574
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:1549
const std::vector< QpMap > & getPRefinementSideMap(const Elem &elem) const
Get the map describing for each side quadrature point (qp) on the refined level which qp on the previ...
Definition: MooseMesh.C:4384
std::unordered_map< dof_id_type, ElemInfo > _elem_to_elem_info
Map connecting elems with their corresponding ElemInfo, we use the element ID as the key...
Definition: MooseMesh.h:1596
libMesh::Node * _node
pointer to the node
Definition: BndNode.h:24
Node * getQuadratureNode(const Elem *elem, const unsigned short int side, const unsigned int qp)
Get a specified quadrature node.
Definition: MooseMesh.C:1662
void printInfo(std::ostream &os=libMesh::out, const unsigned int verbosity=0) const
Calls print_info() on the underlying Mesh.
Definition: MooseMesh.C:3462
std::string & sideset_name(boundary_id_type id)
virtual dof_id_type nNodes() const
Calls n_nodes/elem() on the underlying libMesh mesh object.
Definition: MooseMesh.C:3075
std::pair< T, U > ResultItem
Definition: KDTree.h:24
static MooseEnum partitioning()
returns MooseMesh partitioning options so other classes can use it
Definition: MooseMesh.C:3898
const std::set< boundary_id_type > & get_boundary_ids() const
const std::vector< std::vector< dof_id_type > > & dofIndices() const
Definition: ElemInfo.h:39
virtual Node * add_node(Node *n)=0
virtual const Elem * elem_ptr(const dof_id_type i) const=0
void buildCoarseningMap(const Elem &elem, libMesh::QBase &qrule, libMesh::QBase &qrule_face, int input_side)
Build the coarsening map for a given element type.
Definition: MooseMesh.C:2551
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was set by the user.
virtual unsigned int n_sides() const=0
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:943
std::unordered_map< dof_id_type, std::set< dof_id_type > > getElemIDMapping(const std::string &from_id_name, const std::string &to_id_name) const
Definition: MooseMesh.C:1100
std::map< std::pair< int, libMesh::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 Map key:
Definition: MooseMesh.h:1779
void setBoundaryName(BoundaryID boundary_id, BoundaryName name)
This method sets the boundary name of the boundary based on the id parameter.
Definition: MooseMesh.C:1780
const Elem * neighbor_ptr(unsigned int i) const
void remove_side(const Elem *elem, const unsigned short int side)
Moose::PatchUpdateType _patch_update_strategy
The patch update strategy.
Definition: MooseMesh.h:1571
infix_ostream_iterator< T, charT, traits > & operator=(T const &item)
Definition: InfixIterator.h:46
Physical unit management class with runtime unit string parsing, unit checking, unit conversion...
Definition: Units.h:32
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:2217
bool absolute_fuzzy_equals(const T &var1, const T2 &var2, const Real tol=TOLERANCE *TOLERANCE)
unsigned int level() const
void setCurrentSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:385
std::set< BoundaryID > _mesh_sideset_ids
Definition: MooseMesh.h:1525
void setGeneralAxisymmetricCoordAxes(const std::vector< SubdomainName > &blocks, const std::vector< std::pair< Point, RealVectorValue >> &axes)
Sets the general coordinate axes for axisymmetric blocks.
Definition: MooseMesh.C:4221
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:3201
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _partitioner_overridden
Definition: MooseMesh.h:1437
bool detectOrthogonalDimRanges(Real tol=1e-6)
This routine determines whether the Mesh is a regular orthogonal mesh (i.e.
Definition: MooseMesh.C:1925
virtual const Elem * query_elem_ptr(const dof_id_type i) const=0
void updateActiveSemiLocalNodeRange(std::set< dof_id_type > &ghosted_elems)
Clears the "semi-local" node list and rebuilds it.
Definition: MooseMesh.C:951
subdomain_id_type subdomain_id() const
std::map< const Elem *, std::vector< const Elem * > > _coarsened_element_children
Map of Parent elements to children elements for elements that were just coarsened.
CoordinateSystemType
Definition: MooseTypes.h:810
std::set< dof_id_type > getElemIDsOnBlocks(unsigned int elem_id_index, const std::set< SubdomainID > &blks) const
Return all the unique element IDs for an extra element integer with its index on a set of subdomains...
Definition: MooseMesh.C:1139
unsigned int getBlocksMaxDimension(const std::vector< SubdomainName > &blocks) const
Returns the maximum element dimension on the given blocks.
Definition: MooseMesh.C:2949
std::unique_ptr< libMesh::MeshBase > _mesh
Pointer to underlying libMesh mesh object.
Definition: MooseMesh.h:1433
void max(const T &r, T &o, Request &req) const
libMesh::NodeRange * getActiveNodeRange()
Definition: MooseMesh.C:1252
virtual unsigned short dim() const=0
infix_ostream_iterator< T, charT, traits > & operator*()
Definition: InfixIterator.h:56
const Node * node_ptr(const unsigned int i) const
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:3235
void freeBndNodes()
Definition: MooseMesh.C:376
The definition of the elem_info_iterator struct.
Definition: MooseMesh.h:1953
Real dimensionWidth(unsigned int component) const
Returns the width of the requested dimension.
Definition: MooseMesh.C:2193
PatchUpdateType
Type of patch update strategy for modeling node-face constraints or contact.
Definition: MooseTypes.h:953
bool is_my_variable(unsigned int var_num) const
bool _skip_deletion_repartition_after_refine
Whether or not skip remote deletion and repartition after uniform refinements.
Definition: MooseMesh.h:1463
std::set< SubdomainID > getBoundaryConnectedSecondaryBlocks(const BoundaryID bid) const
Get the list of subdomains associated with the given boundary of its secondary side.
Definition: MooseMesh.C:3532
void add_side(const dof_id_type elem, const unsigned short int side, const boundary_id_type id)
bool _need_delete
Whether we need to delete remote elements after init&#39;ing the EquationSystems.
Definition: MooseMesh.h:1825
const std::string & get_sideset_name(boundary_id_type id) const
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_refinement_map
Definition: MooseMesh.h:1756
void setAxisymmetricCoordAxis(const MooseEnum &rz_coord_axis)
For axisymmetric simulations, set the symmetry coordinate axis.
Definition: MooseMesh.C:4213
std::set< BoundaryID > getSubdomainInterfaceBoundaryIds(const SubdomainID subdomain_id) const
Get the list of boundaries that contact the given subdomain.
Definition: MooseMesh.C:3508
std::vector< const Elem * > _refined_elements
The elements that were just refined.
void mooseWarning(Args &&... args) const
Emits a warning prefixed with object name and type.
Definition: MooseBase.h:295
IntRange< T > make_range(T beg, T end)
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
Definition: MooseMesh.C:820
infix_ostream_iterator< T, charT, traits > & operator++()
Definition: InfixIterator.h:57
libMesh::BoundingBox getInflatedProcessorBoundingBox(Real inflation_multiplier=0.01) const
Get a (slightly inflated) processor bounding box.
Definition: MooseMesh.C:3419
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:1626
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type and optionally a file path to the top-level block p...
Definition: MooseBase.h:267
unsigned int mesh_dimension() const
void setMeshBase(std::unique_ptr< MeshBase > mesh_base)
Method to set the mesh_base object.
Definition: MooseMesh.C:2868
void buildBndElemList()
Definition: MooseMesh.C:1155
std::set< SubdomainID > getInterfaceConnectedBlocks(const BoundaryID bid) const
Get the list of subdomains contacting the given boundary.
Definition: MooseMesh.C:3543
SolverSystem & getSolverSystem(unsigned int sys_num)
Get non-constant reference to a solver system.
std::vector< Real > _ghosted_boundaries_inflation
Definition: MooseMesh.h:1559
virtual SimpleRange< element_iterator > active_subdomain_set_elements_ptr_range(std::set< subdomain_id_type > ss)=0
std::vector< std::unordered_map< SubdomainID, std::set< dof_id_type > > > _block_id_mapping
Unique element integer IDs for each subdomain and each extra element integers.
Definition: MooseMesh.h:1840
IntRange< unsigned short > node_index_range() const
void findAdaptivityQpMaps(const Elem *template_elem, libMesh::QBase &qrule, libMesh::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:2621
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
unsigned int _patch_size
The number of nodes to consider in the NearestNode neighborhood.
Definition: MooseMesh.h:1562
elem_info_iterator ownedElemInfoEnd()
Definition: MooseMesh.C:1536
void buildPeriodicNodeSets(std::map< BoundaryID, std::set< dof_id_type >> &periodic_node_sets, unsigned int var_number, libMesh::PeriodicBoundaries *pbs) const
This routine builds a datastructure of node ids organized by periodic boundary ids.
Definition: MooseMesh.C:1898
virtual std::unique_ptr< libMesh::PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libMesh mesh (default)...
Definition: MooseMesh.C:3733
void addGhostedBoundary(BoundaryID boundary_id)
This will add the boundary ids to be ghosted to this processor.
Definition: MooseMesh.C:3229
virtual dof_id_type maxNodeId() const
Calls max_node/elem_id() on the underlying libMesh mesh object.
Definition: MooseMesh.C:3087
virtual Elem * elem(const dof_id_type i)
Various accessors (pointers/references) for Elem "i".
Definition: MooseMesh.C:3099
libMesh::StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > * getBoundaryElementRange()
Definition: MooseMesh.C:1303
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseBase.h:195
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
face_info_iterator ownedFaceInfoEnd()
Definition: MooseMesh.C:1519
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type...
std::unordered_map< SubdomainID, std::pair< Point, RealVectorValue > > _subdomain_id_to_rz_coord_axis
Map of subdomain ID to general axisymmetric axis.
Definition: MooseMesh.h:1864
Class used for caching additional information for elements such as the volume and centroid...
Definition: ElemInfo.h:25
void paramWarning(const std::string &param, Args... args) const
Emits a warning prefixed with the file and line number of the given param (from the input file) along...
Definition: MooseBase.h:442
MeshBase::node_iterator localNodesEnd()
Definition: MooseMesh.C:3033
const RealVectorValue & getNormalByBoundaryID(BoundaryID id) const
Returns the normal vector associated with a given BoundaryID.
Definition: MooseMesh.C:2812
std::set< SubdomainID > _mesh_subdomains
A set of subdomain IDs currently present in the mesh.
Definition: MooseMesh.h:1516
ConstElemPointerRange * coarsenedElementRange() const
Return a range that is suitable for threaded execution over elements that were just coarsened...
Definition: MooseMesh.C:937
static InputParameters validParams()
Definition: MooseObject.C:25
const Moose::PatchUpdateType & getPatchUpdateStrategy() const
Get the current patch update strategy.
Definition: MooseMesh.C:3413
bool doingPRefinement() const
Query whether we have p-refinement.
Definition: MooseMesh.h:1352
std::vector< std::pair< BoundaryID, BoundaryID > > _paired_boundary
A vector holding the paired boundaries for a regular orthogonal mesh.
Definition: MooseMesh.h:1583
const Point & max() const
void ghostGhostedBoundaries()
Actually do the ghosting of boundaries that need to be ghosted to this processor. ...
Definition: MooseMesh.C:3314
std::set< unsigned int > _ghosted_boundaries
Definition: MooseMesh.h:1558
MeshBase::node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:3027
unsigned int _rz_coord_axis
Storage for RZ axis selection.
Definition: MooseMesh.h:1861
void computeFiniteVolumeCoords() const
Compute the face coordinate value for all FaceInfo and ElemInfo objects.
Definition: MooseMesh.C:3877
void buildNodeList()
Calls BoundaryInfo::build_node_list()/build_side_list() and makes separate copies of Nodes/Elems in t...
Definition: MooseMesh.C:1025
virtual dof_id_type max_node_id() const=0
std::unordered_map< SubdomainID, SubdomainData > _sub_to_data
Holds a map from subdomain ids to associated data.
Definition: MooseMesh.h:1800
std::unique_ptr< libMesh::Partitioner > _custom_partitioner
The custom partitioner.
Definition: MooseMesh.h:1440
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:3591
virtual dof_id_type n_elem() const=0
virtual const Node * node_ptr(const dof_id_type i) const=0
std::map< const Elem *, std::vector< const Elem * > > _coarsened_element_children
Map of Parent elements to child elements for elements that were just coarsened.
Definition: MooseMesh.h:1485
processor_id_type processor_id() const
const std::vector< QpMap > & getPCoarseningMap(const Elem &elem) const
Get the map describing for each volumetric quadrature point (qp) on the coarse level which qp on the ...
Definition: MooseMesh.C:4390
bool isRecovering() const
Whether or not this is a "recover" calculation.
Definition: MooseApp.C:1801
auto min(const L &left, const R &right)
virtual std::size_t numSolverSystems() const override
SearchParams SearchParameters
bool active() const
std::vector< const Elem * > _coarsened_elements
The elements that were just coarsened.
void buildPeriodicNodeMap(std::multimap< dof_id_type, dof_id_type > &periodic_node_map, unsigned int var_number, libMesh::PeriodicBoundaries *pbs) const
This routine builds a multimap of boundary ids to matching boundary ids across all periodic boundarie...
Definition: MooseMesh.C:1813
processor_id_type processor_id() const
std::string getRestartRecoverFileBase() const
The file_base for the recovery file.
Definition: MooseApp.h:487
virtual ElemType type() const=0
bool isParamSetByUser(const std::string &name) const
Test if the supplied parameter is set by a user, as opposed to not set or set to default.
Definition: MooseBase.h:201
dof_id_type node_id(const unsigned int i) const
const Point & point(const unsigned int i) const
uint8_t unique_id_type
const std::unordered_map< boundary_id_type, std::unordered_set< dof_id_type > > & getBoundariesToActiveSemiLocalElemIds() const
Returns a map of boundaries to ids of elements on the boundary.
Definition: MooseMesh.C:1325
const MeshBase::element_iterator activeLocalElementsEnd()
Definition: MooseMesh.C:3057
bool relative_fuzzy_equals(const TypeVector< Real > &rhs, Real tol=TOLERANCE) const
std::unique_ptr< libMesh::ConstNodeRange > _local_node_range
Definition: MooseMesh.h:1498
libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
Definition: MooseMesh.C:1289
void ErrorVector unsigned int
auto index_range(const T &sizable)
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
Definition: MooseMesh.h:1577
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:1509
virtual dof_id_type n_nodes() const=0
bool prepare(const MeshBase *mesh_to_clone)
Calls prepare_for_use() if the underlying MeshBase object isn&#39;t prepared, then communicates various b...
Definition: MooseMesh.C:409
const ElemInfo & elemInfo(const dof_id_type id) const
Accessor for the elemInfo object for a given element ID.
Definition: MooseMesh.C:3871
void setCustomPartitioner(libMesh::Partitioner *partitioner)
Setter for custom partitioner.
Definition: MooseMesh.C:3699
Real _distance
The distance between them.
Definition: MooseMesh.h:81
SubdomainID subdomain_id() const
We return the subdomain ID of the corresponding libmesh element.
Definition: ElemInfo.h:43
dof_id_type get_extra_integer(const unsigned int index) const
const Elem * child_ptr(unsigned int i) const
static MooseEnum elemTypes()
returns MooseMesh element type options
Definition: MooseMesh.C:3906
void meshChanged()
Declares that the MooseMesh has changed, invalidates cached data and rebuilds caches.
Definition: MooseMesh.C:882
void buildPRefinementAndCoarseningMaps(Assembly *assembly)
Definition: MooseMesh.C:2368
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
Definition: MooseMesh.C:1692
std::unique_ptr< libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > > _bnd_node_range
Definition: MooseMesh.h:1500
uint8_t dof_id_type
virtual dof_id_type nElem() const
Definition: MooseMesh.C:3081
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:1178
libMesh::QBase *const & writeableQRuleFace()
Returns the reference to the current quadrature being used on a current face.
Definition: Assembly.h:300
bool isBoundaryFullyExternalToSubdomains(BoundaryID bid, const std::set< SubdomainID > &blk_group) const
Returns whether a boundary (given by its id) is not crossing through a group of blocks, by which we mean that elements on both sides of the boundary are in those blocks.
Definition: MooseMesh.C:1373
std::unique_ptr< SemiLocalNodeRange > _active_semilocal_node_range
Definition: MooseMesh.h:1496
const std::set< SubdomainID > & meshSubdomains() const
Returns a read-only reference to the set of subdomains currently present in the Mesh.
Definition: MooseMesh.C:3171
SubdomainID getSubdomainID(const SubdomainName &subdomain_name) const
Get the associated subdomain ID for the subdomain name.
Definition: MooseMesh.C:1731
unsigned int _max_p_level
Maximum p-refinement level of all elements.
Definition: MooseMesh.h:1879
void setupFiniteVolumeMeshData() const
Sets up the additional data needed for finite volume computations.
Definition: MooseMesh.C:4072
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.
void set_union(T &data, const unsigned int root_id) const
const RemoteElem * remote_elem
virtual unsigned int effectiveSpatialDimension() const
Returns the effective spatial dimension determined by the coordinates actually used by the mesh...
Definition: MooseMesh.C:2934