GCC Code Coverage Report


Directory: ./
File: lib/geogram/mesh/mesh_tetrahedralize.cpp
Date: 2026-09-07 02:36:43
Exec Total Coverage
Lines: 0 98 0.0%
Functions: 0 1 0.0%
Branches: 0 244 0.0%

Line Branch Exec Source
1 /*
2 * Copyright (c) 2000-2022 Inria
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * * Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 * * Redistributions in binary form must reproduce the above copyright notice,
11 * this list of conditions and the following disclaimer in the documentation
12 * and/or other materials provided with the distribution.
13 * * Neither the name of the ALICE Project-Team nor the names of its
14 * contributors may be used to endorse or promote products derived from this
15 * software without specific prior written permission.
16 *
17 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 * POSSIBILITY OF SUCH DAMAGE.
28 *
29 * Contact: Bruno Levy
30 *
31 * https://www.inria.fr/fr/bruno-levy
32 *
33 * Inria,
34 * Domaine de Voluceau,
35 * 78150 Le Chesnay - Rocquencourt
36 * FRANCE
37 *
38 */
39
40 #include <geogram/mesh/mesh_tetrahedralize.h>
41 #include <geogram/mesh/mesh_surface_intersection.h>
42 #include <geogram/mesh/mesh_repair.h>
43 #include <geogram/mesh/mesh_geometry.h>
44 #include <geogram/mesh/mesh_fill_holes.h>
45 #include <geogram/mesh/mesh_io.h>
46 #include <geogram/mesh/mesh.h>
47 #include <geogram/delaunay/delaunay.h>
48 #include <geogram/delaunay/delaunay_tetgen.h>
49 #include <geogram/basic/logger.h>
50 #include <geogram/basic/command_line.h>
51
52 namespace GEO {
53
54 bool mesh_tetrahedralize(
55 Mesh& M, const MeshTetrahedralizeParameters& parameters
56 ) {
57
58 bool verbose = parameters.verbose;
59 bool preprocess = parameters.preprocess;
60 double epsilon = parameters.preprocess_merge_vertices_epsilon;
61 double max_hole_area = parameters.preprocess_fill_hole_max_area;
62 bool refine = parameters.refine;
63 double quality = parameters.refine_quality;
64 bool keep_regions = parameters.keep_regions;
65
66 if(!DelaunayFactory::has_creator("tetgen")) {
67 Logger::err("TetMeshing")
68 << "Not supported in this version" << std::endl;
69 Logger::err("TetMeshing")
70 << "(need to recompile with tetgen support)"
71 << std::endl;
72 return false;
73 }
74
75 if(!M.facets.are_simplices()) {
76 Logger::warn("TetMeshing")
77 << "Mesh is not triangulated (triangulating it)"
78 << std::endl;
79 tessellate_facets(M,3);
80 }
81
82 // in percent of bbox diagonal
83 epsilon *= (0.01 * bbox_diagonal(M));
84
85 max_hole_area *= (0.01 * Geom::mesh_area(M));
86
87 bool ok = true;
88 Delaunay_var delaunay;
89
90 for(index_t iter=0; iter<5; ++iter) {
91 if(iter != 0 && verbose) {
92 Logger::warn("Tetrahedralize")
93 << "Retrying, because tetgen may have moved some vertices"
94 << std::endl;
95 }
96
97 if(preprocess) {
98
99 // Snap to floating-point coords
100 if(epsilon > 0) {
101 index_t N = M.vertices.nb() * M.vertices.dimension();
102 double* coords = M.vertices.point_ptr(0);
103 for(index_t i=0; i<N; ++i) {
104 coords[i] = double(float(coords[i]));
105 }
106 }
107
108 mesh_repair(M, MESH_REPAIR_DEFAULT, epsilon);
109 fill_holes(M, max_hole_area);
110 MeshSurfaceIntersection intersection(M);
111 intersection.set_verbose(verbose);
112 intersection.intersect();
113 intersection.remove_internal_shells();
114 if(iter == 0 && parameters.preprocess_merge_coplanar_facets) {
115 intersection.simplify_coplanar_facets();
116 }
117 mesh_repair(M, MESH_REPAIR_DEFAULT, epsilon);
118 }
119
120 if(verbose) {
121 Logger::out("TetMeshing") << "Tetrahedralizing..." << std::endl;
122 }
123
124 delaunay = Delaunay::create(3,"tetgen");
125 delaunay->set_refine(refine);
126 delaunay->set_quality(quality);
127 delaunay->set_constraints(&M);
128 delaunay->set_keep_regions(keep_regions);
129
130 try {
131 ok = true;
132 delaunay->set_vertices(0,nullptr); // No additional vertex
133 ok = ok &&
134 delaunay->nb_vertices() != 0 &&
135 delaunay->nb_cells() != 0;
136 if(ok) {
137 break;
138 }
139 } catch(const Delaunay::InvalidInput& error_report) {
140 geo_argused(error_report);
141 if(verbose) {
142 Logger::warn("Tetrahedralize") << "Encountered error"
143 << std::endl;
144 }
145 ok = false;
146 if(!preprocess) {
147 break;
148 }
149 }
150 }
151
152 if(!ok) {
153 Logger::err("Tetrahedralize") << "failed" << std::endl;
154 return false;
155 }
156
157 vector<double> pts(delaunay->nb_vertices() * 3);
158 vector<index_t> tet2v(delaunay->nb_cells() * 4);
159 for(index_t v = 0; v < delaunay->nb_vertices(); ++v) {
160 pts[3 * v] = delaunay->vertex_ptr(v)[0];
161 pts[3 * v + 1] = delaunay->vertex_ptr(v)[1];
162 pts[3 * v + 2] = delaunay->vertex_ptr(v)[2];
163 }
164 for(index_t t = 0; t < delaunay->nb_cells(); ++t) {
165 tet2v[4 * t] = index_t(delaunay->cell_vertex(t, 0));
166 tet2v[4 * t + 1] = index_t(delaunay->cell_vertex(t, 1));
167 tet2v[4 * t + 2] = index_t(delaunay->cell_vertex(t, 2));
168 tet2v[4 * t + 3] = index_t(delaunay->cell_vertex(t, 3));
169 }
170
171 if(pts.size() == 0 || tet2v.size() == 0) {
172 Logger::err("Tetrahedralize")
173 << "Did not generate any tetrahedron"
174 << std::endl;
175 return false;
176 }
177
178 M.cells.assign_tet_mesh(3, pts, tet2v, true);
179
180 if(keep_regions) {
181 Attribute<index_t> region(M.cells.attributes(), "region");
182 for(index_t t: M.cells) {
183 region[t] = delaunay->region(t);
184 }
185 }
186
187 M.cells.connect();
188 if(verbose) {
189 M.show_stats("TetMeshing");
190 }
191 return true;
192 }
193
194 }
195