libMesh::TetGenMeshInterface Class Reference

#include <mesh_tetgen_interface.h>

Public Member Functions

 TetGenMeshInterface (UnstructuredMesh &mesh)
 
 ~TetGenMeshInterface ()
 
void set_switches (const std::string &)
 
void triangulate_pointset ()
 
void pointset_convexhull ()
 
void triangulate_conformingDelaunayMesh (double quality_constraint=0., double volume_constraint=0.)
 
void triangulate_conformingDelaunayMesh_carvehole (const std::vector< Point > &holes, double quality_constraint=0., double volume_constraint=0.)
 

Protected Member Functions

void fill_pointlist (TetGenWrapper &wrapper)
 
void assign_nodes_to_elem (unsigned *node_labels, Elem *elem)
 
unsigned check_hull_integrity ()
 
void process_hull_integrity_result (unsigned result)
 
void delete_2D_hull_elements ()
 

Protected Attributes

UnstructuredMesh_mesh
 
std::vector< unsigned > _sequential_to_libmesh_node_map
 
MeshSerializer _serializer
 
std::string _switches
 

Detailed Description

Class TetGenMeshInterface provides an interface for tetrahedrization 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

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

Constructor. Takes a reference to the mesh.

Definition at line 38 of file mesh_tetgen_interface.C.

38  :
39  _mesh(mesh),
41  _switches("Q")
42 {
43 }
MeshBase & mesh

Member Function Documentation

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 382 of file mesh_tetgen_interface.C.

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

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

383 {
384  for (unsigned int j=0; j<elem->n_nodes(); ++j)
385  {
386  // Get the mapped node index to ask the Mesh for
387  unsigned mapped_node_id = _sequential_to_libmesh_node_map[ node_labels[j] ];
388 
389  Node * current_node = this->_mesh.node_ptr( mapped_node_id );
390 
391  elem->set_node(j) = current_node;
392  }
393 }
std::vector< unsigned > _sequential_to_libmesh_node_map
virtual const Node * node_ptr(const dof_id_type i) const =0
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-NULL 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 NULL-neighbor is found

Definition at line 399 of file mesh_tetgen_interface.C.

References _mesh, libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), end, libmesh_nullptr, libMesh::MeshBase::n_elem(), libMesh::Elem::n_neighbors(), libMesh::Elem::neighbor(), libMesh::TRI3, and libMesh::Elem::type().

Referenced by triangulate_conformingDelaunayMesh_carvehole(), and ~TetGenMeshInterface().

400 {
401  // Check for easy return: if the Mesh is empty (i.e. if
402  // somebody called triangulate_conformingDelaunayMesh on
403  // a Mesh with no elements, then hull integrity check must
404  // fail...
405  if (_mesh.n_elem() == 0)
406  return 3;
407 
408  MeshBase::element_iterator it = this->_mesh.elements_begin();
409  const MeshBase::element_iterator end = this->_mesh.elements_end();
410 
411  for (; it != end ; ++it)
412  {
413  Elem * elem = *it;
414 
415  // Check for proper element type
416  if (elem->type() != TRI3)
417  {
418  //libmesh_error_msg("ERROR: Some of the elements in the original mesh were not TRI3!");
419  return 1;
420  }
421 
422  for (unsigned int i=0; i<elem->n_neighbors(); ++i)
423  {
424  if (elem->neighbor(i) == libmesh_nullptr)
425  {
426  // libmesh_error_msg("ERROR: Non-convex hull, cannot be tetrahedralized.");
427  return 2;
428  }
429  }
430  }
431 
432  // If we made it here, return success!
433  return 0;
434 }
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
virtual element_iterator elements_begin()=0
virtual element_iterator elements_end()=0
virtual dof_id_type n_elem() const =0
void libMesh::TetGenMeshInterface::delete_2D_hull_elements ( )
protected

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

Definition at line 459 of file mesh_tetgen_interface.C.

References _mesh, libMesh::MeshBase::delete_elem(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), end, libMesh::MeshBase::get_boundary_info(), libMesh::BoundaryInfo::regenerate_id_sets(), libMesh::TRI3, and libMesh::Elem::type().

