CGAL 6.1 - 3D Constrained Triangulations
Loading...
Searching...
No Matches
User Manual

Author
Laurent Rineau and Jane Tournois


Constrained Triangulations in 3D

3D triangulations partition space and are useful in many applications. In some cases, it is important to ensure that specific faces, such as those representing the sharp features of an object, appear in the output. When a triangulation exactly respects these constraints, it is called a constrained triangulation. However, it is sometimes only possible to preserve the geometry of the constraints, but not their exact combinatorics. In such cases, additional points, called Steiner points, must be inserted. This process results in a conforming triangulation.

This package implements an algorithm for constructing conforming triangulations of 3D polygonal constraints. Specifically, it requires that these piecewise linear constraints are provided as a piecewise linear complex (PLC). The resulting triangulations are of type Triangulation_3, as described in the chapter 3D Triangulations.

The article by Cohen-Steiner et al. [2] discusses the problem of constructing conforming Delaunay triangulations and proposes an algorithm to address it. Si et al.'s work [7], [8], [4], presents an algorithm for computing conforming constrained Delaunay triangulations in 3D.

Definitions

This section introduces the key concepts necessary to understand and use this package effectively.

Piecewise Linear Complex

A piecewise linear complex (PLC) is the three-dimensional generalization of a planar straight-line graph. It consists of a finite set of vertices, edges, and polygonal faces that satisfy the following properties:

  • The vertices and edges of the PLC form a simplicial complex: two edges may intersect only at a shared vertex.
  • The boundary of each polygonal face in the PLC is an ordered list of vertices from the PLC, forming one closed loop.
  • Each polygonal face must be a simple polygon, i.e., its edges don't intersect, except consecutive edges, which intersect at their common vertex.
  • Each polygonal face must be planar, meaning all its vertices lie on the same plane.
  • Each polygonal face may have one or more holes, each of them also represented by an ordered list of vertices from the PLC, forming a closed loop.
  • If two polygons in the PLC intersect, their intersection is a union of edges and vertices from the PLC. In particular, the interiors of two polygons cannot overlap.

Polygons in a PLC may be non-convex and may have holes.

Figure 50.2 A piecewise linear complex, composed of planar faces connected by edges and vertices.


Conforming Constrained Delaunay Triangulation

The algorithms developed in this package are designed to compute a constrained Delaunay triangulation that contains a given set of polygonal constraints in 3D as a subcomplex.

A triangulation is a Delaunay triangulation if the circumscribing sphere of any simplex in the triangulation contains no vertex in its interior (see chapter 3D Triangulations for more details on Delaunay triangulations).

A constrained Delaunay triangulation of a PLC is a constrained triangulation that is as close as possible to being Delaunay, given that some faces are marked as constrained. More precisely, a triangulation is constrained Delaunay if, for any simplex \(s\) of the triangulation, the interior of its circumscribing sphere contains no vertex of the triangulation that is visible from any point in the interior of the simplex \(s\). Two points are visible if the open line segment joining them does not intersect any polygonal face of the PLC, except for polygons that are coplanar with the segment.

In 3D, constrained triangulations do not always exist. This can be demonstrated using the example of Schönhardt polyhedra [5] (see Figure 50.3), [1]. Shewchuk [6] demonstrated that for any PLC, there exists a refined version of the original PLC that admits a constrained Delaunay triangulation. This refinement is achieved by adding Steiner vertices to the input edges and polygons. The constrained triangulation built on this refined PLC is known as a conforming constrained Delaunay triangulation (CCDT for short). Figure 50.4 illustrates an example of a conforming constrained Delaunay triangulation constructed from a PLC.

Figure 50.3 A Schönhardt polyhedron.


Figure 50.4 Left: PLC (360 vertices); Right: CCDT (2452 vertices).


The algorithm implemented in this package is based on the work of Hang Si et al., who developed particular algorithms for constructing conforming constrained Delaunay triangulations from PLCs. The corresponding implementation takes with floating point numbers as coordinates [7], [8], [4].

Software Design

Representation of Piecewise Linear Complexes

There is no universal or canonical way to represent all possible PLCs in CGAL.

Any polyhedral surface is a PLC, so any model of FaceListGraph, such as CGAL::Surface_mesh, can be used to represent such a PLC. In this representation, the geometric structure of the PLC is directly mapped to the elements of the CGAL::Surface_mesh:

  • vertices of the PLC geometrically correspond to vertices of the surface mesh,
  • edges of the PLC correspond to edges of the surface mesh,
  • and polygonal faces of the PLC correspond to faces of the surface mesh, covering the surface of a geometric object. However, PLCs represented in this way are limited to being manifold (that is, each edge belongs to exactly two faces), and their faces cannot have holes.

