10#include <range/v3/algorithm/equal.hpp>
11#include <range/v3/range/conversion.hpp>
12#include <unordered_map>
31 auto const& nodes = subdomain_mesh.
getNodes();
32 std::vector<MathLib::Point3dWithID*> subdomain_points{begin(nodes),
35 auto const& bulk_node_ids =
41 "Expected to find exactly one node in the bulk mesh for each node "
42 "of the subdomain; Found {:d} nodes in the bulk mesh out of {:d} "
43 "nodes in the subdomain.",
54 std::vector<std::size_t>
const& node_ids,
55 std::vector<std::vector<std::size_t>>
const& connected_element_ids_per_node)
60 std::unordered_map<std::size_t, int> element_counts(8);
61 for (
auto const node_id : node_ids)
63 for (
auto const element_id : connected_element_ids_per_node[node_id])
65 element_counts[element_id]++;
73 auto const nnodes = node_ids.size();
74 std::vector<std::size_t> element_ids;
75 for (
auto const& pair : element_counts)
77 if (pair.second ==
static_cast<int>(nnodes))
79 element_ids.push_back(pair.first);
95 std::vector<std::vector<std::size_t>> bulk_element_ids_map(
99 std::vector<std::vector<std::size_t>> connected_element_ids_per_node(
103 connected_element_ids_per_node[node_id] =
108 auto const& elements = subdomain_mesh.
getElements();
109#pragma omp parallel for
110 for (std::ptrdiff_t j = 0; j < std::ssize(elements); ++j)
112 auto*
const e = elements[j];
113 std::vector<std::size_t> element_node_ids(e->getNumberOfBaseNodes());
114 for (
unsigned n = 0; n < e->getNumberOfBaseNodes(); ++n)
118 std::vector<std::size_t> element_node_ids_bulk(
119 e->getNumberOfBaseNodes());
120 std::transform(begin(element_node_ids), end(element_node_ids),
121 begin(element_node_ids_bulk),
122 [&bulk_node_ids](std::size_t
const id)
123 {
return bulk_node_ids[id]; });
126 element_node_ids_bulk, connected_element_ids_per_node);
128 if (bulk_element_ids.empty())
130 ERR(
"No element could be found for the subdomain element {:d}. "
131 "Corresponding bulk mesh node ids are:",
133 for (
auto const i : element_node_ids_bulk)
138 "Expect at least one element to be found in the bulk mesh.");
141 bulk_element_ids_map[e->getID()] = std::move(bulk_element_ids);
144 return bulk_element_ids_map;
150 std::vector<std::size_t>
const& values,
154 if (!properties.existsPropertyVector<std::size_t>(property_name))
156 addPropertyToMesh(mesh, property_name, mesh_item_type, 1, values);
163 auto& original_property =
165 if (ranges::equal(original_property, values))
168 "There is already a '{:s}' property present in the subdomain mesh "
169 "'{:s}' and it is equal to the newly computed values.",
170 property_name, mesh.
getName());
178 "There is already a '{:s}' property present in the subdomain mesh "
179 "'{:s}' and it is not equal to the newly computed values.",
183 if (!force_overwrite)
185 OGS_FATAL(
"The force overwrite flag was not specified, exiting.");
188 INFO(
"Overwriting '{:s}' property.", property_name);
189 original_property.assign(values.begin(), values.end());
198 bool const force_overwrite =
false)
202 auto const& bulk_node_ids =
203 identifySubdomainMeshNodes(subdomain_mesh, mesh_node_searcher);
204 INFO(
"identifySubdomainMesh(): identifySubdomainMeshNodes took {:g} s",
207 updateOrCheckExistingSubdomainProperty(
212 auto const& bulk_element_ids =
213 identifySubdomainMeshElements(subdomain_mesh, bulk_mesh);
214 INFO(
"identifySubdomainMesh(): identifySubdomainMeshElements took {:g} s",
222 if (all_of(begin(bulk_element_ids), end(bulk_element_ids),
223 [](std::vector<std::size_t>
const& v) {
return v.size() == 1; }))
227 std::vector<std::size_t> unique_bulk_element_ids;
228 unique_bulk_element_ids.reserve(bulk_element_ids.size());
229 transform(begin(bulk_element_ids), end(bulk_element_ids),
230 back_inserter(unique_bulk_element_ids),
231 [](std::vector<std::size_t>
const& v) {
return v[0]; });
233 updateOrCheckExistingSubdomainProperty(
244 std::vector<std::size_t> flat_bulk_element_ids;
245 flat_bulk_element_ids.reserve(2 * bulk_element_ids.size());
246 std::vector<std::size_t> number_of_bulk_element_ids;
247 number_of_bulk_element_ids.reserve(bulk_element_ids.size());
249 for (std::vector<std::size_t>
const& v : bulk_element_ids)
251 number_of_bulk_element_ids.push_back(v.size());
252 flat_bulk_element_ids.insert(end(flat_bulk_element_ids), begin(v),
256 updateOrCheckExistingSubdomainProperty(
257 subdomain_mesh,
"number_bulk_elements", number_of_bulk_element_ids,
259 updateOrCheckExistingSubdomainProperty(
Definition of the Element class.
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition of the Mesh class.
Definition of the Node class.
Definition of the RunTime class.
double elapsed() const
Get the elapsed time in seconds.
void start()
Start the timer.
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Properties & getProperties()
const std::string getName() const
Get name of the mesh.
std::size_t getNumberOfNodes() const
Get the number of nodes.
std::vector< Element const * > const & getElementsConnectedToNode(std::size_t node_id) const
std::size_t getNumberOfElements() const
Get the number of elements.
PropertyVector< T > const * getPropertyVector(std::string_view name) const
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
constexpr std::string_view getBulkIDString(MeshItemType mesh_item_type)
std::size_t getNodeIndex(Element const &element, unsigned const idx)
PropertyVector< std::size_t > const * bulkNodeIDs(Mesh const &mesh)
std::vector< std::size_t > identifySubdomainMeshNodes(MeshLib::Mesh const &subdomain_mesh, MeshGeoToolsLib::MeshNodeSearcher const &mesh_node_searcher)
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.
std::vector< std::vector< std::size_t > > identifySubdomainMeshElements(MeshLib::Mesh const &subdomain_mesh, MeshLib::Mesh const &bulk_mesh)
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)