Referenced by triangulate_conformingDelaunayMesh_carvehole(), and ~TetGenMeshInterface().

460 {
461  MeshBase::element_iterator it = this->_mesh.elements_begin();
462  const MeshBase::element_iterator end = this->_mesh.elements_end();
463 
464  for (; it != end ; ++it)
465  {
466  Elem * elem = *it;
467 
468  // Check for proper element type. Yes, we legally delete elements while
469  // iterating over them because no entries from the underlying container
470  // are actually erased.
471  if (elem->type() == TRI3)
472  _mesh.delete_elem(elem);
473  }
474 
475  // We just removed any boundary info associated with hull element
476  // edges, so let's update the boundary id caches.
478 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:111
IterBase * end
virtual element_iterator elements_begin()=0
virtual void delete_elem(Elem *e)=0
virtual element_iterator elements_end()=0
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 355 of file mesh_tetgen_interface.C.

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

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

356 {
357  // fill input structure with point set data:
358  wrapper.allocate_pointlist( this->_mesh.n_nodes() );
359 
360  // Make enough space to store a mapping between the implied sequential
361  // node numbering used in tetgen and libmesh's (possibly) non-sequential
362  // numbering scheme.
365 
366  {
367  unsigned index = 0;
368  MeshBase::node_iterator it = this->_mesh.nodes_begin();
369  const MeshBase::node_iterator end = this->_mesh.nodes_end();
370  for ( ; it != end; ++it)
371  {
372  _sequential_to_libmesh_node_map[index] = (*it)->id();
373  wrapper.set_node(index++, (**it)(0), (**it)(1), (**it)(2));
374  }
375  }
376 }
std::vector< unsigned > _sequential_to_libmesh_node_map
IterBase * end
virtual node_iterator nodes_begin()=0
virtual node_iterator nodes_end()=0
virtual dof_id_type n_nodes() const =0
void libMesh::TetGenMeshInterface::pointset_convexhull ( )

Method invokes TetGen library to compute a Delaunay tetrahedrization from the nodes point set. Stores only 2D hull surface elements.

Definition at line 113 of file mesh_tetgen_interface.C.

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

Referenced by ~TetGenMeshInterface().

114 {
115  // class tetgen_wrapper allows library access on a basic level
116  TetGenWrapper tetgen_wrapper;
117 
118  // Copy Mesh's node points into TetGen data structure
119  this->fill_pointlist(tetgen_wrapper);
120 
121  // Run TetGen triangulation method:
122  // Q = quiet, no terminal output
123  // Note: if no switch is used, the input must be a list of 3D points
124  // (.node file) and the Delaunay tetrahedralization of this point set
125  // will be generated. In this particular function, we are throwing
126  // away the tetrahedra generated by TetGen, and keeping only the
127  // convex hull...
128  tetgen_wrapper.set_switches(_switches);
129  tetgen_wrapper.run_tetgen();
130  unsigned int num_elements = tetgen_wrapper.get_numberoftrifaces();
131 
132  // Delete *all* old elements. Yes, we legally delete elements while
133  // iterating over them because no entries from the underlying container
134  // are actually erased.
135  {
136  MeshBase::element_iterator it = this->_mesh.elements_begin();
137  const MeshBase::element_iterator end = this->_mesh.elements_end();
138  for ( ; it != end; ++it)
139  this->_mesh.delete_elem (*it);
140  }
141 
142  // We just removed any boundary info associated with element faces
143  // or edges, so let's update the boundary id caches.
145 
146  // Add the 2D elements which comprise the convex hull back to the mesh.
147  // Vector that temporarily holds the node labels defining element.
148  unsigned int node_labels[3];
149 
150  for (unsigned int i=0; i<num_elements; ++i)
151  {
152  Elem * elem = new Tri3;
153 
154  // Get node labels associated with this element
155  for (unsigned int j=0; j<elem->n_nodes(); ++j)
156  node_labels[j] = tetgen_wrapper.get_triface_node(i,j);
157 
158  this->assign_nodes_to_elem(node_labels, elem);
159 
160  // Finally, add this element to the mesh.
161  this->_mesh.add_elem(elem);
162  }
163 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:111
IterBase * end
virtual element_iterator elements_begin()=0
virtual void delete_elem(Elem *e)=0
virtual Elem * add_elem(Elem *e)=0
virtual element_iterator elements_end()=0
void assign_nodes_to_elem(unsigned *node_labels, Elem *elem)
void fill_pointlist(TetGenWrapper &wrapper)
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 440 of file mesh_tetgen_interface.C.

References libMesh::err.

Referenced by triangulate_conformingDelaunayMesh_carvehole(), and ~TetGenMeshInterface().

441 {
442  if (result != 0)
443  {
444  libMesh::err << "Error! Conforming Delaunay mesh tetrahedralization requires a convex hull." << std::endl;
445 
446  if (result==1)
447  {
448  libMesh::err << "Non-TRI3 elements were found in the input Mesh. ";
449  libMesh::err << "A constrained Delaunay triangulation requires a convex hull of TRI3 elements." << std::endl;
450  }
451 
452  libmesh_error_msg("Consider calling TetGenMeshInterface::pointset_convexhull() followed by Mesh::find_neighbors() first.");
453  }
454 }
OStreamProxy err(std::cerr)
void libMesh::TetGenMeshInterface::set_switches ( const std::string &  switches)

Method to set swithes to tetgen, allowing for different behaviours

Definition at line 45 of file mesh_tetgen_interface.C.

References _switches.

Referenced by ~TetGenMeshInterface().

46 {
47  // set the tetgen switch manually:
48  // p = tetrahedralizes a piecewise linear complex (see definition in user manual)
49  // Q = quiet, no terminal output
50  // q = specify a minimum radius/edge ratio
51  // a = tetrahedron volume constraint
52  // V = verbose output
53  // for full list of options and their meaning: see the tetgen manual
54  // (http://wias-berlin.de/software/tetgen/1.5/doc/manual/manual005.html)
55  _switches = switches;
56 }
void libMesh::TetGenMeshInterface::triangulate_conformingDelaunayMesh ( double  quality_constraint = 0.,
double  volume_constraint = 0. 
)

Method invokes TetGen library to compute a Delaunay tetrahedrization from the nodes point set. Boundary constraints are taken from elements array.

Definition at line 169 of file mesh_tetgen_interface.C.

References triangulate_conformingDelaunayMesh_carvehole().

Referenced by ~TetGenMeshInterface().

171 {
172  // start triangulation method with empty holes list:
173  std::vector<Point> noholes;
174  triangulate_conformingDelaunayMesh_carvehole(noholes, quality_constraint, volume_constraint);
175 }
void triangulate_conformingDelaunayMesh_carvehole(const std::vector< Point > &holes, double quality_constraint=0., double volume_constraint=0.)
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 tetrahedrization from the nodes point set. Boundary constraints are taken from elements array. Include carve-out functionality.

Definition at line 179 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(), check_hull_integrity(), delete_2D_hull_elements(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), end, 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(), libMesh::Elem::n_nodes(), libMesh::Elem::node_id(), process_hull_integrity_result(), libMesh::TetGenWrapper::run_tetgen(), libMesh::TetGenWrapper::set_hole(), libMesh::TetGenWrapper::set_switches(), libMesh::TetGenWrapper::set_vertex(), and libMesh::x.

Referenced by triangulate_conformingDelaunayMesh(), and ~TetGenMeshInterface().

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

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

Definition at line 59 of file mesh_tetgen_interface.C.

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

Referenced by ~TetGenMeshInterface().

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

Member Data Documentation

UnstructuredMesh& libMesh::TetGenMeshInterface::_mesh
protected
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().

MeshSerializer libMesh::TetGenMeshInterface::_serializer
protected

Tetgen only operates on serial meshes.

Definition at line 165 of file mesh_tetgen_interface.h.

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

Parameter controling 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: