anonymous_namespace{IdentifySubdomainMesh.cpp} Namespace Reference

Functions
std::vector< std::size_t >	identifySubdomainMeshNodes (MeshLib::Mesh const &subdomain_mesh, MeshGeoToolsLib::MeshNodeSearcher const &mesh_node_searcher)
std::vector< std::size_t >	findElementsInMesh (std::vector< std::size_t > const &node_ids, std::vector< std::vector< std::size_t > > const &connected_element_ids_per_node)
std::vector< std::vector< std::size_t > >	identifySubdomainMeshElements (MeshLib::Mesh const &subdomain_mesh, MeshLib::Mesh const &bulk_mesh)
void	updateOrCheckExistingSubdomainProperty (MeshLib::Mesh &mesh, std::string_view property_name, std::vector< std::size_t > const &values, MeshLib::MeshItemType const mesh_item_type, bool const force_overwrite)
	Updates or checks the existing mesh's property with the given values.

Function Documentation

findElementsInMesh()

std::vector< std::size_t > anonymous_namespace{IdentifySubdomainMesh.cpp}::findElementsInMesh	(	std::vector< std::size_t > const &	node_ids,
		std::vector< std::vector< std::size_t > > const &	connected_element_ids_per_node )

Find all elements which consist of all of the given node ids.

Note

It is recommended to include only base nodes of elements into the search node_ids.

Definition at line 49 of file IdentifySubdomainMesh.cpp.

{
    if (node_ids.empty())
    {
        return {};
    }
 
    // An element belongs to the subdomain element when it is connected to all
    // of the given nodes. Use the node with the fewest connected elements as
    // the candidate set and keep only those candidates that are connected to
    // every other node. This avoids any allocation that scales with the bulk
    // mesh size, so the work per call stays cheap inside the parallel loop
    // over subdomain elements. Only local state is used, hence thread-safe.
    auto const& smallest_elements =
        connected_element_ids_per_node[*std::min_element(
            begin(node_ids), end(node_ids),
            [&connected_element_ids_per_node](std::size_t const a,
                                              std::size_t const b)
            {
                return connected_element_ids_per_node[a].size() <
                       connected_element_ids_per_node[b].size();
            })];
 
    std::vector<std::size_t> element_ids;
    for (auto const element_id : smallest_elements)
    {
        bool const shared_by_all = std::all_of(
            begin(node_ids), end(node_ids),
            [&connected_element_ids_per_node,
             element_id](std::size_t const node_id)
            {
                auto const& elements = connected_element_ids_per_node[node_id];
                return std::find(begin(elements), end(elements), element_id) !=
                       end(elements);
            });
        if (shared_by_all)
        {
            element_ids.push_back(element_id);
        }
    }
 
    // Sort descending so the result is deterministic across compilers and
    // parallel execution ordering.
    std::sort(element_ids.begin(), element_ids.end(),
              std::greater<std::size_t>());
 
    return element_ids;
}

Referenced by identifySubdomainMeshElements().

identifySubdomainMeshElements()

std::vector< std::vector< std::size_t > > anonymous_namespace{IdentifySubdomainMesh.cpp}::identifySubdomainMeshElements	(	MeshLib::Mesh const &	subdomain_mesh,
		MeshLib::Mesh const &	bulk_mesh )

Tries to find all elements' ids in the bulk mesh. The subdomain elements' nodes are mapped to the bulk mesh and are then used to identify the bulk mesh elements.

Definition at line 103 of file IdentifySubdomainMesh.cpp.

{
    auto const& bulk_node_ids = *MeshLib::bulkNodeIDs(subdomain_mesh);
 
    // Allocate space for all elements for random insertion.
    std::vector<std::vector<std::size_t>> bulk_element_ids_map(
        subdomain_mesh.getNumberOfElements());
 
    // For each node a vector of connected element ids of that node.
    std::vector<std::vector<std::size_t>> connected_element_ids_per_node(
        bulk_mesh.getNumberOfNodes());
    for (auto const node_id : bulk_mesh.getNodes() | MeshLib::views::ids)
    {
        connected_element_ids_per_node[node_id] =
            bulk_mesh.getElementsConnectedToNode(node_id) |
            MeshLib::views::ids | ranges::to<std::vector>;
    }
 
    auto const& elements = subdomain_mesh.getElements();
#pragma omp parallel for
    for (std::ptrdiff_t j = 0; j < std::ssize(elements); ++j)
    {
        auto* const e = elements[j];
        std::vector<std::size_t> element_node_ids(e->getNumberOfBaseNodes());
        for (unsigned n = 0; n < e->getNumberOfBaseNodes(); ++n)
        {
            element_node_ids[n] = MeshLib::getNodeIndex(*e, n);
        }
        std::vector<std::size_t> element_node_ids_bulk(
            e->getNumberOfBaseNodes());
        std::transform(begin(element_node_ids), end(element_node_ids),
                       begin(element_node_ids_bulk),
                       [&bulk_node_ids](std::size_t const id)
                       { return bulk_node_ids[id]; });
 
        std::vector<std::size_t> bulk_element_ids = findElementsInMesh(
            element_node_ids_bulk, connected_element_ids_per_node);
 
        if (bulk_element_ids.empty())
        {
            ERR("No element could be found for the subdomain element {:d}. "
                "Corresponding bulk mesh node ids are:",
                e->getID());
            for (auto const i : element_node_ids_bulk)
            {
                ERR("\t{:d}", i);
            }
            OGS_FATAL(
                "Expect at least one element to be found in the bulk mesh.");
        }
 
        bulk_element_ids_map[e->getID()] = std::move(bulk_element_ids);
    }
 
    return bulk_element_ids_map;
}

