24 #include <unordered_map> 25 #include <unordered_set> 28 #include "libmesh/elem_range.h" 29 #include "libmesh/mesh_base.h" 30 #include "libmesh/replicated_mesh.h" 31 #include "libmesh/distributed_mesh.h" 32 #include "libmesh/node_range.h" 33 #include "libmesh/nanoflann.hpp" 34 #include "libmesh/vector_value.h" 35 #include "libmesh/point.h" 36 #include "libmesh/partitioner.h" 49 class PeriodicBoundaries;
51 class GhostingFunctor;
59 "QUAD QUAD4 QUAD8 QUAD9 " 61 "HEX HEX8 HEX20 HEX27 " 62 "TET TET4 TET10 TET14 " 63 "PRISM PRISM6 PRISM15 PRISM18 " 64 "PYRAMID PYRAMID5 PYRAMID13 PYRAMID14";
121 virtual std::unique_ptr<MooseMesh>
safeClone()
const = 0;
142 template <
typename T>
149 void setMeshBase(std::unique_ptr<MeshBase> mesh_base);
201 const unsigned short int side)
const;
209 const Elem *
getLowerDElem(
const Elem *,
unsigned short int)
const;
238 const std::map<dof_id_type, std::vector<dof_id_type>> &
nodeToElemMap();
253 struct bnd_node_iterator;
254 struct const_bnd_node_iterator;
256 struct bnd_elem_iterator;
257 struct const_bnd_elem_iterator;
283 std::vector<unsigned short int> & sl,
284 std::vector<boundary_id_type> & il);
288 std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>>
buildSideList();
294 std::vector<std::tuple<dof_id_type, unsigned short int, boundary_id_type>>
455 const std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>> &
461 const std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>> &
490 const std::set<SubdomainID> & blk_group)
const;
553 bool prepare(
const MeshBase * mesh_to_clone);
662 const MeshBase &
getMesh()
const;
663 const MeshBase &
getMesh(
const std::string &
name)
const;
687 const Node *
addUniqueNode(
const Point & p, Real tol = 1e-6);
703 const unsigned short int side,
704 const unsigned int qp,
706 const Point & point);
737 std::vector<BoundaryID>
getBoundaryIDs(
const std::vector<BoundaryName> & boundary_name,
738 bool generate_unknown =
false)
const;
754 std::vector<SubdomainID>
755 getSubdomainIDs(
const std::vector<SubdomainName> & subdomain_names)
const;
756 std::set<SubdomainID>
getSubdomainIDs(
const std::set<SubdomainName> & subdomain_names)
const;
780 std::vector<SubdomainName>
799 unsigned int var_number,
806 unsigned int var_number,
890 const std::vector<std::vector<QpMap>> &
1038 bool use_distributed_mesh,
1065 virtual std::unique_ptr<libMesh::PointLocatorBase>
getPointLocator()
const;
1129 bool areElemIDsIdentical(
const std::string & id_name1,
const std::string & id_name2)
const;
1134 std::set<dof_id_type>
getAllElemIDs(
unsigned int elem_id_index)
const;
1141 const std::set<SubdomainID> & blks)
const;
1143 std::unordered_map<dof_id_type, std::set<dof_id_type>>
1144 getElemIDMapping(
const std::string & from_id_name,
const std::string & to_id_name)
const;
1150 const std::vector<const FaceInfo *> &
faceInfo()
const;
1153 struct face_info_iterator;
1154 struct const_face_info_iterator;
1164 struct elem_info_iterator;
1165 struct const_elem_info_iterator;
1184 const std::vector<FaceInfo> &
allFaceInfo()
const;
1218 const std::map<boundary_id_type, std::vector<dof_id_type>> &
nodeSetNodes()
const;
1235 const std::map<SubdomainID, Moose::CoordinateSystemType> &
getCoordSystem()
const;
1263 const std::vector<std::pair<Point, RealVectorValue>> & axes);
1270 const std::pair<Point, RealVectorValue> &
1321 const std::unordered_map<std::pair<const Elem *, unsigned short int>,
const Elem *> &
1433 std::unique_ptr<libMesh::MeshBase>
_mesh;
1499 std::unique_ptr<libMesh::StoredRange<MooseMesh::const_bnd_node_iterator, const BndNode *>>
1501 std::unique_ptr<libMesh::StoredRange<MooseMesh::const_bnd_elem_iterator, const BndElement *>>
1545 std::unordered_map<boundary_id_type, std::unordered_set<dof_id_type>>
_bnd_elem_ids;
1548 std::map<dof_id_type, std::map<unsigned int, std::map<dof_id_type, Node *>>>
1611 mutable std::unordered_map<std::pair<const Elem *, unsigned int>,
FaceInfo *>
1683 void mapPoints(
const std::vector<Point> & from,
1684 const std::vector<Point> & to,
1685 std::vector<QpMap> & qp_map);
1712 std::vector<std::vector<QpMap>> & refinement_map,
1713 std::vector<std::pair<unsigned int, QpMap>> & coarsen_map,
1722 const std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>> &)
const;
1725 const std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>> &)
const;
1752 std::map<std::pair<int, libMesh::ElemType>, std::vector<std::vector<QpMap>>>
1755 std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>>
1757 std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>>
1761 std::map<libMesh::ElemType, std::map<std::pair<int, int>, std::vector<std::vector<QpMap>>>>
1778 std::map<std::pair<int, libMesh::ElemType>, std::vector<std::pair<unsigned int, QpMap>>>
1781 std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>>
1783 std::map<std::pair<libMesh::ElemType, unsigned int>, std::vector<QpMap>>
1810 std::unordered_map<std::pair<const Elem *, unsigned short int>,
const Elem *>
1883 template <
typename T>
1890 mooseAssert(
_coord_transform,
"The coordinate transformation object is null.");
1913 template <
typename PredType,
typename IterType>
1925 const FaceInfo * const,
1926 const FaceInfo * const &,
1927 const FaceInfo * const *>
1930 template <
typename PredType,
typename IterType>
1957 template <
typename PredType,
typename IterType>
1969 const ElemInfo * const,
1970 const ElemInfo * const &,
1971 const ElemInfo * const *>
1974 template <
typename PredType,
typename IterType>
2000 template <
typename PredType,
typename IterType>
2017 template <
typename PredType,
typename IterType>
2043 template <
typename PredType,
typename IterType>
2056 BndElement * const &,
2057 BndElement * const *>
2060 template <
typename PredType,
typename IterType>
2088 template <
typename T>
2098 mooseError(
"A MooseMesh object is being asked to build a libMesh mesh that is a different " 2099 "parallel type than the libMesh mesh that it wraps. This is not allowed. Please " 2100 "create another MooseMesh object to wrap the new libMesh mesh");
2107 dim = getParam<MooseEnum>(
"dim");
2115 if (!getParam<bool>(
"allow_renumbering"))
2116 mesh->allow_renumbering(
false);
2124 if (!
parameters().isParamSetByAddParam(
"partitioner"))
2125 mooseError(
"If partitioner block is provided, partitioner keyword cannot be used!");
2128 mooseError(
"Custom partitioner requested but not set!");
2153 return getMesh().has_elem_integer(id_name);
2160 mooseError(
"Mesh does not have element ID for ", id_name);
2161 return getMesh().get_elem_integer_index(id_name);
2172 inline const std::vector<const FaceInfo *> &
2178 inline const std::vector<FaceInfo> &
2184 inline const std::map<boundary_id_type, std::vector<dof_id_type>> &
2190 inline const std::unordered_map<std::pair<const Elem *, unsigned short int>,
const Elem *> &
2199 return libmesh_map_find(
_sub_to_data, subdomain_id).is_lower_d;
ParallelType _parallel_type
Can be set to DISTRIBUTED, REPLICATED, or DEFAULT.
static InputParameters validParams()
Typical "Moose-style" constructor and copy constructor.
virtual std::string getFileName() const
Returns the name of the mesh file read to produce this mesh if any or an empty string otherwise...
virtual bnd_node_iterator bndNodesEnd()
virtual bnd_elem_iterator bndElemsEnd()
std::vector< std::vector< Real > > _bounds
The bounds in each dimension of the mesh for regular orthogonal meshes.
const std::map< boundary_id_type, std::vector< dof_id_type > > & nodeSetNodes() const
std::set< Node * > _semilocal_node_list
Used for generating the semilocal node range.
std::map< dof_id_type, Node * > _quadrature_nodes
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...
virtual Real getMaxInDimension(unsigned int component) const
virtual unsigned int spatialDimension() const
Returns MeshBase::spatial_dimension.
bool _node_to_elem_map_built
std::unique_ptr< libMesh::NodeRange > _active_node_range
std::vector< Node * > _extreme_nodes
A vector containing the nodes at the corners of a regular orthogonal mesh.
Node * addQuadratureNode(const Elem *elem, const unsigned short int side, const unsigned int qp, BoundaryID bid, const Point &point)
Adds a fictitious "QuadratureNode".
std::vector< std::tuple< dof_id_type, unsigned short int, boundary_id_type > > buildSideList()
As above, but uses the non-deprecated std::tuple interface.
const_bnd_elem_iterator(const bnd_elem_iterator &rhs)
bool isFiniteVolumeInfoDirty() const
const std::set< BoundaryID > & meshNodesetIds() const
Returns a read-only reference to the set of nodesets currently present in the Mesh.
void buildElemIDInfo()
Build extra data for faster access to the information of extra element integers.
std::vector< FaceInfo > _all_face_info
FaceInfo object storing information for face based loops.
std::vector< const FaceInfo * > _face_info
Holds only those FaceInfo objects that have processor_id equal to this process's id, e.g.
bool allowRemoteElementRemoval() const
Whether we are allow remote element removal.
bool getConstructNodeListFromSideList()
Return construct node list from side list boolean.
virtual Real getMinInDimension(unsigned int component) const
Returns the min or max of the requested dimension respectively.
libMesh::ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
virtual void onMeshChanged()
Declares a callback function that is executed at the conclusion of meshChanged(). ...
bool prepared() const
Setter/getter for whether the mesh is prepared.
void needsPrepareForUse()
If this method is called, we will call libMesh's prepare_for_use method when we call Moose's prepare ...
const std::set< BoundaryID > & getBoundaryIDs() const
Returns a const reference to a set of all user-specified boundary IDs.
bool _is_nemesis
True if a Nemesis Mesh was read in.
std::vector< SubdomainName > _provided_coord_blocks
Set for holding user-provided coordinate system type block names.
virtual MooseMesh & clone() const
Clone method.
void allowRecovery(bool allow)
Set whether or not this mesh is allowed to read a recovery file.
const std::vector< QpMap > & getPCoarseningMapHelper(const Elem &elem, const std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap >> &) const
const std::string LIST_GEOM_ELEM
bool isCustomPartitionerRequested() const
Setter and getter for _custom_partitioner_requested.
bool _need_ghost_ghosted_boundaries
A parallel mesh generator such as DistributedRectilinearMeshGenerator already make everything ready...
A class for creating restricted objects.
bool isDisplaced() const
whether this mesh is a displaced mesh
unsigned int _uniform_refine_level
The level of uniform refinement requested (set to zero if AMR is disabled)
const std::vector< QpMap > & getPRefinementMapHelper(const Elem &elem, const std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap >> &) const
std::vector< dof_id_type > _min_ids
Minimum integer ID for each extra element integer.
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.
const std::set< SubdomainID > & interiorLowerDBlocks() const
const MooseUnits & lengthUnit() const
void checkDuplicateSubdomainNames()
Loop through all subdomain IDs and check if there is name duplication used for the subdomains with sa...
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.
unsigned int _max_leaf_size
The definition of the const_bnd_elem_iterator struct.
RealVectorValue _half_range
A convenience vector used to hold values in each dimension representing half of the range...
Keeps track of stuff related to assembling.
const_bnd_elem_iterator(const IterType &d, const IterType &e, const PredType &p)
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).
void setCoordData(const MooseMesh &other_mesh)
Set the coordinate system data to that of other_mesh.
const std::unordered_map< std::pair< const Elem *, unsigned short int >, const Elem * > & getLowerDElemMap() const
This function attempts to return the map from a high-order element side to its corresponding lower-d ...
bool _finite_volume_info_dirty
virtual Elem * elemPtr(const dof_id_type i)
elem_info_iterator(const IterType &d, const IterType &e, const PredType &p)
const_bnd_node_iterator(const MooseMesh::bnd_node_iterator &rhs)
The definition of the bnd_elem_iterator struct.
std::map< SubdomainID, Moose::CoordinateSystemType > & _coord_sys
Type of coordinate system per subdomain.
const std::vector< const ElemInfo * > & elemInfoVector() const
Accessor for the element info objects owned by this process.
bool isBoundaryNode(dof_id_type node_id) const
Returns true if the requested node is in the list of boundary nodes, false otherwise.
face_info_iterator ownedFaceInfoBegin()
Iterators to owned faceInfo objects.
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 ...
std::vector< std::tuple< dof_id_type, unsigned short int, boundary_id_type > > buildActiveSideList() const
Calls BoundaryInfo::build_active_side_list.
static void setPartitioner(MeshBase &mesh_base, MooseEnum &partitioner, bool use_distributed_mesh, const InputParameters ¶ms, MooseObject &context_obj)
Method for setting the partitioner on the passed in mesh_base object.
void buildLowerDMesh()
Build lower-d mesh for all sides.
const std::set< BoundaryID > & meshSidesetIds() const
Returns a read-only reference to the set of sidesets currently present in the Mesh.
void setParallelType(ParallelType parallel_type)
Allow to change parallel type.
std::unordered_map< const Elem *, unsigned short int > _lower_d_elem_to_higher_d_elem_side
void cacheFVElementalDoFs() const
Cache the DoF indices for FV variables on each element.
void cacheFaceInfoVariableOwnership() const
Cache if variables live on the elements connected by the FaceInfo objects.
bool _custom_partitioner_requested
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.
Moose::CoordinateSystemType getUniqueCoordSystem() const
Get the coordinate system from the mesh, it must be the same in all subdomains otherwise this will er...
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.
const_face_info_iterator(const MooseMesh::face_info_iterator &rhs)
const_bnd_node_iterator(const IterType &d, const IterType &e, const PredType &p)
bool isSemiLocal(Node *const node) const
Returns true if the node is semi-local.
libMesh::StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > ConstBndElemRange
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...
std::unique_ptr< libMesh::StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > > _bnd_elem_range
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...
const std::vector< std::pair< unsigned int, QpMap > > & getCoarseningMap(const Elem &elem, int input_side)
Get the coarsening map for a given element type.
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.
bool _is_changed
true if mesh is changed (i.e. after adaptivity step)
bool _doing_p_refinement
Whether we have p-refinement (as opposed to h-refinement)
void determineUseDistributedMesh()
Determine whether to use a distributed mesh.
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.
const std::string & getBoundaryName(BoundaryID boundary_id)
Return the name of the boundary given the id.
void cacheChangedLists()
Cache information about what elements were refined and coarsened in the previous step.
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...
std::set< SubdomainID > neighbor_subs
Neighboring subdomain ids.
unsigned int getGhostingPatchSize() const
Getter for the ghosting_patch_size parameter.
std::vector< BndNode * >::iterator bnd_node_iterator_imp
const_face_info_iterator(const IterType &d, const IterType &e, const PredType &p)
The definition of the bnd_node_iterator struct.
face_info_iterator(const IterType &d, const IterType &e, const PredType &p)
std::vector< const ElemInfo * > _elem_info
Holds only those ElemInfo objects that have processor_id equal to this process's id, e.g.
void buildHRefinementAndCoarseningMaps(Assembly *assembly)
bool usingGeneralAxisymmetricCoordAxes() const
Returns true if general axisymmetric coordinate axes are being used.
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:
const std::set< SubdomainID > & getBlockConnectedBlocks(const SubdomainID subdomain_id) const
Get the list of subdomains neighboring a given subdomain.
virtual const Node * queryNodePtr(const dof_id_type i) const
libMesh::StoredRange< std::set< Node * >::iterator, Node * > SemiLocalNodeRange
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.
Helper object for holding qp mapping info.
bool _has_lower_d
Whether there are any lower-dimensional blocks that are manifolds of higher-dimensional block faces...
const std::vector< Real > & getGhostedBoundaryInflation() const
Return a writable reference to the _ghosted_boundaries_inflation vector.
std::unique_ptr< ConstElemPointerRange > _refined_elements
The elements that were just refined.
dof_id_type maxElementID(unsigned int elem_id_index) const
Return the maximum element ID for an extra element integer with its accessing index.
std::vector< std::vector< bool > > _id_identical_flag
Flags to indicate whether or not any two extra element integers are the same.
virtual dof_id_type maxElemId() const
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
virtual bnd_elem_iterator bndElemsBegin()
Return iterators to the beginning/end of the boundary elements list.
The definition of the const_bnd_node_iterator struct.
void setUniformRefineLevel(unsigned int, bool deletion=true)
Set uniform refinement level.
MooseEnum _partitioner_name
The partitioner used on this mesh.
std::unique_ptr< libMesh::ConstElemRange > _active_local_elem_range
A range for use with threading.
std::map< dof_id_type, std::set< SubdomainID > > _block_node_list
list of nodes that belongs to a specified block (domain)
bool areElemIDsIdentical(const std::string &id_name1, const std::string &id_name2) const
Whether or not two extra element integers are identical.
bool _node_to_active_semilocal_elem_map_built
std::map< boundary_id_type, std::set< dof_id_type > > _bnd_node_ids
Map of sets of node IDs in each boundary.
ConstElemPointerRange * refinedElementRange() const
Return a range that is suitable for threaded execution over elements that were just refined...
virtual void init()
Initialize the Mesh object.
The definition of the const_face_info_iterator struct.
bool _skip_refine_when_use_split
Whether or not to skip uniform refinements when using a pre-split mesh.
unsigned int getHigherDSide(const Elem *elem) const
Returns the local side ID of the interior parent aligned with the lower dimensional element...
std::unordered_map< std::pair< const Elem *, unsigned int >, FaceInfo * > _elem_side_to_face_info
Map from elem-side pair to FaceInfo.
void detectPairedSidesets()
This routine detects paired sidesets of a regular orthogonal mesh (.i.e.
const std::set< SubdomainID > & getNodeBlockIds(const Node &node) const
Return list of blocks to which the given node belongs.
const Parallel::Communicator & _communicator
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_coarsening_side_map
void setMeshBoundaryIDs(std::set< BoundaryID > boundary_IDs)
Sets the set of BoundaryIDs Is called by AddAllSideSetsByNormals.
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] ...
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.
bool _linear_finite_volume_dofs_cached
void markFiniteVolumeInfoDirty()
Mark the finite volume information as dirty.
std::basic_ostream< charT, traits > * os
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.
std::vector< std::unique_ptr< libMesh::GhostingFunctor > > _ghosting_functors
Deprecated (DO NOT USE)
std::set< Elem * > _ghost_elems_from_ghost_boundaries
Set of elements ghosted by ghostGhostedBoundaries.
virtual void buildMesh()=0
Must be overridden by child classes.
void setPartitionerHelper(MeshBase *mesh=nullptr)
void deleteRemoteElements()
Delete remote elements.
unsigned int _to
The qp to map to.
libMesh::ConstNodeRange * getLocalNodeRange()
virtual const Node & nodeRef(const dof_id_type i) const
bool _allow_recovery
Whether or not this Mesh is allowed to read a recovery file.
virtual const std::string & name() const
Get the name of the class.
std::vector< SubdomainID > getSubdomainIDs(const std::vector< SubdomainName > &subdomain_names) const
Get the associated subdomainIDs for the subdomain names that are passed in.
void buildNodeListFromSideList()
Calls BoundaryInfo::build_node_list_from_side_list().
const std::string & getSubdomainName(SubdomainID subdomain_id) const
Return the name of a block given an id.
void setPatchUpdateStrategy(Moose::PatchUpdateType patch_update_strategy)
Set the patch size update strategy.
void buildFiniteVolumeInfo() const
Builds the face and elem info vectors that store meta-data needed for looping over and doing calculat...
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.
const std::pair< Point, RealVectorValue > & getGeneralAxisymmetricCoordAxis(SubdomainID subdomain_id) const
Gets the general axisymmetric coordinate axis for a block.
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...
std::set< SubdomainID > _lower_d_interior_blocks
Mesh blocks for interior lower-d elements in different types.
virtual Elem * queryElemPtr(const dof_id_type i)
bool isLowerD(const SubdomainID subdomain_id) const
elem_info_iterator ownedElemInfoBegin()
Iterators to owned faceInfo objects.
void setIsCustomPartitionerRequested(bool cpr)
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
std::set< BoundaryID > boundary_ids
The boundary ids that are attached.
unsigned int _max_h_level
Maximum h-refinement level of all elements.
virtual bnd_node_iterator bndNodesBegin()
Return iterators to the beginning/end of the boundary nodes list.
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".
bnd_elem_iterator(const IterType &d, const IterType &e, const PredType &p)
bool is_lower_d
Whether this subdomain is a lower-dimensional manifold of a higher-dimensional subdomain.
void attachRelationshipManagers(Moose::RelationshipManagerType rm_type, bool attach_geometric_rm_final=false)
Attach the relationship managers of the given type Note: Geometric relationship managers that are sup...
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
const MeshBase * getMeshPtr() const
This data structure is used to store geometric and variable related metadata about each cell face in ...
void updateCoordTransform()
Update the coordinate transformation object based on our coordinate system data.
const std::vector< const FaceInfo * > & faceInfo() const
Accessor for local FaceInfo objects.
dof_id_type minElementID(unsigned int elem_id_index) const
Return the minimum element ID for an extra element integer with its accessing index.
const std::vector< FaceInfo > & allFaceInfo() const
Accessor for all FaceInfo objects.
bool _use_distributed_mesh
False by default.
std::unique_ptr< std::map< BoundaryID, RealVectorValue > > _boundary_to_normal_map
The boundary to normal map - valid only when AddAllSideSetsByNormals is active.
std::vector< BndElement * >::const_iterator const_bnd_elem_iterator_imp
bool hasElementID(const std::string &id_name) const
Whether mesh has an extra element integer with a given name.
bool _built_from_other_mesh
Whether or not this mesh was built from another mesh.
virtual dof_id_type nLocalNodes() const
bool _allow_remote_element_removal
Whether to allow removal of remote elements.
bnd_node_iterator(const IterType &d, const IterType &e, const PredType &p)
virtual bool skipNoncriticalPartitioning() const
SemiLocalNodeRange * getActiveSemiLocalNodeRange() const
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.
void setSubdomainName(SubdomainID subdomain_id, const SubdomainName &name)
This method sets the name for subdomain_id to name.
void clearQuadratureNodes()
Clear out any existing quadrature nodes.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
const std::map< SubdomainID, Moose::CoordinateSystemType > & getCoordSystem() const
Get the map from subdomain ID to coordinate system type, e.g.
std::vector< BndNode * > _bnd_nodes
array of boundary nodes
std::unique_ptr< T > buildTypedMesh(unsigned int dim=libMesh::invalid_uint)
Shortcut method to construct a unique pointer to a libMesh mesh instance.
Every object that can be built by the factory should be derived from this class.
std::vector< std::shared_ptr< RelationshipManager > > _relationship_managers
The list of active geometric relationship managers (bound to the underlying MeshBase object)...
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh me...
unsigned int _from
The qp to map from.
std::pair< const Node *, BoundaryID > PeriodicNodeInfo
Helper type for building periodic node maps.
virtual SubdomainID nSubdomains() const
const_elem_info_iterator(const IterType &d, const IterType &e, const PredType &p)
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 ...
unsigned int getAxisymmetricRadialCoord() const
Returns the desired radial direction for RZ coordinate transformation.
bool isPartitionerForced() const
Tell the user if the partitioner was overriden for any reason.
const std::set< unsigned int > & getGhostedBoundaries() const
Return a writable reference to the set of ghosted boundary IDs.
bool _construct_node_list_from_side_list
Whether or not to allow generation of nodesets from sidesets.
boundary_id_type BoundaryID
unsigned int sideWithBoundaryID(const Elem *const elem, const BoundaryID boundary_id) const
Calls BoundaryInfo::side_with_boundary_id().
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.
unsigned int getMaxLeafSize() const
Getter for the maximum leaf size parameter.
virtual const Node * nodePtr(const dof_id_type i) const
std::vector< dof_id_type > _max_ids
Maximum integer ID for each extra element integer.
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
void update()
Calls buildNodeListFromSideList(), buildNodeList(), and buildBndElemList().
const std::pair< BoundaryID, BoundaryID > * getPairedBoundaryMapping(unsigned int component)
This function attempts to return the paired boundary ids for the given component. ...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
std::set< BoundaryID > _mesh_nodeset_ids
std::unique_ptr< MooseAppCoordTransform > _coord_transform
A coordinate transformation object that describes how to transform this problem's coordinate system i...
std::set< SubdomainID > getBoundaryConnectedBlocks(const BoundaryID bid) const
Get the list of subdomains associated with the given boundary.
void checkCoordinateSystems()
Performs a sanity check for every element in the mesh.
ParallelType getParallelType() const
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_refinement_side_map
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_coarsening_map
const bool _is_split
Whether or not we are using a (pre-)split mesh (automatically DistributedMesh)
void setCoordSystem(const std::vector< SubdomainName > &blocks, const MultiMooseEnum &coord_sys)
Set the coordinate system for the provided blocks to coord_sys.
std::unique_ptr< ConstElemPointerRange > _coarsened_elements
The elements that were just coarsened.
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.
std::set< BoundaryID > _mesh_boundary_ids
A set of boundary IDs currently present in the mesh.
const Node * addUniqueNode(const Point &p, Real tol=1e-6)
Add a new node to the mesh.
std::vector< BndNode > _extra_bnd_nodes
MooseApp & _app
The MOOSE application this is associated with.
bool _moose_mesh_prepared
True if prepare has been called on the mesh.
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.
The definition of the face_info_iterator struct.
unsigned int uniformRefineLevel() const
Returns the level of uniform refinement requested (zero if AMR is disabled).
std::vector< BndElement * > _bnd_elems
array of boundary elems
bool _coord_system_set
Whether the coordinate system has been set.
const std::set< SubdomainID > & boundaryLowerDBlocks() const
std::vector< SubdomainName > getSubdomainNames(const std::vector< SubdomainID > &subdomain_ids) const
Get the associated subdomainNames for the subdomain ids that are passed in.
MooseMesh & operator=(const MooseMesh &other_mesh)=delete
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...
virtual dof_id_type nActiveLocalElem() const
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...
bool hasSecondOrderElements()
check if the mesh has SECOND order elements
unsigned int getPatchSize() const
Getter for the patch_size parameter.
void buildRefinementAndCoarseningMaps(Assembly *assembly)
Create the refinement and coarsening maps necessary for projection of stateful material properties wh...
bool _parallel_type_overridden
bool isRegularOrthogonal()
Getter to query if the mesh was detected to be regular and orthogonal.
Interface for objects interacting with the PerfGraph.
MeshBase::element_iterator activeLocalElementsBegin()
Calls active_local_nodes_begin/end() on the underlying libMesh mesh object.
std::vector< Node * > _node_map
Vector of all the Nodes in the mesh for determining when to add a new point.
std::map< dof_id_type, std::map< unsigned int, std::map< dof_id_type, Node * > > > _elem_to_side_to_qp_to_quadrature_nodes
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...
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...
bool skipRefineWhenUseSplit() const
Whether or not skip uniform refinements when using a pre-split mesh.
std::vector< BndElement * >::iterator bnd_elem_iterator_imp
Node * getQuadratureNode(const Elem *elem, const unsigned short int side, const unsigned int qp)
Get a specified quadrature node.
void printInfo(std::ostream &os=libMesh::out, const unsigned int verbosity=0) const
Calls print_info() on the underlying Mesh.
libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > ConstBndNodeRange
Some useful StoredRange typedefs.
virtual dof_id_type nNodes() const
Calls n_nodes/elem() on the underlying libMesh mesh object.
static MooseEnum partitioning()
returns MooseMesh partitioning options so other classes can use it
MooseAppCoordTransform & coordTransform()
void buildCoarseningMap(const Elem &elem, libMesh::QBase &qrule, libMesh::QBase &qrule_face, int input_side)
Build the coarsening map for a given element type.
const std::vector< const Elem * > & coarsenedElementChildren(const Elem *elem) const
Get the newly removed children element ids for an element that was just coarsened.
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
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:
void setBoundaryName(BoundaryID boundary_id, BoundaryName name)
This method sets the boundary name of the boundary based on the id parameter.
Moose::PatchUpdateType _patch_update_strategy
The patch update strategy.
Physical unit management class with runtime unit string parsing, unit checking, unit conversion...
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...
std::set< BoundaryID > _mesh_sideset_ids
void setGeneralAxisymmetricCoordAxes(const std::vector< SubdomainName > &blocks, const std::vector< std::pair< Point, RealVectorValue >> &axes)
Sets the general coordinate axes for axisymmetric blocks.
void setBoundaryToNormalMap(std::unique_ptr< std::map< BoundaryID, RealVectorValue >> boundary_map)
Sets the mapping between BoundaryID and normal vector Is called by AddAllSideSetsByNormals.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _partitioner_overridden
bool detectOrthogonalDimRanges(Real tol=1e-6)
This routine determines whether the Mesh is a regular orthogonal mesh (i.e.
RelationshipManagers are used for describing what kinds of non-local resources are needed for an obje...
void updateActiveSemiLocalNodeRange(std::set< dof_id_type > &ghosted_elems)
Clears the "semi-local" node list and rebuilds it.
void isDisplaced(bool is_displaced)
Set whether this mesh is a displaced mesh.
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...
unsigned int getBlocksMaxDimension(const std::vector< SubdomainName > &blocks) const
Returns the maximum element dimension on the given blocks.
virtual std::unique_ptr< MooseMesh > safeClone() const =0
A safer version of the clone() method that hands back an allocated object wrapped in a smart pointer...
bool _is_displaced
Whether this mesh is displaced.
bool hasMeshBase() const
Whether mesh base object was constructed or not.
std::unique_ptr< libMesh::MeshBase > _mesh
Pointer to underlying libMesh mesh object.
libMesh::NodeRange * getActiveNodeRange()
virtual unsigned int nPartitions() 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...
The definition of the elem_info_iterator struct.
Real dimensionWidth(unsigned int component) const
Returns the width of the requested dimension.
PatchUpdateType
Type of patch update strategy for modeling node-face constraints or contact.
bool _skip_deletion_repartition_after_refine
Whether or not skip remote deletion and repartition after uniform refinements.
std::set< SubdomainID > getBoundaryConnectedSecondaryBlocks(const BoundaryID bid) const
Get the list of subdomains associated with the given boundary of its secondary side.
bool _need_delete
Whether we need to delete remote elements after init'ing the EquationSystems.
std::map< std::pair< libMesh::ElemType, unsigned int >, std::vector< QpMap > > _elem_type_to_p_refinement_map
void setAxisymmetricCoordAxis(const MooseEnum &rz_coord_axis)
For axisymmetric simulations, set the symmetry coordinate axis.
std::set< BoundaryID > getSubdomainInterfaceBoundaryIds(const SubdomainID subdomain_id) const
Get the list of boundaries that contact the given subdomain.
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
unsigned int maxPLevel() const
Returns the maximum p-refinement level of all elements.
libMesh::BoundingBox getInflatedProcessorBoundingBox(Real inflation_multiplier=0.01) const
Get a (slightly inflated) processor bounding box.
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...
void setMeshBase(std::unique_ptr< MeshBase > mesh_base)
Method to set the mesh_base object.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
std::set< SubdomainID > getInterfaceConnectedBlocks(const BoundaryID bid) const
Get the list of subdomains contacting the given boundary.
bool isParallelTypeForced() const
Tell the user if the distribution was overriden for any reason.
std::vector< Real > _ghosted_boundaries_inflation
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.
void needGhostGhostedBoundaries(bool needghost)
Whether or not we want to ghost ghosted boundaries.
const InputParameters & parameters() const
Get the parameters of the object.
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's are closest to each other between a parent and it...
bool needsRemoteElemDeletion() const
Whether we need to delete remote elements.
std::vector< BndNode * >::const_iterator const_bnd_node_iterator_imp
unsigned int _patch_size
The number of nodes to consider in the NearestNode neighborhood.
bool skipDeletionRepartitionAfterRefine() const
Return a flag indicating whether or not we should skip remote deletion and repartition after uniform ...
elem_info_iterator ownedElemInfoEnd()
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.
virtual std::unique_ptr< libMesh::PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libMesh mesh (default)...
void addGhostedBoundary(BoundaryID boundary_id)
This will add the boundary ids to be ghosted to this processor.
virtual dof_id_type maxNodeId() const
Calls max_node/elem_id() on the underlying libMesh mesh object.
virtual Elem * elem(const dof_id_type i)
Various accessors (pointers/references) for Elem "i".
libMesh::StoredRange< MooseMesh::const_bnd_elem_iterator, const BndElement * > * getBoundaryElementRange()
virtual dof_id_type nActiveElem() const
const_elem_info_iterator(const MooseMesh::elem_info_iterator &rhs)
face_info_iterator ownedFaceInfoEnd()
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.
const MooseEnum & partitionerName() const
Class used for caching additional information for elements such as the volume and centroid...
MeshBase::node_iterator localNodesEnd()
const RealVectorValue & getNormalByBoundaryID(BoundaryID id) const
Returns the normal vector associated with a given BoundaryID.
void doingPRefinement(bool doing_p_refinement)
Indicate whether the kind of adaptivity we're doing is p-refinement.
std::set< SubdomainID > _mesh_subdomains
A set of subdomain IDs currently present in the mesh.
ConstElemPointerRange * coarsenedElementRange() const
Return a range that is suitable for threaded execution over elements that were just coarsened...
const Moose::PatchUpdateType & getPatchUpdateStrategy() const
Get the current patch update strategy.
bool doingPRefinement() const
Query whether we have p-refinement.
std::vector< std::pair< BoundaryID, BoundaryID > > _paired_boundary
A vector holding the paired boundaries for a regular orthogonal mesh.
virtual bool skipPartitioning() const
unsigned int maxHLevel() const
Returns the maximum h-refinement level of all elements.
void ghostGhostedBoundaries()
Actually do the ghosting of boundaries that need to be ghosted to this processor. ...
std::set< unsigned int > _ghosted_boundaries
MeshBase::node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
unsigned int _rz_coord_axis
Storage for RZ axis selection.
void computeFiniteVolumeCoords() const
Compute the face coordinate value for all FaceInfo and ElemInfo objects.
bool _distribution_overridden
void buildNodeList()
Calls BoundaryInfo::build_node_list()/build_side_list() and makes separate copies of Nodes/Elems in t...
std::unordered_map< SubdomainID, SubdomainData > _sub_to_data
Holds a map from subdomain ids to associated data.
std::unique_ptr< libMesh::Partitioner > _custom_partitioner
The custom partitioner.
bool isBoundaryElem(dof_id_type elem_id) const
Returns true if the requested element is in the list of boundary elements, false otherwise.
unsigned int getElementIDIndex(const std::string &id_name) const
Return the accessing integer for an extra element integer with its name.
void needsRemoteElemDeletion(bool need_delete)
Set whether we need to delete remote elements.
std::map< const Elem *, std::vector< const Elem * > > _coarsened_element_children
Map of Parent elements to child elements for elements that were just coarsened.
virtual bool isDistributedMesh() const
Returns the final Mesh distribution type.
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 ...
OStreamProxy out(std::cout)
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...
unsigned int _ghosting_patch_size
The number of nearest neighbors to consider for ghosting purposes when iteration patch update strateg...
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.
const MeshBase::element_iterator activeLocalElementsEnd()
std::unique_ptr< libMesh::ConstNodeRange > _local_node_range
libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > * getBoundaryNodeRange()
The definition of the const_elem_info_iterator struct.
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
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.
bool prepare(const MeshBase *mesh_to_clone)
Calls prepare_for_use() if the underlying MeshBase object isn't prepared, then communicates various b...
const ElemInfo & elemInfo(const dof_id_type id) const
Accessor for the elemInfo object for a given element ID.
void setCustomPartitioner(libMesh::Partitioner *partitioner)
Setter for custom partitioner.
Real _distance
The distance between them.
unsigned int nFace() const
accessors for the FaceInfo objects
static MooseEnum elemTypes()
returns MooseMesh element type options
void meshChanged()
Declares that the MooseMesh has changed, invalidates cached data and rebuilds caches.
void buildPRefinementAndCoarseningMaps(Assembly *assembly)
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
std::unique_ptr< libMesh::StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode * > > _bnd_node_range
virtual dof_id_type nElem() const
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...
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.
std::unique_ptr< SemiLocalNodeRange > _active_semilocal_node_range
const std::set< SubdomainID > & meshSubdomains() const
Returns a read-only reference to the set of subdomains currently present in the Mesh.
SubdomainID getSubdomainID(const SubdomainName &subdomain_name) const
Get the associated subdomain ID for the subdomain name.
unsigned int _max_p_level
Maximum p-refinement level of all elements.
void setupFiniteVolumeMeshData() const
Sets up the additional data needed for finite volume computations.
virtual unsigned int effectiveSpatialDimension() const
Returns the effective spatial dimension determined by the coordinates actually used by the mesh...