Skip to content

Commit

Permalink
remap_subdomain_ids capability in stitch_meshes
Browse files Browse the repository at this point in the history
I'm going to try to default it to "true", on the theory that if we don't
see a lot of CI failures then this will be a cheap feature to enable by
default whereas if we do see a lot of CI failures then this will be an
important feature to enable by default.
  • Loading branch information
roystgnr committed Aug 7, 2023
1 parent 5dd7c84 commit 81303ec
Show file tree
Hide file tree
Showing 2 changed files with 154 additions and 9 deletions.
23 changes: 20 additions & 3 deletions include/mesh/unstructured_mesh.h
Original file line number Diff line number Diff line change
Expand Up @@ -199,6 +199,15 @@ class UnstructuredMesh : public MeshBase
* and other_mesh and node IDs in the stitched mesh because the number of nodes (and hence
* the node IDs) in the stitched mesh depend on how many nodes are stitched.
*
* If \p remap_subdomain_ids is true then we assume that some
* subdomain ids might have been autogenerated, so we remap them as
* necessary, treating subdomain names as the important thing for
* consistency; if we have missing names and cannot infer a
* consistent resolution to an id conflict then we exit with an
* error. If \p remap_subdomain_ids is false then we revert to the
* older libMesh behavior: leave all subdomain ids alone and woe
* unto you if you weren't keeping track of them.
*
* \returns the count of how many nodes were merged between the two meshes.
* This can be zero in the case of no matching nodes or if
* \p merge_boundary_nodes_all_or_nothing was active and relevant.
Expand All @@ -211,7 +220,8 @@ class UnstructuredMesh : public MeshBase
bool verbose=true,
bool use_binary_search=true,
bool enforce_all_nodes_match_on_boundaries=false,
bool merge_boundary_nodes_all_or_nothing=false);
bool merge_boundary_nodes_all_or_nothing=false,
bool remap_subdomain_ids=true);

/**
* Similar to stitch_meshes, except that we stitch two adjacent surfaces within this mesh.
Expand All @@ -234,12 +244,18 @@ class UnstructuredMesh : public MeshBase
* \p other_mesh has element or node extra_integer data, any names
* for that data which do not already exist on \p this mesh are
* added so that all such data can be copied.
*
* If an \p id_remapping map is provided, then element subdomain ids
* in \p other_mesh will be converted using it before adding them to
* \p this mesh.
*/
virtual void copy_nodes_and_elements (const MeshBase & other_mesh,
const bool skip_find_neighbors = false,
dof_id_type element_id_offset = 0,
dof_id_type node_id_offset = 0,
unique_id_type unique_id_offset = 0);
unique_id_type unique_id_offset = 0,
std::unordered_map<subdomain_id_type, subdomain_id_type> *
id_remapping = nullptr);