A PLC can also be represented as a polygon soup: a collection of vertices and a set of polygons, where each polygon is defined by an ordered list of vertices, without explicit connectivity information between polygons. For a polygon soup to represent a valid PLC, its polygons must satisfy the properties described in the previous section. This approach allows for the representation of non-manifold geometries; however, polygons in a polygon soup cannot have holes.

This package also provides a way to group polygons into distinct surface patches using a property map, named plc_face_id. Each polygon can be assigned a patch identifier, allowing multiple polygons to form a continuous surface patch, which may include holes. Some necessary geometric conditions must be satisfied for these patches to be used in the conforming constrained Delaunay triangulation construction:

  • Each patch must be planar, meaning all polygonal faces in the patch lie on the same plane;
  • The polygonal faces of the patch must not intersect except at their shared edges.

When this property map is provided, the input PLC is interpreted in terms of its polygonal faces, edges and vertices as follows:

  • Each polygonal face of the PLC is defined as the union of input polygons sharing the same patch identifier;
  • The edges of the PLC are those from the surface mesh or polygon soup that satisfy one of the following conditions: – they are adjacent to only one polygonal face; – they are adjacent to two polygonal faces with different patch identifiers; – they are adjacent to more than two polygonal faces with differing patch identifiers, indicating non-manifold features of the PLC.
  • The vertices of the PLC are the ones lying on the boundaries of surface patches in the original surface mesh or polygon soup.

API

This package provides a primary class, CGAL::Conforming_constrained_Delaunay_triangulation_3. This class is templated by a geometric traits class and an underlying triangulation class, allowing for flexibility and customization.

In addition to the main class, the package includes several auxiliary classes that define the types of vertices, cells, and associated metadata used within the triangulation. These supporting classes enable users to extend or adapt the triangulation data structure to their specific needs.

Two overloads of the constructor function CGAL::make_conforming_constrained_Delaunay_triangulation_3() are provided to facilitate the creation of a CGAL::Conforming_constrained_Delaunay_triangulation_3 object from either a surface mesh or a polygon soup.

Traits and Kernel Choice

The requirements for geometric objects and operations are specified by the traits class concept ConformingConstrainedDelaunayTriangulationTraits_3. Any CGAL kernel is a model of this concept. However, because this package builds upon the 3D Triangulation package, it inherits the requirement that the traits class must provide exact predicates.

A key aspect of this algorithm is the creation of new points, known as Steiner points, which are inserted on the segments and polygons of the input PLC. If a traits class with inexact constructions is used, it cannot be guaranteed that these points will lie exactly on the intended segments or polygons. As a result, the output will only approximate the input, with the accuracy limited by the rounding of the computed Steiner points.

Furthermore, when using inexact constructions, the algorithm may fail if the input PLC contains non-adjacent simplices that are too close to each other. In such cases, the triangulation process will emit an error if the distance between simplices falls below an internally computed threshold. An error message describing the involved simplices will be displayed on the standard output. If the issue is caused by poorly shaped triangles, functions such as CGAL::Polygon_mesh_processing::remove_almost_degenerate_faces() may help resolve the problem.

Examples

Build a Conforming Constrained Delaunay Triangulation

The following example illustrates how to use the helper function CGAL::make_conforming_constrained_Delaunay_triangulation_3() to construct a conforming constrained Delaunay triangulation from a given PLC.

The triangulation is saved in the MEDIT file format, using the CGAL::IO::write_MEDIT() function.


File Constrained_triangulation_3/conforming_constrained_Delaunay_triangulation_3.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/polygon_mesh_io.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
#include <CGAL/IO/write_MEDIT.h>
#include <cassert>
int main(int argc, char* argv[])
{
auto filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/mpi.off");
if(!CGAL::IO::read_polygon_mesh(filename, mesh)) {
std::cerr << "Error: cannot read file " << filename << std::endl;
return EXIT_FAILURE;
}
std::cout << "Read " << mesh.number_of_vertices() << " vertices and "
<< mesh.number_of_faces() << " faces" << std::endl;
std::cout << "Number of vertices in the CDT: "
<< ccdt.triangulation().number_of_vertices() << '\n';
std::cout << "Number of constrained facets in the CDT: "
<< ccdt.number_of_constrained_facets() << '\n';
std::ofstream ofs(argc > 2 ? argv[2] : "out.mesh");
ofs.precision(17);
auto tr = std::move(ccdt).triangulation();
// Now `tr` is a valid `CGAL::Triangulation_3` object that can be used for further processing.
// and the triangulation of `ccdt` is empty.
std::cout << "Number of vertices in the triangulation `tr`: "
<< tr.number_of_vertices() << '\n';
std::cout << "Number of vertices in `ccdt`: "
<< ccdt.triangulation().number_of_vertices() << '\n';
assert(ccdt.triangulation().number_of_vertices() == 0);
return EXIT_SUCCESS;
}
size_type number_of_vertices() const
size_type number_of_faces() const
bool read_polygon_mesh(const std::string &fname, Graph &g, const NamedParameters &np=parameters::default_values())
auto make_conforming_constrained_Delaunay_triangulation_3(const PolygonMesh &mesh, const NamedParameters &np=parameters::default_values())
creates a 3D constrained Delaunay triangulation conforming to the faces of a polygon mesh.
Definition: make_conforming_constrained_Delaunay_triangulation_3.h:146
void write_MEDIT(std::ostream &os, const T3 &t3, const NamedParameters &np=parameters::default_values())

Build a Conforming Constrained Delaunay Triangulation from a Polygon Soup

You can also construct a conforming constrained Delaunay triangulation from a polygon soup. The following example demonstrates how to create such a triangulation from a collection of polygons without explicit connectivity information.


File Constrained_triangulation_3/ccdt_3_from_soup.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/polygon_soup_io.h>
#include <CGAL/IO/write_MEDIT.h>
#include <CGAL/draw_constrained_triangulation_3.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
#include <vector>
int main(int argc, char* argv[])
{
auto filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/cubes.off");
std::vector<K::Point_3> points;
std::vector<std::vector<std::size_t>> polygons;
if(!CGAL::IO::read_polygon_soup(filename, points, polygons)) {
std::cerr << "Error: cannot read file " << filename << std::endl;
return EXIT_FAILURE;
}
std::cout << "Read " << points.size() << " vertices and "
<< polygons.size() << " polygons" << std::endl;
std::cout << "Number of vertices in the CDT: "
<< ccdt.triangulation().number_of_vertices() << '\n'
<< "Number of constrained facets in the CDT: "
<< ccdt.number_of_constrained_facets() << '\n';
std::ofstream ofs(argc > 2 ? argv[2] : "out.mesh");
ofs.precision(17);
CGAL::draw(ccdt);
return EXIT_SUCCESS;
}
bool read_polygon_soup(const std::string &fname, PointRange &points, PolygonRange &polygons, const NamedParameters &np=parameters::default_values())
void draw(const T3 &at3, const GSOptions &gso)

Build a Conforming Constrained Delaunay Triangulation with Known Polygon Identifiers

If the user already knows the set of polygonal face identifiers to associate with each PLC face, this information can be provided and preserved throughout the construction of the conforming constrained Delaunay triangulation.

The following example demonstrates how to detect planar surface patches and remesh them as coarsely as possible using CGAL::Polygon_mesh_processing::remesh_planar_patches() from the Polygon Mesh Processing package. The resulting patches and segmentation are then used to build a conforming constrained Delaunay triangulation.

When the named parameter plc_face_id is specified, each constrained facet in the 3D triangulation is assigned to the corresponding input PLC face, as identified by the provided property map. If this parameter is not specified, each input polygonal face is assigned a unique face index.

Figure Figure 50.5 shows the benefit of using the plc_face_id property map. On the last line of the figure, the input PLC is enriched with a segmentation of the planar faces, provided via the plc_face_id property map. In the resulting conforming constrained Delaunay triangulation, only the boundary edges of the PLC faces are constrained, while the other edges never get inserted as edges of the 3D triangulation.

Without the plc_face_id property map, all edges of the PLC faces are constrained, each PLC face is considered as a constraint, possibly resulting in a 3D triangulation with surfaces that are more refined than necessary.


