libMesh
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
libMesh::TriangleInterface Class Reference

A C++ interface between LibMesh and the Triangle library written by J.R. More...

#include <mesh_triangle_interface.h>

Inheritance diagram for libMesh::TriangleInterface:
[legend]

Public Types

enum  TriangulationType { GENERATE_CONVEX_HULL = 0, PSLG = 1, INVALID_TRIANGULATION_TYPE }
 The TriangulationType is used with the general triangulate function defined below. More...
 

Public Member Functions

 TriangleInterface (UnstructuredMesh &mesh)
 The constructor. More...
 
virtual ~TriangleInterface ()=default
 Empty destructor. More...
 
virtual void triangulate () override
 Internally, this calls Triangle's triangulate routine. More...
 
std::string & extra_flags ()
 Sets and/or gets additional flags to be passed to triangle. More...
 
ElemTypeelem_type ()
 Sets and/or gets the desired element type. More...
 
Realdesired_area ()
 Sets and/or gets the desired triangle area. More...
 
virtual void set_desired_area_function (FunctionBase< Real > *)
 Set a function giving desired triangle area as a function of position. More...
 
virtual FunctionBase< Real > * get_desired_area_function ()
 Get the function giving desired triangle area as a function of position, or nullptr if no such function has been set. More...
 
Realminimum_angle ()
 Sets and/or gets the minimum desired angle. More...
 
TriangulationTypetriangulation_type ()
 Sets and/or gets the desired triangulation type. More...
 
bool & insert_extra_points ()
 Sets and/or gets the flag for inserting add'l points. More...
 
void set_interpolate_boundary_points (int n_points)
 Complicated setter, for compatibility with insert_extra_points() More...
 
int get_interpolate_boundary_points () const
 Complicated getter, for compatibility with insert_extra_points() More...
 
virtual void set_refine_boundary_allowed (bool)
 Set whether or not the triangulation is allowed to refine the mesh boundary when refining the interior. More...
 
virtual bool refine_boundary_allowed () const
 Get whether or not the triangulation is allowed to refine the mesh boundary when refining the interior. More...
 
bool & smooth_after_generating ()
 Sets/gets flag which tells whether to do two steps of Laplace mesh smoothing after generating the grid. More...
 
bool & quiet ()
 Whether not to silence internal messages to stdout. More...
 
void attach_hole_list (const std::vector< Hole *> *holes)
 Attaches a vector of Hole* pointers which will be meshed around. More...
 
void set_verify_hole_boundaries (bool v)
 Verifying that hole boundaries don't cross the outer boundary or each other is something like O(N_bdys^2*N_points_per_bdy^2), so we only do it if requested. More...
 
bool get_verify_hole_boundaries () const
 
void attach_boundary_marker (const std::vector< int > *markers)
 Attaches boundary markers. More...
 
void attach_region_list (const std::vector< Region *> *regions)
 Attaches regions for using attribute to set subdomain IDs and better controlling the triangle sizes within the regions. More...
 
void set_outer_boundary_ids (std::set< std::size_t > bdy_ids)
 A set of ids to allow on the outer boundary loop: interpreted as boundary ids of 2D elements and/or subdomain ids of 1D edges. More...
 
const std::set< std::size_t > & get_outer_boundary_ids () const
 

Public Attributes

std::vector< std::pair< unsigned int, unsigned int > > segments
 When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly through the ordering of the points, you can use the segments vector to specify the segments explicitly, Ex: unit square numbered counter-clockwise starting from origin segments[0] = (0,1) segments[1] = (1,2) segments[2] = (2,3) segments[3] = (3,0) (For the above case you could actually use the implicit ordering!) More...
 

Protected Member Functions

void elems_to_segments ()
 Helper function to create PSLG segments from our other boundary-defining options (1D mesh edges, 2D mesh boundary sides), if no segments already exist. More...
 
void nodes_to_segments (dof_id_type max_node_id)
 Helper function to create PSLG segments from our node ordering, up to the maximum node id, if no segments already exist. More...
 
void insert_any_extra_boundary_points ()
 Helper function to add extra points (midpoints of initial segments) to a PSLG triangulation. More...
 
void verify_holes (const Hole &outer_bdy)
 Helper function to check holes for intersections if requested. More...
 
unsigned int total_hole_points ()
 Helper function to count points in and verify holes. More...
 

Protected Attributes

UnstructuredMesh_mesh
 Reference to the mesh which is to be created by triangle. More...
 
const std::vector< Hole * > * _holes
 A pointer to a vector of Hole*s. More...
 
const std::vector< int > * _markers
 Boundary markers. More...
 
const std::vector< Region * > * _regions
 A pointer to a vector of Regions*s. More...
 
std::set< std::size_t > _bdy_ids
 A set of ids to allow on the outer boundary loop. More...
 
ElemType _elem_type
 The type of elements to generate. More...
 
Real _desired_area
 The desired area for the elements in the resulting mesh. More...
 
Real _minimum_angle
 Minimum angle in triangles. More...
 
TriangulationType _triangulation_type
 The type of triangulation to perform: choices are: convex hull PSLG. More...
 
bool _insert_extra_points
 Flag which tells whether or not to insert additional nodes before triangulation. More...
 
int _interpolate_boundary_points
 Flag which tells how many additional nodes should be inserted between each pair of original mesh points. More...
 
bool _smooth_after_generating
 Flag which tells whether we should smooth the mesh after it is generated. More...
 
bool _quiet
 Flag which tells if we want to suppress stdout outputs. More...
 
bool _verify_hole_boundaries
 Flag which tells if we want to check hole geometry. More...
 

Private Attributes

std::string _extra_flags
 Additional flags to be passed to triangle. More...
 
MeshSerializer _serializer
 Triangle only operates on serial meshes. More...
 

Detailed Description

A C++ interface between LibMesh and the Triangle library written by J.R.

Shewchuk.

Author
John W. Peterson
Date
2011

Definition at line 41 of file mesh_triangle_interface.h.

Member Enumeration Documentation

◆ TriangulationType

The TriangulationType is used with the general triangulate function defined below.

Enumerator
GENERATE_CONVEX_HULL 

