OGS 6.2.1-76-gbb689931b
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 
23 #include "Elements/Element.h"
24 #include "Elements/Tri.h"
25 #include "Elements/Quad.h"
26 #include "Elements/Tet.h"
27 #include "Elements/Hex.h"
28 #include "Elements/Pyramid.h"
29 #include "Elements/Prism.h"
30 
31 namespace MeshLib
32 {
33 Mesh::Mesh(std::string name,
34  std::vector<Node*>
35  nodes,
36  std::vector<Element*>
37  elements,
38  Properties const& properties,
39  const std::size_t n_base_nodes)
40  : _id(_counter_value - 1),
41  _mesh_dimension(0),
42  _edge_length(std::numeric_limits<double>::max(), 0),
43  _node_distance(std::numeric_limits<double>::max(), 0),
44  _name(std::move(name)),
45  _nodes(std::move(nodes)),
46  _elements(std::move(elements)),
47  _n_base_nodes(n_base_nodes),
48  _properties(properties)
49 {
50  assert(_n_base_nodes <= _nodes.size());
51  this->resetNodeIDs();
52  this->resetElementIDs();
53  if (_n_base_nodes == 0)
54  {
56  }
57  if ((_n_base_nodes == 0 && hasNonlinearElement()) || isNonlinear())
58  {
59  this->checkNonlinearNodeIDs();
60  }
61  this->setDimension();
63  //this->setNodesConnectedByEdges();
65  this->setElementNeighbors();
66 
67  this->calcEdgeLengthRange();
68 }
69 
70 Mesh::Mesh(const Mesh &mesh)
72  _edge_length(mesh._edge_length.first, mesh._edge_length.second),
73  _node_distance(mesh._node_distance.first, mesh._node_distance.second),
77 {
78  const std::vector<Node*>& nodes (mesh.getNodes());
79  const std::size_t nNodes (nodes.size());
80  for (unsigned i = 0; i < nNodes; ++i)
81  {
82  _nodes[i] = new Node(*nodes[i]);
83  }
84 
85  const std::vector<Element*>& elements (mesh.getElements());
86  const std::size_t nElements (elements.size());
87  for (unsigned i=0; i<nElements; ++i)
88  {
89  const std::size_t nElemNodes = elements[i]->getNumberOfNodes();
90  _elements[i] = elements[i]->clone();
91  for (unsigned j = 0; j < nElemNodes; ++j)
92  {
93  _elements[i]->_nodes[j] = _nodes[elements[i]->getNode(j)->getID()];
94  }
95  }
96 
97  if (_mesh_dimension == 0)
98  {
99  this->setDimension();
100  }
102  //this->setNodesConnectedByEdges();
103  //this->setNodesConnectedByElements();
104  this->setElementNeighbors();
105 }
106 
108 {
109  const std::size_t nElements (_elements.size());
110  for (std::size_t i = 0; i < nElements; ++i)
111  {
112  delete _elements[i];
113  }
114 
115  const std::size_t nNodes (_nodes.size());
116  for (std::size_t i = 0; i < nNodes; ++i)
117  {
118  delete _nodes[i];
119  }
120 }
121 
122 void Mesh::addNode(Node* node)
123 {
124  _nodes.push_back(node);
125 }
126 
128 {
129  _elements.push_back(elem);
130 
131  // add element information to nodes
132  unsigned nNodes (elem->getNumberOfNodes());
133  for (unsigned i = 0; i < nNodes; ++i)
134  {
135  elem->_nodes[i]->addElement(elem);
136  }
137 }
138 
140 {
141  const std::size_t nNodes(_nodes.size());
142  for (std::size_t i = 0; i < nNodes; ++i)
143  {
144  _nodes[i]->setID(i);
145  }
146 }
147 
149 {
150  std::size_t max_basenode_ID = 0;
151  for (Element const* e : _elements)
152  {
153  for (std::size_t i = 0; i < e->getNumberOfBaseNodes(); i++)
154  {
155  max_basenode_ID = std::max(max_basenode_ID, e->getNodeIndex(i));
156  }
157  }
158  _n_base_nodes = max_basenode_ID + 1;
159 }
160 
162 {
163  const std::size_t nElements (this->_elements.size());
164  for (unsigned i = 0; i < nElements; ++i)
165  {
166  _elements[i]->setID(i);
167  }
168 }
169 
171 {
172  const std::size_t nElements (_elements.size());
173  for (unsigned i = 0; i < nElements; ++i)
174  {
176  {
177  _mesh_dimension = _elements[i]->getDimension();
178  }
179  }
180 }
181 
183 {
184  for (auto& element : _elements)
185  {
186  const unsigned nNodes(element->getNumberOfNodes());
187  for (unsigned j = 0; j < nNodes; ++j)
188  {
189  element->_nodes[j]->addElement(element);
190  }
191  }
192 }
193 
195 {
196  for (auto& node : _nodes)
197  {
198  if (node)
199  {
200  node->clearElements();
201  }
202  }
204 }
205 
207 {
208  this->_edge_length.first = std::numeric_limits<double>::max();
209  this->_edge_length.second = 0;
210  double min_length(0);
211  double max_length(0);
212  const std::size_t nElems (this->getNumberOfElements());
213  for (std::size_t i=0; i<nElems; ++i)
214  {
215  _elements[i]->computeSqrEdgeLengthRange(min_length, max_length);
216  this->_edge_length.first = std::min(this->_edge_length.first, min_length);
217  this->_edge_length.second = std::max(this->_edge_length.second, max_length);
218  }
219  this->_edge_length.first = sqrt(this->_edge_length.first);
220  this->_edge_length.second = sqrt(this->_edge_length.second);
221 }
222 
224 {
225  std::vector<Element*> neighbors;
226  for (auto element : _elements)
227  {
228  // create vector with all elements connected to current element (includes lots of doubles!)
229  const std::size_t nNodes (element->getNumberOfBaseNodes());
230  for (unsigned n(0); n<nNodes; ++n)
231  {
232  std::vector<Element*> const& conn_elems ((element->getNode(n)->getElements()));
233  neighbors.insert(neighbors.end(), conn_elems.begin(), conn_elems.end());
234  }
235  std::sort(neighbors.begin(), neighbors.end());
236  auto const neighbors_new_end = std::unique(neighbors.begin(), neighbors.end());
237 
238  for (auto neighbor = neighbors.begin(); neighbor != neighbors_new_end; ++neighbor)
239  {
240  boost::optional<unsigned> const opposite_face_id = element->addNeighbor(*neighbor);
241  if (opposite_face_id)
242  {
243  (*neighbor)->setNeighbor(element, *opposite_face_id);
244  }
245  }
246  neighbors.clear();
247  }
248 }
249 
251 {
252  const std::size_t nNodes (this->_nodes.size());
253  for (unsigned i=0; i<nNodes; ++i)
254  {
255  MeshLib::Node* node (_nodes[i]);
256  std::vector<MeshLib::Node*> conn_set;
257  const std::vector<MeshLib::Element*> &conn_elems (node->getElements());
258  const std::size_t nConnElems (conn_elems.size());
259  for (unsigned j=0; j<nConnElems; ++j)
260  {
261  MeshLib::Element* conn_ele = conn_elems[j];
262  const unsigned idx (conn_ele->getNodeIDinElement(node));
263  const unsigned nElemNodes (conn_ele->getNumberOfBaseNodes());
264  for (unsigned k(0); k<nElemNodes; ++k)
265  {
266  MeshLib::Node const* node_k = conn_ele->getNode(k);
267  bool is_in_vector (false);
268  const std::size_t nConnNodes (conn_set.size());
269  for (unsigned l(0); l < nConnNodes; ++l)
270  {
271  if (node_k == conn_set[l])
272  {
273  is_in_vector = true;
274  }
275  }
276  if (is_in_vector)
277  {
278  continue;
279  }
280 
281  if (conn_ele->getNumberOfBaseNodes() == conn_ele->getNumberOfNodes())
282  {
283  if (conn_ele->isEdge(idx, k))
284  {
285  conn_set.push_back(const_cast<MeshLib::Node*>(node_k));
286  }
287  }
288  else
289  {
290  for (unsigned l=0; l<conn_ele->getNumberOfEdges(); l++)
291  {
292  std::unique_ptr<Element const> edge(conn_ele->getEdge(l));
293  unsigned match = 0;
294  for (unsigned m=0; m<edge->getNumberOfBaseNodes(); m++)
295  {
296  auto edge_node = edge->getNode(m);
297  if (edge_node == node || edge_node == node_k)
298  {
299  match++;
300  }
301  }
302  if (match != 2)
303  {
304  continue;
305  }
306  conn_set.push_back(const_cast<MeshLib::Node*>(node_k));
307  for (unsigned m = edge->getNumberOfBaseNodes();
308  m < edge->getNumberOfNodes();
309  m++)
310  {
311  conn_set.push_back(
312  const_cast<MeshLib::Node*>(edge->getNode(m)));
313  }
314  break;
315  }
316  }
317 
318  }
319  }
320  node->setConnectedNodes(conn_set);
321  }
322 }
323 
325 {
326  // Allocate temporary space for adjacent nodes.
327  std::vector<Node*> adjacent_nodes;
328  for (Node* const node : _nodes)
329  {
330  adjacent_nodes.clear();
331 
332  // Get all elements, to which this node is connected.
333  std::vector<Element*> const& conn_elems = node->getElements();
334 
335  // And collect all elements' nodes.
336  for (Element const* const element : conn_elems)
337  {
338  Node* const* const single_elem_nodes = element->getNodes();
339  std::size_t const nnodes = element->getNumberOfNodes();
340  for (std::size_t n = 0; n < nnodes; n++)
341  {
342  adjacent_nodes.push_back(single_elem_nodes[n]);
343  }
344  }
345 
346  // Make nodes unique and sorted by their ids.
347  // This relies on the node's id being equivalent to it's address.
348  std::sort(adjacent_nodes.begin(), adjacent_nodes.end(),
349  [](Node* a, Node* b) { return a->getID() < b->getID(); });
350  auto const last = std::unique(adjacent_nodes.begin(), adjacent_nodes.end());
351  adjacent_nodes.erase(last, adjacent_nodes.end());
352 
353  node->setConnectedNodes(adjacent_nodes);
354  }
355 }
356 
358 {
359  for (MeshLib::Element const* e : _elements)
360  {
361  for (unsigned i=e->getNumberOfBaseNodes(); i<e->getNumberOfNodes(); i++)
362  {
363  if (e->getNodeIndex(i) >= getNumberOfBaseNodes())
364  {
365  continue;
366  }
367 
368  WARN(
369  "Found a nonlinear node whose ID (%d) is smaller than the "
370  "number of base node IDs (%d). Some functions may not work "
371  "properly.",
372  e->getNodeIndex(i), getNumberOfBaseNodes());
373  return;
374  }
375  }
376 }
377 
379 {
380  return std::any_of(std::begin(_elements), std::end(_elements),
381  [](Element const* const e) {
382  return e->getNumberOfNodes() != e->getNumberOfBaseNodes();
383  });
384 }
385 
387  std::string const& property_name,
388  double factor)
389 {
390  if (!mesh.getProperties().existsPropertyVector<double>(property_name))
391  {
392  WARN("Did not find PropertyVector '%s' for scaling.",
393  property_name.c_str());
394  return;
395  }
396  for (auto& v :
397  *mesh.getProperties().getPropertyVector<double>(property_name))
398  {
399  v *= factor;
400  }
401 }
402 
404 {
405  auto const& properties = mesh.getProperties();
406  return properties.existsPropertyVector<int>("MaterialIDs",
408  ? properties.getPropertyVector<int>(
409  "MaterialIDs", MeshLib::MeshItemType::Cell, 1)
410  : nullptr;
411 }
412 
413 std::unique_ptr<MeshLib::Mesh> createMeshFromElementSelection(
414  std::string mesh_name, std::vector<MeshLib::Element*> const& elements)
415 {
416  DBUG("Found %d elements in the mesh", elements.size());
417 
418  // Store bulk element ids for each of the new elements.
419  std::vector<std::size_t> bulk_element_ids;
420  bulk_element_ids.reserve(elements.size());
421  std::transform(begin(elements), end(elements),
422  std::back_inserter(bulk_element_ids),
423  [&](auto const& e) { return e->getID(); });
424 
425  // original node pointers to newly created nodes.
426  std::unordered_map<const MeshLib::Node*, MeshLib::Node*> nodes_map;
427  nodes_map.reserve(
428  elements.size()); // There will be at least one node per element.
429 
430  for (auto& e : elements)
431  {
432  // For each node find a cloned node in map or create if there is none.
433  unsigned const n_nodes = e->getNumberOfNodes();
434  for (unsigned i = 0; i < n_nodes; ++i)
435  {
436  const MeshLib::Node* n = e->getNode(i);
437  auto const it = nodes_map.find(n);
438  if (it == nodes_map.end())
439  {
440  auto new_node_in_map = nodes_map[n] = new MeshLib::Node(*n);
441  e->setNode(i, new_node_in_map);
442  }
443  else
444  {
445  e->setNode(i, it->second);
446  }
447  }
448  }
449 
450  // Copy the unique nodes pointers.
451  std::vector<MeshLib::Node*> element_nodes;
452  element_nodes.reserve(nodes_map.size());
453  std::transform(begin(nodes_map), end(nodes_map),
454  std::back_inserter(element_nodes),
455  [](auto const& pair) { return pair.second; });
456 
457  // Store bulk node ids for each of the new nodes.
458  std::vector<std::size_t> bulk_node_ids;
459  bulk_node_ids.reserve(nodes_map.size());
460  std::transform(begin(nodes_map), end(nodes_map),
461  std::back_inserter(bulk_node_ids),
462  [](auto const& pair) { return pair.first->getID(); });
463 
464  auto mesh = std::make_unique<MeshLib::Mesh>(
465  std::move(mesh_name), std::move(element_nodes), std::move(elements));
466  assert(mesh != nullptr);
467 
468  addPropertyToMesh(*mesh, "bulk_element_ids", MeshLib::MeshItemType::Cell, 1,
469  bulk_element_ids);
470  addPropertyToMesh(*mesh, "bulk_node_ids", MeshLib::MeshItemType::Node, 1,
471  bulk_node_ids);
472 
473  return mesh;
474 }
475 } // namespace MeshLib
virtual unsigned getNumberOfBaseNodes() const =0
virtual ~Mesh()
Destructor.
Definition: Mesh.cpp:107
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:99
virtual bool isEdge(unsigned i, unsigned j) const =0
Returns true if these two indeces form an edge and false otherwise.
void setNodesConnectedByEdges()
Definition: Mesh.cpp:250
std::pair< double, double > _edge_length
The minimal and maximal edge length over all elements in the mesh.
Definition: Mesh.h:178
std::string _name
Definition: Mesh.h:181
const std::vector< Element * > & getElements() const
Get all elements the node is part of.
Definition: Node.h:67
std::vector< Element * > _elements
Definition: Mesh.h:183
std::size_t getNumberOfBaseNodes() const
Get the number of base nodes.
Definition: Mesh.h:126
void setConnectedNodes(std::vector< Node *> &connected_nodes)
Definition: Node.h:95
virtual unsigned getNodeIDinElement(const MeshLib::Node *node) const
Returns the position of the given node in the node array of this element.
Definition: Element.cpp:145
Definition of the Tet class.
std::vector< Node * > _nodes
Definition: Mesh.h:182
void setNodesConnectedByElements()
Definition: Mesh.cpp:324
void resetElementsConnectedToNodes()
Definition: Mesh.cpp:194
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:105
MeshLib::Properties & getProperties()
Definition: Mesh.h:134
Definition of the Mesh class.
void resetNodeIDs()
Resets the IDs of all mesh-nodes to their position in the node vector.
Definition: Mesh.cpp:139
std::size_t getID() const
Definition: Point3dWithID.h:62
void scaleMeshPropertyVector(MeshLib::Mesh &mesh, std::string const &property_name, double factor)
Definition: Mesh.cpp:386
void checkNonlinearNodeIDs() const
Check if all the nonlinear nodes are stored at the end of the node vector.
Definition: Mesh.cpp:357
Definition of the Tri class.
Node ** _nodes
Definition: Element.h:211
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties...
Definition: Properties.h:37
Definition of the RunTime class.
const std::string getName() const
Get name of the mesh.
Definition: Mesh.h:102
std::pair< double, double > _node_distance
The minimal and maximal distance of nodes within an element over all elements in the mesh...
Definition: Mesh.h:180
std::unique_ptr< MeshLib::Mesh > createMeshFromElementSelection(std::string mesh_name, std::vector< MeshLib::Element *> const &elements)
Definition: Mesh.cpp:413
void recalculateMaxBaseNodeId()
Finds the maximum id among all of the base nodes.
Definition: Mesh.cpp:148
void addElement(Element *elem)
Definition: Node.h:87
Definition of the Pyramid class.
std::size_t _n_base_nodes
Definition: Mesh.h:184
std::size_t const _id
Definition: Mesh.h:175
unsigned _mesh_dimension
Definition: Mesh.h:176
Interface for heuristic search length strategy.
Definition: ProjectData.h:30
Properties _properties
Definition: Mesh.h:185
bool hasNonlinearElement() const
Check if the mesh contains any nonlinear element.
Definition: Mesh.cpp:378
bool isNonlinear() const
Return true if the mesh has any nonlinear nodes.
Definition: Mesh.h:132
virtual unsigned getNumberOfEdges() const =0
Get the number of edges for this element.
void calcEdgeLengthRange()
Set the minimum and maximum length over the edges of the mesh.
Definition: Mesh.cpp:206
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition: Mesh.cpp:403
virtual const Element * getEdge(unsigned i) const =0
Returns the i-th edge of the element.
PropertyVector< T > const * getPropertyVector(std::string const &name) const
Definition: Properties.h:119
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
bool existsPropertyVector(std::string const &name) const
Definition: Properties.h:79
Definition of the Quad class.
void addPropertyToMesh(MeshLib::Mesh &mesh, std::string const &name, MeshLib::MeshItemType item_type, std::size_t number_of_components, std::vector< T > const &values)
Definition: Mesh.h:219
void setElementsConnectedToNodes()
Fills in the neighbor-information for nodes (i.e. which element each node belongs to)...
Definition: Mesh.cpp:182
void setDimension()
Sets the dimension of the mesh.
Definition: Mesh.cpp:170
Mesh(std::string name, std::vector< Node *> nodes, std::vector< Element *> elements, Properties const &properties=Properties(), const std::size_t n_base_nodes=0)
Definition: Mesh.cpp:33
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:81
std::size_t getNumberOfElements() const
Get the number of elements.
Definition: Mesh.h:96
virtual unsigned getNumberOfNodes() const =0
Returns the number of all nodes including both linear and nonlinear nodes.
void setElementNeighbors()
Definition: Mesh.cpp:223
Definition of the Element class.
void addNode(Node *node)
Add a node to the mesh.
Definition: Mesh.cpp:122
static std::size_t _counter_value
Definition: Counter.h:23
Definition of the Prism class.
Definition of the Hex class.
const Node * getNode(unsigned i) const
Definition: Element.cpp:158
void addElement(Element *elem)
Add an element to the mesh.
Definition: Mesh.cpp:127
void resetElementIDs()
Resets the IDs of all mesh-elements to their position in the element vector.
Definition: Mesh.cpp:161