OGS
Mesh.cpp
Go to the documentation of this file.
1 
15 #include "Mesh.h"
16 
17 #include <memory>
18 #include <unordered_map>
19 #include <utility>
20 
21 #include "BaseLib/RunTime.h"
22 #include "Elements/Element.h"
23 #include "Elements/Hex.h"
24 #include "Elements/Prism.h"
25 #include "Elements/Pyramid.h"
26 #include "Elements/Quad.h"
27 #include "Elements/Tet.h"
28 #include "Elements/Tri.h"
29 
31 static std::size_t global_mesh_counter = 0;
32 
33 namespace MeshLib
34 {
35 std::vector<std::vector<Element const*>> findElementsConnectedToNodes(
36  Mesh const& mesh)
37 {
38  std::vector<std::vector<Element const*>> elements_connected_to_nodes;
39  auto const& nodes = mesh.getNodes();
40  elements_connected_to_nodes.resize(nodes.size());
41 
42  for (auto const* element : mesh.getElements())
43  {
44  unsigned const number_nodes(element->getNumberOfNodes());
45  for (unsigned j = 0; j < number_nodes; ++j)
46  {
47  auto const node_id = element->getNode(j)->getID();
48  elements_connected_to_nodes[node_id].push_back(element);
49  }
50  }
51  return elements_connected_to_nodes;
52 }
53 
54 Mesh::Mesh(std::string name,
55  std::vector<Node*>
56  nodes,
57  std::vector<Element*>
58  elements,
59  Properties const& properties)
60  : _id(global_mesh_counter++),
61  _mesh_dimension(0),
62  _edge_length(std::numeric_limits<double>::max(), 0),
63  _node_distance(std::numeric_limits<double>::max(), 0),
64  _name(std::move(name)),
65  _nodes(std::move(nodes)),
66  _elements(std::move(elements)),
67  _properties(properties)
68 {
69  this->resetNodeIDs();
70  this->resetElementIDs();
71  this->setDimension();
72 
74 
75  this->setElementNeighbors();
76  this->calcEdgeLengthRange();
77 }
78 
79 Mesh::Mesh(const Mesh& mesh)
80  : _id(global_mesh_counter++),
81  _mesh_dimension(mesh.getDimension()),
82  _edge_length(mesh._edge_length.first, mesh._edge_length.second),
83  _node_distance(mesh._node_distance.first, mesh._node_distance.second),
84  _name(mesh.getName()),
85  _nodes(mesh.getNumberOfNodes()),
86  _elements(mesh.getNumberOfElements()),
87  _properties(mesh._properties)
88 {
89  const std::vector<Node*>& nodes(mesh.getNodes());
90  const std::size_t nNodes(nodes.size());
91  for (unsigned i = 0; i < nNodes; ++i)
92  {
93  _nodes[i] = new Node(*nodes[i]);
94  }
95 
96  const std::vector<Element*>& elements(mesh.getElements());
97  const std::size_t nElements(elements.size());
98  for (unsigned i = 0; i < nElements; ++i)
99  {
100  const std::size_t nElemNodes = elements[i]->getNumberOfNodes();
101  _elements[i] = elements[i]->clone();
102  for (unsigned j = 0; j < nElemNodes; ++j)
103  {
104  _elements[i]->setNode(j, _nodes[elements[i]->getNode(j)->getID()]);
105  }
106  }
107 
108  if (_mesh_dimension == 0)
109  {
110  this->setDimension();
111  }
113  this->setElementNeighbors();
114 }
115 
117 {
118  const std::size_t nElements(_elements.size());
119  for (std::size_t i = 0; i < nElements; ++i)
120  {
121  delete _elements[i];
122  }
123 
124  const std::size_t nNodes(_nodes.size());
125  for (std::size_t i = 0; i < nNodes; ++i)
126  {
127  delete _nodes[i];
128  }
129 }
130 
132 {
133  _elements.push_back(elem);
134 }
135 
137 {
138  const std::size_t nNodes(_nodes.size());
139  for (std::size_t i = 0; i < nNodes; ++i)
140  {
141  _nodes[i]->setID(i);
142  }
143 }
144 
146 {
147  const std::size_t nElements(this->_elements.size());
148  for (unsigned i = 0; i < nElements; ++i)
149  {
150  _elements[i]->setID(i);
151  }
152 }
153 
155 {
156  const std::size_t nElements(_elements.size());
157  for (unsigned i = 0; i < nElements; ++i)
158  {
160  {
161  _mesh_dimension = _elements[i]->getDimension();
162  }
163  }
164 }
165 
167 {
168  const std::size_t nElems(getNumberOfElements());
169  for (std::size_t i = 0; i < nElems; ++i)
170  {
171  auto const& [min_length, max_length] =
173  _edge_length.first = std::min(_edge_length.first, min_length);
174  _edge_length.second = std::max(_edge_length.second, max_length);
175  }
176  _edge_length.first = std::sqrt(_edge_length.first);
177  _edge_length.second = std::sqrt(_edge_length.second);
178 }
179 
181 {
182  std::vector<Element const*> neighbors;
183  for (auto element : _elements)
184  {
185  // create vector with all elements connected to current element
186  // (includes lots of doubles!)
187  const std::size_t nNodes(element->getNumberOfBaseNodes());
188  for (unsigned n(0); n < nNodes; ++n)
189  {
190  auto const& conn_elems(
191  _elements_connected_to_nodes[element->getNode(n)->getID()]);
192  neighbors.insert(neighbors.end(), conn_elems.begin(),
193  conn_elems.end());
194  }
195  std::sort(neighbors.begin(), neighbors.end());
196  auto const neighbors_new_end =
197  std::unique(neighbors.begin(), neighbors.end());
198 
199  for (auto neighbor = neighbors.begin(); neighbor != neighbors_new_end;
200  ++neighbor)
201  {
202  std::optional<unsigned> const opposite_face_id =
203  element->addNeighbor(const_cast<Element*>(*neighbor));
204  if (opposite_face_id)
205  {
206  const_cast<Element*>(*neighbor)->setNeighbor(element,
207  *opposite_face_id);
208  }
209  }
210  neighbors.clear();
211  }
212 }
213 
214 std::size_t Mesh::getNumberOfBaseNodes() const
215 {
216  return std::count_if(begin(_nodes), end(_nodes),
217  [this](auto const* const node) {
218  return isBaseNode(
219  *node,
220  _elements_connected_to_nodes[node->getID()]);
221  });
222 }
223 
225 {
226  return std::any_of(
227  std::begin(_elements), std::end(_elements),
228  [](Element const* const e)
229  { return e->getNumberOfNodes() != e->getNumberOfBaseNodes(); });
230 }
231 
232 std::vector<MeshLib::Element const*> const& Mesh::getElementsConnectedToNode(
233  std::size_t const node_id) const
234 {
235  return _elements_connected_to_nodes[node_id];
236 }
237 
238 std::vector<MeshLib::Element const*> const& Mesh::getElementsConnectedToNode(
239  Node const& node) const
240 {
241  return _elements_connected_to_nodes[node.getID()];
242 }
243 
245  std::string const& property_name,
246  double factor)
247 {
248  if (!mesh.getProperties().existsPropertyVector<double>(property_name))
249  {
250  WARN("Did not find PropertyVector '{:s}' for scaling.", property_name);
251  return;
252  }
253  auto& pv = *mesh.getProperties().getPropertyVector<double>(property_name);
254  std::transform(pv.begin(), pv.end(), pv.begin(),
255  [factor](auto const& v) { return v * factor; });
256 }
257 
259 {
260  auto const& properties = mesh.getProperties();
261  return properties.existsPropertyVector<int>("MaterialIDs",
263  ? properties.getPropertyVector<int>(
264  "MaterialIDs", MeshLib::MeshItemType::Cell, 1)
265  : nullptr;
266 }
267 
268 std::unique_ptr<MeshLib::Mesh> createMeshFromElementSelection(
269  std::string mesh_name, std::vector<MeshLib::Element*> const& elements)
270 {
271  DBUG("Found {:d} elements in the mesh", elements.size());
272 
273  // Store bulk element ids for each of the new elements.
274  std::vector<std::size_t> bulk_element_ids;
275  bulk_element_ids.reserve(elements.size());
276  std::transform(begin(elements), end(elements),
277  std::back_inserter(bulk_element_ids),
278  [&](auto const& e) { return e->getID(); });
279 
280  // original node pointers to newly created nodes.
281  std::unordered_map<const MeshLib::Node*, MeshLib::Node*> nodes_map;
282  nodes_map.reserve(
283  elements.size()); // There will be at least one node per element.
284 
285  for (auto& e : elements)
286  {
287  // For each node find a cloned node in map or create if there is none.
288  unsigned const n_nodes = e->getNumberOfNodes();
289  for (unsigned i = 0; i < n_nodes; ++i)
290  {
291  const MeshLib::Node* n = e->getNode(i);
292  auto const it = nodes_map.find(n);
293  if (it == nodes_map.end())
294  {
295  auto new_node_in_map = nodes_map[n] = new MeshLib::Node(*n);
296  e->setNode(i, new_node_in_map);
297  }
298  else
299  {
300  e->setNode(i, it->second);
301  }
302  }
303  }
304 
305  // Copy the unique nodes pointers.
306  std::vector<MeshLib::Node*> element_nodes;
307  element_nodes.reserve(nodes_map.size());
308  std::transform(begin(nodes_map), end(nodes_map),
309  std::back_inserter(element_nodes),
310  [](auto const& pair) { return pair.second; });
311 
312  // Store bulk node ids for each of the new nodes.
313  std::vector<std::size_t> bulk_node_ids;
314  bulk_node_ids.reserve(nodes_map.size());
315  std::transform(begin(nodes_map), end(nodes_map),
316  std::back_inserter(bulk_node_ids),
317  [](auto const& pair) { return pair.first->getID(); });
318 
319  auto mesh = std::make_unique<MeshLib::Mesh>(
320  std::move(mesh_name), std::move(element_nodes), std::move(elements));
321  assert(mesh != nullptr);
322 
323  addPropertyToMesh(*mesh, "bulk_element_ids", MeshLib::MeshItemType::Cell, 1,
324  bulk_element_ids);
325  addPropertyToMesh(*mesh, "bulk_node_ids", MeshLib::MeshItemType::Node, 1,
326  bulk_node_ids);
327 
328  return mesh;
329 }
330 
331 std::vector<std::vector<Node*>> calculateNodesConnectedByElements(
332  Mesh const& mesh)
333 {
334  auto const elements_connected_to_nodes = findElementsConnectedToNodes(mesh);
335 
336  std::vector<std::vector<Node*>> nodes_connected_by_elements;
337  auto const& nodes = mesh.getNodes();
338  nodes_connected_by_elements.resize(nodes.size());
339  for (std::size_t i = 0; i < nodes.size(); ++i)
340  {
341  auto& adjacent_nodes = nodes_connected_by_elements[i];
342  auto const* node = nodes[i];
343 
344  // Get all elements, to which this node is connected.
345  auto const& connected_elements =
346  elements_connected_to_nodes[node->getID()];
347 
348  // And collect all elements' nodes.
349  for (Element const* const element : connected_elements)
350  {
351  Node* const* const single_elem_nodes = element->getNodes();
352  std::size_t const nnodes = element->getNumberOfNodes();
353  for (std::size_t n = 0; n < nnodes; n++)
354  {
355  adjacent_nodes.push_back(single_elem_nodes[n]);
356  }
357  }
358 
359  // Make nodes unique and sorted by their ids.
360  // This relies on the node's id being equivalent to it's address.
361  std::sort(adjacent_nodes.begin(), adjacent_nodes.end(),
362  [](Node* a, Node* b) { return a->getID() < b->getID(); });
363  auto const last =
364  std::unique(adjacent_nodes.begin(), adjacent_nodes.end());
365  adjacent_nodes.erase(last, adjacent_nodes.end());
366  }
367  return nodes_connected_by_elements;
368 }
369 
370 bool isBaseNode(Node const& node,
371  std::vector<Element const*> const& elements_connected_to_node)
372 {
373  // Check if node is connected.
374  if (elements_connected_to_node.empty())
375  {
376  return true;
377  }
378 
379  // In a mesh a node always belongs to at least one element.
380  auto const e = elements_connected_to_node[0];
381 
382  auto const n_base_nodes = e->getNumberOfBaseNodes();
383  auto const local_index = getNodeIDinElement(*e, &node);
384  return local_index < n_base_nodes;
385 }
386 } // namespace MeshLib
Definition of the Element class.
Definition of the Hex class.
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:27
void WARN(char const *fmt, Args const &... args)
Definition: Logging.h:37
static std::size_t global_mesh_counter
Mesh counter used to uniquely identify meshes by id.
Definition: Mesh.cpp:31
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::size_t getNumberOfBaseNodes() const
Get the number of base nodes.
Definition: Mesh.cpp:214
void calcEdgeLengthRange()
Set the minimum and maximum length over the edges of the mesh.
Definition: Mesh.cpp:166
std::vector< std::vector< Element const * > > _elements_connected_to_nodes
Definition: Mesh.h:155
Mesh(std::string name, std::vector< Node * > nodes, std::vector< Element * > elements, Properties const &properties=Properties())
Definition: Mesh.cpp:54
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:95
void addElement(Element *elem)
Add an element to the mesh.
Definition: Mesh.cpp:131
std::pair< double, double > _edge_length
The minimal and maximal edge length over all elements in the mesh.
Definition: Mesh.h:145
unsigned _mesh_dimension
Definition: Mesh.h:143
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:71
std::size_t getID() const
Get id of the mesh.
Definition: Mesh.h:110
void resetNodeIDs()
Resets the IDs of all mesh-nodes to their position in the node vector.
Definition: Mesh.cpp:136
std::vector< Element * > _elements
Definition: Mesh.h:152
void setDimension()
Sets the dimension of the mesh.
Definition: Mesh.cpp:154
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:98
void resetElementIDs()
Resets the IDs of all mesh-elements to their position in the element vector.
Definition: Mesh.cpp:145
void setElementNeighbors()
Definition: Mesh.cpp:180
virtual ~Mesh()
Destructor.
Definition: Mesh.cpp:116
Properties & getProperties()
Definition: Mesh.h:123
std::vector< Element const * > const & getElementsConnectedToNode(std::size_t node_id) const
Definition: Mesh.cpp:232
bool hasNonlinearElement() const
Check if the mesh contains any nonlinear element.
Definition: Mesh.cpp:224
const Node * getNode(std::size_t idx) const
Get the node with the given index.
Definition: Mesh.h:74
std::size_t getNumberOfElements() const
Get the number of elements.
Definition: Mesh.h:86
std::vector< Node * > _nodes
Definition: Mesh.h:151
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties....
Definition: Properties.h:36
PropertyVector< T > const * getPropertyVector(std::string const &name) const
bool existsPropertyVector(std::string const &name) const
constexpr bool any_of(List const &values)
Checks if any of the elements in the given list is true.
Definition: Algorithm.h:325
static const double v
std::vector< std::vector< Node * > > calculateNodesConnectedByElements(Mesh const &mesh)
Definition: Mesh.cpp:331
std::vector< std::vector< Element const * > > findElementsConnectedToNodes(Mesh const &mesh)
Definition: Mesh.cpp:35
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition: Mesh.cpp:258
std::pair< double, double > computeSqrEdgeLengthRange(Element const &element)
Compute the minimum and maximum squared edge length for this element.
Definition: Element.cpp:162
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:212
std::unique_ptr< MeshLib::Mesh > createMeshFromElementSelection(std::string mesh_name, std::vector< MeshLib::Element * > const &elements)
Definition: Mesh.cpp:268
void addPropertyToMesh(Mesh &mesh, std::string const &name, MeshItemType item_type, std::size_t number_of_components, std::vector< T > const &values)
Definition: Mesh.h:193
bool isBaseNode(Node const &node, std::vector< Element const * > const &elements_connected_to_node)
Definition: Mesh.cpp:370
void scaleMeshPropertyVector(MeshLib::Mesh &mesh, std::string const &property_name, double factor)
Definition: Mesh.cpp:244
unsigned getDimension(MeshLib::MeshElemType eleType)