libMesh
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
libMesh::TetGenMeshInterface Class Reference

Class TetGenMeshInterface provides an interface for tetrahedralization of meshes using the TetGen library. More...

#include <mesh_tetgen_interface.h>

Public Member Functions

 TetGenMeshInterface (UnstructuredMesh &mesh)
 Constructor. More...
 
 ~TetGenMeshInterface ()=default
 Empty destructor. More...
 
void set_switches (std::string new_switches)
 Method to set switches to tetgen, allowing for different behaviours. More...
 
void triangulate_pointset ()
 Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set. More...
 
void pointset_convexhull ()
 Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set. More...
 
void triangulate_conformingDelaunayMesh (double quality_constraint=0., double volume_constraint=0.)
 Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set. More...
 
void triangulate_conformingDelaunayMesh_carvehole (const std::vector< Point > &holes, double quality_constraint=0., double volume_constraint=0.)
 Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set. More...
 

Protected Member Functions

void fill_pointlist (TetGenWrapper &wrapper)
 This function copies nodes from the _mesh into TetGen's pointlist. More...
 
void assign_nodes_to_elem (unsigned *node_labels, Elem *elem)
 Assigns the node IDs contained in the 'node_labels' array to 'elem'. More...
 
unsigned check_hull_integrity ()
 This function checks the integrity of the current set of elements in the Mesh to see if they comprise a convex hull, that is: More...
 
void process_hull_integrity_result (unsigned result)
 This function prints an informative message and crashes based on the output of the check_hull_integrity() function. More...
 
void delete_2D_hull_elements ()
 Delete original convex hull elements from the Mesh after performing a Delaunay tetrahedralization. More...
 

Protected Attributes

UnstructuredMesh_mesh
 Local reference to the mesh we are working with. More...
 
std::vector< unsigned > _sequential_to_libmesh_node_map
 We should not assume libmesh nodes are numbered sequentially... More...
 
MeshSerializer _serializer
 Tetgen only operates on serial meshes. More...
 
std::string _switches
 Parameter controlling the behaviour of tetgen. More...
 

Detailed Description

Class TetGenMeshInterface provides an interface for tetrahedralization of meshes using the TetGen library.

For information about TetGen cf. TetGen home page.

Author
Steffen Petersen
Date
2004
Author
John W. Peterson
Date
2011

Definition at line 54 of file mesh_tetgen_interface.h.

Constructor & Destructor Documentation

◆ TetGenMeshInterface()

libMesh::TetGenMeshInterface::TetGenMeshInterface ( UnstructuredMesh mesh)
explicit

Constructor.

Takes a reference to the mesh.

Definition at line 40 of file mesh_tetgen_interface.C.

40  :
41  _mesh(mesh),
43  _switches("Q")
44 {
45 }
MeshSerializer _serializer
Tetgen only operates on serial meshes.
std::string _switches
Parameter controlling the behaviour of tetgen.
MeshBase & mesh
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.

◆ ~TetGenMeshInterface()

libMesh::TetGenMeshInterface::~TetGenMeshInterface ( )
default

Empty destructor.

Member Function Documentation

◆ assign_nodes_to_elem()

void libMesh::TetGenMeshInterface::assign_nodes_to_elem ( unsigned *  node_labels,
Elem elem 
)
protected

Assigns the node IDs contained in the 'node_labels' array to 'elem'.

Definition at line 381 of file mesh_tetgen_interface.C.

References _mesh, _sequential_to_libmesh_node_map, libMesh::Elem::node_index_range(), libMesh::MeshBase::node_ptr(), and libMesh::Elem::set_node().

Referenced by pointset_convexhull(), triangulate_conformingDelaunayMesh_carvehole(), and triangulate_pointset().

382 {
383  for (auto j : elem->node_index_range())
384  {
385  // Get the mapped node index to ask the Mesh for
386  unsigned mapped_node_id = _sequential_to_libmesh_node_map[ node_labels[j] ];
387 
388  Node * current_node = this->_mesh.node_ptr( mapped_node_id );
389 
390  elem->set_node(j) = current_node;
391  }
392 }
std::vector< unsigned > _sequential_to_libmesh_node_map
We should not assume libmesh nodes are numbered sequentially...
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual const Node * node_ptr(const dof_id_type i) const =0