First generate a convex hull from the set of points passed in, and then triangulate this set of points.

This is probably the most common type of usage.

PSLG 

Triangulate the interior of a Planar Straight Line Graph, which is defined implicitly by the order of the "points" vector: a straight line is assumed to lie between each successive pair of points, with an additional line joining the final and first points.

Explicitly telling the triangulator to add additional points may be important for this option.

INVALID_TRIANGULATION_TYPE 

Does nothing, used as a "null" value.

Definition at line 61 of file triangulator_interface.h.

62  {
69 
80  PSLG = 1,
81 
86  };
First generate a convex hull from the set of points passed in, and then triangulate this set of point...
Triangulate the interior of a Planar Straight Line Graph, which is defined implicitly by the order of...

Constructor & Destructor Documentation

◆ TriangleInterface()

libMesh::TriangleInterface::TriangleInterface ( UnstructuredMesh mesh)
explicit

The constructor.

A reference to the mesh containing the points which are to be triangulated must be provided. Unless otherwise specified, a convex hull will be computed for the set of input points and the convex hull will be meshed.

Definition at line 48 of file mesh_triangle_interface.C.

50  _extra_flags(""),
52 {}
std::string _extra_flags
Additional flags to be passed to triangle.
MeshBase & mesh
UnstructuredMesh & _mesh
Reference to the mesh which is to be created by triangle.
TriangulatorInterface(UnstructuredMesh &mesh)
The constructor.
MeshSerializer _serializer
Triangle only operates on serial meshes.

◆ ~TriangleInterface()

virtual libMesh::TriangleInterface::~TriangleInterface ( )
virtualdefault

Empty destructor.

Member Function Documentation

◆ attach_boundary_marker()

void libMesh::TriangulatorInterface::attach_boundary_marker ( const std::vector< int > *  markers)
inlineinherited

Attaches boundary markers.

If segments is set, the number of markers must be equal to the size of segments, otherwise, it is equal to the number of points.

Definition at line 236 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_markers.

236 { _markers = markers; }
const std::vector< int > * _markers
Boundary markers.

◆ attach_hole_list()

void libMesh::TriangulatorInterface::attach_hole_list ( const std::vector< Hole *> *  holes)
inlineinherited

◆ attach_region_list()

void libMesh::TriangulatorInterface::attach_region_list ( const std::vector< Region *> *  regions)
inlineinherited

Attaches regions for using attribute to set subdomain IDs and better controlling the triangle sizes within the regions.

Definition at line 242 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_regions.

242 { _regions = regions; }
const std::vector< Region * > * _regions
A pointer to a vector of Regions*s.

◆ desired_area()

Real& libMesh::TriangulatorInterface::desired_area ( )
inlineinherited

Sets and/or gets the desired triangle area.

Set to zero to disable area constraint.

If a desired_area_function is set, then desired_area() should be used to set a minimum desired area; this will reduce "false negatives" by suggesting how finely to sample desired_area_function inside large triangles, where ideally the desired_area_function will be satisfied in the triangle interior and not just at the triangle vertices.

Definition at line 127 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_desired_area.

Referenced by libMesh::MeshTools::Generation::build_delaunay_square(), MeshTriangulationTest::commonSettings(), libMesh::Poly2TriTriangulator::insert_refinement_points(), libMesh::Poly2TriTriangulator::should_refine_elem(), MeshTriangulationTest::testPoly2TriRefinementBase(), MeshTriangulationTest::testTriangulatorInterp(), and triangulate_domain().

127 {return _desired_area;}
Real _desired_area
The desired area for the elements in the resulting mesh.

◆ elem_type()

ElemType& libMesh::TriangulatorInterface::elem_type ( )
inlineinherited

Sets and/or gets the desired element type.

Definition at line 114 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_elem_type.

Referenced by libMesh::MeshTools::Generation::build_delaunay_square().

114 {return _elem_type;}
ElemType _elem_type
The type of elements to generate.

◆ elems_to_segments()

void libMesh::TriangulatorInterface::elems_to_segments ( )
protectedinherited

Helper function to create PSLG segments from our other boundary-defining options (1D mesh edges, 2D mesh boundary sides), if no segments already exist.

