21 #include "libmesh/libmesh_config.h" 24 #include "libmesh/boundary_info.h" 25 #include "libmesh/libmesh_logging.h" 26 #include "libmesh/elem.h" 27 #include "libmesh/ghost_point_neighbors.h" 28 #include "libmesh/mesh_base.h" 29 #include "libmesh/mesh_communication.h" 30 #include "libmesh/mesh_serializer.h" 31 #include "libmesh/mesh_tools.h" 32 #include "libmesh/parallel.h" 33 #include "libmesh/parallel_algebra.h" 34 #include "libmesh/parallel_fe_type.h" 35 #include "libmesh/partitioner.h" 36 #include "libmesh/point_locator_base.h" 37 #include "libmesh/sparse_matrix.h" 38 #include "libmesh/threads.h" 39 #include "libmesh/enum_elem_type.h" 40 #include "libmesh/enum_point_locator_type.h" 41 #include "libmesh/enum_to_string.h" 42 #include "libmesh/point_locator_nanoflann.h" 43 #include "libmesh/elem_side_builder.h" 50 #include <unordered_map> 65 _default_mapping_data(0),
68 _count_lower_dim_elems_in_point_locator(true),
70 #ifdef LIBMESH_ENABLE_UNIQUE_ID
71 _next_unique_id(
DofObject::invalid_unique_id),
74 _skip_noncritical_partitioning(false),
76 _skip_renumber_nodes_and_elements(false),
77 _skip_find_neighbors(false),
78 _allow_remote_element_removal(true),
79 _spatial_dimension(d),
81 _point_locator_close_to_point_tol(0.)
85 libmesh_assert_less_equal (LIBMESH_DIM, 3);
86 libmesh_assert_greater_equal (LIBMESH_DIM, d);
95 _n_parts (other_mesh._n_parts),
96 _default_mapping_type(other_mesh._default_mapping_type),
97 _default_mapping_data(other_mesh._default_mapping_data),
98 _is_prepared (other_mesh._is_prepared),
100 _count_lower_dim_elems_in_point_locator(other_mesh._count_lower_dim_elems_in_point_locator),
102 #ifdef LIBMESH_ENABLE_UNIQUE_ID
103 _next_unique_id(other_mesh._next_unique_id),
107 _interior_mesh((other_mesh._interior_mesh == &other_mesh) ?
108 this : other_mesh._interior_mesh),
109 _skip_noncritical_partitioning(other_mesh._skip_noncritical_partitioning),
110 _skip_all_partitioning(other_mesh._skip_all_partitioning),
111 _skip_renumber_nodes_and_elements(other_mesh._skip_renumber_nodes_and_elements),
112 _skip_find_neighbors(other_mesh._skip_find_neighbors),
113 _allow_remote_element_removal(other_mesh._allow_remote_element_removal),
114 _elem_dims(other_mesh._elem_dims),
115 _elem_default_orders(other_mesh._elem_default_orders),
116 _supported_nodal_order(other_mesh._supported_nodal_order),
117 _elemset_codes_inverse_map(other_mesh._elemset_codes_inverse_map),
118 _all_elemset_ids(other_mesh._all_elemset_ids),
119 _spatial_dimension(other_mesh._spatial_dimension),
121 _point_locator_close_to_point_tol(other_mesh._point_locator_close_to_point_tol)
129 if (gf == other_default_ghosting)
135 std::shared_ptr<GhostingFunctor> clone_gf = gf->
clone();
144 clone_gf->set_mesh(
this);
149 libmesh_deprecated();
165 LOG_SCOPE(
"operator=(&&)",
"MeshBase");
170 _n_parts = other_mesh.n_partitions();
176 #ifdef LIBMESH_ENABLE_UNIQUE_ID 189 _elem_dims = std::move(other_mesh.elem_dimensions());
213 "Elemset code " << code <<
" not found in _elmset_codes container.");
214 libmesh_assert_equal_to(it->second, &
set);
228 LOG_SCOPE(
"operator==()",
"MeshBase");
318 auto it = other_rows.find(other_node);
319 if (it == other_rows.end())
322 const auto & other_row = it->second;
323 if (row.size() != other_row.size())
328 const auto & [elem_pair, coef] = row[i];
329 const auto & [other_elem_pair, other_coef] = other_row[i];
332 if (elem_pair.first->id() !=
333 other_elem_pair.first->id() ||
335 other_elem_pair.second ||
341 for (
const auto & [elemset_code, elemset_ptr] : this->
_elemset_codes)
343 it == other_mesh.
_elemset_codes.end() || *elemset_ptr != *it->second)
375 return cast_int<unsigned int>(*
_elem_dims.rbegin());
386 parallel_object_only();
390 "Specified element dimensions does not match true element dimensions!");
404 const auto & inserted_id_set = it1->first;
414 auto [it2, inserted2] =
_elemset_codes.emplace(code, &inserted_id_set);
418 libmesh_error_msg_if(!inserted2 && it2->second != &inserted_id_set,
419 "The elemset code " << code <<
" already exists with a different id_set.");
432 id_set_to_fill.clear();
436 id_set_to_fill.insert(it->second->begin(), it->second->end());
447 std::vector<dof_id_type> ret;
450 ret.push_back(pr.first);
473 "Expected _elemset_codes_inverse_map entry for elemset code " << old_code);
492 for (
auto & elem : this->element_ptr_range())
495 elem->get_extra_integer(elemset_index);
497 if (elemset_code == old_code)
498 elem->set_extra_integer(elemset_index, new_code);
510 "Cannot change elemset id " << old_id <<
511 " to " << new_id <<
", " << new_id <<
" already exists.");
515 std::map<MeshBase::elemset_type, dof_id_type> new_elemset_codes_inverse_map;
518 auto id_set_copy = id_set;
519 if (id_set_copy.count(old_id))
522 id_set_copy.erase(old_id);
523 id_set_copy.insert(new_id);
527 new_elemset_codes_inverse_map.emplace(id_set_copy, elemset_code);
585 const std::vector<dof_id_type> * default_values)
587 libmesh_assert(!default_values || default_values->size() == names.size());
590 std::unordered_map<std::string, std::size_t> name_indices;
594 std::vector<unsigned int> returnval(names.size());
596 bool added_an_integer =
false;
599 const std::string &
name = names[i];
600 if (
const auto it = name_indices.find(
name);
601 it != name_indices.end())
603 returnval[i] = it->second;
610 name_indices[
name] = returnval[i];
614 added_an_integer =
true;
618 if (allocate_data && added_an_integer)
632 libmesh_error_msg(
"Unknown elem integer " <<
name);
674 const std::vector<dof_id_type> * default_values)
676 libmesh_assert(!default_values || default_values->size() == names.size());
679 std::unordered_map<std::string, std::size_t> name_indices;
683 std::vector<unsigned int> returnval(names.size());
685 bool added_an_integer =
false;
688 const std::string &
name = names[i];
689 if (
const auto it = name_indices.find(
name);
690 it != name_indices.end())
692 returnval[i] = it->second;
699 name_indices[
name] = returnval[i];
703 added_an_integer =
true;
707 if (allocate_data && added_an_integer)
721 libmesh_error_msg(
"Unknown node integer " <<
name);
740 LOG_SCOPE(
"remove_orphaned_nodes()",
"MeshBase");
743 std::unordered_set<Node *> connected_nodes;
747 for (
const auto & element : this->element_ptr_range())
748 for (
auto & n : element->node_ref_range())
749 connected_nodes.insert(&n);
751 for (
const auto & node : this->node_ptr_range())
752 if (!connected_nodes.count(node))
758 #ifdef LIBMESH_ENABLE_DEPRECATED 761 libmesh_deprecated();
766 if (skip_renumber_nodes_and_elements)
780 libmesh_deprecated();
785 if (skip_renumber_nodes_and_elements)
790 #endif // LIBMESH_ENABLE_DEPRECATED 796 LOG_SCOPE(
"prepare_for_use()",
"MeshBase");
798 parallel_object_only();
864 #if defined(DEBUG) && defined(LIBMESH_ENABLE_UNIQUE_ID) 904 #ifdef LIBMESH_ENABLE_UNIQUE_ID 963 parallel_object_only();
967 for (
const auto & elem : this->active_local_element_ptr_range())
968 ids.insert(elem->subdomain_id());
975 for (
const auto & elem : this->active_unpartitioned_element_ptr_range())
976 ids.insert(elem->subdomain_id());
1000 parallel_object_only();
1002 std::set<subdomain_id_type> ids;
1006 return cast_int<subdomain_id_type>(ids.size());
1013 std::set<subdomain_id_type> ids;
1017 return cast_int<subdomain_id_type>(ids.size());
1031 this->pid_nodes_end (proc_id)));
1044 this->pid_elements_end (proc_id)));
1053 this->active_pid_elements_end (proc_id)));
1062 for (
const auto & elem : this->element_ptr_range())
1063 ne += elem->n_sub_elem();
1074 for (
const auto & elem : this->active_element_ptr_range())
1075 ne += elem->n_sub_elem();
1084 std::ostringstream oss;
1086 oss <<
" Mesh Information:" <<
'\n';
1090 oss <<
" elem_dimensions()={";
1093 std::ostream_iterator<unsigned int>(oss,
", "));
1094 oss << cast_int<unsigned int>(*
_elem_dims.rbegin());
1100 oss <<
" elem_default_orders()={";
1103 std::ostream_iterator<std::string>(oss,
", "),
1105 {
return Utility::enum_to_string<Order>(o); });
1112 <<
" n_nodes()=" << this->
n_nodes() <<
'\n' 1114 <<
" n_elem()=" << this->
n_elem() <<
'\n' 1116 #ifdef LIBMESH_ENABLE_AMR 1117 oss <<
" n_active_elem()=" << this->
n_active_elem() <<
'\n';
1120 oss <<
" n_subdomains()=" <<
static_cast<std::size_t
>(this->
n_subdomains()) <<
'\n';
1122 oss <<
" n_local_subdomains()= " <<
static_cast<std::size_t
>(this->
n_local_subdomains()) <<
'\n';
1123 oss <<
" n_elemsets()=" <<
static_cast<std::size_t
>(this->
n_elemsets()) <<
'\n';
1126 oss <<
" n_partitions()=" <<
static_cast<std::size_t
>(this->
n_partitions()) <<
'\n' 1127 <<
" n_processors()=" << static_cast<std::size_t>(this->
n_processors()) <<
'\n' 1129 <<
" processor_id()=" << static_cast<std::size_t>(this->
processor_id()) <<
'\n' 1130 <<
" is_prepared()=" << (this->
is_prepared() ?
"true" :
"false") <<
'\n' 1131 <<
" is_replicated()=" << (this->
is_replicated() ?
"true" :
"false") <<
'\n';
1137 libmesh_parallel_only(this->
comm());
1139 oss <<
"\n Detailed global get_info() (verbosity > 0) is reduced and output to only rank 0.";
1143 const auto elem_type_helper = [](
const std::set<int> &
elem_types) {
1144 std::stringstream ss;
1156 const auto include_object = [
this, &global](
const DofObject & dof_object) {
1157 return this->
processor_id() == dof_object.processor_id() ||
1168 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Bounding Box:\n" 1169 <<
" Minimum: " << bbox.min() <<
"\n" 1170 <<
" Maximum: " << bbox.max() <<
"\n" 1171 <<
" Delta: " << (bbox.max() - bbox.min()) <<
"\n";
1175 for (
const Elem * elem : this->active_local_element_ptr_range())
1181 for (
const Elem * elem : this->active_unpartitioned_element_ptr_range())
1189 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Element Type(s):\n " 1200 std::size_t num_nodes = 0;
1203 std::map<boundary_id_type, NodesetInfo> nodeset_info_map;
1206 if (!include_object(*node))
1209 NodesetInfo &
info = nodeset_info_map[id];
1220 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Nodesets:\n";
1221 if (nodeset_ids.empty())
1226 for (
const auto id : nodeset_ids)
1228 NodesetInfo &
info = nodeset_info_map[id];
1233 this->
comm().
sum(info.num_nodes);
1236 this->
comm().
min(info.bbox.min());
1237 this->
comm().
max(info.bbox.max());
1241 const bool has_name = nodeset_name_map.count(
id) && nodeset_name_map.at(
id).size();
1242 const std::string
name = has_name ? nodeset_name_map.at(
id) :
"";
1246 if (global ? this->
processor_id() == 0 : info.num_nodes > 0)
1248 oss <<
" Nodeset " << id;
1250 oss <<
" (" <<
name <<
")";
1251 oss <<
", " <<
info.num_nodes <<
" " << (global ?
"" :
"local ") <<
"nodes\n";
1255 oss <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box minimum: " 1256 <<
info.bbox.min() <<
"\n" 1257 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box maximum: " 1258 <<
info.bbox.max() <<
"\n" 1259 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box delta: " 1260 << (
info.bbox.max() -
info.bbox.min()) <<
"\n";
1273 std::size_t num_sides = 0;
1275 std::set<int> side_elem_types;
1277 std::set<dof_id_type> elem_ids;
1278 std::set<dof_id_type> node_ids;
1282 std::map<boundary_id_type, SidesetInfo> sideset_info_map;
1285 const Elem * elem = pair.first;
1286 if (!include_object(*elem))
1289 const auto id = pair.second.second;
1290 SidesetInfo &
info = sideset_info_map[id];
1292 const auto s = pair.second.first;
1293 const Elem & side = side_builder(*elem, s);
1296 info.side_elem_types.insert(side.
type());
1297 info.elem_types.insert(elem->
type());
1298 info.elem_ids.insert(elem->
id());
1301 if (include_object(node))
1302 info.node_ids.insert(node.id());
1314 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Sidesets:\n";
1315 if (sideset_ids.empty())
1319 for (
const auto id : sideset_ids)
1321 SidesetInfo &
info = sideset_info_map[id];
1323 auto num_elems =
info.elem_ids.size();
1324 auto num_nodes =
info.node_ids.size();
1329 this->
comm().
sum(info.num_sides);
1336 this->
comm().
sum(info.volume);
1337 this->
comm().
min(info.bbox.min());
1338 this->
comm().
max(info.bbox.max());
1342 const bool has_name = sideset_name_map.count(
id) && sideset_name_map.at(
id).size();
1343 const std::string
name = has_name ? sideset_name_map.at(
id) :
"";
1347 if (global ? this->
processor_id() == 0 : info.num_sides > 0)
1349 oss <<
" Sideset " << id;
1351 oss <<
" (" <<
name <<
")";
1352 oss <<
", " <<
info.num_sides <<
" sides (" << elem_type_helper(
info.side_elem_types) <<
")" 1353 <<
", " << num_elems <<
" " << (global ?
"" :
"local ") <<
"elems (" << elem_type_helper(
info.elem_types) <<
")" 1354 <<
", " << num_nodes <<
" " << (global ?
"" :
"local ") <<
"nodes\n";
1358 oss <<
" " << (global ?
"Side" :
"Local side") <<
" volume: " <<
info.volume <<
"\n" 1359 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box minimum: " 1360 <<
info.bbox.min() <<
"\n" 1361 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box maximum: " 1362 <<
info.bbox.max() <<
"\n" 1363 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box delta: " 1364 << (
info.bbox.max() -
info.bbox.min()) <<
"\n";
1377 std::size_t num_edges = 0;
1378 std::set<int> edge_elem_types;
1381 std::map<boundary_id_type, EdgesetInfo> edgeset_info_map;
1382 std::unique_ptr<const Elem> edge;
1386 const Elem * elem = pair.first;
1387 if (!include_object(*elem))
1390 const auto id = pair.second.second;
1391 EdgesetInfo &
info = edgeset_info_map[id];
1396 info.edge_elem_types.insert(edge->type());
1399 info.bbox.union_with(edge->loose_bounding_box());
1405 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Edgesets:\n";
1406 if (edgeset_ids.empty())
1411 for (
const auto id : edgeset_ids)
1413 EdgesetInfo &
info = edgeset_info_map[id];
1418 this->
comm().
sum(info.num_edges);
1422 this->
comm().
min(info.bbox.min());
1423 this->
comm().
min(info.bbox.max());
1427 const bool has_name = edgeset_name_map.count(
id) && edgeset_name_map.at(
id).size();
1428 const std::string
name = has_name ? edgeset_name_map.at(
id) :
"";
1432 if (global ? this->
processor_id() == 0 : info.num_edges > 0)
1434 oss <<
" Edgeset " << id;
1436 oss <<
" (" <<
name <<
")";
1437 oss <<
", " <<
info.num_edges <<
" " << (global ?
"" :
"local ") <<
"edges (" 1438 << elem_type_helper(
info.edge_elem_types) <<
")\n";
1442 oss <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box minimum: " 1443 <<
info.bbox.min() <<
"\n" 1444 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box maximum: " 1445 <<
info.bbox.max() <<
"\n" 1446 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box delta: " 1447 << (
info.bbox.max() -
info.bbox.min()) <<
"\n";
1453 std::set<subdomain_id_type> subdomains;
1454 for (
const Elem * elem : this->active_element_ptr_range())
1455 if (include_object(*elem))
1456 subdomains.insert(elem->subdomain_id());
1461 struct SubdomainInfo
1463 std::size_t num_elems = 0;
1466 std::set<dof_id_type> active_node_ids;
1467 #ifdef LIBMESH_ENABLE_AMR 1468 std::size_t num_active_elems = 0;
1472 std::map<subdomain_id_type, SubdomainInfo> subdomain_info_map;
1473 for (
const Elem * elem : this->element_ptr_range())
1474 if (include_object(*elem))
1476 SubdomainInfo &
info = subdomain_info_map[elem->subdomain_id()];
1479 info.elem_types.insert(elem->type());
1481 #ifdef LIBMESH_ENABLE_AMR 1483 ++
info.num_active_elems;
1486 for (
const Node & node : elem->node_ref_range())
1487 if (include_object(node) && node.active())
1488 info.active_node_ids.insert(node.id());
1490 if (verbosity > 1 && elem->active())
1492 info.volume += elem->volume();
1498 oss <<
"\n " << (global ?
"" :
"Local ") <<
"Mesh Subdomains:\n";
1500 for (
const auto id : subdomains)
1502 SubdomainInfo &
info = subdomain_info_map[id];
1504 auto num_active_nodes =
info.active_node_ids.size();
1509 this->
comm().
sum(info.num_elems);
1510 #ifdef LIBMESH_ENABLE_AMR 1511 this->
comm().
sum(info.num_active_elems);
1513 this->
comm().
sum(num_active_nodes);
1517 this->
comm().
min(info.bbox.min());
1518 this->
comm().
max(info.bbox.max());
1519 this->
comm().
sum(info.volume);
1525 const bool has_name = subdomain_name_map.count(
id);
1526 const std::string
name = has_name ? subdomain_name_map.at(
id) :
"";
1532 oss <<
" Subdomain " << id;
1534 oss <<
" (" <<
name <<
")";
1535 oss <<
": " <<
info.num_elems <<
" " << (global ?
"" :
"local ") <<
"elems " 1536 <<
"(" << elem_type_helper(
info.elem_types);
1537 #ifdef LIBMESH_ENABLE_AMR 1538 oss <<
", " <<
info.num_active_elems <<
" active";
1540 oss <<
"), " << num_active_nodes <<
" " << (global ?
"" :
"local ") <<
"active nodes\n";
1543 oss <<
" " << (global ?
"Volume" :
"Local volume") <<
": " <<
info.volume <<
"\n";
1544 oss <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box minimum: " 1545 <<
info.bbox.min() <<
"\n" 1546 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box maximum: " 1547 <<
info.bbox.max() <<
"\n" 1548 <<
" " << (global ?
"Bounding" :
"Local bounding") <<
" box delta: " 1549 << (
info.bbox.max() -
info.bbox.min()) <<
"\n";
1554 oss <<
" " << (global ?
"Global" :
"Local") <<
" mesh volume = " <<
volume <<
"\n";
1564 os << this->
get_info(verbosity, global)
1621 parallel_object_only();
1623 unsigned int max_proc_id=0;
1625 for (
const auto & elem : this->active_local_element_ptr_range())
1626 max_proc_id = std::max(max_proc_id, static_cast<unsigned int>(elem->processor_id()));
1647 parallel_object_only();
1649 #ifdef LIBMESH_ENABLE_NANOFLANN_POINTLOCATOR 1662 #ifdef LIBMESH_ENABLE_NANOFLANN_POINTLOCATOR 1700 static const std::string empty;
1706 return iter->second;
1716 if (sbd_name ==
name)
1724 #ifdef LIBMESH_ENABLE_DEPRECATED 1727 libmesh_deprecated();
1731 #endif // LIBMESH_ENABLE_DEPRECATED 1736 parallel_object_only();
1745 for (
const auto & elem : this->active_element_ptr_range())
1747 _elem_dims.insert(cast_int<unsigned char>(elem->dim()));
1753 static_cast<int>(elem->supported_nodal_order())));
1780 for (
const auto & node : this->node_ptr_range())
1790 #if LIBMESH_DIM == 2 1799 if ((*node)(2) != 0.)
1814 parallel_object_only();
1824 const bool separate_interior_mesh = (&(this->
interior_mesh()) !=
this);
1827 std::unordered_map<dof_id_type, std::vector<dof_id_type>> node_to_elem;
1829 for (
const auto & elem : this->element_ptr_range())
1832 for (
auto n :
make_range(elem->n_vertices()))
1834 libmesh_assert_less (elem->id(), this->
max_elem_id());
1836 node_to_elem[elem->node_id(n)].push_back(elem->id());
1841 for (
const auto & element : this->element_ptr_range())
1844 if (element->dim()>=LIBMESH_DIM || element->interior_parent())
1851 std::vector<std::set<dof_id_type>> neighbors( element->n_vertices() );
1853 bool found_interior_parents =
false;
1855 for (
auto n :
make_range(element->n_vertices()))
1857 std::vector<dof_id_type> & element_ids = node_to_elem[element->node_id(n)];
1858 for (
const auto & eid : element_ids)
1859 if (this->
elem_ref(eid).
dim() == element->dim()+1)
1860 neighbors[n].insert(eid);
1862 if (neighbors[n].size()>0)
1864 found_interior_parents =
true;
1871 found_interior_parents =
false;
1880 if (found_interior_parents)
1882 std::set<dof_id_type> & neighbors_0 = neighbors[0];
1883 for (
const auto & interior_parent_id : neighbors_0)
1885 found_interior_parents =
false;
1886 for (
auto n :
make_range(1u, element->n_vertices()))
1888 if (neighbors[n].count(interior_parent_id))
1890 found_interior_parents =
true;
1894 found_interior_parents =
false;
1899 if (found_interior_parents)
1901 element->set_interior_parent(this->
elem_ptr(interior_parent_id));
1910 if (separate_interior_mesh)
1911 libmesh_not_implemented_msg
1912 (
"interior_parent() values in multiple meshes are unsupported.");
1943 for (
auto elem : this->element_ptr_range())
1952 for (
auto node : this->node_ptr_range())
1957 std::pair<std::vector<unsigned int>, std::vector<unsigned int>>
1960 std::pair<std::vector<unsigned int>, std::vector<unsigned int>> returnval;
1992 (sf.second)->set_mesh(
this);
1998 if (other_mesh.partitioner())
2016 for (
const auto & other_node : other_mesh.node_ptr_range())
2021 if (*other_node != *node)
2024 for (
const auto & node : this->node_ptr_range())
2028 for (
const auto & other_elem : other_mesh.element_ptr_range())
2033 if (!other_elem->topologically_equal(*elem))
2036 for (
const auto & elem : this->element_ptr_range())
2046 dof_id_type n_local_rows=0, n_unpartitioned_rows=0;
2051 n_unpartitioned_rows++;
2056 this->
comm().
sum(n_local_rows);
2058 return n_unpartitioned_rows + n_local_rows;
2065 LOG_SCOPE(
"copy_constraint_rows(mesh)",
"MeshBase");
2070 for (
const auto & [other_node, other_node_constraints] : other_constraint_rows)
2072 const Node *
const our_node = this->
node_ptr(other_node->id());
2074 for (
const auto & [other_inner_key_pair, constraint_value] : other_node_constraints)
2076 const auto & [other_elem, local_node_id] = other_inner_key_pair;
2077 const Elem *
const our_elem = this->
elem_ptr(other_elem->id());
2078 our_node_constraints.emplace_back(std::make_pair(our_elem, local_node_id), constraint_value);
2085 template <
typename T>
2088 bool precondition_constraint_operator)
2090 LOG_SCOPE(
"copy_constraint_rows(mat)",
"MeshBase");
2100 libmesh_error_msg_if(this->
n_nodes() != constraint_operator.
m(),
2101 "Constraint operator matrix with " <<
2102 constraint_operator.
m() <<
2103 "rows does not match this mesh with " <<
2113 std::map<dof_id_type, dof_id_type> existing_unconstrained_columns;
2114 std::set<dof_id_type> existing_unconstrained_nodes;
2120 std::vector<std::pair<dof_id_type, Real>>>
2122 columns_type columns(constraint_operator.
n());
2130 std::unordered_map<dof_id_type, Real> column_sums;
2136 std::vector<numeric_index_type> indices;
2137 std::vector<T> values;
2139 constraint_operator.
get_row(i, indices, values);
2140 libmesh_assert_equal_to(indices.size(), values.size());
2142 if (indices.size() == 1 &&
2148 if (existing_unconstrained_columns.find(indices[0]) !=
2149 existing_unconstrained_columns.end())
2151 const auto j = indices[0];
2152 columns[j].emplace_back(i, 1);
2156 existing_unconstrained_nodes.insert(i);
2157 existing_unconstrained_columns.emplace(indices[0],i);
2163 const auto j = indices[jj];
2165 libmesh_assert_equal_to(coef, values[jj]);
2166 columns[j].emplace_back(i, coef);
2174 std::vector<columns_type> all_columns;
2179 for (
auto & [j, subcol] : all_columns[p])
2180 for (
auto [i, v] : subcol)
2181 columns[j].emplace_back(i,v);
2186 std::unordered_map<const Node *, std::pair<dof_id_type, unsigned int>> node_to_elem_ptrs;
2193 for (
const Elem * elem : this->element_ptr_range())
2198 const Node * node = elem->node_ptr(n);
2199 if (existing_unconstrained_nodes.count(node->
id()))
2200 node_to_elem_ptrs.emplace(node, std::make_pair(elem->id(), n));
2206 for (
auto j :
make_range(constraint_operator.
n()))
2209 if (existing_unconstrained_columns.count(j))
2218 Real total_scaling = 0;
2219 unsigned int total_entries = 0;
2223 std::map<processor_id_type, int> pids;
2225 auto & column = columns[j];
2226 for (
auto [i, r] : column)
2229 const Point constrained_pt = constrained_node;
2230 newpt += r*constrained_pt;
2236 if (precondition_constraint_operator)
2237 column_sums[j] = total_scaling;
2239 libmesh_error_msg_if
2241 "Empty column " << j <<
2242 " found in constraint operator matrix");
2248 newpt /= total_scaling;
2250 newpt /= total_entries;
2254 elem->set_node(0, n);
2255 elem->subdomain_id() = new_sbd_id;
2265 node_to_elem_ptrs.emplace(n, std::make_pair(added_elem->
id(), 0));
2266 existing_unconstrained_columns.emplace(j,n->
id());
2272 for (
auto [pid, count] : pids)
2273 if (count >= n_pids)
2284 std::vector<std::pair<std::pair<dof_id_type, unsigned int>,
Real>>>
2285 indexed_constraint_rows;
2290 if (existing_unconstrained_nodes.count(i))
2293 std::vector<numeric_index_type> indices;
2294 std::vector<T> values;
2296 constraint_operator.
get_row(i, indices, values);
2298 std::vector<std::pair<std::pair<dof_id_type, unsigned int>,
Real>> constraint_row;
2303 existing_unconstrained_columns[indices[jj]];
2309 auto p = node_to_elem_ptrs[&constraining_node];
2312 libmesh_assert_equal_to(coef, values[jj]);
2318 if (precondition_constraint_operator)
2319 if (
auto sum_it = column_sums.find(indices[jj]);
2320 sum_it != column_sums.end())
2322 const Real scaling = sum_it->second;
2328 constraint_row.emplace_back(std::make_pair(p, coef));
2331 indexed_constraint_rows.emplace(i, std::move(constraint_row));
2337 for (
auto & [node_id, indexed_row] : indexed_constraint_rows)
2343 for (
auto [p, coef] : indexed_row)
2346 constraint_row.emplace_back
2347 (std::make_pair(std::make_pair(elem, p.second), coef));
2351 std::move(constraint_row));
2357 bool print_nonlocal)
const 2359 parallel_object_only();
2361 std::string local_constraints =
2366 this->
comm().
send(0, local_constraints);
2370 os <<
"Processor 0:\n";
2371 os << local_constraints;
2376 os <<
"Processor " << p <<
":\n";
2377 os << local_constraints;
2386 std::ostringstream os;
2393 os <<
"Mesh Constraint Rows:" 2398 const bool local = (node->processor_id() == this->
processor_id());
2401 if (!print_nonlocal && !local)
2404 os <<
"Constraints for " << (local ?
"Local" :
"Ghost") <<
" Node " << node->id()
2407 for (
const auto & [elem_and_node, coef] : row)
2408 os <<
" ((" << elem_and_node.first->id() <<
',' << elem_and_node.second <<
"), " << coef <<
")\t";
2420 template LIBMESH_EXPORT
void 2422 bool precondition_constraint_operator);
2424 #ifdef LIBMESH_USE_COMPLEX_NUMBERS 2425 template LIBMESH_EXPORT
void 2427 bool precondition_constraint_operator);
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
bool operator==(const MeshBase &other_mesh) const
This tests for exactly-equal data in all the senses that a mathematician would care about (element co...
std::set< subdomain_id_type > _mesh_subdomains
We cache the subdomain ids of the elements present in the mesh.
bool has_node_integer(std::string_view name) const
std::vector< unsigned int > add_elem_integers(const std::vector< std::string > &names, bool allocate_data=true, const std::vector< dof_id_type > *default_values=nullptr)
Register integer data (of type dof_id_type) to be added to each element in the mesh, one string name for each new integer.
bool closed()
Checks that the library has been closed.
ElemType
Defines an enum for geometric element types.
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
const std::set< boundary_id_type > & get_side_boundary_ids() const
virtual bool subclass_locally_equals(const MeshBase &other_mesh) const =0
Shim to allow operator == (&) to behave like a virtual function without having to be one...
bool _skip_renumber_nodes_and_elements
If this is true then renumbering will be kept to a minimum.
Order
defines an enum for polynomial orders.
constraint_rows_type & get_constraint_rows()
Constraint rows accessors.
virtual dof_id_type n_active_elem() const =0
A Node is like a Point, but with more information.
This abstract base class defines the interface by which library code and user code can report associa...
const MeshBase & interior_mesh() const
Real get_point_locator_close_to_point_tol() const
dof_id_type n_elem_on_proc(const processor_id_type proc) const
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
std::vector< std::string > _elem_integer_names
The array of names for integer data associated with each element in the mesh.
std::unique_ptr< PointLocatorBase > sub_point_locator() const
void remove_orphaned_nodes()
Removes any orphaned nodes, nodes not connected to any elements.
std::vector< dof_id_type > get_elemset_codes() const
Return a vector of all elemset codes defined on the mesh.
virtual void all_second_order_range(const SimpleRange< element_iterator > &range, const bool full_ordered=true)=0
Converts a set of this Mesh's elements defined by range from FIRST order to SECOND order...
The IntRange templated class is intended to make it easy to loop over integers which are indices of a...
void detect_interior_parents()
Search the mesh for elements that have a neighboring element of dim+1 and set that element as the int...
dof_id_type n_active_elem_on_proc(const processor_id_type proc) const
static constexpr Real TOLERANCE
MeshBase * _interior_mesh
Defaulting to this, a pointer to the mesh used to generate boundary elements on this.
std::vector< std::pair< std::pair< const Elem *, unsigned int >, Real > > constraint_rows_mapped_type
subdomain_id_type n_local_subdomains() const
void copy_cached_data(const MeshBase &other_mesh)
Helper class to copy cached data, to synchronize with a possibly unprepared other_mesh.
void set_spatial_dimension(unsigned char d)
Sets the "spatial dimension" of the Mesh.
std::map< dof_id_type, const MeshBase::elemset_type * > _elemset_codes
Map from "element set code" to list of set ids to which that element belongs (and vice-versa)...
bool skip_noncritical_partitioning() const
unsigned int add_elem_integer(std::string name, bool allocate_data=true, dof_id_type default_value=DofObject::invalid_id)
Register an integer datum (of type dof_id_type) to be added to each element in the mesh...
void add_elemset_code(dof_id_type code, MeshBase::elemset_type id_set)
Tabulate a user-defined "code" for elements which belong to the element sets specified in id_set...
const std::map< boundary_id_type, std::string > & get_sideset_name_map() const
bool has_elem_integer(std::string_view name) const
bool get_count_lower_dim_elems_in_point_locator() const
Get the current value of _count_lower_dim_elems_in_point_locator.
void remove_ghosting_functor(GhostingFunctor &ghosting_functor)
Removes a functor which was previously added to the set of ghosting functors.
This is the base class from which all geometric element types are derived.
virtual BoundingBox loose_bounding_box() const
dof_id_type n_local_nodes() const
virtual std::unique_ptr< Partitioner > & partitioner()
A partitioner to use at each prepare_for_use()
constraint_rows_type _constraint_rows
std::ostream & operator<<(std::ostream &os, const OrderWrapper &order)
Overload stream operators.
void copy_constraint_rows(const MeshBase &other_mesh)
Copy the constraints from the other mesh to this mesh.
const Parallel::Communicator & comm() const
std::vector< unsigned int > add_node_integers(const std::vector< std::string > &names, bool allocate_data=true, const std::vector< dof_id_type > *default_values=nullptr)
Register integer data (of type dof_id_type) to be added to each node in the mesh. ...
bool _skip_noncritical_partitioning
If this is true then no partitioning should be done with the possible exception of orphaned nodes...
unsigned char _spatial_dimension
The "spatial dimension" of the Mesh.
std::unique_ptr< BoundaryInfo > boundary_info
This class holds the boundary information.
bool _allow_remote_element_removal
If this is false then even on DistributedMesh remote elements will not be deleted during mesh prepara...
The libMesh namespace provides an interface to certain functionality in the library.
std::map< MeshBase::elemset_type, dof_id_type > _elemset_codes_inverse_map
dof_id_type get_elemset_code(const MeshBase::elemset_type &id_set) const
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
dof_id_type n_unpartitioned_elem() const
std::vector< std::string > _node_integer_names
The array of names for integer data associated with each node in the mesh.
bool in_threads
A boolean which is true iff we are in a Threads:: function It may be useful to assert(!Threadsin_thre...
Real distance(const Point &p)
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
Add a new Node at Point p to the end of the vertex array, with processor_id procid.
dof_id_type n_local_elem() const
MeshBase::elemset_type _all_elemset_ids
dof_id_type n_constraint_rows() const
bool nodes_and_elements_equal(const MeshBase &other_mesh) const
Tests for equality of all elements and nodes in the mesh.
uint8_t processor_id_type
This is the MeshBase class.
virtual numeric_index_type row_stop() const =0
subdomain_id_type get_id_by_name(std::string_view name) const
void get_elemsets(dof_id_type elemset_code, MeshBase::elemset_type &id_set_to_fill) const
Look up the element sets for a given elemset code and vice-versa.
void change_elemset_code(dof_id_type old_code, dof_id_type new_code)
Replace elemset code "old_code" with "new_code".
virtual void all_complete_order_range(const SimpleRange< element_iterator > &range)=0
Converts a set of elements in this (conforming, non-refined) mesh into "complete" order elements...
const std::set< boundary_id_type > & get_node_boundary_ids() const
ParallelObject & operator=(const ParallelObject &libmesh_dbg_var(other))
"Assignment" operator.
std::set< GhostingFunctor * > _ghosting_functors
The list of all GhostingFunctor objects to be used when distributing a DistributedMesh.
Order supported_nodal_order() const
unique_id_type _next_unique_id
The next available unique id for assigning ids to DOF objects.
unsigned int _n_parts
The number of partitions the mesh has.
unsigned int get_elem_integer_index(std::string_view name) const
processor_id_type n_processors() const
virtual bool is_serial() const
static const subdomain_id_type invalid_subdomain_id
A static integral constant representing an invalid subdomain id.
virtual void find_neighbors(const bool reset_remote_elements=false, const bool reset_current_list=true)=0
Locate element face (edge in 2D) neighbors.
Status receive(const unsigned int dest_processor_id, T &buf, const MessageTag &tag=any_tag) const
const std::map< boundary_id_type, std::string > & get_nodeset_name_map() const
unsigned char _default_mapping_data
The default mapping data (unused with Lagrange, used for nodal weight lookup index with rational base...
virtual Order supported_nodal_order() const
void min(const T &r, T &o, Request &req) const
const std::map< subdomain_id_type, std::string > & get_subdomain_name_map() const
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
dof_id_type n_sub_elem() const
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static std::unique_ptr< Elem > build(const ElemType type, Elem *p=nullptr)
ElemMappingType _default_mapping_type
The default mapping type (typically Lagrange) between master and physical space to assign to newly ad...
virtual void update_parallel_id_counts()=0
Updates parallel caches so that methods like n_elem() accurately reflect changes on other processors...
void print_info(std::ostream &os=libMesh::out, const unsigned int verbosity=0, const bool global=true) const
Prints relevant information about the mesh.
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
Helper function that allows us to treat a homogenous pair as a range.
virtual numeric_index_type m() const =0
std::unique_ptr< Partitioner > _partitioner
A partitioner to use at each prepare_for_use().
virtual const Node * query_node_ptr(const dof_id_type i) const =0
const std::map< boundary_id_type, std::string > & get_edgeset_name_map() const
virtual dof_id_type max_elem_id() const =0
void clear_point_locator()
Releases the current PointLocator object.
void subdomain_ids(std::set< subdomain_id_type > &ids, const bool global=true) const
Constructs a list of all subdomain identifiers in the local mesh if global == false, and in the global mesh if global == true (default).
std::set< Order > _elem_default_orders
We cache the (default) order of the geometric elements present in the mesh.
bool allow_find_neighbors() const
The BoundaryInfo class contains information relevant to boundary conditions including storing faces...
static const dof_id_type invalid_id
An invalid id to distinguish an uninitialized DofObject.
bool _skip_find_neighbors
If this is true then we will skip find_neighbors in prepare_for_use.
virtual void delete_node(Node *n)=0
Removes the Node n from the mesh.
std::pair< std::vector< unsigned int >, std::vector< unsigned int > > merge_extra_integer_names(const MeshBase &other)
Merge extra-integer arrays from an other mesh.
unsigned int n_elemsets() const
Returns the number of unique elemset ids which have been added via add_elemset_code(), which is the size of the _all_elemset_ids set.
Order _supported_nodal_order
We cache the maximum nodal order supported by all the mesh's elements (the minimum supported_nodal_or...
std::set< unsigned char > _elem_dims
We cache the dimension of the elements present in the mesh.
unsigned int get_node_integer_index(std::string_view name) const
std::string & subdomain_name(subdomain_id_type id)
const std::set< unsigned char > & elem_dimensions() const
Helper for building element sides that minimizes the construction of new elements.
dof_id_type n_nodes_on_proc(const processor_id_type proc) const
std::string get_local_constraints(bool print_nonlocal=false) const
Gets a string reporting all mesh constraint rows local to this processor.
void size_node_extra_integers()
Size extra-integer arrays of all nodes in the mesh.
An object whose state is distributed along a set of processors.
void regenerate_id_sets()
Clears and regenerates the cached sets of ids.
virtual void clear()
Deletes all the element and node data that is currently stored.
std::string enum_to_string(const T e)
bool skip_partitioning() const
std::map< GhostingFunctor *, std::shared_ptr< GhostingFunctor > > _shared_functors
Hang on to references to any GhostingFunctor objects we were passed in shared_ptr form...
SimpleRange< NodeRefIter > node_ref_range()
Returns a range with all nodes of an element, usable in range-based for loops.
unsigned int n_partitions() const
Defines a Cartesian bounding box by the two corner extremum.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
std::vector< dof_id_type > _node_integer_default_values
The array of default initialization values for integer data associated with each node in the mesh...
unsigned int recalculate_n_partitions()
In a few (very rare) cases, the user may have manually tagged the elements with specific processor ID...
std::set< GhostingFunctor * >::const_iterator ghosting_functors_begin() const
Beginning of range of ghosting functors.
virtual void get_row(numeric_index_type i, std::vector< numeric_index_type > &indices, std::vector< T > &values) const =0
Get a row from the matrix.
std::set< elemset_id_type > elemset_type
Typedef for the "set" container used to store elemset ids.
std::vector< dof_id_type > _elem_integer_default_values
The array of default initialization values for integer data associated with each element in the mesh...
subdomain_id_type n_subdomains() const
std::string get_info(const unsigned int verbosity=0, const bool global=true) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual const Elem * query_elem_ptr(const dof_id_type i) const =0
void max(const T &r, T &o, Request &req) const
void post_dofobject_moves(MeshBase &&other_mesh)
Moves any superclass data (e.g.
Temporarily serialize a DistributedMesh for non-distributed-mesh capable code paths.
virtual void all_complete_order()
Calls the range-based version of this function with a range consisting of all elements in the mesh...
virtual numeric_index_type row_start() const =0
unsigned int spatial_dimension() const
void send(const unsigned int dest_processor_id, const T &buf, const MessageTag &tag=no_tag) const
std::map< const Node *, constraint_rows_mapped_type > constraint_rows_type
void print_constraint_rows(std::ostream &os=libMesh::out, bool print_nonlocal=false) const
Prints (from processor 0) all mesh constraint rows.
const std::set< boundary_id_type > & get_edge_boundary_ids() const
virtual bool is_replicated() const
bool _is_prepared
Flag indicating if the mesh has been prepared for use.
virtual const Elem & elem_ref(const dof_id_type i) const
void set_count_lower_dim_elems_in_point_locator(bool count_lower_dim_elems)
In the point locator, do we count lower dimensional elements when we refine point locator regions...
virtual Real volume() const
MeshBase(const Parallel::Communicator &comm_in, unsigned char dim=1)
Constructor.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
void all_second_order(const bool full_ordered=true)
Calls the range-based version of this function with a range consisting of all elements in the mesh...
dof_id_type n_active_sub_elem() const
Same as n_sub_elem(), but only counts active elements.
The DofObject defines an abstract base class for objects that have degrees of freedom associated with...
virtual ~MeshBase()
Destructor.
void set_elem_dimensions(std::set< unsigned char > elem_dims)
Most of the time you should not need to call this, as the element dimensions will be set automaticall...
unsigned int mesh_dimension() const
std::unique_ptr< GhostingFunctor > _default_ghosting
The default geometric GhostingFunctor, used to implement standard libMesh element ghosting behavior...
void change_elemset_id(elemset_id_type old_id, elemset_id_type new_id)
Replace elemset id "old_id" with "new_id".
bool initialized()
Checks that library initialization has been done.
unsigned int add_node_integer(std::string name, bool allocate_data=true, dof_id_type default_value=DofObject::invalid_id)
Register an integer datum (of type dof_id_type) to be added to each node in the mesh.
virtual std::unique_ptr< GhostingFunctor > clone() const
A clone() is needed because GhostingFunctor can not be shared between different meshes.
virtual const Node & node_ref(const dof_id_type i) const
bool _skip_all_partitioning
If this is true then no partitioning should be done.
bool on_command_line(std::string arg)
bool allow_renumbering() const
virtual void delete_remote_elements()
When supported, deletes all nonlocal elements of the mesh except for "ghosts" which touch a local ele...
void union_with(const Point &p)
Enlarges this bounding box to include the given point.
virtual void update_post_partitioning()
Recalculate any cached data after elements and nodes have been repartitioned.
void reinit_ghosting_functors()
Loops over ghosting functors and calls mesh_reinit()
virtual dof_id_type n_elem() const =0
virtual std::unique_ptr< Elem > build_edge_ptr(const unsigned int i)=0
virtual const Node * node_ptr(const dof_id_type i) const =0
processor_id_type processor_id() const
virtual Order default_order() const =0
std::unique_ptr< PointLocatorBase > _point_locator
A PointLocator class for this mesh.
virtual void redistribute()
Redistribute elements between processors.
processor_id_type processor_id() const
Real _point_locator_close_to_point_tol
If nonzero, we will call PointLocatorBase::set_close_to_point_tol() on any PointLocators that we crea...
virtual ElemType type() const =0
static std::unique_ptr< PointLocatorBase > build(PointLocatorType t, const MeshBase &mesh, const PointLocatorBase *master=nullptr)
Builds an PointLocator for the mesh mesh.
bool _count_lower_dim_elems_in_point_locator
Do we count lower dimensional elements in point locator refinement? This is relevant in tree-based po...
A Point defines a location in LIBMESH_DIM dimensional Real space.
bool locally_equals(const MeshBase &other_mesh) const
This behaves the same as operator==, but only for the local and ghosted aspects of the mesh; i...
std::map< subdomain_id_type, std::string > _block_id_to_name
This structure maintains the mapping of named blocks for file formats that support named blocks...
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
virtual dof_id_type n_nodes() const =0
virtual void renumber_nodes_and_elements()=0
After partitioning a mesh it is useful to renumber the nodes and elements so that they lie in contigu...
This class implements the original default geometry ghosting requirements in libMesh: point neighbors...
void size_elem_extra_integers()
Size extra-integer arrays of all elements in the mesh.
virtual numeric_index_type n() const =0
std::set< GhostingFunctor * >::const_iterator ghosting_functors_end() const
End of range of ghosting functors.
void add_ghosting_functor(GhostingFunctor &ghosting_functor)
Adds a functor which can specify ghosting requirements for use on distributed meshes.
MeshBase & operator=(const MeshBase &)=delete
Copy and move assignment are not allowed because MeshBase subclasses manually manage memory (Elems an...
void set_point_locator_close_to_point_tol(Real val)
Set value used by PointLocatorBase::close_to_point_tol().
void set_union(T &data, const unsigned int root_id) const