◆ check_hull_integrity()

unsigned libMesh::TetGenMeshInterface::check_hull_integrity ( )
protected

This function checks the integrity of the current set of elements in the Mesh to see if they comprise a convex hull, that is:

  • If they are all TRI3 elements
  • They all have non-nullptr neighbors
Returns
  • 0 if the mesh forms a valid convex hull
  • 1 if a non-TRI3 element is found
  • 2 if an element with a nullptr-neighbor is found

Definition at line 398 of file mesh_tetgen_interface.C.

References _mesh, libMesh::MeshBase::n_elem(), and libMesh::TRI3.

Referenced by triangulate_conformingDelaunayMesh_carvehole().

399 {
400  // Check for easy return: if the Mesh is empty (i.e. if
401  // somebody called triangulate_conformingDelaunayMesh on
402  // a Mesh with no elements, then hull integrity check must
403  // fail...
404  if (_mesh.n_elem() == 0)
405  return 3;
406 
407  for (auto & elem : this->_mesh.element_ptr_range())
408  {
409  // Check for proper element type
410  if (elem->type() != TRI3)
411  {
412  //libmesh_error_msg("ERROR: Some of the elements in the original mesh were not TRI3!");
413  return 1;
414  }
415 
416  for (auto neigh : elem->neighbor_ptr_range())
417  {
418  if (neigh == nullptr)
419  {
420  // libmesh_error_msg("ERROR: Non-convex hull, cannot be tetrahedralized.");
421  return 2;
422  }
423  }
424  }
425 
426  // If we made it here, return success!
427  return 0;
428 }
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual dof_id_type n_elem() const =0

◆ delete_2D_hull_elements()

void libMesh::TetGenMeshInterface::delete_2D_hull_elements ( )
protected

Delete original convex hull elements from the Mesh after performing a Delaunay tetrahedralization.

Definition at line 453 of file mesh_tetgen_interface.C.

References _mesh, libMesh::MeshBase::delete_elem(), libMesh::MeshBase::get_boundary_info(), libMesh::BoundaryInfo::regenerate_id_sets(), and libMesh::TRI3.

Referenced by triangulate_conformingDelaunayMesh_carvehole().

454 {
455  for (auto & elem : this->_mesh.element_ptr_range())
456  {
457  // Check for proper element type. Yes, we legally delete elements while
458  // iterating over them because no entries from the underlying container
459  // are actually erased.
460  if (elem->type() == TRI3)
461  _mesh.delete_elem(elem);
462  }
463 
464  // We just removed any boundary info associated with hull element
465  // edges, so let's update the boundary id caches.
467 }
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:159
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual void delete_elem(Elem *e)=0
Removes element e from the mesh.
void regenerate_id_sets()
Clears and regenerates the cached sets of ids.

◆ fill_pointlist()

void libMesh::TetGenMeshInterface::fill_pointlist ( TetGenWrapper wrapper)
protected

This function copies nodes from the _mesh into TetGen's pointlist.

Takes some pains to ensure that non-sequential node numberings (which can happen with e.g. DistributedMesh) are handled.

Definition at line 353 of file mesh_tetgen_interface.C.

References _mesh, _sequential_to_libmesh_node_map, libMesh::TetGenWrapper::allocate_pointlist(), libMesh::MeshBase::n_nodes(), and libMesh::TetGenWrapper::set_node().

Referenced by pointset_convexhull(), triangulate_conformingDelaunayMesh_carvehole(), and triangulate_pointset().

