13 #include "libmesh/elem.h" 14 #include "libmesh/boundary_info.h" 15 #include "libmesh/parallel.h" 16 #include "libmesh/parallel_algebra.h" 17 #include "libmesh/utility.h" 31 std::map<boundary_id_type, boundary_id_type> same_name_ids;
33 auto populate_map = [](
const std::map<boundary_id_type, std::string> & map,
34 std::map<boundary_id_type, boundary_id_type> & same_ids)
36 for (
const auto & pair_outer : map)
37 for (
const auto & pair_inner : map)
39 if (pair_outer.second == pair_inner.second && pair_outer.first != pair_inner.first &&
40 same_ids.find(pair_inner.first) == same_ids.end())
41 same_ids[pair_outer.first] = pair_inner.first;
44 populate_map(side_bd_name_map, same_name_ids);
45 populate_map(node_bd_name_map, same_name_ids);
47 for (
const auto & [id1, id2] : same_name_ids)
61 std::vector<boundary_id_type> old_ids;
64 for (
auto & elem :
as_range(
mesh.level_elements_begin(0),
mesh.level_elements_end(0)))
66 unsigned int n_sides = elem->n_sides();
70 if (std::find(old_ids.begin(), old_ids.end(), old_id) != old_ids.end())
72 std::vector<boundary_id_type> new_ids(old_ids);
73 std::replace(new_ids.begin(), new_ids.end(), old_id, new_id);
77 boundary_info.
add_side(elem, s, new_ids);
80 boundary_info.
add_side(elem, s, new_ids);
95 std::vector<boundary_id_type>
97 const std::vector<BoundaryName> & boundary_name,
98 bool generate_unknown)
104 std::vector<boundary_id_type>
106 const std::vector<BoundaryName> & boundary_name,
107 bool generate_unknown,
108 const std::set<BoundaryID> & mesh_boundary_ids)
120 if (generate_unknown)
124 max_boundary_local_id = bids.empty() ? 0 : *(bids.rbegin());
130 BoundaryID max_boundary_id = mesh_boundary_ids.empty() ? 0 : *(mesh_boundary_ids.rbegin());
133 max_boundary_id > max_boundary_local_id ? max_boundary_id : max_boundary_local_id;
135 std::vector<BoundaryID> ids(boundary_name.size());
138 if (boundary_name[i] ==
"ANY_BOUNDARY_ID")
140 ids.assign(mesh_boundary_ids.begin(), mesh_boundary_ids.end());
142 mooseWarning(
"You passed \"ANY_BOUNDARY_ID\" in addition to other boundary_names. This " 143 "may be a logic error.");
147 if (boundary_name[i].empty() && !generate_unknown)
148 mooseError(
"Incoming boundary name is empty and we are not generating unknown boundary IDs. " 160 if (generate_unknown &&
163 id = ++max_boundary_id;
168 id = getIDFromName<BoundaryName, BoundaryID>(boundary_name[i]);
178 const std::vector<BoundaryName> & boundary_name,
179 bool generate_unknown)
182 return std::set<BoundaryID>(boundaries.begin(), boundaries.end());
185 std::vector<subdomain_id_type>
188 std::vector<SubdomainID> ids(subdomain_names.size());
196 std::set<subdomain_id_type>
199 std::set<SubdomainID> ids;
200 for (
const auto &
name : subdomain_names)
209 if (boundary_name.empty())
215 id = getIDFromName<BoundaryName, BoundaryID>(boundary_name);
223 if (subdomain_name ==
"ANY_BLOCK_ID")
224 mooseError(
"getSubdomainID() does not work with \"ANY_BLOCK_ID\"");
227 if (subdomain_name.empty())
233 id = getIDFromName<SubdomainName, SubdomainID>(subdomain_name);
241 for (
const auto & elem :
mesh.element_ptr_range())
242 if (elem->subdomain_id() == old_id)
243 elem->subdomain_id() = new_id;
254 for (
const auto & elem :
255 as_range(
mesh.active_local_elements_begin(),
mesh.active_local_elements_end()))
257 Real elem_vol = elem->volume();
258 centroid_pt += (elem->true_centroid()) * elem_vol;
263 centroid_pt /= vol_tmp;
267 std::unordered_map<dof_id_type, dof_id_type>
269 const std::set<SubdomainID> & block_ids,
270 std::vector<ExtraElementIDName> extra_ids)
273 const bool block_restricted = !block_ids.empty();
275 ExtraElementIDName id_name = extra_ids.back();
276 extra_ids.pop_back();
280 if (extra_ids.empty())
283 std::vector<dof_id_type> ids;
285 std::set<dof_id_type> ids_set;
286 for (
const auto & elem :
mesh.active_element_ptr_range())
288 if (block_restricted && block_ids.find(elem->subdomain_id()) == block_ids.end())
290 const auto id = elem->get_extra_integer(id_index);
294 ids.assign(ids_set.begin(), ids_set.end());
298 std::unordered_map<dof_id_type, dof_id_type> parsed_ids;
299 for (
auto & elem :
mesh.active_element_ptr_range())
301 if (block_restricted && block_ids.find(elem->subdomain_id()) == block_ids.end())
303 parsed_ids[elem->id()] = std::distance(
304 ids.begin(), std::lower_bound(ids.begin(), ids.end(), elem->get_extra_integer(id_index)));
310 const auto base_parsed_ids =
313 std::vector<std::pair<dof_id_type, dof_id_type>> unique_ids;
315 std::set<std::pair<dof_id_type, dof_id_type>> unique_ids_set;
316 for (
const auto & elem :
mesh.active_element_ptr_range())
318 if (block_restricted && block_ids.find(elem->subdomain_id()) == block_ids.end())
320 const dof_id_type id1 = libmesh_map_find(base_parsed_ids, elem->id());
321 const dof_id_type id2 = elem->get_extra_integer(id_index);
322 const std::pair<dof_id_type, dof_id_type> ids = std::make_pair(id1, id2);
323 unique_ids_set.insert(ids);
326 unique_ids.assign(unique_ids_set.begin(), unique_ids_set.end());
329 std::unordered_map<dof_id_type, dof_id_type> parsed_ids;
331 for (
const auto & elem :
mesh.active_element_ptr_range())
333 if (block_restricted && block_ids.find(elem->subdomain_id()) == block_ids.end())
335 const dof_id_type id1 = libmesh_map_find(base_parsed_ids, elem->id());
336 const dof_id_type id2 = elem->get_extra_integer(id_index);
339 std::lower_bound(unique_ids.begin(), unique_ids.end(), std::make_pair(id1, id2)));
340 parsed_ids[elem->id()] = new_id;
349 for (
const auto & pt : vec_pts)
358 return isCoPlanar(vec_pts, plane_nvec, vec_pts.front());
366 std::vector<Point> vec_pts_nonzero{vec_pts[0]};
369 vec_pts_nonzero.push_back(vec_pts[i]);
371 if (vec_pts_nonzero.size() <= 3)
375 for (
const auto i :
make_range(vec_pts_nonzero.size() - 1))
377 const Point tmp_pt = (vec_pts_nonzero[i] - vec_pts_nonzero[0])
378 .cross(vec_pts_nonzero[i + 1] - vec_pts_nonzero[0]);
394 std::set<SubdomainID> preexisting_subdomain_ids;
396 if (preexisting_subdomain_ids.empty())
400 const auto highest_subdomain_id =
401 *std::max_element(preexisting_subdomain_ids.begin(), preexisting_subdomain_ids.end());
403 "A SubdomainID with max possible value was found");
404 return highest_subdomain_id + 1;
412 if (boundary_ids.empty())
414 return (*boundary_ids.rbegin() + 1);
420 std::set<SubdomainID> mesh_blocks;
461 std::vector<dof_id_type> & elem_id_list,
462 std::vector<dof_id_type> & midpoint_node_list,
463 std::vector<dof_id_type> & ordered_node_list,
464 std::vector<dof_id_type> & ordered_elem_id_list)
467 bool is_flipped =
false;
469 mooseAssert(node_assm.size(),
"Node list must not be empty");
470 ordered_node_list.push_back(node_assm.front().first);
472 ordered_node_list.push_back(midpoint_node_list.front());
473 ordered_node_list.push_back(node_assm.front().second);
474 ordered_elem_id_list.push_back(elem_id_list.front());
476 node_assm.erase(node_assm.begin());
477 midpoint_node_list.erase(midpoint_node_list.begin());
478 elem_id_list.erase(elem_id_list.begin());
479 const unsigned int node_assm_size_0 = node_assm.size();
480 for (
unsigned int i = 0; i < node_assm_size_0; i++)
483 dof_id_type end_node_id = ordered_node_list.back();
484 auto isMatch1 = [end_node_id](std::pair<dof_id_type, dof_id_type> old_id_pair)
485 {
return old_id_pair.first == end_node_id; };
486 auto isMatch2 = [end_node_id](std::pair<dof_id_type, dof_id_type> old_id_pair)
487 {
return old_id_pair.second == end_node_id; };
488 auto result = std::find_if(node_assm.begin(), node_assm.end(), isMatch1);
490 if (result == node_assm.end())
493 result = std::find_if(node_assm.begin(), node_assm.end(), isMatch2);
500 if (result != node_assm.end())
502 const auto elem_index = std::distance(node_assm.begin(), result);
504 ordered_node_list.push_back(midpoint_node_list[elem_index]);
505 ordered_node_list.push_back(match_first ? (*result).second : (*result).first);
506 node_assm.erase(result);
507 midpoint_node_list.erase(midpoint_node_list.begin() + elem_index);
508 ordered_elem_id_list.push_back(elem_id_list[elem_index]);
509 elem_id_list.erase(elem_id_list.begin() + elem_index);
519 throw MooseException(
"The node list provided has more than one segments.");
523 std::reverse(ordered_node_list.begin(), ordered_node_list.end());
524 std::reverse(midpoint_node_list.begin(), midpoint_node_list.end());
525 std::reverse(ordered_elem_id_list.begin(), ordered_elem_id_list.end());
534 std::vector<dof_id_type> & elem_id_list,
535 std::vector<dof_id_type> & ordered_node_list,
536 std::vector<dof_id_type> & ordered_elem_id_list)
540 node_assm, elem_id_list, dummy_midpoint_node_list, ordered_node_list, ordered_elem_id_list);
553 const std::string & class_name,
554 const unsigned int num_sections,
555 const unsigned int num_integers,
556 const std::vector<std::vector<std::vector<dof_id_type>>> & elem_integers_swaps,
557 std::vector<std::unordered_map<dof_id_type, dof_id_type>> & elem_integers_swap_pairs)
559 elem_integers_swap_pairs.reserve(num_sections * num_integers);
562 const auto & elem_integer_swaps = elem_integers_swaps[i];
563 std::vector<std::unordered_map<dof_id_type, dof_id_type>> elem_integer_swap_pairs;
567 "elem_integers_swaps",
569 elem_integer_swap_pairs,
577 elem_integers_swap_pairs.insert(elem_integers_swap_pairs.end(),
578 elem_integer_swap_pairs.begin(),
579 elem_integer_swap_pairs.end());
583 std::unique_ptr<ReplicatedMesh>
586 auto poly_mesh = std::make_unique<ReplicatedMesh>(input_mesh.
comm());
590 std::unordered_map<dof_id_type, dof_id_type> old_new_node_map;
591 for (
const auto & bside : side_list)
593 if (std::get<2>(bside) != boundary_id)
596 const Elem * elem = input_mesh.
elem_ptr(std::get<0>(bside));
597 const auto side = std::get<1>(bside);
599 auto copy = side_elem->build(side_elem->type());
601 for (
const auto i : side_elem->node_index_range())
603 auto & n = side_elem->node_ref(i);
605 if (old_new_node_map.count(n.id()))
606 copy->set_node(i, poly_mesh->node_ptr(old_new_node_map[n.id()]));
609 Node * node = poly_mesh->add_point(side_elem->point(i));
610 copy->set_node(i, node);
611 old_new_node_map[n.id()] = node->
id();
614 poly_mesh->add_elem(copy.release());
616 poly_mesh->skip_partitioning(
true);
617 poly_mesh->prepare_for_use();
618 if (poly_mesh->n_elem() == 0)
619 mooseError(
"The input mesh does not have a boundary with id ", boundary_id);
void swapNodesInElem(Elem &elem, const unsigned int nd1, const unsigned int nd2)
std::string name(const ElemQuality q)
void remove_id(boundary_id_type id, bool global=false)
void makeOrderedNodeList(std::vector< std::pair< dof_id_type, dof_id_type >> &node_assm, std::vector< dof_id_type > &elem_id_list, std::vector< dof_id_type > &midpoint_node_list, std::vector< dof_id_type > &ordered_node_list, std::vector< dof_id_type > &ordered_elem_id_list)
Convert a list of sides in the form of a vector of pairs of node ids into a list of ordered nodes bas...
std::unordered_map< dof_id_type, dof_id_type > getExtraIDUniqueCombinationMap(const MeshBase &mesh, const std::set< SubdomainID > &block_ids, std::vector< ExtraElementIDName > extra_ids)
virtual Node *& set_node(const unsigned int i)
std::set< BoundaryID > getBoundaryIDSet(const libMesh::MeshBase &mesh, const std::vector< BoundaryName > &boundary_name, bool generate_unknown)
Gets the boundary IDs into a set with their names.
virtual const char * what() const
Get out the error message.
std::unique_ptr< ReplicatedMesh > buildBoundaryMesh(const ReplicatedMesh &input_mesh, const boundary_id_type boundary_id)
auto norm() const -> decltype(std::norm(Real()))
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
bool hasBoundaryName(const MeshBase &input_mesh, const BoundaryName &name)
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i)=0
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
const BoundaryID INVALID_BOUNDARY_ID
void mooseWarning(Args &&... args)
Emit a warning message with the given stringified, concatenated args.
bool doesMapContainValue(const std::map< T1, T2 > &the_map, const T2 &value)
This routine is a simple helper function for searching a map by values instead of keys...
const std::map< boundary_id_type, std::string > & get_sideset_name_map() const
const Parallel::Communicator & comm() const
bool hasBoundaryID(const MeshBase &input_mesh, const BoundaryID id)
void boundary_ids(const Node *node, std::vector< boundary_id_type > &vec_to_fill) const
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...
const BoundaryInfo & get_boundary_info() const
void renumber_id(boundary_id_type old_id, boundary_id_type new_id)
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)
subdomain_id_type get_id_by_name(std::string_view name) const
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
bool hasSubdomainID(const MeshBase &input_mesh, const SubdomainID &id)
const SubdomainID INVALID_BLOCK_ID
boundary_id_type get_id_by_name(std::string_view name) const
BoundaryID getBoundaryID(const BoundaryName &boundary_name, const MeshBase &mesh)
virtual const Elem * elem_ptr(const dof_id_type i) const override final
unsigned int get_elem_integer_index(std::string_view name) const
virtual bool is_serial() const
TypeVector< Real > unit() const
const std::map< boundary_id_type, std::string > & get_nodeset_name_map() const
boundary_id_type BoundaryID
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
void subdomain_ids(std::set< subdomain_id_type > &ids, const bool global=true) const
static const dof_id_type invalid_id
void idSwapParametersProcessor(const std::string &class_name, const std::string &id_name, const std::vector< std::vector< T >> &id_swaps, std::vector< std::unordered_map< T, T >> &id_swap_pairs, const unsigned int row_index_shift=0)
Reprocess the swap related input parameters to make pairs out of them to ease further processing...
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.
void changeBoundaryId(MeshBase &mesh, const boundary_id_type old_id, const boundary_id_type new_id, bool delete_prev)
bool hasSubdomainName(const MeshBase &input_mesh, const SubdomainName &name)
bool isCoPlanar(const std::vector< Point > vec_pts)
const std::set< boundary_id_type > & get_boundary_ids() const
void changeSubdomainId(MeshBase &mesh, const subdomain_id_type old_id, const subdomain_id_type new_id)
void remove_side(const Elem *elem, const unsigned short int side)
Provides a way for users to bail out of the current solve.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void max(const T &r, T &o, Request &req) const
const Node * node_ptr(const unsigned int i) const
virtual bool is_replicated() const
void add_side(const dof_id_type elem, const unsigned short int side, const boundary_id_type id)
const std::set< boundary_id_type > & get_global_boundary_ids() const
IntRange< T > make_range(T beg, T end)
void extraElemIntegerSwapParametersProcessor(const std::string &class_name, const unsigned int num_sections, const unsigned int num_integers, const std::vector< std::vector< std::vector< dof_id_type >>> &elem_integers_swaps, std::vector< std::unordered_map< dof_id_type, dof_id_type >> &elem_integers_swap_pairs)
Reprocess the elem_integers_swaps into maps so they are easier to use.
Point meshCentroidCalculator(const MeshBase &mesh)
void mergeBoundaryIDsWithSameName(MeshBase &mesh)
SubdomainID getNextFreeSubdomainID(MeshBase &input_mesh)
bool isDigits(const std::string &str)
Courtesy https://stackoverflow.com/a/8889045 and https://en.cppreference.com/w/cpp/string/byte/isdigi...
BoundaryID getNextFreeBoundaryID(MeshBase &input_mesh)
auto index_range(const T &sizable)
void set_union(T &data, const unsigned int root_id) const