OGS
Mesh.cpp
Go to the documentation of this file.
1
15#include "Mesh.h"
16
17#include <memory>
18#include <range/v3/algorithm/contains.hpp>
19#include <range/v3/numeric.hpp>
20#include <range/v3/range/conversion.hpp>
21#include <range/v3/view/enumerate.hpp>
22#include <range/v3/view/indirect.hpp>
23#include <range/v3/view/map.hpp>
24#include <unordered_map>
25#include <utility>
26
27#include "BaseLib/RunTime.h"
28#include "Elements/Element.h"
29#include "Elements/Hex.h"
30#include "Elements/Prism.h"
31#include "Elements/Pyramid.h"
32#include "Elements/Quad.h"
33#include "Elements/Tet.h"
34#include "Elements/Tri.h"
37
39static std::size_t global_mesh_counter = 0;
40
41namespace MeshLib
42{
43using namespace ranges;
44
45std::vector<std::vector<Element const*>> findElementsConnectedToNodes(
46 Mesh const& mesh)
47{
48 std::vector<std::vector<Element const*>> elements_connected_to_nodes;
49 auto const& nodes = mesh.getNodes();
50 elements_connected_to_nodes.resize(nodes.size());
51
52 for (auto const* element : mesh.getElements())
53 {
54 for (auto const node_id : element->nodes() | views::ids)
55 {
56 elements_connected_to_nodes[node_id].push_back(element);
57 }
58 }
59 return elements_connected_to_nodes;
60}
61
62Mesh::Mesh(std::string name,
63 std::vector<Node*>
64 nodes,
65 std::vector<Element*>
66 elements,
67 Properties const& properties)
68 : _id(global_mesh_counter++),
69 _mesh_dimension(0),
70 _node_distance(std::numeric_limits<double>::max(), 0),
71 _name(std::move(name)),
72 _nodes(std::move(nodes)),
73 _elements(std::move(elements)),
74 _properties(properties)
75{
76 this->resetNodeIDs();
77 this->resetElementIDs();
78 this->setDimension();
79
81
82 this->setElementNeighbors();
83}
84
85Mesh::Mesh(const Mesh& mesh)
86 : _id(global_mesh_counter++),
87 _mesh_dimension(mesh.getDimension()),
88 _node_distance(mesh._node_distance.first, mesh._node_distance.second),
89 _name(mesh.getName()),
90 _nodes(mesh.getNumberOfNodes()),
91 _elements(mesh.getNumberOfElements()),
92 _properties(mesh._properties)
93{
94 const std::vector<Node*>& nodes(mesh.getNodes());
95 const std::size_t nNodes(nodes.size());
96 for (unsigned i = 0; i < nNodes; ++i)
97 {
98 _nodes[i] = new Node(*nodes[i]);
99 }
100
101 const std::vector<Element*>& elements(mesh.getElements());
102 const std::size_t nElements(elements.size());
103 for (unsigned i = 0; i < nElements; ++i)
104 {
105 _elements[i] = elements[i]->clone();
106 for (auto const& [j, node_id] :
107 elements[i]->nodes() | views::ids | ranges::views::enumerate)
108 {
109 _elements[i]->setNode(static_cast<unsigned>(j), _nodes[node_id]);
110 }
111 }
112
113 if (_mesh_dimension == 0)
114 {
115 this->setDimension();
116 }
118 this->setElementNeighbors();
119}
120
121Mesh::Mesh(Mesh&& mesh) = default;
122
124{
125 _elements.clear();
126 _nodes.clear();
127}
128
130{
131 const std::size_t nElements(_elements.size());
132 for (std::size_t i = 0; i < nElements; ++i)
133 {
134 delete _elements[i];
135 }
136
137 const std::size_t nNodes(_nodes.size());
138 for (std::size_t i = 0; i < nNodes; ++i)
139 {
140 delete _nodes[i];
141 }
142}
143
145{
146 _elements.push_back(elem);
147}
148
150{
151 const std::size_t nNodes(_nodes.size());
152 for (std::size_t i = 0; i < nNodes; ++i)
153 {
154 _nodes[i]->setID(i);
155 }
156}
157
159{
160 const std::size_t nElements(this->_elements.size());
161 for (unsigned i = 0; i < nElements; ++i)
162 {
163 _elements[i]->setID(i);
164 }
165}
166
168{
169 const std::size_t nElements(_elements.size());
170 for (unsigned i = 0; i < nElements; ++i)
171 {
173 {
174 _mesh_dimension = _elements[i]->getDimension();
175 }
176 }
177}
178
179std::pair<double, double> minMaxEdgeLength(
180 std::vector<Element*> const& elements)
181{
182 auto min_max = [](auto const a, auto const b) -> std::pair<double, double> {
183 return {std::min(a.first, b.first), std::max(a.second, b.second)};
184 };
185
186 using limits = std::numeric_limits<double>;
187 auto const bounds = ranges::accumulate(
188 elements, std::pair{limits::infinity(), -limits::infinity()}, min_max,
189 [](Element const* const e) { return computeSqrEdgeLengthRange(*e); });
190
191 return {std::sqrt(bounds.first), std::sqrt(bounds.second)};
192}
193
195{
196 std::vector<Element const*> neighbors;
197 for (auto element : _elements)
198 {
199 // create vector with all elements connected to current element
200 // (includes lots of doubles!)
201 const std::size_t nNodes(element->getNumberOfBaseNodes());
202 for (unsigned n(0); n < nNodes; ++n)
203 {
204 auto const& conn_elems(
205 _elements_connected_to_nodes[element->getNode(n)->getID()]);
206 neighbors.insert(neighbors.end(), conn_elems.begin(),
207 conn_elems.end());
208 }
209 std::sort(neighbors.begin(), neighbors.end());
210 auto const neighbors_new_end =
211 std::unique(neighbors.begin(), neighbors.end());
212
213 for (auto neighbor = neighbors.begin(); neighbor != neighbors_new_end;
214 ++neighbor)
215 {
216 std::optional<unsigned> const opposite_face_id =
217 element->addNeighbor(const_cast<Element*>(*neighbor));
218 if (opposite_face_id)
219 {
220 const_cast<Element*>(*neighbor)->setNeighbor(element,
221 *opposite_face_id);
222 }
223 }
224 neighbors.clear();
225 }
226}
227
229{
230 return std::count_if(begin(_nodes), end(_nodes),
231 [this](auto const* const node) {
232 return isBaseNode(
233 *node,
234 _elements_connected_to_nodes[node->getID()]);
235 });
236}
237
239{
240 return std::any_of(
241 std::begin(_elements), std::end(_elements),
242 [](Element const* const e)
243 { return e->getNumberOfNodes() != e->getNumberOfBaseNodes(); });
244}
245
246std::vector<MeshLib::Element const*> const& Mesh::getElementsConnectedToNode(
247 std::size_t const node_id) const
248{
249 return _elements_connected_to_nodes[node_id];
250}
251
252std::vector<MeshLib::Element const*> const& Mesh::getElementsConnectedToNode(
253 Node const& node) const
254{
255 return _elements_connected_to_nodes[node.getID()];
256}
257
259{
260 auto const& properties = mesh.getProperties();
261 if (properties.existsPropertyVector<int>("MaterialIDs",
263 {
264 return properties.getPropertyVector<int>(
265 "MaterialIDs", MeshLib::MeshItemType::Cell, 1);
266 }
267 if (properties.hasPropertyVector("MaterialIDs"))
268 {
269 WARN(
270 "The 'MaterialIDs' mesh property exists but is either of wrong "
271 "type (must be int), or it is not defined on element / cell data.");
272 }
273 return nullptr;
274}
275
277{
278 return const_cast<PropertyVector<int>*>(
279 MeshLib::materialIDs(std::as_const(mesh)));
280}
281
283{
284 auto const& properties = mesh.getProperties();
285 return properties.getPropertyVector<std::size_t>(
288}
289
291{
292 auto const& properties = mesh.getProperties();
293 return properties.getPropertyVector<std::size_t>(
296}
297
298std::vector<std::vector<Node*>> calculateNodesConnectedByElements(
299 Mesh const& mesh)
300{
301 auto const elements_connected_to_nodes = findElementsConnectedToNodes(mesh);
302
303 std::vector<std::vector<Node*>> nodes_connected_by_elements;
304 auto const& nodes = mesh.getNodes();
305 nodes_connected_by_elements.resize(nodes.size());
306 for (std::size_t i = 0; i < nodes.size(); ++i)
307 {
308 auto& adjacent_nodes = nodes_connected_by_elements[i];
309 auto const node_id = nodes[i]->getID();
310
311 // Get all elements, to which this node is connected.
312 auto const& connected_elements = elements_connected_to_nodes[node_id];
313
314 // And collect all elements' nodes.
315 for (Element const* const element : connected_elements)
316 {
317 Node* const* const single_elem_nodes = element->getNodes();
318 std::size_t const nnodes = element->getNumberOfNodes();
319 for (std::size_t n = 0; n < nnodes; n++)
320 {
321 adjacent_nodes.push_back(single_elem_nodes[n]);
322 }
323 }
324
325 // Make nodes unique and sorted by their ids.
326 // This relies on the node's id being equivalent to it's address.
327 std::sort(adjacent_nodes.begin(), adjacent_nodes.end(),
328 idsComparator<Node*>);
329 auto const last =
330 std::unique(adjacent_nodes.begin(), adjacent_nodes.end());
331 adjacent_nodes.erase(last, adjacent_nodes.end());
332 }
333 return nodes_connected_by_elements;
334}
335
336bool isBaseNode(Node const& node,
337 std::vector<Element const*> const& elements_connected_to_node)
338{
339 // Check if node is connected.
340 if (elements_connected_to_node.empty())
341 {
342 return true;
343 }
344
345 // In a mesh a node always belongs to at least one element.
346 auto const e = elements_connected_to_node[0];
347
348 auto const n_base_nodes = e->getNumberOfBaseNodes();
349 auto const local_index = getNodeIDinElement(*e, &node);
350 return local_index < n_base_nodes;
351}
352
353} // namespace MeshLib
Definition of Duplicate functions.
Definition of the Element class.
Definition of the Hex class.
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:40
static std::size_t global_mesh_counter
Mesh counter used to uniquely identify meshes by id.
Definition: Mesh.cpp:39
Definition of the Mesh class.
Definition of the Prism class.
Definition of the Pyramid class.
Definition of the Quad class.
Definition of the RunTime class.
std::string getName(std::string const &line)
Returns the name/title from the "Zone"-description.
Definition of the Tet class.
Definition of the Tri class.
std::size_t getID() const
Definition: Point3dWithID.h:62
void setNeighbor(Element *neighbor, unsigned const face_id)
Definition: Element.cpp:32
virtual unsigned getNumberOfNodes() const =0
virtual unsigned getNumberOfBaseNodes() const =0
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:102
std::vector< std::vector< Element const * > > _elements_connected_to_nodes
Definition: Mesh.h:158
Mesh(std::string name, std::vector< Node * > nodes, std::vector< Element * > elements, Properties const &properties=Properties())
Definition: Mesh.cpp:62
void addElement(Element *elem)
Add an element to the mesh.
Definition: Mesh.cpp:144
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:105
unsigned _mesh_dimension
Definition: Mesh.h:150
std::size_t computeNumberOfBaseNodes() const
Get the number of base nodes.
Definition: Mesh.cpp:228
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:84
void resetNodeIDs()
Resets the IDs of all mesh-nodes to their position in the node vector.
Definition: Mesh.cpp:149
Properties & getProperties()
Definition: Mesh.h:130
std::vector< Element * > _elements
Definition: Mesh.h:155
void setDimension()
Sets the dimension of the mesh.
Definition: Mesh.cpp:167
void resetElementIDs()
Resets the IDs of all mesh-elements to their position in the element vector.
Definition: Mesh.cpp:158
void setElementNeighbors()
Definition: Mesh.cpp:194
virtual ~Mesh()
Destructor.
Definition: Mesh.cpp:129
std::vector< Element const * > const & getElementsConnectedToNode(std::size_t node_id) const
Definition: Mesh.cpp:246
void shallowClean()
Definition: Mesh.cpp:123
bool hasNonlinearElement() const
Check if the mesh contains any nonlinear element.
Definition: Mesh.cpp:238
std::vector< Node * > _nodes
Definition: Mesh.h:154
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties....
Definition: Properties.h:36
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.
Definition: Mesh.h:217
std::vector< std::vector< Element const * > > findElementsConnectedToNodes(Mesh const &mesh)
Definition: Mesh.cpp:45
std::vector< std::vector< Node * > > calculateNodesConnectedByElements(Mesh const &mesh)
Definition: Mesh.cpp:298
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition: Mesh.cpp:258
constexpr std::string_view getBulkIDString(MeshItemType mesh_item_type)
Definition: Properties.h:188
PropertyVector< std::size_t > const * bulkElementIDs(Mesh const &mesh)
Definition: Mesh.cpp:290
std::pair< double, double > computeSqrEdgeLengthRange(Element const &element)
Compute the minimum and maximum squared edge length for this element.
Definition: Element.cpp:156
std::pair< double, double > minMaxEdgeLength(std::vector< Element * > const &elements)
Returns the minimum and maximum edge length for given elements.
Definition: Mesh.cpp:179
unsigned getNodeIDinElement(Element const &element, const MeshLib::Node *node)
Returns the position of the given node in the node array of this element.
Definition: Element.cpp:206
bool isBaseNode(Node const &node, std::vector< Element const * > const &elements_connected_to_node)
Definition: Mesh.cpp:336
PropertyVector< std::size_t > const * bulkNodeIDs(Mesh const &mesh)
Definition: Mesh.cpp:282