00001 00002 // // 00003 // TetGen // 00004 // // 00005 // A Quality Tetrahedral Mesh Generator and 3D Delaunay Triangulator // 00006 // // 00007 // Version 1.4 // 00008 // January 14, 2006 // 00009 // // 00010 // Copyright 2002, 2004, 2005, 2006 // 00011 // Hang Si // 00012 // Rathausstr. 9, 10178 Berlin, Germany // 00013 // si@wias-berlin.de // 00014 // // 00015 // You can obtain TetGen via internet: http://tetgen.berlios.de. It may be // 00016 // freely copied, modified, and redistributed under the copyright notices // 00017 // given in the file LICENSE. // 00018 // // 00019 // TetGen computes Delaunay tetrahedralizations, constrained Delaunay tetra- // 00020 // hedralizations, and quality Delaunay tetrahedral meshes. The latter are // 00021 // nicely graded and whose tetrahedra have radius-edge ratio bounded. Such // 00022 // meshes are suitable for finite element and finite volume methods. // 00023 // // 00024 // TetGen incorporates a suit of geometrical and mesh generation algorithms. // 00025 // A brief description of algorithms used in TetGen is found in the first // 00026 // section of the user's manual. References are given for users who are // 00027 // interesting in these approaches. The main references are given below: // 00028 // // 00029 // The efficient Delaunay tetrahedralization algorithm is: H. Edelsbrunner // 00030 // and N. R. Shah, "Incremental Topological Flipping Works for Regular // 00031 // Triangulations". Algorithmica 15: 223-241, 1996. // 00032 // // 00033 // The constrained Delaunay tetrahedralization algorithm is described in: // 00034 // H. Si and K. Gaertner, "Meshing Piecewise Linear Complexs by Constrain- // 00035 // ed Delaunay Tetrahedralizations". Proceedings of the 14th International // 00036 // Mesh Roundtable. September 2005. // 00037 // // 00038 // The Delaunay refinement algorithm is from: J. R. Shewchuk, "Tetrahedral // 00039 // Mesh Generation by Delaunay Refinement". Proceedings of the 14th Annual // 00040 // Symposium on Computational Geometry, pages 86-95, 1998. // 00041 // // 00042 // The mesh data structure of TetGen is a combination of two types of mesh // 00043 // data structures. The tetrahedron-based mesh data structure introduced // 00044 // by Shewchuk is eligible to implement tetrahedralization algorithms. The // 00045 // triangle-edge data structure provided by Muecke is used for represent- // 00046 // ing boundary elements: subfaces and subsegments. These data types are // 00047 // handled through a set of fast mesh manipulation primitives. // 00048 // // 00049 // J. R. Shewchuk, "Delaunay Refinement Mesh Generation". PhD thesis, // 00050 // Carnegie Mellon University, 1997. // 00051 // // 00052 // E. P. Muecke, "Shapes and Implementations in Three-Dimensional // 00053 // Geometry". PhD thesis, Univ. of Illinois, Urbana, Illinois, 1993. // 00054 // // 00055 // The research of mesh generation is definitly on the move. Many State-of- // 00056 // the-art algorithms need implementing and evaluating. I heartily welcome // 00057 // any new algorithm especially for generating quality conforming Delaunay // 00058 // meshes and anisotropic conforming Delaunay meshes. // 00059 // // 00060 // TetGen is supported by the "pdelib" project of Weierstrass Institute for // 00061 // Applied Analysis and Stochastics (WIAS) in Berlin. It is a collection // 00062 // of software components for solving non-linear partial differential // 00063 // equations including 2D and 3D mesh generators, sparse matrix solvers, // 00064 // and scientific visualization tools, etc. For more information please // 00065 // see: http://www.wias-berlin.de/software/pdelib. // 00066 // // 00068 00070 // // 00071 // tetgen.h // 00072 // // 00073 // Header file of the TetGen library. Also is the user-level header file. // 00074 // // 00076 00078 // // 00079 // TetGen Library Overview // 00080 // // 00081 // TetGen library is comprised by several data types and global functions. // 00082 // // 00083 // There are three main data types: tetgenio, tetgenbehavior, and tetgenmesh.// 00084 // Tetgenio is used to pass data into and out of TetGen library; tetgenbeha- // 00085 // vior keeps the runtime options and thus controls the behaviors of TetGen; // 00086 // tetgenmesh, the biggest data type I've ever defined, contains mesh data // 00087 // structures and mesh traversing and transformation operators. These data // 00088 // types are defined as C++ classes. // 00089 // // 00090 // There are few global functions. tetrahedralize() is provided for calling // 00091 // TetGen from another program. Two functions: orient3d() and insphere() are // 00092 // incorporated from a public C code provided by Shewchuk. They performing // 00093 // exact geometrical tests. // 00094 // // 00096 00097 #ifndef tetgenH 00098 #define tetgenH 00099 00100 // To compile TetGen as a library instead of an executable program, define 00101 // the TETLIBRARY symbol. 00102 00103 // #define TETLIBRARY 00104 00105 // Uncomment the following line to disable assert macros. These macros are 00106 // inserted in places where I hope to catch bugs. 00107 00108 // #define NDEBUG 00109 00110 // To insert lots of self-checks for internal errors, define the SELF_CHECK 00111 // symbol. This will slow down the program significantly. 00112 00113 // #define SELF_CHECK 00114 00115 // For single precision ( which will save some memory and reduce paging ), 00116 // define the symbol SINGLE by using the -DSINGLE compiler switch or by 00117 // writing "#define SINGLE" below. 00118 // 00119 // For double precision ( which will allow you to refine meshes to a smaller 00120 // edge length), leave SINGLE undefined. 00121 00122 // #define SINGLE 00123 00124 #ifdef SINGLE 00125 #define REAL float 00126 #else 00127 #define REAL double 00128 #endif // not defined SINGLE 00129 00130 // Here are the most general used head files for C/C++ programs. 00131 00132 #include <stdio.h> // Standard IO: FILE, NULL, EOF, printf(), ... 00133 #include <stdlib.h> // Standard lib: abort(), system(), getenv(), ... 00134 #include <string.h> // String lib: strcpy(), strcat(), strcmp(), ... 00135 #include <math.h> // Math lib: sin(), sqrt(), pow(), ... 00136 #include <assert.h> 00137 00139 // // 00140 // The tetgenio data type // 00141 // // 00142 // Used to pass data into and out of the library of TetGen. // 00143 // // 00144 // If you want to program with the library of TetGen, it's necessary for you // 00145 // to understand the tetgenio data type, while the other two data types can // 00146 // be hidden through calling the global function "tetrahedralize()". As you // 00147 // will see, that tetgenio is just a collection of arrays to storing points, // 00148 // tetrahedra, (triangular) faces, boundary markers, and so forth. They are // 00149 // used to store data from input & output files of TetGen. If you understand // 00150 // TetGen's file formats, then it is easy to understand these arrays. // 00151 // // 00152 // Once an object of tetgenio is declared (or created), all arrays of it are // 00153 // automatically initialized to NULLs (by routine initialize()). Before they // 00154 // can be used, one has to first allocate enough memory for them, i.e., use // 00155 // either 'malloc()' or 'new' operator. On deletion of the object, one needs // 00156 // to free the memory occupied by these arrays. Routine deinitialize() will // 00157 // be automatically called. It will deallocate the memory for an array if it // 00158 // is not a NULL. However, it assumes that the memory is allocated by 'new' // 00159 // (C++ operator). If you use malloc(), you should free() them and set the // 00160 // pointers to NULLs before reaching deinitialize(). // 00161 // // 00162 // In all cases, the first item in any array is stored starting at index [0].// 00163 // However, that item is item number `firstnumber' which may be '0' or '1'. // 00164 // Be sure to set the 'firstnumber' be '1' if your indices pointing into the // 00165 // pointlist is starting from '1'. Default, it is initialized be '0'. // 00166 // // 00167 // Tetgenio also contains routines for reading and writing TetGen's files as // 00168 // well. Both the library of TetGen and TetView use these routines to parse // 00169 // input files, i.e., .node, .poly, .smesh, .ele, .face, and .edge files. // 00170 // Other routines are provided mainly for debugging purpose. // 00171 // // 00173 00174 class tetgenio { 00175 00176 public: 00177 00178 // Maximum number of characters in a file name (including the null). 00179 enum {FILENAMESIZE = 1024}; 00180 00181 // Maxi. numbers of chars in a line read from a file (incl. the null). 00182 enum {INPUTLINESIZE = 1024}; 00183 00184 // The polygon data structure. A "polygon" is a planar polygon. It can 00185 // be arbitrary shaped (convex or non-convex) and bounded by non- 00186 // crossing segments, i.e., the number of vertices it has indictes the 00187 // same number of edges. 00188 // 'vertexlist' is a list of vertex indices (integers), its length is 00189 // indicated by 'numberofvertices'. The vertex indices are odered in 00190 // either counterclockwise or clockwise way. 00191 typedef struct { 00192 int *vertexlist; 00193 int numberofvertices; 00194 } polygon; 00195 00196 static void init(polygon* p) { 00197 p->vertexlist = (int *) NULL; 00198 p->numberofvertices = 0; 00199 } 00200 00201 // The facet data structure. A "facet" is a planar facet. It is used 00202 // to represent a planar straight line graph (PSLG) in two dimension. 00203 // A PSLG contains a list of polygons. It also may conatin holes in it, 00204 // indicated by a list of hole points (their coordinates). 00205 typedef struct { 00206 polygon *polygonlist; 00207 int numberofpolygons; 00208 REAL *holelist; 00209 int numberofholes; 00210 } facet; 00211 00212 static void init(facet* f) { 00213 f->polygonlist = (polygon *) NULL; 00214 f->numberofpolygons = 0; 00215 f->holelist = (REAL *) NULL; 00216 f->numberofholes = 0; 00217 } 00218 00219 // The periodic boundary condition group data structure. A "pbcgroup" 00220 // contains the definition of a pbc and the list of pbc point pairs. 00221 // 'fmark1' and 'fmark2' are the facetmarkers of the two pbc facets f1 00222 // and f2, respectively. 'transmat' is the transformation matrix which 00223 // maps a point in f1 into f2. An array of pbc point pairs are saved 00224 // in 'pointpairlist'. The first point pair is at indices [0] and [1], 00225 // followed by remaining pairs. Two integers per pair. 00226 typedef struct { 00227 int fmark1, fmark2; 00228 REAL transmat[4][4]; 00229 int numberofpointpairs; 00230 int *pointpairlist; 00231 } pbcgroup; 00232 00233 public: 00234 00235 // Items are numbered starting from 'firstnumber' (0 or 1), default is 0. 00236 int firstnumber; 00237 00238 // Dimension of the mesh (2 or 3), default is 3. 00239 int mesh_dim; 00240 00241 // `pointlist': An array of point coordinates. The first point's x 00242 // coordinate is at index [0] and its y coordinate at index [1], its 00243 // z coordinate is at index [2], followed by the coordinates of the 00244 // remaining points. Each point occupies three REALs. 00245 // `pointattributelist': An array of point attributes. Each point's 00246 // attributes occupy `numberofpointattributes' REALs. 00247 // 'addpointlist': An array of additional point coordinates. 00248 // 'addpointattributelist': An array of attributes for addition points. 00249 // `pointmarkerlist': An array of point markers; one int per point. 00250 REAL *pointlist; 00251 REAL *pointattributelist; 00252 REAL *addpointlist; 00253 REAL *addpointattributelist; 00254 int *pointmarkerlist; 00255 int numberofpoints; 00256 int numberofpointattributes; 00257 int numberofaddpoints; 00258 00259 // `elementlist': An array of element (triangle or tetrahedron) corners. 00260 // The first element's first corner is at index [0], followed by its 00261 // other corners in counterclockwise order, followed by any other 00262 // nodes if the element represents a nonlinear element. Each element 00263 // occupies `numberofcorners' ints. 00264 // `elementattributelist': An array of element attributes. Each 00265 // element's attributes occupy `numberofelementattributes' REALs. 00266 // `elementconstraintlist': An array of constraints, i.e. triangle's 00267 // area or tetrahedron's volume; one REAL per element. Input only. 00268 // `neighborlist': An array of element neighbors; 3 or 4 ints per 00269 // element. Output only. 00270 int *tetrahedronlist; 00271 REAL *tetrahedronattributelist; 00272 REAL *tetrahedronvolumelist; 00273 int *neighborlist; 00274 int numberoftetrahedra; 00275 int numberofcorners; 00276 int numberoftetrahedronattributes; 00277 00278 // `facetlist': An array of facets. Each entry is a structure of facet. 00279 // `facetmarkerlist': An array of facet markers; one int per facet. 00280 facet *facetlist; 00281 int *facetmarkerlist; 00282 int numberoffacets; 00283 00284 // `holelist': An array of holes. The first hole's x, y and z 00285 // coordinates are at indices [0], [1] and [2], followed by the 00286 // remaining holes. Three REALs per hole. 00287 REAL *holelist; 00288 int numberofholes; 00289 00290 // `regionlist': An array of regional attributes and volume constraints. 00291 // The first constraint's x, y and z coordinates are at indices [0], 00292 // [1] and [2], followed by the regional attribute at index [3], foll- 00293 // owed by the maximum volume at index [4]. Five REALs per constraint. 00294 // Note that each regional attribute is used only if you select the `A' 00295 // switch, and each volume constraint is used only if you select the 00296 // `a' switch (with no number following). 00297 REAL *regionlist; 00298 int numberofregions; 00299 00300 // `facetconstraintlist': An array of facet maximum area constraints. 00301 // Two REALs per constraint. The first one is the facet marker (cast 00302 // it to int), the second is its maximum area bound. 00303 // Note the 'facetconstraintlist' is used only for the 'q' switch. 00304 REAL *facetconstraintlist; 00305 int numberoffacetconstraints; 00306 00307 // `segmentconstraintlist': An array of seg maximum length constraints. 00308 // Three REALs per constraint. The first two are the indices (pointing 00309 // into 'pointlist') of the endpoints of the segment, the third is its 00310 // maximum length bound. 00311 // Note the 'segmentconstraintlist' is used only for the 'q' switch. 00312 REAL *segmentconstraintlist; 00313 int numberofsegmentconstraints; 00314 00315 // `nodeconstraintlist': An array of segment length constraints. Two 00316 // REALs per constraint. The first one is the index (pointing into 00317 // 'pointlist') of the node, the second is its edge length bound. 00318 // Note the 'nodeconstraintlist' is used only for the 'q' switch. 00319 REAL *nodeconstraintlist; 00320 int numberofnodeconstraints; 00321 00322 // 'pbcgrouplist': An array of periodic boundary condition groups. 00323 pbcgroup *pbcgrouplist; 00324 int numberofpbcgroups; 00325 00326 // `trifacelist': An array of triangular face endpoints. The first 00327 // face's endpoints are at indices [0], [1] and [2], followed by the 00328 // remaining faces. Three ints per face. 00329 // `trifacemarkerlist': An array of face markers; one int per face. 00330 int *trifacelist; 00331 int *trifacemarkerlist; 00332 int numberoftrifaces; 00333 00334 // `edgelist': An array of edge endpoints. The first edge's endpoints 00335 // are at indices [0] and [1], followed by the remaining edges. Two 00336 // ints per edge. 00337 // `edgemarkerlist': An array of edge markers; one int per edge. 00338 int *edgelist; 00339 int *edgemarkerlist; 00340 int numberofedges; 00341 00342 public: 00343 00344 // Initialize routine. 00345 void initialize(); 00346 void deinitialize(); 00347 00348 // Input & output routines. 00349 bool load_node_call(FILE* infile, int markers, char* nodefilename); 00350 bool load_node(char* filename); 00351 bool load_addnodes(char* filename); 00352 bool load_pbc(char* filename); 00353 bool load_var(char* filename); 00354 bool load_poly(char* filename); 00355 bool load_off(char* filename); 00356 bool load_ply(char* filename); 00357 bool load_stl(char* filename); 00358 bool load_medit(char* filename); 00359 bool load_plc(char* filename, int object); 00360 bool load_tetmesh(char* filename); 00361 void save_nodes(char* filename); 00362 void save_elements(char* filename); 00363 void save_faces(char* filename); 00364 void save_edges(char* filename); 00365 void save_neighbors(char* filename); 00366 void save_poly(char* filename); 00367 00368 // Read line and parse string functions. 00369 char *readline(char* string, FILE* infile, int *linenumber); 00370 char *findnextfield(char* string); 00371 char *readnumberline(char* string, FILE* infile, char* infilename); 00372 char *findnextnumber(char* string); 00373 00374 // Constructor and destructor. 00375 tetgenio() {initialize();} 00376 ~tetgenio() {deinitialize();} 00377 }; 00378 00380 // // 00381 // The tetgenbehavior data type // 00382 // // 00383 // Used to parse command line switches and file names. // 00384 // // 00385 // It includes a list of variables corresponding to the commandline switches // 00386 // for control the behavior of TetGen. These varibales are all initialized // 00387 // to their default values. // 00388 // // 00389 // parse_commandline() provides an simple interface to set the vaules of the // 00390 // variables. It accepts the standard parameters (e.g., 'argc' and 'argv') // 00391 // that pass to C/C++ main() function. Alternatively a string which contains // 00392 // the command line options can be used as its parameter. // 00393 // // 00394 // You don't need to understand this data type. It can be implicitly called // 00395 // by the global function "tetrahedralize()" defined below. The necessary // 00396 // thing you need to know is the meaning of command line switches of TetGen. // 00397 // They are described in the third section of the user's manual. // 00398 // // 00400 00401 class tetgenbehavior { 00402 00403 public: 00404 00405 // Labels define the objects which are acceptable by TetGen. They are 00406 // recognized by the file extensions. 00407 // - NODES, a list of nodes (.node); 00408 // - POLY, a piecewise linear complex (.poly or .smesh); 00409 // - OFF, a polyhedron (.off, Geomview's file format); 00410 // - PLY, a polyhedron (.ply, file format from gatech); 00411 // - STL, a surface mesh (.stl, stereolithography format); 00412 // - MEDIT, a surface mesh (.mesh, Medit's file format); 00413 // - MESH, a tetrahedral mesh (.ele). 00414 // If no extension is available, the imposed commandline switch 00415 // (-p or -r) implies the object. 00416 00417 enum objecttype {NONE, NODES, POLY, OFF, PLY, STL, MEDIT, MESH}; 00418 00419 // Variables of command line switches. Each variable is corresponding 00420 // to a specific switch and will be properly initialized. Read the 00421 // user's manul to find out the meaning of these switches. 00422 00423 int plc; // '-p' switch, 0. 00424 int refine; // '-r' switch, 0. 00425 int quality; // '-q' switch, 0. 00426 int varvolume; // '-a' switch without number, 0. 00427 int fixedvolume; // '-a' switch with number, 0. 00428 int insertaddpoints; // '-i' switch, 0. 00429 int regionattrib; // '-A' switch, 0. 00430 int detectinter; // '-d' switch, 0. 00431 int offcenter; // '-R' switch, 0. 00432 int zeroindex; // '-z' switch, 0. 00433 int order; // element order, specified after '-o' switch, 1. 00434 int facesout; // '-f' switch, 0. 00435 int edgesout; // '-e' switch, 0. 00436 int neighout; // '-n' switch, 0. 00437 int meditview; // '-g' switch, 0. 00438 int gidview; // '-G' switch, 0. 00439 int geomview; // '-O' switch, 0. 00440 int nobound; // '-B' switch, 0. 00441 int nonodewritten; // '-N' switch, 0. 00442 int noelewritten; // '-E' switch, 0. 00443 int nofacewritten; // '-F' switch, 0. 00444 int noiterationnum; // '-I' switch, 0. 00445 int nomerge; // count of how often '-M' switch is selected, 0. 00446 int nobisect; // count of how often '-Y' switch is selected, 0. 00447 int noflip; // do not perform flips. '-Y' switch. 0. 00448 int nofreevert; // Remove free vertices. '-X' switch. 0. 00449 int nojettison; // do not jettison redundants nodes. '-J' switch. 0. 00450 int steiner; // number after '-S' switch. 0. 00451 int dofullperturb; // do full permutation. '-P' switch, 0. 00452 int docheck; // '-C' switch, 0. 00453 int quiet; // '-Q' switch, 0. 00454 int verbose; // count of how often '-V' switch is selected, 0. 00455 int useshelles; // '-p', '-r', '-q', '-d', or '-c' switch, 0. 00456 REAL minratio; // number after '-q' switch, 2.0. 00457 REAL goodratio; // number calculated from 'minratio', 0.0. 00458 REAL minangle; // minimum angle bound, 20.0. 00459 REAL goodangle; // cosine squared of minangle, 0.0. 00460 REAL maxvolume; // number after '-a' switch, -1.0. 00461 REAL maxdihedral; // number after '-s' switch, 175.0. 00462 REAL alpha3; // number after '-R' switch, 0.6. 00463 REAL epsilon; // number after '-T' switch, 1.0e-8. 00464 enum objecttype object; // determined by -p, or -r switch. NONE. 00465 00466 // Variables used to save command line switches and in/out file names. 00467 char commandline[1024]; 00468 char infilename[1024]; 00469 char outfilename[1024]; 00470 00471 tetgenbehavior(); 00472 ~tetgenbehavior() {} 00473 00474 void versioninfo(); 00475 void syntax(); 00476 void usage(); 00477 00478 // Command line parse routine. 00479 bool parse_commandline(int argc, char **argv); 00480 bool parse_commandline(char *switches) { 00481 return parse_commandline(0, &switches); 00482 } 00483 }; 00484 00486 // // 00487 // Geometric predicates // 00488 // // 00489 // Return one of the values +1, 0, and -1 on basic geometric questions such // 00490 // as the orientation of point sets, in-circle, and in-sphere tests. They // 00491 // are basic units for the composition of geometric algorithms. TetGen uses // 00492 // two 3D geometric predicates, which are the orientation test and the in- // 00493 // sphere test (e.g. the locally Deklaunay test). // 00494 // // 00495 // Orientation test: let a, b, c be a sequence of 3 non-collinear points in // 00496 // R^3. They defines a unique hypeplane H. Let H+ and H- be the two spaces // 00497 // separated by H, which are defined as follows (using the left-hand rule): // 00498 // make a fist using your left hand in such a way that your fingers follow // 00499 // the order of a, b and c, then your thumb is pointing to H+. Given any // 00500 // point d in R^3, the orientation test returns +1 if d lies in H+, -1 if d // 00501 // lies in H-, or 0 if d lies on H. // 00502 // // 00503 // In-sphere test: let a, b, c, d be 4 non-coplanar points in R^3. They // 00504 // defines a unique circumsphere S. Given any point e in R^3, the in-sphere // 00505 // test returns +1 if e lies inside S, or -1 if e lies outside S, or 0 if e // 00506 // lies on S. // 00507 // // 00508 // The correctness of geometric predicates is crucial for the control flow // 00509 // and hence for the correctness and robustness of an implementation of a // 00510 // geometric algorithm. The following routines use arbitrary precision // 00511 // floating-point arithmetic. They are fast and robust. It is provided by J. // 00512 // Schewchuk in public domain (http://www.cs.cmu.edu/~quake/robust.html). // 00513 // The source code are found in a separate file "predicates.cxx". // 00514 // // 00516 00517 REAL exactinit(); 00518 REAL orient3d(REAL *pa, REAL *pb, REAL *pc, REAL *pd); 00519 REAL insphere(REAL *pa, REAL *pb, REAL *pc, REAL *pd, REAL *pe); 00520 00522 // // 00523 // The tetgenmesh data type // 00524 // // 00525 // Includes data types and mesh routines for creating tetrahedral meshes and // 00526 // Delaunay tetrahedralizations, mesh input & output, and so on. // 00527 // // 00528 // An object of tetgenmesh can be used to store a triangular or tetrahedral // 00529 // mesh and its settings. TetGen's functions operates on one mesh each time. // 00530 // This type allows reusing of the same function for different meshes. // 00531 // // 00532 // The mesh data structure (tetrahedron-based and triangle-edge data struct- // 00533 // ures) are declared. There are other accessary data type defined as well, // 00534 // for efficient memory management and link list operations, etc. // 00535 // // 00536 // All algorithms TetGen used are implemented in this data type as member // 00537 // functions. References of these algorithms can be found in user's manual. // 00538 // // 00539 // It's not necessary to understand this type. There is a global function // 00540 // "tetrahedralize()" (defined at the end of this file) implicitly creates // 00541 // the object and calls its member functions according to the command line // 00542 // switches you specified. // 00543 // // 00545 00546 class tetgenmesh { 00547 00548 public: 00549 00550 // Maximum number of characters in a file name (including the null). 00551 enum {FILENAMESIZE = 1024}; 00552 00553 // For efficiency, a variety of data structures are allocated in bulk. 00554 // The following constants determine how many of each structure is 00555 // allocated at once. 00556 enum {VERPERBLOCK = 4092, SUBPERBLOCK = 4092, ELEPERBLOCK = 8188}; 00557 00558 // Used for the point location scheme of Mucke, Saias, and Zhu, to 00559 // decide how large a random sample of tetrahedra to inspect. 00560 enum {SAMPLEFACTOR = 11}; 00561 00562 // Labels that signify two edge rings of a triangle defined in Muecke's 00563 // triangle-edge data structure, one (CCW) traversing edges in count- 00564 // erclockwise direction and one (CW) in clockwise direction. 00565 enum {CCW = 0, CW = 1}; 00566 00567 // Labels that signify whether a record consists primarily of pointers 00568 // or of floating-point words. Used to make decisions about data 00569 // alignment. 00570 enum wordtype {POINTER, FLOATINGPOINT}; 00571 00572 // Labels that signify the type of a vertex. An UNUSEDVERTEX is a vertex 00573 // read from input (.node file or tetgenio structure) or an isolated 00574 // vertex (outside the mesh). It is the default type for a newpoint. 00575 enum verttype {UNUSEDVERTEX, DUPLICATEDVERTEX, NACUTEVERTEX, ACUTEVERTEX, 00576 FREESEGVERTEX, FACETVERTEX, FREESUBVERTEX, FREEVOLVERTEX, 00577 DEADVERTEX = -32768}; 00578 00579 // Labels that signify the type of a subface. A subface is SHARPSUB if 00580 // it is in a facet which forms an acute dihedral angle (less than 90 00581 // degree with other facets). It is SKINNYSUB if it has two segments 00582 // form a small angle (less than 10 degree). 00583 enum shestype {NSHARPNSKINNY, SHARPSUB, SKINNYSUB, SHARPSKINNYSUB}; 00584 00585 // Labels that signify the type of flips can be applied on a face. 00586 // A flipable face has the one of the types T23, T32, T22, and T44. 00587 // Types UNFLIPABLE, NONCONVEX are unflipable. 00588 enum fliptype {T23, T32, T22, T44, UNFLIPABLE, FORBIDDENFACE, 00589 FORBIDDENEDGE, NONCONVEX}; 00590 00591 // Labels that signify the result of triangle-triangle intersection test. 00592 // Two triangles are DISJOINT, or adjoint at a vertex SHAREVERTEX, or 00593 // adjoint at an edge SHAREEDGE, or coincident SHAREFACE or INTERSECT. 00594 enum interresult {DISJOINT, SHAREVERTEX, SHAREEDGE, SHAREFACE, INTERSECT}; 00595 00596 // Labels that signify the result of point location. The result of a 00597 // search indicates that the point falls inside a tetrahedron, inside 00598 // a triangle, on an edge, on a vertex, or outside the mesh. 00599 enum locateresult {INTETRAHEDRON, ONFACE, ONEDGE, ONVERTEX, OUTSIDE}; 00600 00601 // Labels that signify the result of vertex insertion. The result 00602 // indicates that the vertex was inserted with complete success, was 00603 // inserted but encroaches upon a subsegment, was not inserted because 00604 // it lies on a segment, or was not inserted because another vertex 00605 // occupies the same location. 00606 enum insertsiteresult {SUCCESSINTET, SUCCESSONFACE, SUCCESSONEDGE, 00607 DUPLICATEPOINT, OUTSIDEPOINT}; 00608 00609 // Labels that signify the result of direction finding. The result 00610 // indicates that a segment connecting the two query points accross 00611 // an edge of the direction triangle/tetrahedron, across a face of 00612 // the direction tetrahedron, along the left edge of the direction 00613 // triangle/tetrahedron, along the right edge of the direction 00614 // triangle/tetrahedron, or along the top edge of the tetrahedron. 00615 enum finddirectionresult {ACROSSEDGE, ACROSSFACE, LEFTCOLLINEAR, 00616 RIGHTCOLLINEAR, TOPCOLLINEAR, BELOWHULL}; 00617 00619 // // 00620 // The basic mesh element data structures // 00621 // // 00622 // There are four types of mesh elements: tetrahedra, subfaces, subsegments, // 00623 // and points, where subfaces and subsegments are triangles and edges which // 00624 // appear on boundaries. A tetrahedralization of a 3D point set comprises // 00625 // tetrahedra and points; a surface mesh of a 3D domain comprises subfaces // 00626 // subsegments and points. The elements of all the four types consist of a // 00627 // tetrahedral mesh of a 3D domain. However, TetGen uses three data types: // 00628 // 'tetrahedron', 'shellface', and 'point'. A 'tetrahedron' is a tetrahedron;// 00629 // while a 'shellface' can be either a subface or a subsegment; and a 'point'// 00630 // is a point. These three data types, linked by pointers comprise a mesh. // 00631 // // 00632 // A tetrahedron primarily consists of a list of 4 pointers to its corners, // 00633 // a list of 4 pointers to its adjoining tetrahedra, a list of 4 pointers to // 00634 // its adjoining subfaces (when subfaces are needed). Optinoally, (depending // 00635 // on the selected switches), it may contain an arbitrary number of user- // 00636 // defined floating-point attributes, an optional maximum volume constraint // 00637 // (for -a switch), and a pointer to a list of high-order nodes (-o2 switch).// 00638 // Since the size of a tetrahedron is not determined until running time, it // 00639 // is not simply declared as a structure. // 00640 // // 00641 // The data structure of tetrahedron also stores the geometrical information.// 00642 // Let t be a tetrahedron, v0, v1, v2, and v3 be the 4 nodes corresponding // 00643 // to the order of their storage in t. v3 always has a negative orientation // 00644 // with respect to v0, v1, v2 (ie,, v3 lies above the oriented plane passes // 00645 // through v0, v1, v2). Let the 4 faces of t be f0, f1, f2, and f3. Vertices // 00646 // of each face are stipulated as follows: f0 (v0, v1, v2), f1 (v0, v3, v1), // 00647 // f2 (v1, v3, v2), f3 (v2, v3, v0). // 00648 // // 00649 // A subface has 3 pointers to vertices, 3 pointers to adjoining subfaces, 3 // 00650 // pointers to adjoining subsegments, 2 pointers to adjoining tetrahedra, a // 00651 // boundary marker(an integer). Like a tetrahedron, the pointers to vertices,// 00652 // subfaces, and subsegments are ordered in a way that indicates their geom- // 00653 // etric relation. Let s be a subface, v0, v1 and v2 be the 3 nodes corres- // 00654 // ponding to the order of their storage in s, e0, e1 and e2 be the 3 edges,// 00655 // then we have: e0 (v0, v1), e1 (v1, v2), e2 (v2, v0). // 00656 // // 00657 // A subsegment has exactly the same data fields as a subface has, but only // 00658 // uses some of them. It has 2 pointers to its endpoints, 2 pointers to its // 00659 // adjoining (and collinear) subsegments, a pointer to a subface containing // 00660 // it (there may exist any number of subfaces having it, choose one of them // 00661 // arbitrarily). The geometric relation between its endpoints and adjoining // 00662 // subsegments is kept with respect to the storing order of its endpoints. // 00663 // // 00664 // The data structure of point is relatively simple. A point is a list of // 00665 // floating-point numbers, starting with the x, y, and z coords, followed by // 00666 // an arbitrary number of optional user-defined floating-point attributes, // 00667 // an integer boundary marker, an integer for the point type, and a pointer // 00668 // to a tetrahedron (used for speeding up point location). // 00669 // // 00670 // For a tetrahedron on a boundary (or a hull) of the mesh, some or all of // 00671 // the adjoining tetrahedra may not be present. For an interior tetrahedron, // 00672 // often no neighboring subfaces are present, Such absent tetrahedra and // 00673 // subfaces are never represented by the NULL pointers; they are represented // 00674 // by two special records: `dummytet', the tetrahedron fills "outer space", // 00675 // and `dummysh', the vacuous subfaces which are omnipresent. // 00676 // // 00677 // Tetrahedra and adjoining subfaces are glued together through the pointers // 00678 // saved in each data fields of them. Subfaces and adjoining subsegments are // 00679 // connected in the same fashion. However, there are no pointers directly // 00680 // gluing tetrahedra and adjoining subsegments. For the purpose of saving // 00681 // space, the connections between tetrahedra and subsegments are entirely // 00682 // mediated through subfaces. The following part explains how subfaces are // 00683 // connected in TetGen. // 00684 // // 00686 00688 // // 00689 // The subface-subface and subface-subsegment connections // 00690 // // 00691 // Adjoining subfaces sharing a common edge are connected in such a way that // 00692 // they form a face ring around the edge. It is indeed a single linked list // 00693 // which is cyclic, e.g., one can start from any subface in it and traverse // 00694 // back. When the edge is not a subsegment, the ring only has two coplanar // 00695 // subfaces which are pointing to each other. Otherwise, the face ring may // 00696 // have any number of subfaces (and are not all coplanar). // 00697 // // 00698 // How is the face ring formed? Let s be a subsegment, f is one of subfaces // 00699 // containing s as an edge. The direction of s is stipulated from its first // 00700 // endpoint to its second (according to their storage in s). Once the dir of // 00701 // s is determined, the other two edges of f are oriented to follow this dir.// 00702 // The "directional normal" N_f is a vector formed from any point in f and a // 00703 // points orthogonally above f. // 00704 // // 00705 // The face ring of s is a cyclic ordered set of subfaces containing s, i.e.,// 00706 // F(s) = {f1, f2, ..., fn}, n >= 1. Where the order is defined as follows: // 00707 // let fi, fj be two faces in F(s), the "normal-angle", NAngle(i,j) (range // 00708 // from 0 to 360 degree) is the angle between the N_fi and N_fj; then fi is // 00709 // in front of fj (or symbolically, fi < fj) if there exists another fk in // 00710 // F(s), and NAangle(k, i) < NAngle(k, j). The face ring of s is: f1 < f2 < // 00711 // ... < fn < f1. // 00712 // // 00713 // The easiest way to imagine how a face ring is formed is to use the right- // 00714 // hand rule. Make a fist using your right hand with the thumb pointing to // 00715 // the direction of the subsegment. The face ring is connected following the // 00716 // direction of your fingers. // 00717 // // 00718 // The subface and subsegment are also connected through pointers stored in // 00719 // their own data fields. Every subface has a pointer to its adjoining sub- // 00720 // segment. However, a subsegment only has one pointer to a subface which is // 00721 // containing it. Such subface can be chosen arbitrarily, other subfaces are // 00722 // found through the face ring. // 00723 // // 00725 00726 // The tetrahedron data structure. Fields of a tetrahedron contains: 00727 // - a list of four adjoining tetrahedra; 00728 // - a list of four vertices; 00729 // - a list of four subfaces (optional, used for -p switch); 00730 // - a list of user-defined floating-point attributes (optional); 00731 // - a volume constraint (optional, used for -a switch); 00732 // - a pointer to a list of high-ordered nodes (optional, -o2 switch); 00733 00734 typedef REAL **tetrahedron; 00735 00736 // The shellface data structure. Fields of a shellface contains: 00737 // - a list of three adjoining subfaces; 00738 // - a list of three vertices; 00739 // - a list of two adjoining tetrahedra; 00740 // - a list of three adjoining subsegments; 00741 // - a pointer to a badface containing it (used for -q); 00742 // - an area constraint (optional, used for -q); 00743 // - an integer for boundary marker; 00744 // - an integer for type: SHARPSEGMENT, NONSHARPSEGMENT, ...; 00745 // - an integer for pbc group (optional, if in->pbcgrouplist exists); 00746 00747 typedef REAL **shellface; 00748 00749 // The point data structure. It is actually an array of REALs: 00750 // - x, y and z coordinates; 00751 // - a list of user-defined point attributes (optional); 00752 // - an edge length constraint (optional, used for -q); 00753 // - a pointer to a simplex (tet, tri, edge, or vertex); 00754 // - a pointer to a parent (or duplicate) point; 00755 // - a pointer to another pbc point (optional); 00756 // - an integer for boundary marker; 00757 // - an integer for verttype: INPUTVERTEX, FREEVERTEX, ...; 00758 00759 typedef REAL *point; 00760 00762 // // 00763 // The mesh handle (triface, face) data types // 00764 // // 00765 // Two special data types, 'triface' and 'face' are defined for maintaining // 00766 // and updating meshes. They are like pointers (or handles), which allow you // 00767 // to hold one particular part of the mesh, i.e., a tetrahedron, a triangle, // 00768 // an edge and a vertex. However, these data types do not themselves store // 00769 // any part of the mesh. The mesh is made of the data types defined above. // 00770 // // 00771 // Muecke's "triangle-edge" data structure is the prototype for these data // 00772 // types. It allows a universal representation for every tetrahedron, // 00773 // triangle, edge and vertex. For understanding the following descriptions // 00774 // of these handle data structures, readers are required to read both the // 00775 // introduction and implementation detail of "triangle-edge" data structure // 00776 // in Muecke's thesis. // 00777 // // 00778 // A 'triface' represents a face of a tetrahedron and an oriented edge of // 00779 // the face simultaneously. It has a pointer 'tet' to a tetrahedron, an // 00780 // integer 'loc' (range from 0 to 3) as the face index, and an integer 'ver' // 00781 // (range from 0 to 5) as the edge version. A face of the tetrahedron can be // 00782 // uniquly determined by the pair (tet, loc), and an oriented edge of this // 00783 // face can be uniquly determined by the triple (tet, loc, ver). Therefore, // 00784 // different usages of one triface are possible. If we only use the pair // 00785 // (tet, loc), it refers to a face, and if we add the 'ver' additionally to // 00786 // the pair, it is an oriented edge of this face. // 00787 // // 00788 // A 'face' represents a subface and an oriented edge of it simultaneously. // 00789 // It has a pointer 'sh' to a subface, an integer 'shver'(range from 0 to 5) // 00790 // as the edge version. The pair (sh, shver) determines a unique oriented // 00791 // edge of this subface. A 'face' is also used to represent a subsegment, // 00792 // in this case, 'sh' points to the subsegment, and 'shver' indicates the // 00793 // one of two orientations of this subsegment, hence, it only can be 0 or 1. // 00794 // // 00795 // Mesh navigation and updating are accomplished through a set of mesh // 00796 // manipulation primitives which operate on trifaces and faces. They are // 00797 // introduced below. // 00798 // // 00800 00801 class triface { 00802 00803 public: 00804 00805 tetrahedron* tet; 00806 int loc, ver; 00807 00808 // Constructors; 00809 triface() : tet(0), loc(0), ver(0) {} 00810 // Operators; 00811 triface& operator=(const triface& t) { 00812 tet = t.tet; loc = t.loc; ver = t.ver; 00813 return *this; 00814 } 00815 bool operator==(triface& t) { 00816 return tet == t.tet && loc == t.loc && ver == t.ver; 00817 } 00818 bool operator!=(triface& t) { 00819 return tet != t.tet || loc != t.loc || ver != t.ver; 00820 } 00821 }; 00822 00823 class face { 00824 00825 public: 00826 00827 shellface *sh; 00828 int shver; 00829 00830 // Constructors; 00831 face() : sh(0), shver(0) {} 00832 // Operators; 00833 face& operator=(const face& s) { 00834 sh = s.sh; shver = s.shver; 00835 return *this; 00836 } 00837 bool operator==(face& s) {return (sh == s.sh) && (shver == s.shver);} 00838 bool operator!=(face& s) {return (sh != s.sh) || (shver != s.shver);} 00839 }; 00840 00842 // // 00843 // The badface structure // 00844 // // 00845 // A multiple usages structure. Despite of its name, a 'badface' can be used // 00846 // to represent the following objects: // 00847 // - a face of a tetrahedron which is (possibly) non-Delaunay; // 00848 // - an encroached subsegment or subface; // 00849 // - a bad-quality tetrahedron, i.e, has too large radius-edge ratio; // 00850 // - a sliver, i.e., has good radius-edge ratio but nearly zero volume; // 00851 // - a degenerate tetrahedron (see routine checkdegetet()). // 00852 // - a recently flipped face (saved for undoing the flip later). // 00853 // // 00854 // It has the following fields: 'tt' holds a tetrahedron; 'ss' holds a sub- // 00855 // segment or subface; 'cent' is the circumcent of 'tt' or 'ss', 'key' is a // 00856 // special value depending on the use, it can be either the square of the // 00857 // radius-edge ratio of 'tt' or the flipped type of 'tt'; 'forg', 'fdest', // 00858 // 'fapex', and 'foppo' are vertices saved for checking the object in 'tt' // 00859 // or 'ss' is still the same when it was stored; 'noppo' is the fifth vertex // 00860 // of a degenerate point set. 'previtem' and 'nextitem' implement a double // 00861 // link for managing many basfaces. // 00862 // // 00864 00865 struct badface { 00866 triface tt; 00867 face ss; 00868 REAL key; 00869 REAL cent[3]; 00870 point forg, fdest, fapex, foppo; 00871 point noppo; 00872 struct badface *previtem, *nextitem; 00873 }; 00874 00876 // // 00877 // The pbcdata structure // 00878 // // 00879 // A pbcdata stores data of a periodic boundary condition defined on a pair // 00880 // of facets or segments. Let f1 and f2 define a pbcgroup. 'fmark' saves the // 00881 // facet markers of f1 and f2; 'ss' contains two subfaces belong to f1 and // 00882 // f2, respectively. Let s1 and s2 define a segment pbcgroup. 'segid' are // 00883 // the segment ids of s1 and s2; 'ss' contains two segments belong to s1 and // 00884 // s2, respectively. 'transmat' are two transformation matrices. transmat[0] // 00885 // transforms a point of f1 (or s1) into a point of f2 (or s2), transmat[1] // 00886 // does the inverse. // 00887 // // 00889 00890 struct pbcdata { 00891 int fmark[2]; 00892 int segid[2]; 00893 face ss[2]; 00894 REAL transmat[2][4][4]; 00895 }; 00896 00898 // // 00899 // The list, link and queue data structures // 00900 // // 00901 // These data types are used to manipulate a set of (same-typed) data items. // 00902 // For a given set S = {a, b, c, ...}, a list stores the elements of S in a // 00903 // piece of continuous memory. It allows quickly accessing each element of S,// 00904 // thus is suitable for storing a fix-sized set. While a link stores its // 00905 // elements incontinuously. It allows quickly inserting or deleting an item, // 00906 // thus is suitable for storing a size-variable set. A queue is basically a // 00907 // special case of a link where one data element joins the link at the end // 00908 // and leaves in an ordered fashion at the other end. // 00909 // // 00911 00912 // The compfunc data type. "compfunc" is a pointer to a linear-order 00913 // function, which takes two 'void*' arguments and returning an 'int'. 00914 // 00915 // A function: int cmp(const T &, const T &), is said to realize a 00916 // linear order on the type T if there is a linear order <= on T such 00917 // that for all x and y in T satisfy the following relation: 00918 // -1 if x < y. 00919 // comp(x, y) = 0 if x is equivalent to y. 00920 // +1 if x > y. 00921 typedef int (*compfunc) (const void *, const void *); 00922 00923 // The predefined compare functions for primitive data types. They 00924 // take two pointers of the corresponding date type, perform the 00925 // comparation, and return -1, 0 or 1 indicating the default linear 00926 // order of them. 00927 00928 // Compare two 'integers'. 00929 static int compare_2_ints(const void* x, const void* y); 00930 // Compare two 'longs'. 00931 static int compare_2_longs(const void* x, const void* y); 00932 // Compare two 'unsigned longs'. 00933 static int compare_2_unsignedlongs(const void* x, const void* y); 00934 00935 // The function used to determine the size of primitive data types and 00936 // set the corresponding predefined linear order functions for them. 00937 static void set_compfunc(char* str, int* itembytes, compfunc* pcomp); 00938 00940 // // 00941 // List data structure. // 00942 // // 00943 // A 'list' is an array of items with automatically reallocation of memory. // 00944 // It behaves like an array. // 00945 // // 00946 // 'base' is the starting address of the array; The memory unit in list is // 00947 // byte, i.e., sizeof(char). 'itembytes' is the size of each item in byte, // 00948 // so that the next item in list will be found at the next 'itembytes' // 00949 // counted from the current position. // 00950 // // 00951 // 'items' is the number of items stored in list. 'maxitems' indicates how // 00952 // many items can be stored in this list. 'expandsize' is the increasing // 00953 // size (items) when the list is full. // 00954 // // 00955 // 'comp' is a pointer pointing to a linear order function for the list. // 00956 // default it is set to 'NULL'. // 00957 // // 00958 // The index of list always starts from zero, i.e., for a list L contains // 00959 // n elements, the first element is L[0], and the last element is L[n-1]. // 00960 // This feature lets lists likes C/C++ arrays. // 00961 // // 00963 00964 class list { 00965 00966 public: 00967 00968 char *base; 00969 int itembytes; 00970 int items, maxitems, expandsize; 00971 compfunc comp; 00972 00973 public: 00974 00975 list(int itbytes, compfunc pcomp, int mitems = 256, int exsize = 128) { 00976 listinit(itbytes, pcomp, mitems, exsize); 00977 } 00978 list(char* str, int mitems = 256, int exsize = 128) { 00979 set_compfunc(str, &itembytes, &comp); 00980 listinit(itembytes, comp, mitems, exsize); 00981 } 00982 ~list() { free(base); } 00983 00984 void *operator[](int i) { return (void *) (base + i * itembytes); } 00985 00986 void listinit(int itbytes, compfunc pcomp, int mitems, int exsize); 00987 void setcomp(compfunc compf) { comp = compf; } 00988 void clear() { items = 0; } 00989 int len() { return items; } 00990 void *append(void* appitem); 00991 void *insert(int pos, void* insitem); 00992 void del(int pos, int order); 00993 int hasitem(void* checkitem); 00994 void sort(); 00995 }; 00996 00998 // // 00999 // Memorypool data structure. // 01000 // // 01001 // A type used to allocate memory. (It is incorporated from Shewchuk's // 01002 // Triangle program) // 01003 // // 01004 // firstblock is the first block of items. nowblock is the block from which // 01005 // items are currently being allocated. nextitem points to the next slab // 01006 // of free memory for an item. deaditemstack is the head of a linked list // 01007 // (stack) of deallocated items that can be recycled. unallocateditems is // 01008 // the number of items that remain to be allocated from nowblock. // 01009 // // 01010 // Traversal is the process of walking through the entire list of items, and // 01011 // is separate from allocation. Note that a traversal will visit items on // 01012 // the "deaditemstack" stack as well as live items. pathblock points to // 01013 // the block currently being traversed. pathitem points to the next item // 01014 // to be traversed. pathitemsleft is the number of items that remain to // 01015 // be traversed in pathblock. // 01016 // // 01017 // itemwordtype is set to POINTER or FLOATINGPOINT, and is used to suggest // 01018 // what sort of word the record is primarily made up of. alignbytes // 01019 // determines how new records should be aligned in memory. itembytes and // 01020 // itemwords are the length of a record in bytes (after rounding up) and // 01021 // words. itemsperblock is the number of items allocated at once in a // 01022 // single block. items is the number of currently allocated items. // 01023 // maxitems is the maximum number of items that have been allocated at // 01024 // once; it is the current number of items plus the number of records kept // 01025 // on deaditemstack. // 01026 // // 01028 01029 class memorypool { 01030 01031 public: 01032 01033 void **firstblock, **nowblock; 01034 void *nextitem; 01035 void *deaditemstack; 01036 void **pathblock; 01037 void *pathitem; 01038 wordtype itemwordtype; 01039 int alignbytes; 01040 int itembytes, itemwords; 01041 int itemsperblock; 01042 long items, maxitems; 01043 int unallocateditems; 01044 int pathitemsleft; 01045 01046 public: 01047 01048 memorypool(); 01049 memorypool(int, int, enum wordtype, int); 01050 ~memorypool(); 01051 01052 void poolinit(int, int, enum wordtype, int); 01053 void restart(); 01054 void *alloc(); 01055 void dealloc(void*); 01056 void traversalinit(); 01057 void *traverse(); 01058 }; 01059 01061 // // 01062 // Link data structure. // 01063 // // 01064 // A 'link' is a double linked nodes. It uses the memorypool data structure // 01065 // for memory management. Following is an image of a link. // 01066 // // 01067 // head-> ____0____ ____1____ ____2____ _________<-tail // 01068 // |__next___|--> |__next___|--> |__next___|--> |__NULL___| // 01069 // |__NULL___|<-- |__prev___|<-- |__prev___|<-- |__prev___| // 01070 // | | |_ _| |_ _| | | // 01071 // | | |_ Data1 _| |_ Data2 _| | | // 01072 // |_________| |_________| |_________| |_________| // 01073 // // 01074 // The unit size for storage is size of pointer, which may be 4-byte (in 32- // 01075 // bit machine) or 8-byte (in 64-bit machine). The real size of an item is // 01076 // stored in 'linkitembytes'. // 01077 // // 01078 // 'head' and 'tail' are pointers pointing to the first and last nodes. They // 01079 // do not conatin data (See above). // 01080 // // 01081 // 'nextlinkitem' is a pointer pointing to a node which is the next one will // 01082 // be traversed. 'curpos' remembers the position (1-based) of the current // 01083 // traversing node. // 01084 // // 01085 // 'linkitems' indicates how many items in link. Note it is different with // 01086 // 'items' of memorypool. // 01087 // // 01088 // The index of link starts from 1, i.e., for a link K contains n elements, // 01089 // the first element of the link is K[1], and the last element is K[n]. // 01090 // See the above figure. // 01091 // // 01093 01094 class link : public memorypool { 01095 01096 public: 01097 01098 void **head, **tail; 01099 void *nextlinkitem; 01100 int linkitembytes; 01101 int linkitems; 01102 int curpos; 01103 compfunc comp; 01104 01105 public: 01106 01107 link(int _itembytes, compfunc _comp, int itemcount) { 01108 linkinit(_itembytes, _comp, itemcount); 01109 } 01110 link(char* str, int itemcount) { 01111 set_compfunc(str, &linkitembytes, &comp); 01112 linkinit(linkitembytes, comp, itemcount); 01113 } 01114 01115 void linkinit(int _itembytes, compfunc _comp, int itemcount); 01116 void setcomp(compfunc compf) { comp = compf; } 01117 void rewind() { nextlinkitem = *head; curpos = 1; } 01118 void goend() { nextlinkitem = *(tail + 1); curpos = linkitems; } 01119 long len() { return linkitems; } 01120 void clear(); 01121 bool move(int numberofnodes); 01122 bool locate(int pos); 01123 void *add(void* newitem); 01124 void *insert(int pos, void* insitem); 01125 void *del(void* delitem); 01126 void *del(int pos); 01127 void *getitem(); 01128 void *getnitem(int pos); 01129 int hasitem(void* checkitem); 01130 }; 01131 01133 // // 01134 // Queue data structure. // 01135 // // 01136 // A 'queue' is a basically a link. Following is an image of a queue. // 01137 // ___________ ___________ ___________ // 01138 // Pop() <-- |_ _|<--|_ _|<--|_ _| <-- Push() // 01139 // |_ Data0 _| |_ Data1 _| |_ Data2 _| // 01140 // |___________| |___________| |___________| // 01141 // queue head queue tail // 01142 // // 01144 01145 class queue : public link { 01146 01147 public: 01148 01149 queue(int bytes, int count = 256) : link(bytes, NULL, count) {} 01150 queue(char* str, int count = 256) : link(str, count) {} 01151 01152 int empty() { return linkitems == 0; } 01153 void *push(void* newitem) { return link::add(newitem); } 01154 void *bot() { return link::getnitem(1); } 01155 void *pop() { return link::del(1); } 01156 }; 01157 01159 // // 01160 // Following are variables used in 'tetgenmesh' for miscellaneous purposes. // 01161 // // 01163 01164 // Pointer to an object of 'tetgenio', which contains input data. 01165 tetgenio *in; 01166 01167 // Pointer to an object of 'tetgenbehavor', which contains the user- 01168 // defined command line swithes and filenames. 01169 tetgenbehavior *b; 01170 01171 // Variables used to allocate and access memory for tetrahedra, subfaces 01172 // subsegments, points, encroached subfaces, encroached subsegments, 01173 // bad-quality tetrahedra, and so on. 01174 memorypool *tetrahedrons; 01175 memorypool *subfaces; 01176 memorypool *subsegs; 01177 memorypool *points; 01178 memorypool *badsubsegs; 01179 memorypool *badsubfaces; 01180 memorypool *badtetrahedrons; 01181 memorypool *flipstackers; 01182 01183 // Pointer to the 'tetrahedron' that occupies all of "outer space". 01184 tetrahedron *dummytet; 01185 tetrahedron *dummytetbase; // Keep base address so we can free it later. 01186 01187 // Pointer to the omnipresent subface. Referenced by any tetrahedron, 01188 // or subface that isn't connected to a subface at that location. 01189 shellface *dummysh; 01190 shellface *dummyshbase; // Keep base address so we can free it later. 01191 01192 // A point above the plane in which the facet currently being used lies. 01193 // It is used as a reference point for orient3d(). 01194 point *facetabovepointarray, abovepoint; 01195 01196 // Array (size = numberoftetrahedra * 6) for storing high-order nodes of 01197 // tetrahedra (only used when -o2 switch is selected). 01198 point *highordertable; 01199 01200 // Arrays for storing and searching pbc data. 'subpbcgrouptable', (size 01201 // is numberofpbcgroups) for pbcgroup of subfaces. 'segpbcgrouptable', 01202 // a list for pbcgroup of segments. Because a segment can have several 01203 // pbcgroup incident on it, its size is unknown on input, it will be 01204 // found in 'createsegpbcgrouptable()'. 01205 pbcdata *subpbcgrouptable; 01206 list *segpbcgrouptable; 01207 // A map for searching the pbcgroups of a given segment. 'idx2segpglist' 01208 // (size = number of input segments + 1), and 'segpglist'. 01209 int *idx2segpglist, *segpglist; 01210 01211 // List used in Delaunay refinement for keeping tetrahedra which need to 01212 // be classified by the quality measure. 01213 list *qualchecktetlist; 01214 01215 // Array used in Delaunay refinement for keeping the protecting sphere 01216 // radii of the input vertices. 01217 REAL *rpsarray; 01218 01219 // Queues that maintain the bad (badly-shaped or too large) tetrahedra. 01220 // The tails are pointers to the pointers that have to be filled in to 01221 // enqueue an item. The queues are ordered from 63 (highest priority) 01222 // to 0 (lowest priority). 01223 badface *subquefront[2], **subquetail[2]; 01224 badface *tetquefront[64], **tetquetail[64]; 01225 01226 // Pointer to a recently visited tetrahedron. Improves point location 01227 // if proximate points are inserted sequentially. 01228 triface recenttet; 01229 01230 REAL xmax, xmin, ymax, ymin, zmax, zmin; // Bounding box of points. 01231 REAL longest; // The longest possible edge length. 01232 long hullsize; // Number of faces of convex hull. 01233 long insegment; // Number of input segments. 01234 int steinerleft; // Number of Steiner points not yet used. 01235 int edgeboundindex; // Index to find edge bound of a point. 01236 int pointmarkindex; // Index to find boundary marker of a point. 01237 int point2simindex; // Index to find a simplex adjacent to a point. 01238 int point2pbcptindex; // Index to find a pbc point to a point. 01239 int highorderindex; // Index to find extra nodes for highorder elements. 01240 int elemattribindex; // Index to find attributes of a tetrahedron. 01241 int volumeboundindex; // Index to find volume bound of a tetrahedron. 01242 int shmarkindex; // Index to find boundary marker of a subface. 01243 int areaboundindex; // Index to find area bound of a subface. 01244 int checksubfaces; // Are there subfaces in the mesh yet? 01245 int checkpbcs; // Are there periodic boundary conditions? 01246 int varconstraint; // Are there variant (node, seg, facet) constraints? 01247 int nonconvex; // Is current mesh non-convex? 01248 int dupverts; // Are there duplicated vertices? 01249 int unuverts; // Are there unused vertices? 01250 int relverts; // The number of relocated vertices. 01251 int jettisoninverts; // The number of jettisoned input vertices. 01252 int symbolic; // Use symbolic insphere test. 01253 long samples; // Number of random samples for point location. 01254 unsigned long randomseed; // Current random number seed. 01255 REAL macheps; // The machine epsilon. 01256 REAL cosmaxdihed, cosmindihed; // The cosine values of max/min dihedral. 01257 int maxcavfaces, maxcavverts; // The size of the largest cavity. 01258 int expcavcount; // The times of expanding cavitys. 01259 long abovecount; // Number of abovepoints calculation. 01260 long rejectsmalledgecount; // Number of rejections of small edges. 01261 long bowatvolcount, bowatsubcount, bowatsegcount; // Bowyer-Watsons. 01262 long r1count, r2count, r3count, r4count; // Number of rules performed. 01263 long flip23s, flip32s, flip22s, flip44s; // Number of flips performed. 01264 01266 // // 01267 // Fast lookup tables for mesh manipulation primitives. // 01268 // // 01269 // Mesh manipulation primitives (given below) are basic operations on mesh // 01270 // data structures. They answer basic queries on mesh handles, such as "what // 01271 // is the origin (or destination, or apex) of the face?", "what is the next // 01272 // (or previous) edge in the edge ring?", and "what is the next face in the // 01273 // face ring?", and so on. // 01274 // // 01275 // The implementation of teste basic queries can take advangtage of the fact // 01276 // that the mesh data structures additionally store geometric informations. // 01277 // For example, we have ordered the 4 vertices (from 0 to 3) and the 4 faces // 01278 // (from 0 to 3) of a tetrahedron, and for each face of the tetrahedron, a // 01279 // sequence of vertices has stipulated, therefore the origin of any face of // 01280 // the tetrahedron can be quickly determined by a table 'locver2org', which // 01281 // takes the index of the face and the edge version as inputs. A list of // 01282 // fast lookup tables are defined below. They're just like global variables. // 01283 // These tables are initialized at the runtime. // 01284 // // 01286 01287 // For enext() primitive, uses 'ver' as the index. 01288 static int ve[6]; 01289 01290 // For org(), dest() and apex() primitives, uses 'ver' as the index. 01291 static int vo[6], vd[6], va[6]; 01292 01293 // For org(), dest() and apex() primitives, uses 'loc' as the first 01294 // index and 'ver' as the second index. 01295 static int locver2org[4][6]; 01296 static int locver2dest[4][6]; 01297 static int locver2apex[4][6]; 01298 01299 // For oppo() primitives, uses 'loc' as the index. 01300 static int loc2oppo[4]; 01301 01302 // For fnext() primitives, uses 'loc' as the first index and 'ver' as 01303 // the second index, returns an array containing a new 'loc' and a 01304 // new 'ver'. Note: Only valid for 'ver' equals one of {0, 2, 4}. 01305 static int locver2nextf[4][6][2]; 01306 01307 // For enumerating three edges of a triangle. 01308 static int plus1mod3[3]; 01309 static int minus1mod3[3]; 01310 01312 // // 01313 // Mesh manipulation primitives // 01314 // // 01315 // A serial of mesh operations such as topological maintenance, navigation, // 01316 // local modification, etc., is accomplished through a set of mesh manipul- // 01317 // ation primitives. These primitives are indeed very simple functions which // 01318 // take one or two handles ('triface's and 'face's) as parameters, perform // 01319 // basic operations such as "glue two tetrahedra at a face", "return the // 01320 // origin of a tetrahedron", "return the subface adjoining at the face of a // 01321 // tetrahedron", and so on. // 01322 // // 01324 01325 // Primitives for tetrahedra. 01326 inline void decode(tetrahedron ptr, triface& t); 01327 inline tetrahedron encode(triface& t); 01328 inline void sym(triface& t1, triface& t2); 01329 inline void symself(triface& t); 01330 inline void bond(triface& t1, triface& t2); 01331 inline void dissolve(triface& t); 01332 inline point org(triface& t); 01333 inline point dest(triface& t); 01334 inline point apex(triface& t); 01335 inline point oppo(triface& t); 01336 inline void setorg(triface& t, point pointptr); 01337 inline void setdest(triface& t, point pointptr); 01338 inline void setapex(triface& t, point pointptr); 01339 inline void setoppo(triface& t, point pointptr); 01340 inline void esym(triface& t1, triface& t2); 01341 inline void esymself(triface& t); 01342 inline void enext(triface& t1, triface& t2); 01343 inline void enextself(triface& t); 01344 inline void enext2(triface& t1, triface& t2); 01345 inline void enext2self(triface& t); 01346 inline bool fnext(triface& t1, triface& t2); 01347 inline bool fnextself(triface& t); 01348 inline void enextfnext(triface& t1, triface& t2); 01349 inline void enextfnextself(triface& t); 01350 inline void enext2fnext(triface& t1, triface& t2); 01351 inline void enext2fnextself(triface& t); 01352 inline void infect(triface& t); 01353 inline void uninfect(triface& t); 01354 inline bool infected(triface& t); 01355 inline REAL elemattribute(tetrahedron* ptr, int attnum); 01356 inline void setelemattribute(tetrahedron* ptr, int attnum, REAL value); 01357 inline REAL volumebound(tetrahedron* ptr); 01358 inline void setvolumebound(tetrahedron* ptr, REAL value); 01359 01360 // Primitives for subfaces and subsegments. 01361 inline void sdecode(shellface sptr, face& s); 01362 inline shellface sencode(face& s); 01363 inline void spivot(face& s1, face& s2); 01364 inline void spivotself(face& s); 01365 inline void sbond(face& s1, face& s2); 01366 inline void sbond1(face& s1, face& s2); 01367 inline void sdissolve(face& s); 01368 inline point sorg(face& s); 01369 inline point sdest(face& s); 01370 inline point sapex(face& s); 01371 inline void setsorg(face& s, point pointptr); 01372 inline void setsdest(face& s, point pointptr); 01373 inline void setsapex(face& s, point pointptr); 01374 inline void sesym(face& s1, face& s2); 01375 inline void sesymself(face& s); 01376 inline void senext(face& s1, face& s2); 01377 inline void senextself(face& s); 01378 inline void senext2(face& s1, face& s2); 01379 inline void senext2self(face& s); 01380 inline void sfnext(face&, face&); 01381 inline void sfnextself(face&); 01382 inline badface* shell2badface(face& s); 01383 inline void setshell2badface(face& s, badface* value); 01384 inline REAL areabound(face& s); 01385 inline void setareabound(face& s, REAL value); 01386 inline int shellmark(face& s); 01387 inline void setshellmark(face& s, int value); 01388 inline enum shestype shelltype(face& s); 01389 inline void setshelltype(face& s, enum shestype value); 01390 inline int shellpbcgroup(face& s); 01391 inline void setshellpbcgroup(face& s, int value); 01392 inline void sinfect(face& s); 01393 inline void suninfect(face& s); 01394 inline bool sinfected(face& s); 01395 01396 // Primitives for interacting tetrahedra and subfaces. 01397 inline void tspivot(triface& t, face& s); 01398 inline void stpivot(face& s, triface& t); 01399 inline void tsbond(triface& t, face& s); 01400 inline void tsdissolve(triface& t); 01401 inline void stdissolve(face& s); 01402 01403 // Primitives for interacting subfaces and subsegs. 01404 inline void sspivot(face& s, face& edge); 01405 inline void ssbond(face& s, face& edge); 01406 inline void ssdissolve(face& s); 01407 01408 // Primitives for points. 01409 inline int pointmark(point pt); 01410 inline void setpointmark(point pt, int value); 01411 inline enum verttype pointtype(point pt); 01412 inline void setpointtype(point pt, enum verttype value); 01413 inline tetrahedron point2tet(point pt); 01414 inline void setpoint2tet(point pt, tetrahedron value); 01415 inline shellface point2sh(point pt); 01416 inline void setpoint2sh(point pt, shellface value); 01417 inline point point2ppt(point pt); 01418 inline void setpoint2ppt(point pt, point value); 01419 inline point point2pbcpt(point pt); 01420 inline void setpoint2pbcpt(point pt, point value); 01421 inline REAL edgebound(point pt); 01422 inline void setedgebound(point pt, REAL value); 01423 01424 // Advanced primitives. 01425 inline void adjustedgering(triface& t, int direction); 01426 inline void adjustedgering(face& s, int direction); 01427 inline bool isdead(triface* t); 01428 inline bool isdead(face* s); 01429 inline bool isfacehaspoint(triface* t, point testpoint); 01430 inline bool isfacehaspoint(face* t, point testpoint); 01431 inline bool isfacehasedge(face* s, point tend1, point tend2); 01432 inline bool issymexist(triface* t); 01433 bool getnextface(triface*, triface*); 01434 void getnextsface(face*, face*); 01435 void tsspivot(triface*, face*); 01436 void sstpivot(face*, triface*); 01437 bool findorg(triface* t, point dorg); 01438 bool findorg(face* s, point dorg); 01439 void findedge(triface* t, point eorg, point edest); 01440 void findedge(face* s, point eorg, point edest); 01441 void findface(triface *fface, point forg, point fdest, point fapex); 01442 void getonextseg(face* s, face* lseg); 01443 void getseghasorg(face* sseg, point dorg); 01444 point getsubsegfarorg(face* sseg); 01445 point getsubsegfardest(face* sseg); 01446 void printtet(triface*); 01447 void printsh(face*); 01448 01450 // // 01451 // Triangle-triangle intersection test // 01452 // // 01453 // The triangle-triangle intersection test is implemented with exact arithm- // 01454 // etic. It exactly tells whether or not two triangles in three dimensions // 01455 // intersect. Before implementing this test myself, I tried two C codes // 01456 // (implemented by Thomas Moeller and Philippe Guigue, respectively), which // 01457 // are all public available. However both of them failed frequently. Another // 01458 // unconvenience is both codes only tell whether or not the two triangles // 01459 // intersect without distinguishing the cases whether they exactly intersect // 01460 // in interior or they just share a vertex or share an edge. The two latter // 01461 // cases are acceptable and should return not intersection in TetGen. // 01462 // // 01464 01465 enum interresult edge_vert_col_inter(REAL*, REAL*, REAL*); 01466 enum interresult edge_edge_cop_inter(REAL*, REAL*, REAL*, REAL*, REAL*); 01467 enum interresult tri_vert_cop_inter(REAL*, REAL*, REAL*, REAL*, REAL*); 01468 enum interresult tri_edge_cop_inter(REAL*, REAL*, REAL*, REAL*, REAL*, 01469 REAL*); 01470 enum interresult tri_edge_inter_tail(REAL*, REAL*, REAL*, REAL*, REAL*, 01471 REAL, REAL); 01472 enum interresult tri_edge_inter(REAL*, REAL*, REAL*, REAL*, REAL*); 01473 enum interresult tri_tri_inter(REAL*, REAL*, REAL*, REAL*, REAL*, REAL*); 01474 01475 // Geometric predicates 01476 REAL insphere_sos(REAL*, REAL*, REAL*, REAL*, REAL*, int, int, int, int, 01477 int); 01478 bool iscollinear(REAL*, REAL*, REAL*, REAL eps); 01479 bool iscoplanar(REAL*, REAL*, REAL*, REAL*, REAL vol6, REAL eps); 01480 bool iscospheric(REAL*, REAL*, REAL*, REAL*, REAL*, REAL vol24, REAL eps); 01481 01482 // Linear algebra functions 01483 inline REAL dot(REAL* v1, REAL* v2); 01484 inline void cross(REAL* v1, REAL* v2, REAL* n); 01485 bool lu_decmp(REAL lu[4][4], int n, int* ps, REAL* d, int N); 01486 void lu_solve(REAL lu[4][4], int n, int* ps, REAL* b, int N); 01487 01488 // Geometric quantities calculators. 01489 inline REAL distance(REAL* p1, REAL* p2); 01490 REAL shortdistance(REAL* p, REAL* e1, REAL* e2); 01491 REAL shortdistance(REAL* p, REAL* e1, REAL* e2, REAL* e3); 01492 REAL interiorangle(REAL* o, REAL* p1, REAL* p2, REAL* n); 01493 void projpt2edge(REAL* p, REAL* e1, REAL* e2, REAL* prj); 01494 void projpt2face(REAL* p, REAL* f1, REAL* f2, REAL* f3, REAL* prj); 01495 void facenormal(REAL* pa, REAL* pb, REAL* pc, REAL* n, REAL* nlen); 01496 void edgeorthonormal(REAL* e1, REAL* e2, REAL* op, REAL* n); 01497 REAL facedihedral(REAL* pa, REAL* pb, REAL* pc1, REAL* pc2); 01498 void tetalldihedral(point, point, point, point, REAL dihed[6]); 01499 void tetallnormal(point, point, point, point, REAL N[4][3], REAL* volume); 01500 bool circumsphere(REAL*, REAL*, REAL*, REAL*, REAL* cent, REAL* radius); 01501 void inscribedsphere(REAL*, REAL*, REAL*, REAL*, REAL* cent, REAL* radius); 01502 void rotatepoint(REAL* p, REAL rotangle, REAL* p1, REAL* p2); 01503 void spherelineint(REAL* p1, REAL* p2, REAL* C, REAL R, REAL p[7]); 01504 void linelineint(REAL *p1,REAL *p2, REAL *p3, REAL *p4, REAL p[7]); 01505 void planelineint(REAL*, REAL*, REAL*, REAL*, REAL*, REAL*, REAL*); 01506 01507 // Memory managment routines. 01508 void dummyinit(int, int); 01509 void initializepools(); 01510 void tetrahedrondealloc(tetrahedron*); 01511 tetrahedron *tetrahedrontraverse(); 01512 void shellfacedealloc(memorypool*, shellface*); 01513 shellface *shellfacetraverse(memorypool*); 01514 void badfacedealloc(memorypool*, badface*); 01515 badface *badfacetraverse(memorypool*); 01516 void pointdealloc(point); 01517 point pointtraverse(); 01518 void maketetrahedron(triface*); 01519 void makeshellface(memorypool*, face*); 01520 void makepoint(point*); 01521 01522 // Mesh items searching routines. 01523 void makepoint2tetmap(); 01524 void makeindex2pointmap(point*& idx2verlist); 01525 void makesegmentmap(int*& idx2seglist, shellface**& segsperverlist); 01526 void makesubfacemap(int*& idx2facelist, shellface**& facesperverlist); 01527 void maketetrahedronmap(int*& idx2tetlist, tetrahedron**& tetsperverlist); 01528 01529 // Point location routines. 01530 unsigned long randomnation(unsigned int choices); 01531 REAL distance2(tetrahedron* tetptr, point p); 01532 enum locateresult preciselocate(point searchpt, triface* searchtet, long); 01533 enum locateresult locate(point searchpt, triface* searchtet); 01534 enum locateresult adjustlocate(point, triface*, enum locateresult, REAL); 01535 enum locateresult locatesub(point searchpt, face* searchsh, int stopatseg); 01536 enum locateresult adjustlocatesub(point, face*, enum locateresult, REAL); 01537 enum locateresult locateseg(point searchpt, face* searchseg); 01538 enum locateresult adjustlocateseg(point, face*, enum locateresult, REAL); 01539 01541 // // 01542 // Mesh Local Transformation Operators // 01543 // // 01544 // These operators (including flips, insert & remove vertices and so on) are // 01545 // used to transform (or replace) a set of mesh elements into another set of // 01546 // mesh elements. // 01547 // // 01549 01550 // Mesh transformation routines. 01551 enum fliptype categorizeface(triface& horiz); 01552 void enqueueflipface(triface& checkface, queue* flipqueue); 01553 void enqueueflipedge(face& checkedge, queue* flipqueue); 01554 void flip23(triface* flipface, queue* flipqueue); 01555 void flip32(triface* flipface, queue* flipqueue); 01556 void flip22(triface* flipface, queue* flipqueue); 01557 void flip22sub(face* flipedge, queue* flipqueue); 01558 long flip(queue* flipqueue, badface **plastflip, bool, bool, bool); 01559 void undoflip(badface *lastflip); 01560 long flipsub(queue* flipqueue); 01561 01562 void splittetrahedron(point newpoint, triface* splittet, queue* flipqueue); 01563 void unsplittetrahedron(triface* splittet); 01564 void splittetface(point newpoint, triface* splittet, queue* flipqueue); 01565 void unsplittetface(triface* splittet); 01566 void splitsubface(point newpoint, face* splitface, queue* flipqueue); 01567 void unsplitsubface(face* splitsh); 01568 void splittetedge(point newpoint, triface* splittet, queue* flipqueue); 01569 void unsplittetedge(triface* splittet); 01570 void splitsubedge(point newpoint, face* splitsh, queue* flipqueue); 01571 void unsplitsubedge(face* splitsh); 01572 enum insertsiteresult insertsite(point newpoint, triface* searchtet, 01573 bool approx, queue* flipqueue); 01574 void undosite(enum insertsiteresult insresult, triface* splittet, 01575 point torg, point tdest, point tapex, point toppo); 01576 01577 // Delaunay tetrahedralization routines. 01578 void formstarpolyhedron(point pt, list* tetlist, list* verlist); 01579 void collectcavtets(point newpoint, list* cavtetlist); 01580 bool unifypoint(point testpt, triface*, enum locateresult, REAL); 01581 void closeopenface(triface* openface, queue* flipque); 01582 void inserthullsite(point inspoint, triface* horiz, queue* flipque); 01583 void incrflipdelaunay(triface*, point*, long, bool, bool, REAL, queue*); 01584 long delaunizevertices(); 01585 01586 // Surface triangulation routines. 01587 void formstarpolygon(point pt, list* trilist, list* verlist); 01588 void getfacetabovepoint(face* facetsh); 01589 void collectcavsubs(point newpoint, list* cavsublist); 01590 void collectvisiblesubs(int shmark, point inspoint, face* horiz, queue*); 01591 void incrflipdelaunaysub(int shmark, REAL eps, list*, int, REAL*, queue*); 01592 enum finddirectionresult finddirectionsub(face* searchsh, point tend); 01593 void insertsubseg(face* tri); 01594 bool scoutsegmentsub(face* searchsh, point tend); 01595 void flipedgerecursive(face* flipedge, queue* flipqueue); 01596 void constrainededge(face* startsh, point tend, queue* flipqueue); 01597 void recoversegment(point tstart, point tend, queue* flipqueue); 01598 void infecthullsub(memorypool* viri); 01599 void plaguesub(memorypool* viri); 01600 void carveholessub(int holes, REAL* holelist, memorypool* viri); 01601 void triangulate(int shmark, REAL eps, list* ptlist, list* conlist, 01602 int holes, REAL* holelist, memorypool* viri, queue*); 01603 void retrievefacetsubs(list* newshlist, list* boundedgelist); 01604 void unifysegments(); 01605 void mergefacets(queue* flipqueue); 01606 void removefacetvertex(face* remsh, list* oldshlist, list* ptlist, 01607 list* conlist, list* newshlist, list* boundedgelist, 01608 memorypool* viri, queue* flipque); 01609 void removesegmentvertex(face* remseg, face* nremseg, face* newseg, 01610 list* spinshlist, list* ptlist, list* conlist, 01611 list* newshlist, list* boundedgelist, 01612 list* newsegshlist, list* newshlistlist, 01613 list* oldshlistlist, memorypool*, queue*); 01614 void removevolumevertex(point rempt, list* ptlist, list* newshlist, 01615 list* oldtetlist, list* newtetlist, 01616 list* frontlist, list* misfrontlist, queue*); 01617 void removefreepoints(list*, list*, memorypool*, queue*); 01618 long meshsurface(); 01619 01620 // Detect intersecting facets of PLC. 01621 void interecursive(shellface** subfacearray, int arraysize, int axis, 01622 REAL bxmin, REAL bxmax, REAL bymin, REAL bymax, 01623 REAL bzmin, REAL bzmax, int* internum); 01624 void detectinterfaces(); 01625 01626 // Periodic boundary condition supporting routines. 01627 void createsubpbcgrouptable(); 01628 void getsubpbcgroup(face* pbcsub, pbcdata** pd, int *f1, int *f2); 01629 enum locateresult getsubpbcsympoint(point, face*, point, face*); 01630 void createsegpbcgrouptable(); 01631 enum locateresult getsegpbcsympoint(point, face*, point, face*, int); 01632 01633 // Vertex perturbation routines. 01634 REAL randgenerator(REAL range); 01635 bool checksub4cocir(face* testsub, REAL eps, bool once, bool enqflag); 01636 bool checktet4cosph(triface* testtet, REAL eps, bool once, bool enqflag); 01637 void tallcocirsubs(REAL eps, bool enqflag); 01638 void tallcosphtets(list* testtetlist, REAL eps, bool enqflag); 01639 bool tallencsegsfsubs(point testpt, list* cavsublist); 01640 bool tallencsubsfsubs(point testpt, list* cavtetlist); 01641 void collectflipedges(point inspoint, face* splitseg, queue* flipqueue); 01642 void perturbrepairencsegs(queue* flipqueue); 01643 void perturbrepairencsubs(REAL eps, list* cavsublist, queue* flipqueue); 01644 void perturbrepairbadtets(REAL eps, list* cavsublist, list*, queue*); 01645 void incrperturbvertices(REAL eps); 01646 01647 // Segment recovery routines. 01648 void markacutevertices(REAL acuteangle); 01649 void initializerpsarray(); 01650 enum finddirectionresult finddirection(triface* searchtet, point, long); 01651 void getsearchtet(point p1, point p2, triface* searchtet, point* tend); 01652 bool isedgeencroached(point p1, point p2, point testpt, bool degflag); 01653 point scoutrefpoint(triface* searchtet, point tend); 01654 point getsegmentorigin(face* splitseg); 01655 point getsplitpoint(face* splitseg, point refpoint); 01656 void delaunizesegments(); 01657 01658 // Facets recovery routines. 01659 bool insertsubface(face* insertsh, triface* searchtet); 01660 bool tritritest(triface* checktet, point p1, point p2, point p3); 01661 void initializecavity(list* floorlist, list* ceillist, list* frontlist); 01662 void delaunizecavvertices(triface*, list*, list*, list*, queue*); 01663 void retrievenewtets(list* newtetlist); 01664 void insertauxsubface(triface* front, triface* idfront); 01665 bool scoutfront(triface* front, triface* idfront, list* newtetlist); 01666 void gluefronts(triface* front, triface* front1); 01667 bool identifyfronts(list* frontlist, list* misfrontlist, list* newtetlist); 01668 void detachauxsubfaces(list* newtetlist); 01669 void expandcavity(list* frontlist, list* misfrontlist, list* newtetlist, 01670 list* crosstetlist, queue* missingshqueue, queue*); 01671 void carvecavity(list* newtetlist, list* outtetlist, queue* flipque); 01672 void delaunizecavity(list* floorlist, list* ceillist, list* ceilptlist, 01673 list* floorptlist, list* frontlist,list* misfrontlist, 01674 list* newtetlist, list* crosstetlist, queue*, queue*); 01675 void formmissingregion(face* missingsh, list* missingshlist, 01676 list* equatptlist, int* worklist); 01677 void formcavity(list* missingshlist, list* crossedgelist, 01678 list* equatptlist, list* crossshlist, list* crosstetlist, 01679 list* belowfacelist, list* abovefacelist, 01680 list* horizptlist, list* belowptlist, list* aboveptlist, 01681 queue* missingshqueue, int* worklist); 01682 bool scoutcrossingedge(list* missingshlist, list* boundedgelist, 01683 list* crossedgelist, int* worklist); 01684 void rearrangesubfaces(list* missingshlist, list* boundedgelist, 01685 list* equatptlist, int* worklist); 01686 void insertallsubfaces(queue* missingshqueue); 01687 void constrainedfacets(); 01688 01689 // Carving out holes and concavities routines. 01690 void infecthull(memorypool *viri); 01691 void plague(memorypool *viri); 01692 void regionplague(memorypool *viri, REAL attribute, REAL volume); 01693 void removeholetets(memorypool *viri); 01694 void assignregionattribs(); 01695 void carveholes(); 01696 01697 // Steiner points removing routines. 01698 void fillsteinerpolygon(list* oldshlist, list* boundedgelist); 01699 void orientnewsubfaces(list* newshlist, face* orientsh, REAL* norm); 01700 bool constrainedflip(triface* flipface, triface* front, queue* flipque); 01701 bool recoverfront(triface* front, list* newtetlist, queue* flipque); 01702 void repairflips(queue* flipque); 01703 bool constrainedcavity(triface* oldtet, list* floorlist, list* ceillist, 01704 list* ptlist, list* frontlist, list* misfrontlist, 01705 list* newtetlist, queue* flipque); 01706 void expandsteinercavity(point steinpt, REAL eps, list* frontlist, list*); 01707 bool relocatesteiner(point sp, point np, REAL* n, list* frontlist, list*); 01708 void bowatinsertsite(point steinpt, triface* oldtet, list*, list*, queue*); 01709 void removesubsteiner(face* steinsh, list* oldshlist, list* ptlist, 01710 list* conlist, list* newshlist, list* boundedgelist, 01711 list* oldtetlist, list* newtetlist, list* frontlist, 01712 list*, memorypool*, queue*); 01713 void removesegsteiner(face* steinseg, list* spinshlist, list* ptlist, 01714 list* conlist, list* newshlist, list* boundedgelist, 01715 list* newsegshlist, list* newshlistlist, 01716 list* oldshlistlist, list* oldtetlist, 01717 list* newtetlist, list* frontlist, 01718 list* misfrontlist, memorypool *viri, queue*); 01719 void removesteiners(); 01720 01721 // Mesh reconstruction rotuines. 01722 void assignvarconstraints(point *idx2verlist); 01723 long reconstructmesh(); 01724 void insertaddpoints(); 01725 01726 // Mesh repair routines. 01727 bool checktet4dege(triface* testtet, REAL eps, list* degetetlist); 01728 bool finddiagonal(triface* akite); 01729 void removetetbypeeloff(triface *badtet, queue* flipqueue); 01730 void removetetbyflip32(triface *badtet, queue* flipqueue); 01731 bool removetetbyedgereduce(triface *badtet, queue* flipqueue); 01732 bool removekite(triface* akite, bool reduce, queue* flipqueue); 01733 void talldegetets(REAL eps, list* newtetlist, list* degetetlist); 01734 void repairdegetets(REAL eps, list* newtetlist, queue* flipqueue); 01735 01736 // Delaunay refinement routines. 01737 void marksharpsubfaces(REAL dihedbound); 01738 void markskinnysubfaces(REAL anglebound); 01739 badface* dequeuebadtet(); 01740 bool checkseg4encroach(face* testseg, point testpt, point*, bool enqflag); 01741 bool checksub4encroach(face* testsub, point testpt, bool enqflag); 01742 bool checkseg4badqual(face* testseg, bool enqflag); 01743 bool checksub4badqual(face* testsub, bool enqflag); 01744 bool checktet4badqual(triface* testtet, bool enqflag); 01745 bool checktet4sliver(triface* testtet, bool enqflag); 01746 bool checkseg4splitting(face* testseg, point* pencpt); 01747 bool checksub4splitting(face* testsub, bool bqual); 01748 bool tallencsegs(point testpt, list *cavtetlist, list* cavsublist); 01749 bool tallencsubs(point testpt, list *cavtetlist, list* cavsublist); 01750 void tallencsubsfseg(face* testseg); 01751 void tallbadsegs(); 01752 void tallbadsubs(); 01753 void tallbadtetrahedrons(); 01754 void tallslivers(); 01755 void doqualchecktetlist(); 01756 void getsplitpoint1(face* splitseg, REAL* splitpt, point encpt); 01757 void collectvolcavtets(point, list*, list*); 01758 void collectsubcavtets(point, face*, list*, list*); 01759 void collectsegcavtets(point, face*, list*, list*, list*); 01760 void collectbowatcavsubs(face*, list*, list*); 01761 void formbowatcavity(point newpt, list* cavtetlist, list* ceillist, 01762 list* cutceillist); 01763 void formbowatsubcavity(point newpt, list* spinshlist, list** cavsublists, 01764 list** subceillists, bool cutchk); 01765 void bowatinsertsegsite(point newpt, face* splitseg, list* spinshlist); 01766 void bowatinsertsubsite(point newpt, face* splitseg, face* splitsegsh, 01767 list* cavsublist, list* subceillist, 01768 list* ceillist, bool chksubs); 01769 void bowatinsertvolsite(point newpt, face* splitseg, list* cavtetlist, 01770 list* ceillist, list* spinshlist, 01771 list** cavsublists, list** subceillists, 01772 bool chksubs, bool chkbqual, bool chksliver); 01773 void repairencsegs(bool chksubs, bool chkbqual); 01774 void repairencsegs(list* cavtetlist, list* ceillist, list* cutceillist, 01775 list* spinshlist, list* quadtetlist, bool chksubs, 01776 bool chkbqual); 01777 void repairencsubs(list* cavtetlist, list* ceillist, list* cutceillist, 01778 list* spinshlist, list* quadtetlist, bool chksubs); 01779 void repairbadtets(list* cavtetlist, list* ceillist, list* cutceillist, 01780 list* spinshlist, list* quadtetlist); 01781 void repairslivers(list* cavtetlist, list* ceillist, list* cutceillist, 01782 list* spinshlist); 01783 void enforcequality(); 01784 01785 // I/O routines 01786 void transfernodes(); 01787 void jettisonnodes(); 01788 void highorder(); 01789 void outnodes(tetgenio* out); 01790 void outelements(tetgenio* out); 01791 void outfaces(tetgenio* out); 01792 void outhullfaces(tetgenio* out); 01793 void outsubfaces(tetgenio* out); 01794 void outsubsegments(tetgenio* out); 01795 void outneighbors(tetgenio* out); 01796 void outpbcnodes(tetgenio* out); 01797 void outsmesh(char* smfilename); 01798 void outmesh2medit(char* mfilename); 01799 void outmesh2gid(char* gfilename); 01800 void outmesh2off(char* ofilename); 01801 01802 // User interaction routines. 01803 void internalerror(); 01804 void checkmesh(); 01805 void checkshells(); 01806 void checkdelaunay(REAL eps, queue* flipqueue); 01807 void checkdegeneracy(REAL eps); 01808 void checkconforming(); 01809 void algorithmicstatistics(); 01810 void qualitystatistics(); 01811 void statistics(); 01812 01813 public: 01814 01815 // Constructor and destructor. 01816 tetgenmesh(); 01817 ~tetgenmesh(); 01818 01819 }; // End of class tetgenmesh. 01820 01822 // // 01823 // tetrahedralize() Interface for using TetGen's library to generate // 01824 // Delaunay tetrahedralizations, constrained Delaunay // 01825 // tetrahedralizations, quality tetrahedral meshes. // 01826 // // 01827 // Two functions (interfaces) are available. The difference is only the way // 01828 // of passing switches. One directly accepts an object of 'tetgenbehavior', // 01829 // while the other accepts a string which is the same as one can used in the // 01830 // command line. The latter may be more convenient for users who don't want // 01831 // to kown the 'tetgenbehavir' structure. // 01832 // // 01833 // 'in' is an object of 'tetgenio' which contains a PLC you want to tetrahed-// 01834 // ralize or a previously generated tetrahedral mesh you want to refine. It // 01835 // must not be a NULL. 'out' is another object of 'tetgenio' for storing the // 01836 // generated tetrahedral mesh. It can be a NULL. If so, the output will be // 01837 // saved to file(s). // 01838 // // 01840 01841 void tetrahedralize(tetgenbehavior *b, tetgenio *in, tetgenio *out); 01842 void tetrahedralize(char *switches, tetgenio *in, tetgenio *out); 01843 01844 #endif // #ifndef tetgenH
1.7.1