Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 : #include "libmesh/libmesh_config.h"
19 : #ifdef LIBMESH_HAVE_TETGEN
20 :
21 : // Local includes
22 : #include "libmesh/libmesh_common.h"
23 : #include "libmesh/mesh_tetgen_wrapper.h"
24 :
25 : // C++ includes
26 : #include <iostream>
27 : #include <memory>
28 :
29 : namespace libMesh
30 : {
31 :
32 1092 : TetGenWrapper::TetGenWrapper() :
33 1638 : tetgen_output(std::make_unique<tetgenio>())
34 : {
35 1092 : this->tetgen_data.mesh_dim = 3;
36 1092 : this->tetgen_data.numberofpointattributes = 0;
37 1092 : this->tetgen_data.firstnumber = 0;
38 1092 : }
39 :
40 :
41 :
42 : TetGenWrapper::~TetGenWrapper() = default;
43 :
44 :
45 :
46 6376 : void TetGenWrapper::set_node(unsigned i, REAL x, REAL y, REAL z)
47 : {
48 6376 : unsigned index = i*3;
49 6376 : tetgen_data.pointlist[index++] = x;
50 6376 : tetgen_data.pointlist[index++] = y;
51 6376 : tetgen_data.pointlist[index++] = z;
52 6376 : }
53 :
54 :
55 :
56 4 : void TetGenWrapper::set_hole(unsigned i, REAL x, REAL y, REAL z)
57 : {
58 4 : unsigned index = i*3;
59 4 : tetgen_data.holelist[index++] = x;
60 4 : tetgen_data.holelist[index++] = y;
61 4 : tetgen_data.holelist[index++] = z;
62 4 : }
63 :
64 :
65 :
66 1092 : void TetGenWrapper::set_numberofpoints(int i)
67 : {
68 : // This is an int in tetgen, so use an int here even though it should be unsigned
69 1092 : tetgen_data.numberofpoints = i;
70 1092 : }
71 :
72 :
73 :
74 2928 : void TetGenWrapper::get_output_node(unsigned i, REAL & x, REAL & y, REAL & z)
75 : {
76 : // Bounds checking...
77 2928 : libmesh_error_msg_if(i >= static_cast<unsigned>(tetgen_output->numberofpoints),
78 : "Error, requested point "
79 : << i
80 : << ", but there are only "
81 : << tetgen_output->numberofpoints
82 : << " points available.");
83 :
84 2928 : x = tetgen_output->pointlist[3*i];
85 2928 : y = tetgen_output->pointlist[3*i+1];
86 2928 : z = tetgen_output->pointlist[3*i+2];
87 2928 : }
88 :
89 :
90 :
91 1084 : int TetGenWrapper::get_numberoftetrahedra()
92 : {
93 1084 : return tetgen_output->numberoftetrahedra;
94 : }
95 :
96 :
97 :
98 8 : int TetGenWrapper::get_numberoftrifaces()
99 : {
100 8 : return tetgen_output->numberoftrifaces;
101 : }
102 :
103 :
104 :
105 4 : int TetGenWrapper::get_numberofpoints()
106 : {
107 4 : return tetgen_output->numberofpoints;
108 : }
109 :
110 :
111 :
112 53800 : int TetGenWrapper::get_element_node(unsigned i, unsigned j)
113 : {
114 53800 : return tetgen_output->tetrahedronlist[i*4+j];
115 : }
116 :
117 :
118 :
119 288 : int TetGenWrapper::get_triface_node(unsigned i, unsigned j)
120 : {
121 288 : return tetgen_output->trifacelist[i*3+j];
122 : }
123 :
124 :
125 :
126 0 : REAL TetGenWrapper::get_element_attribute(unsigned i)
127 : {
128 0 : libmesh_assert(tetgen_output->numberoftetrahedronattributes>0);
129 0 : return tetgen_output->tetrahedronattributelist[tetgen_output->numberoftetrahedronattributes*i];
130 : }
131 :
132 :
133 :
134 1092 : void TetGenWrapper::allocate_pointlist(int numofpoints)
135 : {
136 : // This is stored as an int in tetgen, so we store it that way as well.
137 1092 : this->set_numberofpoints(numofpoints);
138 :
139 : // Don't try to allocate an array of size zero, this is not portable...
140 1092 : if (this->tetgen_data.numberofpoints > 0)
141 : {
142 : // Is there previously-allocated memory here?
143 1092 : libmesh_error_msg_if(this->tetgen_data.pointlist != nullptr,
144 : "Cannot allocate on top of previously allocated memory!");
145 :
146 : // We allocate memory here, the tetgenio destructor will delete it.
147 1092 : this->tetgen_data.pointlist = new REAL[this->tetgen_data.numberofpoints * 3];
148 : }
149 1092 : }
150 :
151 :
152 :
153 1092 : void TetGenWrapper::set_switches(std::string_view s)
154 : {
155 : // A temporary buffer for passing to the C API, it requires
156 : // a char *, not a const char *...
157 : char buffer[256];
158 :
159 : // Make sure char buffer has enough room
160 1092 : libmesh_error_msg_if(s.size() >= sizeof(buffer)-1,
161 : "Fixed size buffer of length "
162 : << sizeof(buffer)
163 : << " not large enough to hold TetGen switches.");
164 :
165 : // Copy the string, don't forget to terminate!
166 1092 : buffer[ s.copy( buffer , sizeof( buffer ) - 1 ) ] = '\0' ;
167 :
168 1092 : if (!tetgen_be.parse_commandline(buffer))
169 0 : libMesh::out << "TetGen replies: Wrong switches!" << std::endl;
170 1092 : }
171 :
172 :
173 :
174 1092 : void TetGenWrapper::run_tetgen()
175 : {
176 : // Call tetrahedralize from the TetGen library.
177 1092 : tetrahedralize(&tetgen_be, &tetgen_data, tetgen_output.get());
178 1092 : }
179 :
180 :
181 :
182 4 : void TetGenWrapper::set_numberoffacets(int i)
183 : {
184 : // This is stored as an int in TetGen
185 4 : this->tetgen_data.numberoffacets = i;
186 4 : }
187 :
188 :
189 :
190 4 : void TetGenWrapper::set_numberofholes(int i)
191 : {
192 : // This is stored as an int in TetGen
193 4 : this->tetgen_data.numberofholes = i;
194 4 : }
195 :
196 :
197 :
198 0 : void TetGenWrapper::set_numberofregions(int i)
199 : {
200 : // This is stored as an int in TetGen
201 0 : this->tetgen_data.numberofregions = i;
202 0 : }
203 :
204 :
205 :
206 4 : void TetGenWrapper::allocate_facetlist(int numoffacets, int numofholes)
207 : {
208 : // These are both stored as ints in TetGen
209 4 : this->set_numberoffacets(numoffacets);
210 4 : this->set_numberofholes(numofholes);
211 :
212 : // Don't try to allocate an array of size zero, this is not portable...
213 4 : if (this->tetgen_data.numberoffacets > 0)
214 : {
215 : // Is there previously-allocated memory here?
216 4 : libmesh_error_msg_if(this->tetgen_data.facetlist != nullptr,
217 : "Cannot allocate on top of previously allocated memory!");
218 :
219 : // We allocate memory here, the tetgenio destructor cleans it up.
220 4 : this->tetgen_data.facetlist = new tetgenio::facet[this->tetgen_data.numberoffacets];
221 :
222 100 : for (int i=0; i<numoffacets; i++)
223 96 : this->tetgen_data.init(&(this->tetgen_data.facetlist[i]));
224 : }
225 :
226 :
227 : // Don't try to allocate an array of size zero, this is not portable...
228 4 : if (this->tetgen_data.numberofholes > 0)
229 : {
230 : // Is there previously-allocated memory here?
231 4 : libmesh_error_msg_if(this->tetgen_data.holelist != nullptr,
232 : "Cannot allocate on top of previously allocated memory!");
233 :
234 4 : this->tetgen_data.holelist = new REAL[this->tetgen_data.numberofholes * 3];
235 : }
236 4 : }
237 :
238 :
239 :
240 0 : void TetGenWrapper::allocate_regionlist(int numofregions)
241 : {
242 0 : this->set_numberofregions(numofregions);
243 :
244 : // Don't try to allocate an array of size zero, this is not portable...
245 0 : if (this->tetgen_data.numberofregions > 0)
246 : {
247 : // Is there previously-allocated memory here?
248 0 : libmesh_error_msg_if(this->tetgen_data.regionlist != nullptr,
249 : "Cannot allocate on top of previously allocated memory!");
250 :
251 : // We allocate memory here, the tetgenio destructor cleans it up.
252 0 : this->tetgen_data.regionlist = new REAL[this->tetgen_data.numberofregions * 5];
253 : }
254 0 : }
255 :
256 :
257 :
258 96 : void TetGenWrapper::set_facet_numberofpolygons(unsigned i, int num)
259 : {
260 : // numberofpolygons is stored as an int in TetGen
261 96 : this->tetgen_data.facetlist[i].numberofpolygons = num;
262 96 : }
263 :
264 :
265 :
266 96 : void TetGenWrapper::set_facet_numberofholes(unsigned i, int num)
267 : {
268 : // numberofholes is stored as an int in TetGen
269 96 : this->tetgen_data.facetlist[i].numberofholes = num;
270 96 : }
271 :
272 :
273 :
274 :
275 96 : void TetGenWrapper::allocate_facet_polygonlist(unsigned i, int numofpolygons)
276 : {
277 96 : this->set_facet_numberofpolygons(i, numofpolygons);
278 96 : this->set_facet_numberofholes(i, 0);
279 :
280 : // Don't try to create an array of size zero, this isn't portable
281 96 : if (numofpolygons > 0)
282 : {
283 : // Is there previously-allocated memory here?
284 96 : libmesh_error_msg_if(this->tetgen_data.facetlist[i].polygonlist != nullptr,
285 : "Cannot allocate on top of previously allocated memory!");
286 :
287 : // We allocate memory here, the tetgenio destructor cleans it up.
288 96 : this->tetgen_data.facetlist[i].polygonlist = new tetgenio::polygon[numofpolygons];
289 :
290 192 : for (int j=0; j<this->tetgen_data.facetlist[i].numberofpolygons; j++)
291 96 : this->tetgen_data.init(&(this->tetgen_data.facetlist[i].polygonlist[j]));
292 : }
293 96 : }
294 :
295 :
296 :
297 96 : void TetGenWrapper::set_polygon_numberofvertices(unsigned i, unsigned j, int num)
298 : {
299 : // numberofvertices is stored as an int in TetGen
300 96 : this->tetgen_data.facetlist[i].polygonlist[j].numberofvertices = num;
301 96 : }
302 :
303 :
304 :
305 96 : void TetGenWrapper::allocate_polygon_vertexlist(unsigned i, unsigned j, int numofvertices)
306 : {
307 96 : this->set_polygon_numberofvertices(i, j, numofvertices);
308 :
309 : // Don't try to create an array of size zero, this isn't portable
310 96 : if (numofvertices > 0)
311 : {
312 : // Is there previously-allocated memory here?
313 96 : libmesh_error_msg_if(this->tetgen_data.facetlist[i].polygonlist[j].vertexlist != nullptr,
314 : "Cannot allocate on top of previously allocated memory!");
315 :
316 : // We allocate memory here, the tetgenio destructor cleans it up.
317 96 : this->tetgen_data.facetlist[i].polygonlist[j].vertexlist = new int[numofvertices];
318 : }
319 96 : }
320 :
321 :
322 :
323 :
324 288 : void TetGenWrapper::set_vertex(unsigned i, unsigned j, unsigned k, int nodeindex)
325 : {
326 : // vertexlist entries are stored as ints in TetGen
327 288 : this->tetgen_data.facetlist[i].polygonlist[j].vertexlist[k] = nodeindex;
328 288 : }
329 :
330 :
331 :
332 0 : void TetGenWrapper::set_region(unsigned i, REAL x, REAL y, REAL z,
333 : REAL attribute, REAL vol_constraint)
334 : {
335 0 : unsigned index = i*5;
336 0 : tetgen_data.regionlist[index++] = x;
337 0 : tetgen_data.regionlist[index++] = y;
338 0 : tetgen_data.regionlist[index++] = z;
339 0 : tetgen_data.regionlist[index++] = attribute;
340 0 : tetgen_data.regionlist[index++] = vol_constraint;
341 0 : }
342 :
343 : } // namespace libMesh
344 :
345 :
346 : #endif // LIBMESH_HAVE_TETGEN
|