/**
* Move node and elements from other_mesh to this mesh.
Expand Down Expand Up @@ -277,7 +293,8 @@ class UnstructuredMesh : public MeshBase
bool use_binary_search,
bool enforce_all_nodes_match_on_boundaries,
bool skip_find_neighbors,
bool merge_boundary_nodes_all_or_nothing);
bool merge_boundary_nodes_all_or_nothing,
bool remap_subdomain_ids);
};


Expand Down
140 changes: 134 additions & 6 deletions src/mesh/unstructured_mesh.C
Original file line number Diff line number Diff line change
Expand Up @@ -316,7 +316,9 @@ void UnstructuredMesh::copy_nodes_and_elements(const MeshBase & other_mesh,
#ifdef LIBMESH_ENABLE_UNIQUE_ID
unique_id_offset
#endif
)
,
std::unordered_map<subdomain_id_type, subdomain_id_type> *
id_remapping)
{
LOG_SCOPE("copy_nodes_and_elements()", "UnstructuredMesh");

Expand Down Expand Up @@ -388,7 +390,14 @@ void UnstructuredMesh::copy_nodes_and_elements(const MeshBase & other_mesh,
nullptr;
auto el = Elem::build(old->type(), newparent);

el->subdomain_id() = old->subdomain_id();
subdomain_id_type sbd_id = old->subdomain_id();
if (id_remapping)
{
auto remapping_it = id_remapping->find(sbd_id);
if (remapping_it != id_remapping->end())
sbd_id = remapping_it->second;
}
el->subdomain_id() = sbd_id;

// Hold off on trying to set the interior parent because we may actually
// add lower dimensional elements before their interior parents
Expand Down Expand Up @@ -1730,7 +1739,8 @@ UnstructuredMesh::stitch_meshes (const MeshBase & other_mesh,
bool verbose,
bool use_binary_search,
bool enforce_all_nodes_match_on_boundaries,
bool merge_boundary_nodes_all_or_nothing)
bool merge_boundary_nodes_all_or_nothing,
bool remap_subdomain_ids)
{
LOG_SCOPE("stitch_meshes()", "UnstructuredMesh");
return stitching_helper(&other_mesh,
Expand All @@ -1742,7 +1752,8 @@ UnstructuredMesh::stitch_meshes (const MeshBase & other_mesh,
use_binary_search,
enforce_all_nodes_match_on_boundaries,
true,
merge_boundary_nodes_all_or_nothing);
merge_boundary_nodes_all_or_nothing,
remap_subdomain_ids);
}


Expand All @@ -1766,7 +1777,8 @@ UnstructuredMesh::stitch_surfaces (boundary_id_type boundary_id_1,
use_binary_search,
enforce_all_nodes_match_on_boundaries,
true,
merge_boundary_nodes_all_or_nothing);
merge_boundary_nodes_all_or_nothing,
false);
}


Expand All @@ -1780,7 +1792,8 @@ UnstructuredMesh::stitching_helper (const MeshBase * other_mesh,
bool use_binary_search,
bool enforce_all_nodes_match_on_boundaries,
bool skip_find_neighbors,
bool merge_boundary_nodes_all_or_nothing)
bool merge_boundary_nodes_all_or_nothing,
bool remap_subdomain_ids)
{
#ifdef DEBUG
// We rely on neighbor links here
Expand Down Expand Up @@ -2172,6 +2185,121 @@ UnstructuredMesh::stitching_helper (const MeshBase * other_mesh,
for (auto & entry : pr.second)
entry += elem_delta;

// We run into problems when the libMesh subdomain standard (the
// id defines the subdomain; the name was an afterthought) and
// the MOOSE standard (the name defines the subdomain; the id
// might be autogenerated) clash.
//
// Subdomain ids with the same name in both meshes are surely
// meant to represent the same subdomain. We can just merge
// them.
//
// Subdomain ids which don't have a name in either mesh are
// almost surely meant to represent the same subdomain. We'll
// just merge them.
//
// Subdomain ids with different names in different meshes, or
// names with different ids in different meshes, are trickier.
// For backwards compatibility we default to the old "just copy
// all the subdomain ids over" behavior, but if requested we'll
// remap any ids that appear to be clear conflicts, and we'll
// scream and die if we see any ids that are ambiguous due to
// being named in one mesh but not the other.
std::unordered_map<subdomain_id_type, subdomain_id_type> id_remapping;
if (remap_subdomain_ids)
{
const auto & this_map = this->get_subdomain_name_map();
const auto & other_map = other_mesh->get_subdomain_name_map();
std::unordered_map<std::string, subdomain_id_type> other_map_reversed;
for (auto & [sid, sname] : other_map)
other_map_reversed.emplace(sname, sid);

std::unordered_map<std::string, subdomain_id_type> this_map_reversed;
for (auto & [sid, sname] : this_map)
this_map_reversed.emplace(sname, sid);

// We don't require either mesh to be prepared, but that
// means we need to check for subdomains manually.
auto get_subdomains = [](const MeshBase & mesh) {
std::set<subdomain_id_type> all_subdomains;
for (auto & el : mesh.element_ptr_range())
all_subdomains.insert(el->subdomain_id());
return all_subdomains;
};

const auto this_subdomains = get_subdomains(*this);
const auto other_subdomains = get_subdomains(*other_mesh);

for (auto & [sid, sname] : this_map)
{
auto other_reverse_it = other_map_reversed.find(sname);

// The same name with the same id means we're fine. The
// same name with another id means we remap their id to
// ours
if (other_reverse_it != other_map_reversed.end())
{
if (other_reverse_it->second != sid)
id_remapping[other_reverse_it->second] = sid;
}

// The same id with a different name, we'll get to
// later. The same id without any name means we don't
// know what the user wants.
if (other_subdomains.count(sid) && !other_map.count(sid))
libmesh_error_msg("Can't safely stitch with a mesh sharing subdomain id "
<< sid << " but not subdomain name " << sname);
}

subdomain_id_type next_free_id = 0;
// We might try to stitch empty meshes ...
if (!this_subdomains.empty())
next_free_id = *this_subdomains.rbegin() + 1;
if (!other_subdomains.empty())
next_free_id =
std::max(next_free_id,
cast_int<subdomain_id_type>
(*other_subdomains.rbegin() + 1));

for (auto & [sid, sname] : other_map)
{
auto reverse_it = this_map_reversed.find(sname);

// At this point we've figured out any remapping
// necessary for an sname that we share. And we don't
// need to remap any sid we don't share.
if (reverse_it == this_map_reversed.end())
{
// But if we don't have this sname and we do have this
// sid then we can't just merge into that.
if (this_subdomains.count(sid))
{
// If we have this sid with no name, we don't
// know what the user wants.
if (this_map.count(sid))
libmesh_error_msg("Can't safely stitch with a mesh sharing subdomain id "
<< sid << " but under subdomain name " << sname);

// We have this sid under a different name, so
// we just need to give the other elements a new
// id.

// Users might have done crazy things with id
// choice so let's make sure they didn't get too
// crazy.
libmesh_error_msg_if ((!this_subdomains.empty() &&
next_free_id < *this_subdomains.rbegin()) ||
(!other_subdomains.empty() &&
next_free_id < *other_subdomains.rbegin()),
"Subdomain id overflow");

id_remapping[sid] = next_free_id++;
this->subdomain_name(next_free_id) = sname;
}
}
}
}

// Copy mesh data. If we skip the call to find_neighbors(), the lists
// of neighbors will be copied verbatim from the other mesh
this->copy_nodes_and_elements(*other_mesh, skip_find_neighbors,
Expand Down

0 comments on commit 81303ec

Please sign in to comment.