References MeshLib::bulkNodeIDs(), ERR(), findElementsInMesh(), MeshLib::Mesh::getElements(), MeshLib::Mesh::getElementsConnectedToNode(), MeshLib::getNodeIndex(), MeshLib::Mesh::getNodes(), MeshLib::Mesh::getNumberOfElements(), MeshLib::Mesh::getNumberOfNodes(), MeshLib::views::ids, and OGS_FATAL.

identifySubdomainMeshNodes()

std::vector< std::size_t > anonymous_namespace{IdentifySubdomainMesh.cpp}::identifySubdomainMeshNodes	(	MeshLib::Mesh const &	subdomain_mesh,
		MeshGeoToolsLib::MeshNodeSearcher const &	mesh_node_searcher )

Find one-to-one mapping of subdomain nodes to the bulk mesh nodes and returns the map of ids.

Definition at line 22 of file IdentifySubdomainMesh.cpp.

{
    // Convert nodes pointers needed for the mesh_node_searcher algorithm.
    auto const& nodes = subdomain_mesh.getNodes();
    std::vector<MathLib::Point3dWithID*> subdomain_points{begin(nodes),
                                                          end(nodes)};
 
    auto const& bulk_node_ids =
        mesh_node_searcher.getMeshNodeIDs(subdomain_points);
 
    if (bulk_node_ids.size() != subdomain_mesh.getNumberOfNodes())
    {
        OGS_FATAL(
            "Expected to find exactly one node in the bulk mesh for each node "
            "of the subdomain; Found {:d} nodes in the bulk mesh out of {:d} "
            "nodes in the subdomain.",
            bulk_node_ids.size(), subdomain_mesh.getNumberOfNodes());
    }
 
    return bulk_node_ids;
}

References MeshGeoToolsLib::MeshNodeSearcher::getMeshNodeIDs(), MeshLib::Mesh::getNodes(), MeshLib::Mesh::getNumberOfNodes(), and OGS_FATAL.

updateOrCheckExistingSubdomainProperty()

void anonymous_namespace{IdentifySubdomainMesh.cpp}::updateOrCheckExistingSubdomainProperty	(	MeshLib::Mesh &	mesh,
		std::string_view	property_name,
		std::vector< std::size_t > const &	values,
		MeshLib::MeshItemType const	mesh_item_type,
		bool const	force_overwrite )

Updates or checks the existing mesh's property with the given values.

Definition at line 162 of file IdentifySubdomainMesh.cpp.

{
    auto& properties = mesh.getProperties();
    const bool property_name_exists =
        properties.hasPropertyVector(property_name);
    if (!properties.existsPropertyVector<std::size_t>(property_name,
                                                      mesh_item_type, 1))
    {
        if (property_name_exists)
        {
            if (!force_overwrite)
            {
                OGS_FATAL(
                    "A property named '{:s}' already exists on mesh '{:s}', "
                    "but it has a different mesh item type or number of "
                    "components. Use force overwrite to replace it.",
                    property_name, mesh.getName());
            }
 
            WARN(
                "A property named '{:s}' exists on mesh '{:s}' with a "
                "different mesh item type or number of components. "
                "Overwriting it.",
                property_name, mesh.getName());
            properties.removePropertyVector(property_name);
        }
 
        addPropertyToMesh<std::size_t>(mesh, property_name, mesh_item_type, 1,
                                       {values});
        return;
    }
 
    //
    // Check the existing property against new values.
    //
    auto& original_property =
        *properties.getPropertyVector<std::size_t>(property_name);
    if (ranges::equal(original_property, values))
    {
        INFO(
            "There is already a '{:s}' property present in the subdomain mesh "
            "'{:s}' and it is equal to the newly computed values.",
            property_name, mesh.getName());
        return;
    }
 
    //
    // Property differs. Notify and update if forced.
    //
    WARN(
        "There is already a '{:s}' property present in the subdomain mesh "
        "'{:s}' and it is not equal to the newly computed values.",
        property_name,
        mesh.getName());
 
    if (!force_overwrite)
    {
        OGS_FATAL("The force overwrite flag was not specified, exiting.");
    }
 
    INFO("Overwriting '{:s}' property.", property_name);
    original_property.assign(values);
}

References MeshLib::Mesh::getName(), MeshLib::Mesh::getProperties(), MeshLib::Properties::hasPropertyVector(), INFO(), OGS_FATAL, and WARN().