354 {
355  // fill input structure with point set data:
356  wrapper.allocate_pointlist( this->_mesh.n_nodes() );
357 
358  // Make enough space to store a mapping between the implied sequential
359  // node numbering used in tetgen and libmesh's (possibly) non-sequential
360  // numbering scheme.
363 
364  {
365  unsigned index = 0;
366  for (auto & node : this->_mesh.node_ptr_range())
367  {
368  _sequential_to_libmesh_node_map[index] = node->id();
369  wrapper.set_node(index++,
370  REAL((*node)(0)),
371  REAL((*node)(1)),
372  REAL((*node)(2)));
373  }
374  }
375 }
std::vector< unsigned > _sequential_to_libmesh_node_map
We should not assume libmesh nodes are numbered sequentially...
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual dof_id_type n_nodes() const =0

◆ pointset_convexhull()

void libMesh::TetGenMeshInterface::pointset_convexhull ( )

Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set.

Stores only 2D hull surface elements.

Definition at line 115 of file mesh_tetgen_interface.C.

References _mesh, _switches, libMesh::MeshBase::add_elem(), assign_nodes_to_elem(), libMesh::Elem::build(), libMesh::MeshBase::delete_elem(), fill_pointlist(), libMesh::MeshBase::get_boundary_info(), libMesh::TetGenWrapper::get_numberoftrifaces(), libMesh::TetGenWrapper::get_triface_node(), libMesh::BoundaryInfo::regenerate_id_sets(), libMesh::TetGenWrapper::run_tetgen(), libMesh::TetGenWrapper::set_switches(), and libMesh::TRI3.

Referenced by add_cube_convex_hull_to_mesh().

116 {
117  // class tetgen_wrapper allows library access on a basic level
118  TetGenWrapper tetgen_wrapper;
119 
120  // Copy Mesh's node points into TetGen data structure
121  this->fill_pointlist(tetgen_wrapper);
122 
123  // Run TetGen triangulation method:
124  // Q = quiet, no terminal output
125  // Note: if no switch is used, the input must be a list of 3D points
126  // (.node file) and the Delaunay tetrahedralization of this point set
127  // will be generated. In this particular function, we are throwing
128  // away the tetrahedra generated by TetGen, and keeping only the
129  // convex hull...
130  tetgen_wrapper.set_switches(_switches);
131  tetgen_wrapper.run_tetgen();
132  unsigned int num_elements = tetgen_wrapper.get_numberoftrifaces();
133 
134  // Delete *all* old elements. Yes, we legally delete elements while
135  // iterating over them because no entries from the underlying container
136  // are actually erased.
137  for (auto & elem : this->_mesh.element_ptr_range())
138  this->_mesh.delete_elem (elem);
139 
140  // We just removed any boundary info associated with element faces
141  // or edges, so let's update the boundary id caches.
143 
144  // Add the 2D elements which comprise the convex hull back to the mesh.
145  // Vector that temporarily holds the node labels defining element.
146  unsigned int node_labels[3];
147 
148  for (unsigned int i=0; i<num_elements; ++i)
149  {
150  auto elem = Elem::build(TRI3);
151 
152  // Get node labels associated with this element
153  for (auto j : elem->node_index_range())
154  node_labels[j] = tetgen_wrapper.get_triface_node(i,j);
155 
156  this->assign_nodes_to_elem(node_labels, elem.get());
157 
158  // Finally, add this element to the mesh.
159  this->_mesh.add_elem(std::move(elem));
160  }
161 }
std::string _switches
Parameter controlling the behaviour of tetgen.
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:159
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual void delete_elem(Elem *e)=0
Removes element e from the mesh.
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static std::unique_ptr< Elem > build(const ElemType type, Elem *p=nullptr)
Definition: elem.C:273
void assign_nodes_to_elem(unsigned *node_labels, Elem *elem)
Assigns the node IDs contained in the &#39;node_labels&#39; array to &#39;elem&#39;.
void regenerate_id_sets()
Clears and regenerates the cached sets of ids.
void fill_pointlist(TetGenWrapper &wrapper)
This function copies nodes from the _mesh into TetGen&#39;s pointlist.

◆ process_hull_integrity_result()

void libMesh::TetGenMeshInterface::process_hull_integrity_result ( unsigned  result)
protected

This function prints an informative message and crashes based on the output of the check_hull_integrity() function.