File Constrained_triangulation_3/conforming_constrained_Delaunay_triangulation_3_fimap.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/Polygon_mesh_processing/remesh_planar_patches.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
#include <CGAL/IO/write_MEDIT.h>
#include <CGAL/Polygon_mesh_processing/IO/polygon_mesh_io.h>
using face_descriptor = boost::graph_traits<Mesh>::face_descriptor;
int main(int argc, char* argv[])
{
std::string filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/cross.off");
Mesh input;
std::cerr << "Invalid input." << std::endl;
return 1;
}
std::cout << "Read " << input.number_of_vertices() << " vertices and "
<< input.number_of_faces() << " facets\n";
Mesh mesh;
auto plc_facet_map = get(CGAL::face_patch_id_t<int>(), mesh);
// Remesh planar patches and segment the mesh into planar patches
CGAL::parameters::face_patch_map(plc_facet_map)
.do_not_triangulate_faces(true));
filename = argc > 2 ? argv[2] : "mesh.ply";
CGAL::parameters::stream_precision(17));
std::cout << "Wrote segmented mesh to " << filename << "\n";
// Create a conforming constrained Delaunay triangulation from the mesh
CGAL::parameters::plc_face_id(plc_facet_map));
std::cout << "Number of vertices in the CDT: "
<< ccdt.triangulation().number_of_vertices() << '\n'
<< "Number of constrained facets in the CDT: "
<< ccdt.number_of_constrained_facets() << '\n';
filename = argc > 3 ? argv[3] : "out.mesh";
std::ofstream out(filename);
out.precision(17);
CGAL::IO::write_MEDIT(out, ccdt, CGAL::parameters::with_plc_face_id(true));
std::cout << "Wrote CDT to " << filename << "\n";
return EXIT_SUCCESS;
}
bool read_polygon_mesh(const std::string &fname, PolygonMesh &g, const NamedParameters &np=parameters::default_values())
void remesh_planar_patches(const TriangleMeshIn &tm_in, PolygonMeshOut &pm_out, const NamedParametersIn &np_in=parameters::default_values(), const NamedParametersOut &np_out=parameters::default_values())
bool write_polygon_mesh(const std::string &fname, Graph &g, const NamedParameters &np=parameters::default_values())
Default_named_parameters default_values()

Figure 50.5 shows the input and output of this triangulation construction example.

Build a Conforming Constrained Delaunay Triangulation with Detected Polygon Identifiers

If the user does not know the set of polygonal face identifiers to associate with each PLC face, this information can be automatically detected using the CGAL::Polygon_mesh_processing::region_growing_of_planes_on_faces() function from the Polygon Mesh Processing package.

The following example demonstrates how to detect planar surface patches and build a conforming constrained Delaunay triangulation using the detected segmentation. The named parameter plc_face_id is used to associate each facet of the triangulation with the corresponding input PLC face.


File Constrained_triangulation_3/ccdt_3_fimap_region_growing.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/write_MEDIT.h>
#include <CGAL/Polygon_mesh_processing/IO/polygon_mesh_io.h>
#include <CGAL/Polygon_mesh_processing/bbox.h>
#include <CGAL/Polygon_mesh_processing/region_growing.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
using face_descriptor = boost::graph_traits<Mesh>::face_descriptor;
int main(int argc, char* argv[])
{
std::string filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/cross_quad.off");
Mesh mesh;
if(!PMP::IO::read_polygon_mesh(filename, mesh)) {
std::cerr << "Invalid input." << std::endl;
return 1;
}
std::cout << "Read " << mesh.number_of_vertices() << " vertices and "
<< mesh.number_of_faces() << " facets\n";
auto [face_patch_map, _] =mesh.add_property_map<face_descriptor, std::size_t>("f:patch_id");
const auto bbox = CGAL::Polygon_mesh_processing::bbox(mesh);
const double bbox_max_span = (std::max)({bbox.x_span(), bbox.y_span(), bbox.z_span()});
std::cout << "Merging facets into coplanar patches..." << std::endl;
mesh,
face_patch_map,
CGAL::parameters::maximum_distance(bbox_max_span * 1.e-6)
.maximum_angle(5.));
for(auto f: faces(mesh))
{
// if region growing did not assign a patch id, assign one
if(get(face_patch_map, f) == static_cast<std::size_t>(-1)) {
put(face_patch_map, f, number_of_patches++);
}
}
std::cout << "Number of patches: " << number_of_patches << std::endl;
filename = argc > 2 ? argv[2] : "mesh.ply";
CGAL::parameters::stream_precision(17)
.use_binary_mode(false)
.face_patch_map(face_patch_map));
std::cout << "-- Wrote segmented mesh to \"" << filename << "\"\n";
std::cout << "Creating a conforming constrained Delaunay triangulation...\n";
CGAL::parameters::plc_face_id(face_patch_map));
std::cout << "Number of vertices in the CDT: "
<< ccdt.triangulation().number_of_vertices() << '\n'
<< "Number of constrained facets in the CDT: "
<< ccdt.number_of_constrained_facets() << '\n';
// Write the CDT to a file, with the PLC face ids
filename = argc > 3 ? argv[3] : "out.mesh";
std::ofstream out(filename);
out.precision(17);
CGAL::IO::write_MEDIT(out, ccdt, CGAL::parameters::with_plc_face_id(true));
std::cout << "-- Wrote CDT to \"" << filename << "\"\n";
return EXIT_SUCCESS;
}
std::size_t region_growing_of_planes_on_faces(const PolygonMesh &mesh, RegionMap region_map, const NamedParameters &np=parameters::default_values())
CGAL::Bbox_3 bbox(const PolygonMesh &pmesh, const NamedParameters &np=parameters::default_values())

