CGAL::Delaunay_mesher_2<CDT, Criteria>

This class implements a 2D mesh generator.

#include <CGAL/Delaunay_mesher_2.h>

Parameters

The template parameter CDT should be a model of the concept ConstrainedDelaunayTriangulation_2, and type CDT::Face should be a model of the concept MeshFaceBase_2.

The geometric traits class of the instance of CDT has to be a model of the concept DelaunayMeshTraits_2.

The template parameter Criteria should be a model of the concept MeshingCriteria_2. This traits class defines the shape and size criteria for the triangles of the mesh. Criteria::Face_handle has to be the same as CDT::Face_handle.

Using this class

The constructor of the class Delaunay_mesher_2<CDT, Criteria> takes a reference to a CDT as an argument. A call to the refinement method refine_mesh() will refine the constrained Delaunay triangulation into a mesh satisfying the size and shape criteria specified in the traits class. Note that if, during the life time of the Delaunay_mesher_2<CDT, Criteria> object, the triangulation is externally modified, any further call to its member methods may crash. Consider constructing a new Delaunay_mesher_2<CDT, Criteria> object if the triangulation has been modified.

Meshing domain

The domain to be mesh is defined by the constrained edges and a set of seed points. The constrained edges divides the plane into several connected components. The mesh domain is either the union of the bounded connected components including at least one seed, or the union of the bounded connected components that do no contain any seed. Note that the unbounded component of the plane is never meshed.

Types

typedef CDT::Geom_traits Geom_traits; the geometric traits class.

Delaunay_mesher_2<CDT, Criteria>::Seeds_iterator
const iterator over defined seeds. Its value type is Geom_traits::Point_2.

Creation

Delaunay_mesher_2<CDT, Criteria> mesher ( CDT& t, Criteria criteria = Criteria());
Create a new mesher, working on t, with meshing criteria criteria.

Seeds functions

The following functions are used to define seeds.

void mesher.clear_seeds () Sets seeds to the empty set. All finite connected components of the constrained triangulation will be refined.

template<class InputIterator>
void mesher.set_seeds ( InputIterator begin, InputIterator end, const bool mark=false)
Sets seeds to the sequence [begin, end]. If mark=true, the mesh domain is the union of the bounded connected components including at least one seed. If mark=false, the domain is the union of the bounded components including no seed. Note that the unbounded component of the plane is never meshed.
Requirement: The value_type of begin and end is Geom_traits::Point_2.

Seeds_const_iterator mesher.seeds_begin () const Start of the seeds sequence.
Seeds_const_iterator mesher.seeds_end () const Past the end of the seeds sequence.

Meshing methods

void mesher.refine_mesh () Refines the constrained Delaunay triangulation into a mesh satisfying the criteria defined by the traits.

Criteria mesher.get_criteria () Returns a const reference to the criteria traits object.

void mesher.set_criteria ( Criteria criteria)
Assigns criteria to the criteria traits object.

The function set_criteria scans all faces to recalculate the list of bad faces, that are faces not conforming to the meshing criteria. This function actually has an optional argument that permits to prevent this recalculation. The filling of the list of bad faces can then be done by a call to set_bad_faces.

void mesher.set_criteria ( Criteria criteria, bool recalculate_bad_faces)
Assigns criteria to the criteria traits object. If recalculate_bad_faces is false, the list of bad faces is let empty and the function set_bad_faces should be called before refine_mesh.

template <class InputIterator>
void mesher.set_bad_faces ( InputIterator begin, InputIterator end)
This method permits to set the list of bad triangles directly, from the sequence [begin, end], so that the algorithm will not scan the whole set of triangles to find bad ones. To use if there is a non-naive way to find bad triangles.
Requirement: The value_type of begin and end is Face_handle.

Step by step operations

The Delaunay_mesher_2<CDT, Criteria> class allows, for debugging or demos, to play the meshing algorithm step by step, using the following methods.

void mesher.init () This method must be called just before the first call to the following step by step refinement method, that is when all vertices and constrained edges have been inserted into the constrained Delaunay triangulation. It must be called again before any subsequent calls to the step by step refinement method if new vertices or constrained edges have been inserted since the last call.

bool mesher.is_refinement_done () Tests if the step by step refinement algorithm is done. If it returns true, the following calls to step_by_step_refine_mesh will not insert any points, until some new constrained segments or points are inserted in the triangulation and init is called again.

bool mesher.step_by_step_refine_mesh ()
Applies one step of the algorithm, by inserting one point, if the algorithm is not done. Returns false iff no point has been inserted because the algorithm is done.