It is a separate function so that you can check hull integrity without crashing if you desire.

Definition at line 434 of file mesh_tetgen_interface.C.

References libMesh::err.

Referenced by triangulate_conformingDelaunayMesh_carvehole().

435 {
436  if (result != 0)
437  {
438  libMesh::err << "Error! Conforming Delaunay mesh tetrahedralization requires a convex hull." << std::endl;
439 
440  if (result==1)
441  {
442  libMesh::err << "Non-TRI3 elements were found in the input Mesh. ";
443  libMesh::err << "A constrained Delaunay triangulation requires a convex hull of TRI3 elements." << std::endl;
444  }
445 
446  libmesh_error_msg("Consider calling TetGenMeshInterface::pointset_convexhull() followed by Mesh::find_neighbors() first.");
447  }
448 }
OStreamProxy err

◆ set_switches()

void libMesh::TetGenMeshInterface::set_switches ( std::string  new_switches)

Method to set switches to tetgen, allowing for different behaviours.

Definition at line 47 of file mesh_tetgen_interface.C.

References _switches.

Referenced by tetrahedralize_domain().

48 {
49  // set the tetgen switch manually:
50  // p = tetrahedralizes a piecewise linear complex (see definition in user manual)
51  // Q = quiet, no terminal output
52  // q = specify a minimum radius/edge ratio
53  // a = tetrahedron volume constraint
54  // V = verbose output
55  // for full list of options and their meaning: see the tetgen manual
56  // (http://wias-berlin.de/software/tetgen/1.5/doc/manual/manual005.html)
57  _switches = std::move(switches);
58 }
std::string _switches
Parameter controlling the behaviour of tetgen.

◆ triangulate_conformingDelaunayMesh()

void libMesh::TetGenMeshInterface::triangulate_conformingDelaunayMesh ( double  quality_constraint = 0.,
double  volume_constraint = 0. 
)

Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set.

Boundary constraints are taken from elements array.

Definition at line 167 of file mesh_tetgen_interface.C.

References triangulate_conformingDelaunayMesh_carvehole().

169 {
170  // start triangulation method with empty holes list:
171  std::vector<Point> noholes;
172  triangulate_conformingDelaunayMesh_carvehole(noholes, quality_constraint, volume_constraint);
173 }
void triangulate_conformingDelaunayMesh_carvehole(const std::vector< Point > &holes, double quality_constraint=0., double volume_constraint=0.)
Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set...

◆ triangulate_conformingDelaunayMesh_carvehole()

void libMesh::TetGenMeshInterface::triangulate_conformingDelaunayMesh_carvehole ( const std::vector< Point > &  holes,
double  quality_constraint = 0.,
double  volume_constraint = 0. 
)

Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set.

Boundary constraints are taken from elements array. Include carve-out functionality.

Definition at line 177 of file mesh_tetgen_interface.C.

References _mesh, _sequential_to_libmesh_node_map, _switches, libMesh::MeshBase::add_elem(), libMesh::MeshBase::add_point(), libMesh::TetGenWrapper::allocate_facet_polygonlist(), libMesh::TetGenWrapper::allocate_facetlist(), libMesh::TetGenWrapper::allocate_polygon_vertexlist(), assign_nodes_to_elem(), libMesh::Utility::binary_find(), libMesh::Elem::build(), check_hull_integrity(), delete_2D_hull_elements(), distance(), fill_pointlist(), libMesh::TetGenWrapper::get_element_node(), libMesh::TetGenWrapper::get_numberofpoints(), libMesh::TetGenWrapper::get_numberoftetrahedra(), libMesh::TetGenWrapper::get_output_node(), libMesh::DofObject::id(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), process_hull_integrity_result(), libMesh::TetGenWrapper::run_tetgen(), libMesh::TetGenWrapper::set_hole(), libMesh::TetGenWrapper::set_switches(), libMesh::TetGenWrapper::set_vertex(), and libMesh::TET4.

Referenced by tetrahedralize_domain(), and triangulate_conformingDelaunayMesh().