Definition at line 87 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_bdy_ids, libMesh::TriangulatorInterface::_holes, libMesh::TriangulatorInterface::_mesh, libMesh::TriangulatorInterface::_verify_hole_boundaries, libMesh::MeshBase::clear_elems(), libMesh::MeshBase::delete_node(), libMesh::make_range(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::node_ptr(), libMesh::TriangulatorInterface::segments, and libMesh::TriangulatorInterface::verify_holes().

Referenced by triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

88 {
89  // Don't try to override manually specified segments
90  if (!this->segments.empty())
91  return;
92 
93  // If we have edges, they should form the polyline with the ordering
94  // we want. Let's turn them into segments for later use, because
95  // we're going to delete the original elements to replace with our
96  // triangulation.
97  if (_mesh.n_elem())
98  {
99  // Mapping from points to node ids, to back those out from
100  // MeshedHole results later
101  std::map<Point, dof_id_type> point_id_map;
102 
103  for (Node * node : _mesh.node_ptr_range())
104  {
105  // We're not going to support overlapping nodes on the boundary
106  libmesh_error_msg_if
107  (point_id_map.count(*node),
108  "TriangulatorInterface does not support overlapping nodes found at "
109  << static_cast<Point&>(*node));
110 
111  point_id_map.emplace(*node, node->id());
112  }
113 
114  // We'll steal the ordering calculation from
115  // the MeshedHole code
116  const TriangulatorInterface::MeshedHole mh { _mesh, this->_bdy_ids };
117 
118  // And now we're done with elements. Delete them lest they have
119  // dangling pointers to nodes we'll be deleting.
120  _mesh.clear_elems();
121 
122  // If we've specified only a subset of the mesh as our outer
123  // boundary, then we may have nodes that don't actually fall
124  // inside that boundary. Triangulator code doesn't like Steiner
125  // points that aren't inside the triangulation domain, so we
126  // need to get rid of them.
127  std::unordered_set<Node *> nodes_to_delete;
128  if (!this->_bdy_ids.empty())
129  {
130  for (auto & node : _mesh.node_ptr_range())
131  if (!mh.contains(*node))
132  nodes_to_delete.insert(node);
133  }
134 
135  // Make segments from boundary nodes; also make sure we don't
136  // delete them.
137  const std::size_t np = mh.n_points();
138  for (auto i : make_range(np))
139  {
140  const Point pt = mh.point(i);
141  const dof_id_type id0 = libmesh_map_find(point_id_map, pt);
142  nodes_to_delete.erase(_mesh.node_ptr(id0));
143  const Point next_pt = mh.point((np+i+1)%np);
144  const dof_id_type id1 = libmesh_map_find(point_id_map, next_pt);
145  this->segments.emplace_back(id0, id1);
146  }
147 
148  for (Node * node : nodes_to_delete)
149  _mesh.delete_node(node);
150 
151  if (this->_verify_hole_boundaries && _holes)
152  this->verify_holes(mh);
153  }
154 }
bool _verify_hole_boundaries
Flag which tells if we want to check hole geometry.
std::set< std::size_t > _bdy_ids
A set of ids to allow on the outer boundary loop.
UnstructuredMesh & _mesh
Reference to the mesh which is to be created by triangle.
const std::vector< Hole * > * _holes
A pointer to a vector of Hole*s.
virtual void delete_node(Node *n)=0
Removes the Node n from the mesh.
virtual void clear_elems()=0
Deletes all the element data that is currently stored.
std::vector< std::pair< unsigned int, unsigned int > > segments
When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly ...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
void verify_holes(const Hole &outer_bdy)
Helper function to check holes for intersections if requested.
virtual dof_id_type n_elem() const =0
virtual const Node * node_ptr(const dof_id_type i) const =0
uint8_t dof_id_type
Definition: id_types.h:67

◆ extra_flags()

std::string& libMesh::TriangleInterface::extra_flags ( )
inline

Sets and/or gets additional flags to be passed to triangle.

Definition at line 66 of file mesh_triangle_interface.h.

References _extra_flags.

66 {return _extra_flags;}
std::string _extra_flags
Additional flags to be passed to triangle.

◆ get_desired_area_function()

virtual FunctionBase<Real>* libMesh::TriangulatorInterface::get_desired_area_function ( )
inlinevirtualinherited

Get the function giving desired triangle area as a function of position, or nullptr if no such function has been set.

Reimplemented in libMesh::Poly2TriTriangulator.

Definition at line 143 of file triangulator_interface.h.

144  { return nullptr; }

◆ get_interpolate_boundary_points()

int libMesh::TriangulatorInterface::get_interpolate_boundary_points ( ) const
inherited

Complicated getter, for compatibility with insert_extra_points()

Definition at line 75 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_insert_extra_points, and libMesh::TriangulatorInterface::_interpolate_boundary_points.

Referenced by libMesh::TriangulatorInterface::insert_any_extra_boundary_points().

76 {
77  // backwards compatibility - someone might have turned this off via
78  // the old API
80  return 0;
81 
83 }
bool _insert_extra_points
Flag which tells whether or not to insert additional nodes before triangulation.
int _interpolate_boundary_points
Flag which tells how many additional nodes should be inserted between each pair of original mesh poin...

◆ get_outer_boundary_ids()

const std::set<std::size_t>& libMesh::TriangulatorInterface::get_outer_boundary_ids ( ) const
inlineinherited

Definition at line 251 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_bdy_ids.

251 { return _bdy_ids; }
std::set< std::size_t > _bdy_ids
A set of ids to allow on the outer boundary loop.

◆ get_verify_hole_boundaries()

bool libMesh::TriangulatorInterface::get_verify_hole_boundaries ( ) const
inlineinherited

Definition at line 214 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_verify_hole_boundaries.

214 {return _verify_hole_boundaries;}
bool _verify_hole_boundaries
Flag which tells if we want to check hole geometry.

◆ insert_any_extra_boundary_points()

void libMesh::TriangulatorInterface::insert_any_extra_boundary_points ( )
protectedinherited

Helper function to add extra points (midpoints of initial segments) to a PSLG triangulation.

Definition at line 197 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_mesh, libMesh::TriangulatorInterface::_triangulation_type, libMesh::MeshBase::add_point(), libMesh::TriangulatorInterface::get_interpolate_boundary_points(), libMesh::DofObject::id(), libMesh::libmesh_assert(), libMesh::make_range(), libMesh::MeshBase::max_node_id(), libMesh::MeshBase::node_ptr(), libMesh::TriangulatorInterface::PSLG, and libMesh::TriangulatorInterface::segments.

Referenced by triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

198 {
199  // If the initial PSLG is really simple, e.g. an L-shaped domain or
200  // a square/rectangle, the resulting triangulation may be very
201  // "structured" looking. Sometimes this is a problem if your
202  // intention is to work with an "unstructured" looking grid. We can
203  // attempt to work around this limitation by inserting midpoints
204  // into the original PSLG. Inserting additional points into a
205  // set of points meant to be a convex hull usually makes less sense.
206 
207  const int n_interpolated = this->get_interpolate_boundary_points();
208  if ((_triangulation_type==PSLG) && n_interpolated)
209  {
210  // If we were lucky enough to start with contiguous node ids,
211  // let's keep them that way.
213 
214  std::vector<std::pair<unsigned int, unsigned int>> old_segments =
215  std::move(this->segments);
216 
217  // We expect to have converted any elems and/or nodes into
218  // segments by now.
219  libmesh_assert(!old_segments.empty());
220 
221  this->segments.clear();
222 
223  // Insert a new point on each segment at evenly spaced locations
224  // between existing boundary points.
225  // np=index into new points vector
226  // n =index into original points vector
227  for (auto old_segment : old_segments)
228  {
229  Node * begin_node = _mesh.node_ptr(old_segment.first);
230  Node * end_node = _mesh.node_ptr(old_segment.second);
231  dof_id_type current_id = begin_node->id();
232  for (auto i : make_range(n_interpolated))
233  {
234  // new points are equispaced along the original segments
235  const Point new_point =
236  ((n_interpolated-i) * *(Point *)(begin_node) +
237  (i+1) * *(Point *)(end_node)) /
238  (n_interpolated + 1);
239  Node * next_node = _mesh.add_point(new_point, nn++);
240  this->segments.emplace_back(current_id,
241  next_node->id());
242  current_id = next_node->id();
243  }
244  this->segments.emplace_back(current_id,
245  end_node->id());
246  }
247  }
248 }
UnstructuredMesh & _mesh
Reference to the mesh which is to be created by triangle.
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.
libmesh_assert(ctx)
Triangulate the interior of a Planar Straight Line Graph, which is defined implicitly by the order of...
std::vector< std::pair< unsigned int, unsigned int > > segments
When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly ...
int get_interpolate_boundary_points() const
Complicated getter, for compatibility with insert_extra_points()
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
virtual dof_id_type max_node_id() const =0
virtual const Node * node_ptr(const dof_id_type i) const =0
TriangulationType _triangulation_type
The type of triangulation to perform: choices are: convex hull PSLG.
uint8_t dof_id_type
Definition: id_types.h:67

◆ insert_extra_points()

bool& libMesh::TriangulatorInterface::insert_extra_points ( )
inlineinherited

Sets and/or gets the flag for inserting add'l points.

Definition at line 160 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_insert_extra_points.

160 {return _insert_extra_points;}
bool _insert_extra_points
Flag which tells whether or not to insert additional nodes before triangulation.

◆ minimum_angle()

Real& libMesh::TriangulatorInterface::minimum_angle ( )
inlineinherited

Sets and/or gets the minimum desired angle.

Set to zero to disable angle constraint.

Definition at line 150 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_minimum_angle.

Referenced by MeshTriangulationTest::commonSettings(), and libMesh::Poly2TriTriangulator::insert_refinement_points().

150 {return _minimum_angle;}
Real _minimum_angle
Minimum angle in triangles.

◆ nodes_to_segments()

void libMesh::TriangulatorInterface::nodes_to_segments ( dof_id_type  max_node_id)
protectedinherited

Helper function to create PSLG segments from our node ordering, up to the maximum node id, if no segments already exist.

Definition at line 158 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_holes, libMesh::TriangulatorInterface::_mesh, libMesh::TriangulatorInterface::_triangulation_type, libMesh::TriangulatorInterface::_verify_hole_boundaries, libMesh::DofObject::id(), libMesh::MeshBase::point(), libMesh::TriangulatorInterface::PSLG, libMesh::TriangulatorInterface::segments, and libMesh::TriangulatorInterface::verify_holes().

Referenced by triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

159 {
160  // Don't try to override manually specified segments, or try to add
161  // segments if we're doing a convex hull
162  if (!this->segments.empty() || _triangulation_type != PSLG)
163  return;
164 
165  for (auto node_it = _mesh.nodes_begin(),
166  node_end = _mesh.nodes_end();
167  node_it != node_end;)
168  {
169  Node * node = *node_it;
170 
171  // If we're out of boundary nodes, the rest are going to be
172  // Steiner points or hole points
173  if (node->id() >= max_node_id)
174  break;
175 
176  ++node_it;
177 
178  Node * next_node = (node_it == node_end) ?
179  *_mesh.nodes_begin() : *node_it;
180 
181  this->segments.emplace_back(node->id(), next_node->id());
182  }
183 
184  if (this->_verify_hole_boundaries && _holes)
185  {
186  std::vector<Point> outer_pts;
187  for (auto segment : this->segments)
188  outer_pts.push_back(_mesh.point(segment.first));
189 
190  ArbitraryHole ah(outer_pts);
191  this->verify_holes(ah);
192  }
193 }
bool _verify_hole_boundaries
Flag which tells if we want to check hole geometry.
UnstructuredMesh & _mesh
Reference to the mesh which is to be created by triangle.
const std::vector< Hole * > * _holes
A pointer to a vector of Hole*s.
Triangulate the interior of a Planar Straight Line Graph, which is defined implicitly by the order of...
std::vector< std::pair< unsigned int, unsigned int > > segments
When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly ...
virtual const Point & point(const dof_id_type i) const =0
void verify_holes(const Hole &outer_bdy)
Helper function to check holes for intersections if requested.
TriangulationType _triangulation_type
The type of triangulation to perform: choices are: convex hull PSLG.

◆ quiet()

bool& libMesh::TriangulatorInterface::quiet ( )
inlineinherited

Whether not to silence internal messages to stdout.

Definition at line 199 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_quiet.

199 {return _quiet;}
bool _quiet
Flag which tells if we want to suppress stdout outputs.

◆ refine_boundary_allowed()

virtual bool libMesh::TriangulatorInterface::refine_boundary_allowed ( ) const
inlinevirtualinherited

Get whether or not the triangulation is allowed to refine the mesh boundary when refining the interior.

True by default.

Reimplemented in libMesh::Poly2TriTriangulator.

Definition at line 187 of file triangulator_interface.h.

188  { return true; }

◆ set_desired_area_function()

virtual void libMesh::TriangulatorInterface::set_desired_area_function ( FunctionBase< Real > *  )
inlinevirtualinherited

Set a function giving desired triangle area as a function of position.

Set this to nullptr to disable position-dependent area constraint (falling back on desired_area()).

This may not be implemented in all subclasses.

Reimplemented in libMesh::Poly2TriTriangulator.

Definition at line 136 of file triangulator_interface.h.

137  { libmesh_not_implemented(); }

◆ set_interpolate_boundary_points()

void libMesh::TriangulatorInterface::set_interpolate_boundary_points ( int  n_points)
inherited

Complicated setter, for compatibility with insert_extra_points()

Definition at line 61 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_insert_extra_points, libMesh::TriangulatorInterface::_interpolate_boundary_points, and libMesh::libmesh_assert().

Referenced by MeshTriangulationTest::testTriangulatorInterp().

62 {
63  // Maybe we'll reserve a meaning for negatives later?
64  libmesh_assert(n_points >= 0);
65 
67 
68  // backwards compatibility - someone (including us) might want to
69  // query this via the old API.
70  _insert_extra_points = n_points;
71 }
bool _insert_extra_points
Flag which tells whether or not to insert additional nodes before triangulation.
int _interpolate_boundary_points
Flag which tells how many additional nodes should be inserted between each pair of original mesh poin...
libmesh_assert(ctx)

◆ set_outer_boundary_ids()

void libMesh::TriangulatorInterface::set_outer_boundary_ids ( std::set< std::size_t >  bdy_ids)
inlineinherited

A set of ids to allow on the outer boundary loop: interpreted as boundary ids of 2D elements and/or subdomain ids of 1D edges.

If this is empty, then the outer boundary may be constructed from boundary edges of any id!

Definition at line 250 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_bdy_ids.

Referenced by MeshTriangulationTest::testHalfDomain().

250 { _bdy_ids = std::move(bdy_ids); }
std::set< std::size_t > _bdy_ids
A set of ids to allow on the outer boundary loop.

◆ set_refine_boundary_allowed()

virtual void libMesh::TriangulatorInterface::set_refine_boundary_allowed ( bool  )
inlinevirtualinherited

Set whether or not the triangulation is allowed to refine the mesh boundary when refining the interior.

This is true by default, but may be set to false to make the mesh boundary more predictable (and so easier to stitch to other meshes) later.

This may not be implemented in all subclasses.

Reimplemented in libMesh::Poly2TriTriangulator.

Definition at line 180 of file triangulator_interface.h.

181  { libmesh_not_implemented(); }

◆ set_verify_hole_boundaries()

void libMesh::TriangulatorInterface::set_verify_hole_boundaries ( bool  v)
inlineinherited

Verifying that hole boundaries don't cross the outer boundary or each other is something like O(N_bdys^2*N_points_per_bdy^2), so we only do it if requested.

Definition at line 212 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_verify_hole_boundaries.

Referenced by MeshTriangulationTest::commonSettings().

bool _verify_hole_boundaries
Flag which tells if we want to check hole geometry.

◆ smooth_after_generating()

bool& libMesh::TriangulatorInterface::smooth_after_generating ( )
inlineinherited

Sets/gets flag which tells whether to do two steps of Laplace mesh smoothing after generating the grid.

Definition at line 194 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_smooth_after_generating.

Referenced by MeshTriangulationTest::commonSettings(), and triangulate_domain().

194 {return _smooth_after_generating;}
bool _smooth_after_generating
Flag which tells whether we should smooth the mesh after it is generated.

◆ total_hole_points()

unsigned int libMesh::TriangulatorInterface::total_hole_points ( )
protectedinherited

Helper function to count points in and verify holes.

Definition at line 276 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_holes.

Referenced by triangulate().

277 {
278  // If the holes vector is non-nullptr (and non-empty) we need to determine
279  // the number of additional points which the holes will add to the
280  // triangulation.
281  // Note that the number of points is always equal to the number of segments
282  // that form the holes.
283  unsigned int n_hole_points = 0;
284 
285  if (_holes)
286  for (const auto & hole : *_holes)
287  {
288  n_hole_points += hole->n_points();
289  // A hole at least has one enclosure.
290  // Points on enclosures are ordered so that we can add segments implicitly.
291  // Elements in segment_indices() indicates the starting points of all enclosures.
292  // The last element in segment_indices() is the number of total points.
293  libmesh_assert_greater(hole->segment_indices().size(), 1);
294  libmesh_assert_equal_to(hole->segment_indices().back(), hole->n_points());
295  }
296 
297  return n_hole_points;
298 }
const std::vector< Hole * > * _holes
A pointer to a vector of Hole*s.

◆ triangulate()

void libMesh::TriangleInterface::triangulate ( )
overridevirtual

Internally, this calls Triangle's triangulate routine.

Implements libMesh::TriangulatorInterface.

Definition at line 58 of file mesh_triangle_interface.C.

References libMesh::TriangulatorInterface::_desired_area, libMesh::TriangulatorInterface::_elem_type, _extra_flags, libMesh::TriangulatorInterface::_holes, libMesh::TriangulatorInterface::_markers, libMesh::TriangulatorInterface::_mesh, libMesh::TriangulatorInterface::_minimum_angle, libMesh::TriangulatorInterface::_quiet, libMesh::TriangulatorInterface::_regions, libMesh::TriangulatorInterface::_smooth_after_generating, libMesh::TriangulatorInterface::_triangulation_type, libMesh::TriangleWrapper::copy_tri_to_mesh(), libMesh::TriangleWrapper::destroy(), libMesh::TriangulatorInterface::elems_to_segments(), libMesh::Utility::enum_to_string(), libMesh::TriangulatorInterface::GENERATE_CONVEX_HULL, libMesh::TriangleWrapper::init(), libMesh::TriangleWrapper::INPUT, libMesh::TriangulatorInterface::insert_any_extra_boundary_points(), libMesh::TriangulatorInterface::INVALID_TRIANGULATION_TYPE, libMesh::MeshBase::max_node_id(), libMesh::MeshBase::n_nodes(), libMesh::TriangulatorInterface::nodes_to_segments(), libMesh::TriangleWrapper::OUTPUT, libMesh::MeshBase::prepare_for_use(), libMesh::TriangulatorInterface::PSLG, libMesh::TriangulatorInterface::segments, libMesh::MeshBase::set_mesh_dimension(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::TOLERANCE, libMesh::TriangulatorInterface::total_hole_points(), libMesh::TRI3, and libMesh::TRI6.

Referenced by libMesh::MeshTools::Generation::build_delaunay_square(), and triangulate_domain().

59 {
60  // Will the triangulation have holes?
61  const bool have_holes = ((_holes != nullptr) && (!_holes->empty()));
62 
63  unsigned int n_hole_points = this->total_hole_points();
64 
65  // If we have no explicit segments defined, we may get them from
66  // mesh elements
67  this->elems_to_segments();
68 
69  // If we're doing PSLG without segments, construct them from all our
70  // mesh nodes
72 
73  // Insert additional new points in between existing boundary points,
74  // if that is requested and reasonable
76 
77  // Regardless of whether we added additional points, the set of points to
78  // triangulate is now sitting in the mesh.
79 
80  // Triangle data structure for the mesh
81  TriangleWrapper::triangulateio initial;
82  TriangleWrapper::triangulateio final;
83  TriangleWrapper::triangulateio voronoi;
84 
85  // Pseudo-Constructor for the triangle io structs
86  TriangleWrapper::init(initial);
87  TriangleWrapper::init(final);
88  TriangleWrapper::init(voronoi);
89 
90  initial.numberofpoints = _mesh.n_nodes() + n_hole_points;
91  initial.pointlist = static_cast<REAL*>(std::malloc(initial.numberofpoints * 2 * sizeof(REAL)));
92 
94  {
95  // Implicit segment ordering: One segment per point, including hole points
96  if (this->segments.empty())
97  initial.numberofsegments = initial.numberofpoints;
98 
99  // User-defined segment ordering: One segment per entry in the segments vector
100  else
101  initial.numberofsegments = this->segments.size() + n_hole_points;
102  }
103 
105  initial.numberofsegments = n_hole_points; // One segment for each hole point
106 
107  // Allocate space for the segments (2 int per segment)
108  if (initial.numberofsegments > 0)
109  {
110  initial.segmentlist = static_cast<int *> (std::malloc(initial.numberofsegments * 2 * sizeof(int)));
111  if (_markers)
112  initial.segmentmarkerlist = static_cast<int *> (std::malloc(initial.numberofsegments * sizeof(int)));
113  }
114 
115 
116  // Copy all the holes' points and segments into the triangle struct.
117 
118  // The hole_offset is a constant offset into the points vector which points
119  // past the end of the last hole point added.
120  unsigned int hole_offset=0;
121 
122  if (have_holes)
123  for (const auto & hole : *_holes)
124  {
125  for (unsigned int ctr=0, h=0, i=0, hsism=hole->segment_indices().size()-1; i<hsism; ++i)
126  {
127  unsigned int begp = hole_offset + hole->segment_indices()[i];
128  unsigned int endp = hole->segment_indices()[i+1];
129 
130  for (; h<endp; ctr+=2, ++h)
131  {
132  Point p = hole->point(h);
133 
134  const unsigned int index0 = 2*hole_offset+ctr;
135  const unsigned int index1 = 2*hole_offset+ctr+1;
136 
137  // Save the x,y locations in the triangle struct.
138  initial.pointlist[index0] = p(0);
139  initial.pointlist[index1] = p(1);
140 
141  // Set the points which define the segments
142  initial.segmentlist[index0] = hole_offset+h;
143  initial.segmentlist[index1] = (h == endp - 1) ? begp : hole_offset + h + 1; // wrap around
144  if (_markers)
145  // 1 is reserved for boundaries of holes
146  initial.segmentmarkerlist[hole_offset+h] = 1;
147  }
148  }
149 
150  // Update the hole_offset for the next hole
151  hole_offset += hole->n_points();
152  }
153 
154 
155  // Copy all the non-hole points and segments into the triangle struct.
156  std::vector<unsigned int> libmesh_id_to_pointlist_index(_mesh.max_node_id());
157  {
158  dof_id_type ctr=0;
159  for (auto & node : _mesh.node_ptr_range())
160  {
161  dof_id_type index = 2*hole_offset + ctr;
162 
163  // Set x,y values in pointlist
164  initial.pointlist[index] = (*node)(0);
165  initial.pointlist[index+1] = (*node)(1);
166  libmesh_id_to_pointlist_index[node->id()] = ctr/2;
167 
168  // If the user requested a PSLG, the non-hole points are also segments
170  {
171  // Use implicit ordering to define segments
172  if (this->segments.empty())
173  {
174  dof_id_type n = ctr/2; // ctr is always even
175  initial.segmentlist[index] = hole_offset+n;
176  initial.segmentlist[index+1] = (n==_mesh.n_nodes()-1) ? hole_offset : hole_offset+n+1; // wrap around
177  if (_markers)
178  initial.segmentmarkerlist[hole_offset + n] = (*_markers)[n];
179  }
180  }
181 
182  ctr +=2;
183  }
184  }
185 
186 
187  // If the user provided it, use his ordering to define the segments
188  for (std::size_t ctr=0, s=0, ss=this->segments.size(); s<ss; ctr+=2, ++s)
189  {
190  const unsigned int index0 = 2*hole_offset+ctr;
191  const unsigned int index1 = 2*hole_offset+ctr+1;
192 
193  initial.segmentlist[index0] = hole_offset +
194  libmesh_id_to_pointlist_index[this->segments[s].first];
195  initial.segmentlist[index1] = hole_offset +
196  libmesh_id_to_pointlist_index[this->segments[s].second];
197  if (_markers)
198  initial.segmentmarkerlist[hole_offset + s] = (*_markers)[s];
199  }
200 
201 
202 
203  // Tell the input struct about the holes
204  if (have_holes)
205  {
206  initial.numberofholes = _holes->size();
207  initial.holelist = static_cast<REAL*>(std::malloc(initial.numberofholes * 2 * sizeof(REAL)));
208  for (std::size_t i=0, ctr=0, hs=_holes->size(); i<hs; ++i, ctr+=2)
209  {
210  Point inside_point = (*_holes)[i]->inside();
211  initial.holelist[ctr] = inside_point(0);
212  initial.holelist[ctr+1] = inside_point(1);
213  }
214  }
215 
216  if (_regions)
217  {
218  initial.numberofregions = _regions->size();
219  initial.regionlist = static_cast<REAL*>(std::malloc(initial.numberofregions * 4 * sizeof(REAL)));
220  for (std::size_t i=0, ctr=0, rs=_regions->size(); i<rs; ++i, ctr+=4)
221  {
222  Point inside_point = (*_regions)[i]->inside();
223  initial.regionlist[ctr] = inside_point(0);
224  initial.regionlist[ctr+1] = inside_point(1);
225  initial.regionlist[ctr+2] = (*_regions)[i]->attribute();
226  initial.regionlist[ctr+3] = (*_regions)[i]->max_area();
227  }
228  }
229 
230  // Set the triangulation flags.
231  // c ~ enclose convex hull with segments
232  // z ~ use zero indexing
233  // B ~ Suppresses boundary markers in the output
234  // Q ~ run in "quiet" mode
235  // p ~ Triangulates a Planar Straight Line Graph
236  // If the `p' switch is used, `segmentlist' must point to a list of
237  // segments, `numberofsegments' must be properly set, and
238  // `segmentmarkerlist' must either be set to nullptr (in which case all
239  // markers default to zero), or must point to a list of markers.
240  // D ~ Conforming Delaunay: use this switch if you want all triangles
241  // in the mesh to be Delaunay, and not just constrained Delaunay
242  // q ~ Quality mesh generation with no angles smaller than 20 degrees.
243  // An alternate minimum angle may be specified after the q
244  // a ~ Imposes a maximum triangle area constraint.
245  // -P Suppresses the output .poly file. Saves disk space, but you lose the ability to maintain
246  // constraining segments on later refinements of the mesh.
247  // -e Outputs (to an .edge file) a list of edges of the triangulation.
248  // -v Outputs the Voronoi diagram associated with the triangulation.
249  // Create the flag strings, depends on element type
250  std::ostringstream flags;
251 
252  // Default flags always used
253  flags << "z";
254 
255  if (_quiet)
256  flags << "QP";
257  else
258  flags << "V";
259 
260  if (_markers)
261  flags << "ev";
262 
263  // Flags which are specific to the type of triangulation
264  switch (_triangulation_type)
265  {
267  {
268  flags << "c";
269  break;
270  }
271 
272  case PSLG:
273  {
274  flags << "p";
275  break;
276  }
277 
279  libmesh_error_msg("ERROR: INVALID_TRIANGULATION_TYPE selected!");
280 
281  default:
282  libmesh_error_msg("Unrecognized _triangulation_type");
283  }
284 
285 
286  // Flags specific to the type of element
287  switch (_elem_type)
288  {
289  case TRI3:
290  {
291  // do nothing.
292  break;
293  }
294 
295  case TRI6:
296  {
297  flags << "o2";
298  break;
299  }
300 
301  default:
302  libmesh_error_msg("ERROR: Unrecognized triangular element type == " << Utility::enum_to_string(_elem_type));
303  }
304 
305 
306  // If we do have holes and the user asked to GENERATE_CONVEX_HULL,
307  // need to add the p flag so the triangulation respects those segments.
308  if ((_triangulation_type==GENERATE_CONVEX_HULL) && (have_holes))
309  flags << "p";
310 
311  // Finally, add the area constraint
312  if (_desired_area > TOLERANCE)
313  flags << "a" << std::fixed << _desired_area;
314 
315  // add minimum angle constraint
317  flags << "q" << std::fixed << _minimum_angle;
318 
319  if (_regions)
320  flags << "Aa";
321 
322  // add user provided extra flags
323  if (_extra_flags.size() > 0)
324  flags << _extra_flags;
325 
326  // Refine the initial output to conform to the area constraint
327  if (_markers)
328  {
329  // need Voronoi to generate boundary information
330  TriangleWrapper::triangulate(const_cast<char *>(flags.str().c_str()),
331  &initial,
332  &final,
333  &voronoi);
334 
335  // Send the information computed by Triangle to the Mesh.
337  _mesh,
338  _elem_type,
339  &voronoi);
340  }
341  else
342  {
343  TriangleWrapper::triangulate(const_cast<char *>(flags.str().c_str()),
344  &initial,
345  &final,
346  nullptr);
347  // Send the information computed by Triangle to the Mesh.
349  _mesh,
350  _elem_type);
351  }
352 
353 
355 
356 
357  // To the naked eye, a few smoothing iterations usually looks better,
358  // so we do this by default unless the user says not to.
359  if (this->_smooth_after_generating)
360  LaplaceMeshSmoother(_mesh).smooth(2);
361 
362 
363  // Clean up.
367 
368  // Prepare the mesh for use before returning. This ensures (among
369  // other things) that it is partitioned and therefore users can
370  // iterate over local elements, etc.
372 }
void insert_any_extra_boundary_points()
Helper function to add extra points (midpoints of initial segments) to a PSLG triangulation.
static constexpr Real TOLERANCE
const std::vector< int > * _markers
Boundary markers.
void prepare_for_use(const bool skip_renumber_nodes_and_elements, const bool skip_find_neighbors)
Prepare a newly ecreated (or read) mesh for use.
Definition: mesh_base.C:710
Real _minimum_angle
Minimum angle in triangles.
const std::vector< Region * > * _regions
A pointer to a vector of Regions*s.
std::string _extra_flags
Additional flags to be passed to triangle.
UnstructuredMesh & _mesh
Reference to the mesh which is to be created by triangle.
Real _desired_area
The desired area for the elements in the resulting mesh.
void elems_to_segments()
Helper function to create PSLG segments from our other boundary-defining options (1D mesh edges...
ElemType _elem_type
The type of elements to generate.
bool _quiet
Flag which tells if we want to suppress stdout outputs.
void copy_tri_to_mesh(const triangulateio &triangle_data_input, UnstructuredMesh &mesh_output, const ElemType type, const triangulateio *voronoi=nullptr)
Copies triangulation data computed by triangle from a triangulateio object to a LibMesh mesh...
const std::vector< Hole * > * _holes
A pointer to a vector of Hole*s.
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
void set_mesh_dimension(unsigned char d)
Resets the logical dimension of the mesh.
Definition: mesh_base.h:269
First generate a convex hull from the set of points passed in, and then triangulate this set of point...
unsigned int total_hole_points()
Helper function to count points in and verify holes.
std::string enum_to_string(const T e)
Triangulate the interior of a Planar Straight Line Graph, which is defined implicitly by the order of...
std::vector< std::pair< unsigned int, unsigned int > > segments
When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly ...
void nodes_to_segments(dof_id_type max_node_id)
Helper function to create PSLG segments from our node ordering, up to the maximum node id...
bool _smooth_after_generating
Flag which tells whether we should smooth the mesh after it is generated.
void destroy(triangulateio &t, IO_Type)
Frees any memory which has been dynamically allocated by Triangle.
virtual dof_id_type max_node_id() const =0
TriangulationType _triangulation_type
The type of triangulation to perform: choices are: convex hull PSLG.
virtual dof_id_type n_nodes() const =0
uint8_t dof_id_type
Definition: id_types.h:67

◆ triangulation_type()

TriangulationType& libMesh::TriangulatorInterface::triangulation_type ( )
inlineinherited

Sets and/or gets the desired triangulation type.

Definition at line 155 of file triangulator_interface.h.

References libMesh::TriangulatorInterface::_triangulation_type.

Referenced by libMesh::MeshTools::Generation::build_delaunay_square(), MeshTriangulationTest::commonSettings(), and triangulate_domain().

155 {return _triangulation_type;}
TriangulationType _triangulation_type
The type of triangulation to perform: choices are: convex hull PSLG.

◆ verify_holes()

void libMesh::TriangulatorInterface::verify_holes ( const Hole outer_bdy)
protectedinherited

Helper function to check holes for intersections if requested.

Definition at line 251 of file triangulator_interface.C.

References libMesh::TriangulatorInterface::_holes, libMesh::TriangulatorInterface::Hole::contains(), and libMesh::make_range().

Referenced by libMesh::TriangulatorInterface::elems_to_segments(), and libMesh::TriangulatorInterface::nodes_to_segments().

252 {
253  for (const Hole * hole : *_holes)
254  {
255  for (const Hole * hole2 : *_holes)
256  {
257  if (hole == hole2)
258  continue;
259 
260  for (auto i : make_range(hole2->n_points()))
261  if (hole->contains(hole2->point(i)))
262  libmesh_error_msg
263  ("Found point " << hole2->point(i) <<
264  " on one hole boundary and another's interior");
265  }
266 
267  for (auto i : make_range(hole->n_points()))
268  if (!outer_bdy.contains(hole->point(i)))
269  libmesh_error_msg
270  ("Found point " << hole->point(i) <<
271  " on hole boundary but outside outer boundary");
272  }
273 }
const std::vector< Hole * > * _holes
A pointer to a vector of Hole*s.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134

Member Data Documentation

◆ _bdy_ids

std::set<std::size_t> libMesh::TriangulatorInterface::_bdy_ids
protectedinherited

◆ _desired_area

Real libMesh::TriangulatorInterface::_desired_area
protectedinherited

The desired area for the elements in the resulting mesh.

Definition at line 319 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::desired_area(), and triangulate().

◆ _elem_type

ElemType libMesh::TriangulatorInterface::_elem_type
protectedinherited

The type of elements to generate.

(Defaults to TRI3).

Definition at line 314 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::elem_type(), triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

◆ _extra_flags

std::string libMesh::TriangleInterface::_extra_flags
private

Additional flags to be passed to triangle.

Definition at line 72 of file mesh_triangle_interface.h.

Referenced by extra_flags(), and triangulate().

◆ _holes

const std::vector<Hole*>* libMesh::TriangulatorInterface::_holes
protectedinherited

◆ _insert_extra_points

bool libMesh::TriangulatorInterface::_insert_extra_points
protectedinherited

Flag which tells whether or not to insert additional nodes before triangulation.

This can sometimes be used to "de-regularize" the resulting triangulation.

This flag is supported for backwards compatibility; setting _interpolate_boundary_points = 1 is equivalent.

Definition at line 341 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::get_interpolate_boundary_points(), libMesh::TriangulatorInterface::insert_extra_points(), and libMesh::TriangulatorInterface::set_interpolate_boundary_points().

◆ _interpolate_boundary_points

int libMesh::TriangulatorInterface::_interpolate_boundary_points
protectedinherited

Flag which tells how many additional nodes should be inserted between each pair of original mesh points.

Definition at line 347 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::get_interpolate_boundary_points(), and libMesh::TriangulatorInterface::set_interpolate_boundary_points().

◆ _markers

const std::vector<int>* libMesh::TriangulatorInterface::_markers
protectedinherited

◆ _mesh

UnstructuredMesh& libMesh::TriangulatorInterface::_mesh
protectedinherited

◆ _minimum_angle

Real libMesh::TriangulatorInterface::_minimum_angle
protectedinherited

Minimum angle in triangles.

Definition at line 324 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::minimum_angle(), and triangulate().

◆ _quiet

bool libMesh::TriangulatorInterface::_quiet
protectedinherited

Flag which tells if we want to suppress stdout outputs.

Definition at line 358 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::quiet(), and triangulate().

◆ _regions

const std::vector<Region*>* libMesh::TriangulatorInterface::_regions
protectedinherited

A pointer to a vector of Regions*s.

If this is nullptr, there are no regions!

Definition at line 303 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::attach_region_list(), triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

◆ _serializer

MeshSerializer libMesh::TriangleInterface::_serializer
private

Triangle only operates on serial meshes.

Definition at line 77 of file mesh_triangle_interface.h.

◆ _smooth_after_generating

bool libMesh::TriangulatorInterface::_smooth_after_generating
protectedinherited

Flag which tells whether we should smooth the mesh after it is generated.

True by default.

Definition at line 353 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::smooth_after_generating(), triangulate(), and libMesh::Poly2TriTriangulator::triangulate().

◆ _triangulation_type

TriangulationType libMesh::TriangulatorInterface::_triangulation_type
protectedinherited

◆ _verify_hole_boundaries

bool libMesh::TriangulatorInterface::_verify_hole_boundaries
protectedinherited

◆ segments

std::vector<std::pair<unsigned int, unsigned int> > libMesh::TriangulatorInterface::segments
inherited

When constructing a PSLG, if the node numbers do not define the desired boundary segments implicitly through the ordering of the points, you can use the segments vector to specify the segments explicitly, Ex: unit square numbered counter-clockwise starting from origin segments[0] = (0,1) segments[1] = (1,2) segments[2] = (2,3) segments[3] = (3,0) (For the above case you could actually use the implicit ordering!)

Definition at line 229 of file triangulator_interface.h.

Referenced by libMesh::TriangulatorInterface::elems_to_segments(), libMesh::TriangulatorInterface::insert_any_extra_boundary_points(), libMesh::Poly2TriTriangulator::insert_refinement_points(), libMesh::TriangulatorInterface::nodes_to_segments(), MeshTriangulationTest::testTriangulatorSegments(), triangulate(), and libMesh::Poly2TriTriangulator::triangulate_current_points().


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