zhangyiss/tetgen_cmake
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

update

Browse Source
This commit is contained in:
张壹
2023-12-24 10:13:04 +08:00
parent 24a2adda64
commit 554439818b
11 changed files with 45855 additions and 54 deletions
Download Patch File Download Diff File Expand all files Collapse all files

355
src/lib/tetgen.cxx
View File

@@ -2503,7 +2503,9 @@ void tetgenio::save_nodes(const char *filebasename)
char outmtrfilename[FILENAMESIZE];
int i, j;
sprintf(outnodefilename, "%s.node", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outnodefilename, "%s.node", filebasename);
snprintf(outnodefilename, FILENAMESIZE, "%s.node", filebasename);
printf("Saving nodes to %s\n", outnodefilename);
fout = fopen(outnodefilename, "w");
fprintf(fout, "%d %d %d %d\n", numberofpoints, mesh_dim,
@@ -2529,7 +2531,9 @@ void tetgenio::save_nodes(const char *filebasename)
// If the point metrics exist, output them to a .mtr file.
if ((numberofpointmtrs > 0) && (pointmtrlist != (REAL *) NULL)) {
sprintf(outmtrfilename, "%s.mtr", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outmtrfilename, "%s.mtr", filebasename);
snprintf(outmtrfilename, FILENAMESIZE, "%s.mtr", filebasename);
printf("Saving metrics to %s\n", outmtrfilename);
fout = fopen(outmtrfilename, "w");
fprintf(fout, "%d %d\n", numberofpoints, numberofpointmtrs);
@@ -2555,7 +2559,9 @@ void tetgenio::save_elements(const char* filebasename)
char outelefilename[FILENAMESIZE];
int i, j;
sprintf(outelefilename, "%s.ele", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outelefilename, "%s.ele", filebasename);
snprintf(outelefilename, FILENAMESIZE, "%s.ele", filebasename);
printf("Saving elements to %s\n", outelefilename);
fout = fopen(outelefilename, "w");
if (mesh_dim == 3) {
@@ -2602,7 +2608,9 @@ void tetgenio::save_faces(const char* filebasename)
char outfacefilename[FILENAMESIZE];
int i;
sprintf(outfacefilename, "%s.face", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outfacefilename, "%s.face", filebasename);
snprintf(outfacefilename, FILENAMESIZE, "%s.face", filebasename);
printf("Saving faces to %s\n", outfacefilename);
fout = fopen(outfacefilename, "w");
fprintf(fout, "%d %d\n", numberoftrifaces,
@@ -2631,7 +2639,9 @@ void tetgenio::save_edges(char* filebasename)
char outedgefilename[FILENAMESIZE];
int i;
sprintf(outedgefilename, "%s.edge", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outedgefilename, "%s.edge", filebasename);
snprintf(outedgefilename, FILENAMESIZE, "%s.edge", filebasename);
printf("Saving edges to %s\n", outedgefilename);
fout = fopen(outedgefilename, "w");
fprintf(fout, "%d %d\n", numberofedges, edgemarkerlist != NULL ? 1 : 0);
@@ -2659,7 +2669,9 @@ void tetgenio::save_neighbors(char* filebasename)
char outneighborfilename[FILENAMESIZE];
int i;
sprintf(outneighborfilename, "%s.neigh", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outneighborfilename, "%s.neigh", filebasename);
snprintf(outneighborfilename, FILENAMESIZE, "%s.neigh", filebasename);
printf("Saving neighbors to %s\n", outneighborfilename);
fout = fopen(outneighborfilename, "w");
fprintf(fout, "%d %d\n", numberoftetrahedra, mesh_dim + 1);
@@ -2694,7 +2706,9 @@ void tetgenio::save_poly(const char *filebasename)
char outpolyfilename[FILENAMESIZE];
int i, j, k;
sprintf(outpolyfilename, "%s.poly", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outpolyfilename, "%s.poly", filebasename);
snprintf(outpolyfilename, FILENAMESIZE, "%s.poly", filebasename);
printf("Saving poly to %s\n", outpolyfilename);
fout = fopen(outpolyfilename, "w");
@@ -2792,7 +2806,9 @@ void tetgenio::save_faces2smesh(char* filebasename)
char outsmeshfilename[FILENAMESIZE];
int i, j;
sprintf(outsmeshfilename, "%s.smesh", filebasename);
// modified by Yi Zhang on 2023-12-22
//sprintf(outsmeshfilename, "%s.smesh", filebasename);
snprintf(outsmeshfilename, FILENAMESIZE, "%s.smesh", filebasename);
printf("Saving faces to %s\n", outsmeshfilename);
fout = fopen(outsmeshfilename, "w");
@@ -2824,6 +2840,93 @@ void tetgenio::save_faces2smesh(char* filebasename)
fclose(fout);
}
//============================================================================//
// //
// save_mesh2gmsh() Save mesh to a .msh file. //
// //
//============================================================================//
void tetgenio::save_mesh2gmsh(const char *filebasename)
{
FILE *fout;
char outnodefilename[FILENAMESIZE];
char outmtrfilename[FILENAMESIZE];
int i, j;
if (numberofcorners != 4)
{
printf(" Write Gmsh not implemented for order 2 elements \n");
return;
}
// modified by Yi Zhang on 2023-12-22
//sprintf(outnodefilename, "%s.node", filebasename);
snprintf(outnodefilename, FILENAMESIZE, "%s.msh", filebasename);
printf("Saving mesh to %s\n", outnodefilename);
fout = fopen(outnodefilename, "w");
fprintf(fout, "$MeshFormat\n2.2 0 8\n$EndMeshFormat\n");
fprintf(fout, "$Nodes\n");
fprintf(fout, "%d\n", numberofpoints);
for (i = 0; i < numberofpoints; i++)
{
fprintf(fout, "%d %.16g %.16g %.16g\n", i + firstnumber,
pointlist[i * 3], pointlist[i * 3 + 1], pointlist[i * 3 + 2]);
}
fprintf(fout, "$EndNodes\n");
fprintf(fout, "$Elements\n");
fprintf(fout, "%d\n", numberoftetrahedra + numberoftrifaces);
int fac_tag = 0;
if (trifacemarkerlist != NULL) fac_tag = 1;
for (i = 0; i < numberoftrifaces; i++) {
fprintf(fout, "%d 2 %d", i + firstnumber, fac_tag);
if (trifacemarkerlist != NULL) {
fprintf(fout, " %d", trifacemarkerlist[i]);
}
fprintf(fout, " %d %d %d\n", trifacelist[i * 3], trifacelist[i * 3 + 1], trifacelist[i * 3 + 2]);
}
int st_number = numberoftrifaces + firstnumber;
for (i = 0; i < numberoftetrahedra; i++) {
fprintf(fout, "%d 4 %d", i + st_number, numberoftetrahedronattributes);
for (j = 0; j < numberoftetrahedronattributes; j++) {
fprintf(fout, " %g", tetrahedronattributelist[i * numberoftetrahedronattributes + j]);
}
for (j = 0; j < numberofcorners; j++) {
fprintf(fout, " %d", tetrahedronlist[i * numberofcorners + j]);
}
fprintf(fout, "\n");
}
fprintf(fout, "$EndElements\n");
if (trifacemarkerlist != NULL) {
fprintf(fout, "$ElementData\n1\n\"Face Tag\"\n1\n0\n3\n0\n1\n");
fprintf(fout, "%d\n", numberoftrifaces);
for (i = 0; i < numberoftrifaces; i++) {
fprintf(fout, "%d %d\n", i + firstnumber, trifacemarkerlist[i]);
}
fprintf(fout, "$EndElementData\n");
}
for (j = 0; j < numberoftetrahedronattributes; j++) {
//fprintf(fout, "$ElementData\n1\n\"Element Tag\"\n1\n0\n3\n0\n1\n");
fprintf(fout, "$ElementData\n1\n\"Element Tag");
fprintf(fout, " %d\"\n1\n0\n3\n0\n1\n", j+1);
fprintf(fout, "%d\n", numberoftetrahedra);
for (i = 0; i < numberoftetrahedra; i++) {
fprintf(fout, "%d %g\n", i + st_number, tetrahedronattributelist[i * numberoftetrahedronattributes + j]);
}
fprintf(fout, "$EndElementData\n");
}
fclose(fout);
}
//============================================================================//
// //
// readline() Read a nonempty line from a file. //
@@ -2964,7 +3067,7 @@ char* tetgenio::findnextnumber(char *string)
void tetgenbehavior::syntax()
{
printf(" tetgen [-pYrq_Aa_miO_S_T_XMwcdzfenvgkJBNEFICQVh] input_file\n");
printf(" tetgen [-pYrq_Aa_miO_S_T_XMwcdzfenvgGkJBNEFICQVh] input_file\n"); // Modified by Yi Zhang
printf(" -p Tetrahedralizes a piecewise linear complex (PLC).\n");
printf(" -Y Preserves the input surface mesh (does not modify it).\n");
printf(" -r Reconstructs a previously generated mesh.\n");
@@ -2987,6 +3090,7 @@ void tetgenbehavior::syntax()
printf(" -e Outputs all edges to .edge file.\n");
printf(" -n Outputs tetrahedra neighbors to .neigh file.\n");
printf(" -g Outputs mesh to .mesh file for viewing by Medit.\n");
printf(" -G Outputs mesh to .msh file for viewing by Gmsh.\n"); // Modified by Yi Zhang
printf(" -k Outputs mesh to .vtk file for viewing by Paraview.\n");
printf(" -J No jettison of unused vertices from output .node file.\n");
printf(" -B Suppresses output of boundary information.\n");
@@ -3011,10 +3115,10 @@ void tetgenbehavior::usage()
printf("TetGen\n");
printf("A Quality Tetrahedral Mesh Generator and 3D Delaunay ");
printf("Triangulator\n");
printf("Version 1.6\n");
printf("August, 2020\n");
printf("Version 1.6.0_1\n");
printf("December, 2023\n");
printf("\n");
printf("Copyright (C) 2002 - 2020\n");
printf("Copyright (C) 2023\n");
printf("\n");
printf("What Can TetGen Do?\n");
printf("\n");
@@ -3056,7 +3160,11 @@ void tetgenbehavior::usage()
printf("have volume\n of 0.1 or less, and writes the mesh to ");
printf("object.1.node, object.1.ele,\n object.1.face, and object.1.edge\n");
printf("\n");
printf("Please send bugs/comments to Hang Si <si@wias-berlin.de>\n");
printf("Please send bugs/comments to Hang Si <si@wias-berlin.de>\n\n");
printf("Updates:\n");
printf("(1) This version is modified to support outputs of Gmsh .msh file.\n");
printf("(2) Fixed a bug for outputing VTK files when the first node index is zero.\n");
printf("Please send corresponding bugs/comments to Yi Zhang <yizhang-geo@zju.edu.cn>\n");
terminatetetgen(NULL, 0);
}
@@ -3423,6 +3531,8 @@ bool tetgenbehavior::parse_commandline(int argc, char **argv)
neighout++;
} else if (argv[i][j] == 'g') {
meditview = 1;
} else if (argv[i][j] == 'G') { // modified by Yi Zhang
gmshview = 1; // modified by Yi Zhang
} else if (argv[i][j] == 'k') {
if (argv[i][j + 1] == '2') { // -k2
vtksurfview = 1;
@@ -3782,7 +3892,9 @@ bool tetgenbehavior::parse_commandline(int argc, char **argv)
workstring[increment] = '%';
workstring[increment + 1] = 'd';
workstring[increment + 2] = '\0';
sprintf(outfilename, workstring, meshnumber + 1);
// modified by Yi Zhang on 2023-12-22
//sprintf(outfilename, workstring, meshnumber + 1);
snprintf(outfilename, 1024, workstring, meshnumber + 1);
}
// Additional input file name has the end ".a".
strcpy(addinfilename, infilename);
@@ -32989,13 +33101,17 @@ void tetgenmesh::qualitystatistics()
shortest, longest);
printf(" Smallest asp.ratio: %13.5g | Largest asp.ratio: %13.5g\n",
smallestratio, biggestratio);
sprintf(sbuf, "%.17g", biggestfaangle);
// modified by Yi Zhang on 2023-12-22
//sprintf(sbuf, "%.17g", biggestfaangle);
snprintf(sbuf, 128, "%.17g", biggestfaangle);
if (strlen(sbuf) > 8) {
sbuf[8] = '\0';
}
printf(" Smallest facangle: %14.5g | Largest facangle: %s\n",
smallestfaangle, sbuf);
sprintf(sbuf, "%.17g", biggestdiangle);
// modified by Yi Zhang on 2023-12-22
//sprintf(sbuf, "%.17g", biggestdiangle);
snprintf(sbuf, 128, "%.17g", biggestdiangle);
if (strlen(sbuf) > 8) {
sbuf[8] = '\0';
}
@@ -35790,6 +35906,181 @@ void tetgenmesh::outmesh2medit(char* mfilename)
//============================================================================//
// //
// outmesh2gmsh() Write mesh to a .msh file, which can be read and //
// rendered by Gmsh (an open source 3D finite element mesh //
// generator with a built-in CAD engine and post-processor).//
// //
// You can specify a filename (without suffix) in 'mfilename'. If you don't //
// supply a filename (let mfilename be NULL), the default name stored in //
// 'tetgenbehavior' will be used. The output file will have the suffix .msh. //
// //
//============================================================================//
void tetgenmesh::outmesh2gmsh(char* mfilename)
{
FILE *outfile;
char mefilename[FILENAMESIZE];
tetrahedron* tetptr;
face faceloop;
point ptloop, p1, p2, p3, p4;
int pointnumber;
int marker;
if (b->order == 2) {
printf(" Write Gmsh not implemented for order 2 elements \n");
return;
}
if (mfilename != (char *) NULL && mfilename[0] != '\0') {
strcpy(mefilename, mfilename);
} else if (b->outfilename[0] != '\0') {
strcpy(mefilename, b->outfilename);
} else {
strcpy(mefilename, "unnamed");
}
strcat(mefilename, ".msh");
if (!b->quiet) {
printf("Writing %s.\n", mefilename);
}
outfile = fopen(mefilename, "w");
if (outfile == (FILE *) NULL) {
printf("File I/O Error: Cannot create file %s.\n", mefilename);
return;
}
fprintf(outfile, "$MeshFormat\n2.2 0 8\n$EndMeshFormat\n");
fprintf(outfile, "$Nodes\n");
fprintf(outfile, "%ld\n", points->items);
points->traversalinit();
ptloop = pointtraverse();
pointnumber = 1;
while (ptloop != (point) NULL) {
// Point coordinates.
fprintf(outfile, "%d %.17g %.17g %.17g\n", pointnumber, ptloop[0], ptloop[1], ptloop[2]);
setpointmark(ptloop, pointnumber);
ptloop = pointtraverse();
pointnumber++;
}
fprintf(outfile, "$EndNodes\n");
long ntets = tetrahedrons->items - hullsize;
long nfaces = subfaces->items;
fprintf(outfile, "$Elements\n");
fprintf(outfile, "%ld\n", ntets + nfaces);
int *face_marks = new int [nfaces];
double *ele_marks = new double [numelemattrib * ntets];
triface abuttingtet;
int t1ver;
subfaces->traversalinit();
faceloop.sh = shellfacetraverse(subfaces);
int facidx = 1;
while (faceloop.sh != (shellface *) NULL) {
stpivot(faceloop, abuttingtet);
// If there is a tetrahedron containing this subface, orient it so
// that the normal of this face points to inside of the volume by
// right-hand rule.
if (abuttingtet.tet != NULL) {
if (ishulltet(abuttingtet)) {
fsymself(abuttingtet);
}
}
if (abuttingtet.tet != NULL) {
p1 = org (abuttingtet);
p2 = dest(abuttingtet);
p3 = apex(abuttingtet);
} else {
// A dangling subfacet.
p1 = sorg(faceloop);
p2 = sdest(faceloop);
p3 = sapex(faceloop);
}
if (!b->nobound) {
marker = shellmark(faceloop);
face_marks[facidx-1] = marker;
}
fprintf(outfile, "%d 2 1 %d %d %d %d\n",
facidx, marker, pointmark(p1), pointmark(p2), pointmark(p3));
facidx++;
faceloop.sh = shellfacetraverse(subfaces);
}
tetrahedrons->traversalinit();
tetptr = tetrahedrontraverse();
int tetnumber = facidx;
int tetcount = 0;
while (tetptr != (tetrahedron *) NULL) {
if (!b->reversetetori) {
p1 = (point) tetptr[4];
p2 = (point) tetptr[5];
} else {
p1 = (point) tetptr[5];
p2 = (point) tetptr[4];
}
p3 = (point) tetptr[6];
p4 = (point) tetptr[7];
fprintf(outfile, "%d 4 %d", tetnumber, numelemattrib);
for (int i = 0; i < numelemattrib; i++) {
ele_marks[tetcount * numelemattrib + i] = elemattribute(tetptr, i);
fprintf(outfile, " %.17g", elemattribute(tetptr, i));
}
fprintf(outfile, " %d %d %d %d\n",
pointmark(p1), pointmark(p2),
pointmark(p3), pointmark(p4));
// Remember the index of this element (for counting edges).
setelemindex(tetptr, tetnumber);
tetptr = tetrahedrontraverse();
tetnumber++;
tetcount++;
}
fprintf(outfile, "$EndElements\n");
for (int i = 0; i < nfaces; i++)
{
fprintf(outfile, "$ElementData\n1\n\"Face Tag\"\n1\n0\n3\n0\n1\n");
fprintf(outfile, "%ld\n", nfaces);
for (i = 0; i < nfaces; i++) {
fprintf(outfile, "%d %d\n", i + 1, face_marks[i]);
}
fprintf(outfile, "$EndElementData\n");
}
tetnumber = nfaces + 1;
for (int j = 0; j < numelemattrib; j++) {
fprintf(outfile, "$ElementData\n1\n\"Element Tag");
fprintf(outfile, " %d\"\n1\n0\n3\n0\n1\n", j+1);
fprintf(outfile, "%ld\n", ntets);
for (int i = 0; i < ntets; i++) {
fprintf(outfile, "%d %g\n", i + tetnumber, ele_marks[i * numelemattrib + j]);
}
fprintf(outfile, "$EndElementData\n");
}
fclose(outfile);
delete[] face_marks;
delete[] ele_marks;
}
//============================================================================//
// //
@@ -35820,7 +36111,9 @@ void tetgenmesh::outmesh2vtk(char* ofilename, int mesh_idx)
if (ofilename != (char *) NULL && ofilename[0] != '\0') {
//strcpy(vtkfilename, ofilename);
sprintf(vtkfilename, "%s.%d.vtk", ofilename, mesh_idx);
// modified by Yi Zhang on 2023-12-22
//sprintf(vtkfilename, "%s.%d.vtk", ofilename, mesh_idx);
snprintf(vtkfilename, FILENAMESIZE, "%s.%d.vtk", ofilename, mesh_idx);
} else if (b->outfilename[0] != '\0') {
strcpy(vtkfilename, b->outfilename);
strcat(vtkfilename, ".vtk");
@@ -35863,6 +36156,7 @@ void tetgenmesh::outmesh2vtk(char* ofilename, int mesh_idx)
tetrahedrons->traversalinit();
tptr = tetrahedrontraverse();
//elementnumber = firstindex; // in->firstnumber;
int first_id = 1;
while (tptr != (tetrahedron *) NULL) {
if (!b->reversetetori) {
p1 = (point) tptr[4];
@@ -35873,10 +36167,12 @@ void tetgenmesh::outmesh2vtk(char* ofilename, int mesh_idx)
}
p3 = (point) tptr[6];
p4 = (point) tptr[7];
n1 = pointmark(p1) - in->firstnumber;
n2 = pointmark(p2) - in->firstnumber;
n3 = pointmark(p3) - in->firstnumber;
n4 = pointmark(p4) - in->firstnumber;
// modified by Yi Zhang
// VTK's node index allways starts with zero
n1 = pointmark(p1) - first_id; //in->firstnumber;
n2 = pointmark(p2) - first_id; //in->firstnumber;
n3 = pointmark(p3) - first_id; //in->firstnumber;
n4 = pointmark(p4) - first_id; //in->firstnumber;
fprintf(outfile, "%d %4d %4d %4d %4d\n", nnodes, n1, n2, n3, n4);
tptr = tetrahedrontraverse();
}
@@ -35929,7 +36225,9 @@ void tetgenmesh::out_surfmesh_vtk(char* ofilename, int mesh_idx)
if (ofilename != (char *) NULL && ofilename[0] != '\0') {
//strcpy(vtkfilename, ofilename);
sprintf(vtkfilename, "%s.%d.vtk", ofilename, mesh_idx);
// modified by Yi Zhang on 2023-12-22
//sprintf(vtkfilename, "%s.%d.vtk", ofilename, mesh_idx);
snprintf(vtkfilename, FILENAMESIZE, "%s.%d.vtk", ofilename, mesh_idx);
} else if (b->outfilename[0] != '\0') {
strcpy(vtkfilename, b->outfilename);
strcat(vtkfilename, ".surf.vtk");
@@ -35972,6 +36270,7 @@ void tetgenmesh::out_surfmesh_vtk(char* ofilename, int mesh_idx)
subfaces->traversalinit();
faceloop.sh = shellfacetraverse(subfaces);
//facenumber = firstindex; // in->firstnumber;
int first_id = 1;
while (faceloop.sh != (shellface *) NULL) {
stpivot(faceloop, abuttingtet);
// If there is a tetrahedron containing this subface, orient it so
@@ -35993,9 +36292,11 @@ void tetgenmesh::out_surfmesh_vtk(char* ofilename, int mesh_idx)
tapex = sapex(faceloop);
}
n1 = pointmark(torg) - in->firstnumber;
n2 = pointmark(tdest) - in->firstnumber;
n3 = pointmark(tapex) - in->firstnumber;
// modified by Yi Zhang
// VTK's node index allways starts with zero
n1 = pointmark(torg) - first_id; //in->firstnumber;
n2 = pointmark(tdest) - first_id; //in->firstnumber;
n3 = pointmark(tapex) - first_id; //in->firstnumber;
fprintf(outfile, "%d %4d %4d %4d\n", nnodes, n1, n2, n3);
//facenumber++;
@@ -36452,6 +36753,10 @@ void tetrahedralize(tetgenbehavior *b, tetgenio *in, tetgenio *out,
m.outmesh2medit(b->outfilename);
}
// modified by Yi Zhang
if (!out && b->gmshview) {
m.outmesh2gmsh(b->outfilename);
}
if (!out && b->vtkview) {
m.outmesh2vtk(NULL, 0); // b->outfilename