180 {
181  // Before calling this function, the Mesh must contain a convex hull
182  // of TRI3 elements which define the boundary.
183  unsigned hull_integrity_check = check_hull_integrity();
184 
185  // Possibly die if hull integrity check failed
186  this->process_hull_integrity_result(hull_integrity_check);
187 
188  // class tetgen_wrapper allows library access on a basic level
189  TetGenWrapper tetgen_wrapper;
190 
191  // Copy Mesh's node points into TetGen data structure
192  this->fill_pointlist(tetgen_wrapper);
193 
194  // >>> fill input structure "tetgenio" with facet data:
195  int facet_num = this->_mesh.n_elem();
196 
197  // allocate memory in "tetgenio" structure:
198  tetgen_wrapper.allocate_facetlist
199  (facet_num, cast_int<int>(holes.size()));
200 
201 
202  // Set up tetgen data structures with existing facet information
203  // from the convex hull.
204  {
205  int insertnum = 0;
206  for (auto & elem : this->_mesh.element_ptr_range())
207  {
208  tetgen_wrapper.allocate_facet_polygonlist(insertnum, 1);
209  tetgen_wrapper.allocate_polygon_vertexlist(insertnum, 0, 3);
210 
211  for (auto j : elem->node_index_range())
212  {
213  // We need to get the sequential index of elem->node_ptr(j), but
214  // it should already be stored in _sequential_to_libmesh_node_map...
215  unsigned libmesh_node_id = elem->node_id(j);
216 
217  // The libmesh node IDs may not be sequential, but can we assume
218  // they are at least in order??? We will do so here.
219  std::vector<unsigned>::iterator node_iter =
222  libmesh_node_id);
223 
224  // Check to see if not found: this could also indicate the sequential
225  // node map is not sorted...
226  libmesh_error_msg_if(node_iter == _sequential_to_libmesh_node_map.end(),
227  "Global node " << libmesh_node_id << " not found in sequential node map!");
228 
229  int sequential_index = cast_int<int>
231  node_iter));
232 
233  // Debugging:
234  // libMesh::out << "libmesh_node_id=" << libmesh_node_id
235  // << ", sequential_index=" << sequential_index
236  // << std::endl;
237 
238  tetgen_wrapper.set_vertex(insertnum, // facet number
239  0, // polygon (always 0)
240  j, // local vertex index in tetgen input
241  sequential_index);
242  }
243 
244  // Go to next facet in polygonlist
245  insertnum++;
246  }
247  }
248 
249 
250 
251  // fill hole list (if there are holes):
252  if (holes.size() > 0)
253  {
254  unsigned hole_index = 0;
255  for (Point ihole : holes)
256  tetgen_wrapper.set_hole(hole_index++,
257  REAL(ihole(0)),
258  REAL(ihole(1)),
259  REAL(ihole(2)));
260  }
261 
262 
263  // Run TetGen triangulation method
264 
265  // Assemble switches: we append the user's switches (if any) to
266  // - 'p' tetrahedralize a piecewise linear complex
267  // - 'C' check consistency of mesh (avoid inverted elements)
268  // (see definition and further options in user manual
269  // http://wias-berlin.de/software/tetgen/1.5/doc/manual/manual005.html )
270  std::ostringstream oss;
271  oss << "pC";
272  oss << _switches;
273 
274  if (quality_constraint != 0)
275  oss << "q" << std::fixed << quality_constraint;
276 
277  if (volume_constraint != 0)
278  oss << "a" << std::fixed << volume_constraint;
279 
280  std::string params = oss.str();
281 
282  tetgen_wrapper.set_switches(params); // TetGen switches: Piecewise linear complex, Quiet mode
283  tetgen_wrapper.run_tetgen();
284 
285  // => nodes:
286  unsigned int old_nodesnum = this->_mesh.n_nodes();
287  REAL x=0., y=0., z=0.;
288  const unsigned int num_nodes = tetgen_wrapper.get_numberofpoints();
289 
290  // Debugging:
291  // libMesh::out << "Original mesh had " << old_nodesnum << " nodes." << std::endl;
292  // libMesh::out << "Reserving space for " << num_nodes << " total nodes." << std::endl;
293 
294  // Reserve space for additional nodes in the node map
295  _sequential_to_libmesh_node_map.reserve(num_nodes);
296 
297  // Add additional nodes to the Mesh.
298  // Original code had i<=num_nodes here (Note: the indexing is:
299  // foo[3*i], [3*i+1], [3*i+2]) But according to the TetGen docs, "In
300  // all cases, the first item in any array is stored starting at
301  // index [0]."
302  for (unsigned int i=old_nodesnum; i<num_nodes; i++)
303  {
304  // Fill in x, y, z values
305  tetgen_wrapper.get_output_node(i, x,y,z);
306 
307  // Catch the node returned by add_point()... this will tell us the ID
308  // assigned by the Mesh.
309  Node * new_node = this->_mesh.add_point ( Point(x,y,z) );
310 
311  // Store this new ID in our sequential-to-libmesh node mapping array
312  _sequential_to_libmesh_node_map.push_back( new_node->id() );
313  }
314 
315  // Debugging:
316  // std::copy(_sequential_to_libmesh_node_map.begin(),
317  // _sequential_to_libmesh_node_map.end(),
318  // std::ostream_iterator<unsigned>(std::cout, " "));
319  // std::cout << std::endl;
320 
321 
322  // => tetrahedra:
323  const unsigned int num_elements = tetgen_wrapper.get_numberoftetrahedra();
324 
325  // Vector that temporarily holds the node labels defining element connectivity.
326  unsigned int node_labels[4];
327 
328  for (unsigned int i=0; i<num_elements; i++)
329  {
330  // TetGen only supports Tet4 elements.
331  auto elem = Elem::build(TET4);
332 
333  // Fill up the the node_labels vector
334  for (auto j : elem->node_index_range())
335  node_labels[j] = tetgen_wrapper.get_element_node(i,j);
336 
337  // Associate nodes with this element
338  this->assign_nodes_to_elem(node_labels, elem.get());
339 
340  // Finally, add this element to the mesh
341  this->_mesh.add_elem(std::move(elem));
342  }
343 
344  // Delete original convex hull elements. Is there ever a case where
345  // we should not do this?
346  this->delete_2D_hull_elements();
347 }
std::string _switches
Parameter controlling the behaviour of tetgen.
std::vector< unsigned > _sequential_to_libmesh_node_map
We should not assume libmesh nodes are numbered sequentially...
void process_hull_integrity_result(unsigned result)
This function prints an informative message and crashes based on the output of the check_hull_integri...
Real distance(const Point &p)
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id)=0
Add a new Node at Point p to the end of the vertex array, with processor_id procid.
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static std::unique_ptr< Elem > build(const ElemType type, Elem *p=nullptr)
Definition: elem.C:273
void assign_nodes_to_elem(unsigned *node_labels, Elem *elem)
Assigns the node IDs contained in the &#39;node_labels&#39; array to &#39;elem&#39;.
ForwardIterator binary_find(ForwardIterator first, ForwardIterator last, const T &value)
The STL provides std::binary_search() which returns true or false depending on whether the searched-f...
Definition: utility.h:265
unsigned check_hull_integrity()
This function checks the integrity of the current set of elements in the Mesh to see if they comprise...
void fill_pointlist(TetGenWrapper &wrapper)
This function copies nodes from the _mesh into TetGen&#39;s pointlist.
virtual dof_id_type n_elem() const =0
void delete_2D_hull_elements()
Delete original convex hull elements from the Mesh after performing a Delaunay tetrahedralization.
virtual dof_id_type n_nodes() const =0