Preprocessing the Input for Conforming Constrained Delaunay Triangulations

Given a PLC, the algorithms in this package can construct a conforming constrained Delaunay triangulation, provided the input surface can be represented as a valid surface mesh or a collection of surface meshes, and does not contain self-intersections. Several preprocessing functions are available in the Polygon Mesh Processing package to help ensure these preconditions are met.

The following example demonstrates how to construct a conforming constrained Delaunay triangulation from an input mesh that is not triangulated and may contain self-intersections, using autorefinement.


File Constrained_triangulation_3/ccdt_3_after_autorefinement.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/polygon_mesh_io.h>
#include <CGAL/IO/write_MEDIT.h>
#include <CGAL/Polygon_mesh_processing/autorefinement.h>
#include <CGAL/Polygon_mesh_processing/polygon_soup_self_intersections.h>
#include <CGAL/Polygon_mesh_processing/self_intersections.h>
#include <CGAL/Surface_mesh/Surface_mesh.h>
#include <CGAL/draw_constrained_triangulation_3.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
using Point = K::Point_3;
using Surface_mesh = CGAL::Surface_mesh<Point>;
int main(int argc, char* argv[])
{
const auto filename = (argc > 1) ? argv[1]
: CGAL::data_file_path("meshes/spheres_intersecting.off");
if(!CGAL::IO::read_polygon_mesh(filename, mesh)) {
std::cerr << "Error: cannot read file " << filename << std::endl;
return EXIT_FAILURE;
}
std::cout << "Number of facets in " << filename << ": "
<< mesh.number_of_faces() << "\n";
if(PMP::does_self_intersect(mesh))
{
std::cout << "Mesh self-intersects, performing autorefine...\n";
// use a polygon soup as container as the output will most likely be non-manifold
std::vector<Point> points;
std::vector<std::vector<std::size_t>> polygons;
PMP::polygon_mesh_to_polygon_soup(mesh, points, polygons);
PMP::autorefine_triangle_soup(points, polygons);
std::cout << "Number of input triangles after autorefine: "
<< polygons.size() << "\n";
if(PMP::does_polygon_soup_self_intersect(points, polygons))
{
std::cerr << "Error: mesh still self-intersects after autorefine\n";
std::cerr << "Run autorefinement again with an exact kernel ";
std::cerr << "such as CGAL::Exact_predicates_exact_constructions_kernel \n";
return EXIT_FAILURE;
}
}
else
{
}
std::cout << "Number of constrained facets in the CDT: "
<< ccdt.number_of_constrained_facets() << '\n';
std::ofstream ofs(argc > 2 ? argv[2] : "out.mesh");
ofs.precision(17);
CGAL::draw(ccdt);
return EXIT_SUCCESS;
}
This class template represents a 3D conforming constrained Delaunay triangulation.
Definition: Conforming_constrained_Delaunay_triangulation_3.h:548
Triangulation::size_type number_of_constrained_facets() const
returns the number of constrained facets in the triangulation.
Definition: Conforming_constrained_Delaunay_triangulation_3.h:948

The function CGAL::Polygon_mesh_processing::does_self_intersect() is used to detect self-intersections, but it requires the input mesh to be triangulated. Therefore, the input mesh must first be triangulated using CGAL::Polygon_mesh_processing::triangulate_faces() before performing the self-intersection check.

If self-intersections are found, the triangulated mesh is converted into a triangle soup, which is then processed with CGAL::Polygon_mesh_processing::autorefine_triangle_soup() to resolve the self-intersections.

Remeshing a Conforming Constrained Delaunay Triangulation

