GCC Code Coverage Report

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