◆ triangulate_pointset()

void libMesh::TetGenMeshInterface::triangulate_pointset ( )

Method invokes TetGen library to compute a Delaunay tetrahedralization from the nodes point set.

Definition at line 61 of file mesh_tetgen_interface.C.

References _mesh, _switches, libMesh::MeshBase::add_elem(), assign_nodes_to_elem(), libMesh::Elem::build(), fill_pointlist(), libMesh::TetGenWrapper::get_element_node(), libMesh::TetGenWrapper::get_numberoftetrahedra(), libMesh::TetGenWrapper::run_tetgen(), libMesh::TetGenWrapper::set_switches(), and libMesh::TET4.

62 {
63  // class tetgen_wrapper allows library access on a basic level:
64  TetGenWrapper tetgen_wrapper;
65 
66  // fill input structure with point set data:
67  this->fill_pointlist(tetgen_wrapper);
68 
69  // Run TetGen triangulation method:
70  // Q = quiet, no terminal output
71  // V = verbose, more terminal output
72  // Note: if no switch is used, the input must be a list of 3D points
73  // (.node file) and the Delaunay tetrahedralization of this point set
74  // will be generated.
75 
76  // Can we apply quality and volume constraints in
77  // triangulate_pointset()?. On at least one test problem,
78  // specifying any quality or volume constraints here causes tetgen
79  // to segfault down in the insphere method: a nullptr is passed
80  // to the routine.
81  std::ostringstream oss;
82  oss << _switches;
83  // oss << "V"; // verbose operation
84  //oss << "q" << std::fixed << 2.0; // quality constraint
85  //oss << "a" << std::fixed << 100.; // volume constraint
86  tetgen_wrapper.set_switches(oss.str());
87 
88  // Run tetgen
89  tetgen_wrapper.run_tetgen();
90 
91  // save elements to mesh structure, nodes will not be changed:
92  const unsigned int num_elements = tetgen_wrapper.get_numberoftetrahedra();
93 
94  // Vector that temporarily holds the node labels defining element.
95  unsigned int node_labels[4];
96 
97  for (unsigned int i=0; i<num_elements; ++i)
98  {
99  auto elem = Elem::build(TET4);
100 
101  // Get the nodes associated with this element
102  for (auto j : elem->node_index_range())
103  node_labels[j] = tetgen_wrapper.get_element_node(i,j);
104 
105  // Associate the nodes with this element
106  this->assign_nodes_to_elem(node_labels, elem.get());
107 
108  // Finally, add this element to the mesh.
109  this->_mesh.add_elem(std::move(elem));
110  }
111 }
std::string _switches
Parameter controlling the behaviour of tetgen.
UnstructuredMesh & _mesh
Local reference to the mesh we are working with.
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static std::unique_ptr< Elem > build(const ElemType type, Elem *p=nullptr)
Definition: elem.C:273
void assign_nodes_to_elem(unsigned *node_labels, Elem *elem)
Assigns the node IDs contained in the &#39;node_labels&#39; array to &#39;elem&#39;.
void fill_pointlist(TetGenWrapper &wrapper)
This function copies nodes from the _mesh into TetGen&#39;s pointlist.