After constructing the triangulation, you can improve its quality or adapt it to a specific sizing field by applying the CGAL::tetrahedral_isotropic_remeshing() function from the Tetrahedral Remeshing package.

The following example demonstrates how to remesh a conforming constrained Delaunay triangulation.


File Constrained_triangulation_3/remesh_constrained_Delaunay_triangulation_3.cpp

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Conforming_constrained_Delaunay_triangulation_cell_base_3.h>
#include <CGAL/Conforming_constrained_Delaunay_triangulation_vertex_base_3.h>
#include <CGAL/Tetrahedral_remeshing/Remeshing_cell_base_3.h>
#include <CGAL/Tetrahedral_remeshing/Remeshing_vertex_base_3.h>
#include <CGAL/make_conforming_constrained_Delaunay_triangulation_3.h>
#include <CGAL/tetrahedral_remeshing.h>
#include <CGAL/IO/File_medit.h>
#include <CGAL/IO/polygon_mesh_io.h>
#include <CGAL/IO/write_MEDIT.h>
#include <CGAL/property_map.h>
#include <unordered_set>
#include <fstream>
#include <string>
// Triangulation for Remeshing
using CCDT_Tr = CCDT::Triangulation;
using Vertex_handle = Triangulation_3::Vertex_handle;
using Vertex_pair = std::pair<Vertex_handle, Vertex_handle>;
using Constraints_set = std::unordered_set<Vertex_pair, boost::hash<Vertex_pair>>;
int main(int argc, char* argv[])
{
std::string filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/mpi.off");
if(!CGAL::IO::read_polygon_mesh(filename, mesh))
{
std::cerr << "Error: cannot read file " << filename << std::endl;
return EXIT_FAILURE;
}
CCDT ccdt = CGAL::make_conforming_constrained_Delaunay_triangulation_3<CCDT>(mesh);
std::ofstream out("ccdt.mesh");
out.close();
Constraints_set constraints;
Constraints_pmap constraints_pmap(constraints);
namespace np = CGAL::parameters;
namespace Tet_remesh = CGAL::Tetrahedral_remeshing;
Tr tr = Tet_remesh::get_remeshing_triangulation(std::move(ccdt),
np::edge_is_constrained_map(constraints_pmap));
std::cout << "Number of vertices in tr: " << tr.number_of_vertices() << std::endl;
1., // target edge length
np::number_of_iterations(3)
.edge_is_constrained_map(constraints_pmap));
std::cout << "Number of vertices in tr: "
<< tr.number_of_vertices() << std::endl;
std::ofstream ofs("tr.mesh");
return EXIT_SUCCESS;
}
Cell base class for the 3D conforming constrained Delaunay triangulation.
Definition: Conforming_constrained_Delaunay_triangulation_cell_base_3.h:42
Vertex base class for the 3D conforming constrained Delaunay triangulation.
Definition: Conforming_constrained_Delaunay_triangulation_vertex_base_3.h:43
Triangulation_data_structure::Vertex_handle Vertex_handle
void tetrahedral_isotropic_remeshing(CGAL::Triangulation_3< Traits, TDS, SLDS > &tr, const SizingFunction &sizing, const NamedParameters &np=parameters::default_values())
Definition: Conforming_constrained_Delaunay_triangulation_3.h:3995

Figures

The following table of figures (Figure 50.5) illustrates some results of the examples provided in this package. The left column shows the input PLC, while the right column displays the resulting conforming constrained Delaunay triangulation.

From top to bottom, the lines show different input PLC, from the same input triangulated surface and, for each of them, the resulting conforming constrained Delaunay triangulation. The input data are:

On the fourth line, the input PLC is remeshed using CGAL::Polygon_mesh_processing::isotropic_remeshing(). The resulting conforming constrained Delaunay triangulation contains fewer vertices than the input remeshed and segmented input PLC. This reduction occurs because only the boundary edges of the PLC faces are marked as constraints in the triangulation; interior edges that do not lie on the boundaries of surface patches (as defined by plc_face_id) are ignored. As a result, these non-boundary edges are omitted from the triangulation, leading to a coarser triangulation.

Figure 50.5 A collection of conforming constrained Delaunay triangulations built from different inputs. The left column shows the input PLC, while the right column displays the resulting 3D triangulation.


Implementation History

The initial version of this package was implemented by Laurent Rineau and released in CGAL 6.1 (2025). Jane Tournois contributed to the documentation and helped improve the API. The package design and algorithms are grounded in the theoretical work of Hang Si et al. on meshing algorithms [7][4].