Member Data Documentation

◆ _mesh

UnstructuredMesh& libMesh::TetGenMeshInterface::_mesh
protected

◆ _sequential_to_libmesh_node_map

std::vector<unsigned> libMesh::TetGenMeshInterface::_sequential_to_libmesh_node_map
protected

We should not assume libmesh nodes are numbered sequentially...

This is not the default behavior of DistributedMesh, for example, unless you specify node IDs explicitly. So this array allows us to keep a mapping between the sequential numbering in tetgen_data.pointlist.

Definition at line 160 of file mesh_tetgen_interface.h.

Referenced by assign_nodes_to_elem(), fill_pointlist(), and triangulate_conformingDelaunayMesh_carvehole().

◆ _serializer

MeshSerializer libMesh::TetGenMeshInterface::_serializer
protected

Tetgen only operates on serial meshes.

Definition at line 165 of file mesh_tetgen_interface.h.

◆ _switches

std::string libMesh::TetGenMeshInterface::_switches
protected

Parameter controlling the behaviour of tetgen.

By default quiet.

Definition at line 171 of file mesh_tetgen_interface.h.

Referenced by pointset_convexhull(), set_switches(), triangulate_conformingDelaunayMesh_carvehole(), and triangulate_pointset().


The documentation for this class was generated from the following files: