Skip to content
Snippets Groups Projects
Commit f842ab46 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

simpler tetgen integration by Jozef Vesely
parent 9e2a5cdd
No related branches found
No related tags found
No related merge requests found
// $Id: CommandLine.cpp,v 1.62 2005-08-09 23:41:12 geuzaine Exp $
// $Id: CommandLine.cpp,v 1.63 2005-08-22 00:29:11 geuzaine Exp $
//
// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
//
......@@ -69,7 +69,7 @@ void Print_Usage(char *name){
Msg(DIRECT, " -1, -2, -3 Perform batch 1D, 2D and 3D mesh generation");
Msg(DIRECT, " -saveall Save all elements (discard physical group definitions)");
Msg(DIRECT, " -o file Specify mesh output file name");
Msg(DIRECT, " -format string Set output mesh format (msh, unv, gref, p3d)");
Msg(DIRECT, " -format string Set output mesh format (msh, unv, gref, stl, p3d)");
Msg(DIRECT, " -algo string Select mesh algorithm (iso, tri, aniso, netgen, tetgen)");
Msg(DIRECT, " -smooth int Set number of mesh smoothing steps");
Msg(DIRECT, " -optimize Optimize quality of tetrahedral elements");
......@@ -386,6 +386,9 @@ void Get_Options(int argc, char *argv[])
!strcmp(argv[i], "GREF") || !strcmp(argv[i], "Gref")) {
CTX.mesh.format = FORMAT_GREF;
}
else if(!strcmp(argv[i], "stl") || !strcmp(argv[i], "STL")) {
CTX.mesh.format = FORMAT_STL;
}
else if(!strcmp(argv[i], "p3d") ||
!strcmp(argv[i], "P3D") || !strcmp(argv[i], "Plot3D")) {
CTX.mesh.format = FORMAT_P3D;
......
......@@ -930,8 +930,8 @@ StringXNumber MeshOptions_Number[] = {
{ F|O, "PointType" , opt_mesh_point_type , 0. ,
"Display mesh vertices as solid color dots (0) or 3D spheres (1)" },
//{ F|O, "Quality" , opt_mesh_quality , 0.0 ,
// "Target quality for tetrahedral elements (not fully functional)" },
{ F|O, "Quality" , opt_mesh_quality , 1.0 ,
"Target quality for tetrahedral elements (currently only used by Tetgen)" },
{ F|O, "QualityInf" , opt_mesh_quality_inf , 0.0 ,
"Only display elements whose quality measure is greater than QualityInf" },
{ F|O, "QualitySup" , opt_mesh_quality_sup , 0.0 ,
......
// $Id: Solvers.cpp,v 1.37 2005-08-09 23:41:13 geuzaine Exp $
// $Id: Solvers.cpp,v 1.38 2005-08-22 00:29:11 geuzaine Exp $
//
// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
//
......@@ -56,6 +56,8 @@ SolverInfo SINFO[MAXSOLVERS];
// each new connection, but that's still a bit of a nightmare to
// maintain in a portable way on all the platforms.)
// FIXME: we should reimplement this using Fl::add_fd()
int WaitForData(int socket, int num, int pollint, double waitint)
{
struct pollfd pfd;
......
// $Id: 3D_Mesh_Tetgen.cpp,v 1.2 2005-07-03 08:02:24 geuzaine Exp $
// $Id: 3D_Mesh_Tetgen.cpp,v 1.3 2005-08-22 00:29:11 geuzaine Exp $
//
// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
//
......@@ -20,11 +20,10 @@
// Please report all bugs and problems to <gmsh@geuz.org>.
//
// Contributor(s):
// Nicolas Tardieu
// Jozef Vesely
//
#include "Gmsh.h"
#include "Geo.h"
#include "Mesh.h"
#include "Create.h"
#include "Numeric.h"
......@@ -41,259 +40,144 @@ int Mesh_Tetgen(Volume * v)
Msg(GERROR, "Tetgen is not compiled in this version of Gmsh");
return 0;
}
#else
#include "tetgen.h"
class Tetgen{
private:
List_T *_surverts, *_volverts;
Volume *_vol;
tetgenio *in, *out;
public:
Tetgen(Volume *vol);
~Tetgen();
void MeshVolume();
void TransferVolumeMesh();
};
Tetgen::Tetgen(Volume *vol)
: _volverts(0), _vol(vol)
{
// ===================================
// TRANSFER GMSH => TETGEN
// ===================================
// creates Tetgen mesh structure
tetgenio *in = new tetgenio;
tetgenio *out = new tetgenio;
tetgenio::facet *f;
tetgenio::polygon *p;
// All indices start from 1.
in->firstnumber = 1;
// Get all surface vertices (the same vertex can belong to several
// surfaces...)
Tree_T *tree = Tree_Create(sizeof(Vertex*), compareVertex);
for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
Surface *s;
List_Read(_vol->Surfaces, i, &s);
Tree_Unit(tree, s->Vertices);
}
_surverts = Tree2List(tree);
List_Sort(_surverts, compareVertex);
Tree_Delete(tree);
// Tetgen points input data
in->numberofpoints = List_Nbr(_surverts);
in->pointlist = new REAL[in->numberofpoints * 3];
// Transfer the vertices
for(int i = 0; i < List_Nbr(_surverts); i++){
Vertex *v;
List_Read(_surverts, i, &v);
in->pointlist[3*i + 0] = v->Pos.X;
in->pointlist[3*i + 1] = v->Pos.Y;
in->pointlist[3*i + 2] = v->Pos.Z;
}
int Mesh_Tetgen(Volume * vol) {
tetgenio in, out;
char opts[128];
// Count number of facets
int NbOfFacets=0;
for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
Surface *s;
List_Read(_vol->Surfaces, i, &s);
List_T *simplist = Tree2List(s->Simplexes);
NbOfFacets = NbOfFacets + List_Nbr(simplist);
}
// Tetgen elements input data
in->numberoffacets = NbOfFacets;
in->facetlist = new tetgenio::facet[in->numberoffacets];
in->facetmarkerlist = new int[in->numberoffacets];
// Transfert all surface simplices
int currentFacet=0;
for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
Surface *s;
List_Read(_vol->Surfaces, i, &s);
int sign;
List_Read(_vol->SurfacesOrientations, i, &sign);
List_T *simplist = Tree2List(s->Simplexes);
for(int j = 0; j < List_Nbr(simplist); j++) {
Simplex *simp;
List_Read(simplist, j, &simp);
int tmp[3], index[3];
if(sign > 0){
index[0] = 0;
index[1] = 1;
index[2] = 2;
}
else{
index[0] = 0;
index[1] = 2;
index[2] = 1;
}
tmp[0] = 1 + List_ISearch(_surverts, &simp->V[index[0]], compareVertex);
tmp[1] = 1 + List_ISearch(_surverts, &simp->V[index[1]], compareVertex);
tmp[2] = 1 + List_ISearch(_surverts, &simp->V[index[2]], compareVertex);
if(CTX.mesh.algo3d != DELAUNAY_TETGEN)
return 0;
f = &in->facetlist[currentFacet];
f->numberofpolygons = 1;
f->polygonlist = new tetgenio::polygon[f->numberofpolygons];
f->numberofholes = 0;
f->holelist = NULL;
p = &f->polygonlist[0];
p->numberofvertices = 3;
p->vertexlist = new int[p->numberofvertices];
p->vertexlist[0] = tmp[0];
p->vertexlist[1] = tmp[1];
p->vertexlist[2] = tmp[2];
in->facetmarkerlist[currentFacet] = 0;
currentFacet=currentFacet+1;
}
List_Delete(simplist);
if(THEM->BGM.Typ == ONFILE){
Msg(GERROR, "Tetgen is not ready to be used with a background mesh");
return 0;
}
// In order to check the input...
in->save_nodes("barin");
in->save_poly("barin");
Msg(STATUS3, "Meshing volume %d with experimental tetgen", vol->Num);
// Get all surface vertices (from all surfaces)
Tree_T *treeVrtx = Tree_Create(sizeof(Vertex*), compareVertex);
Tree_T *treeSimp = Tree_Create(sizeof(Simplex*), compareSimplex);
for(int i = 0; i < List_Nbr(vol->Surfaces); i++) {
Surface *sur;
List_Read(vol->Surfaces, i, &sur);
// ===================================
// MESHING...
// ===================================
tetrahedralize("pqV", in, out);
Tree_Unit(treeVrtx, sur->Vertices);
Tree_Unit(treeSimp, sur->Simplexes);
// DELETE all triangles and vertices,
// tetgen will refine also surface mesh
Tree_Delete(sur->Simplexes);
Tree_Delete(sur->Vertices);
sur->Simplexes = Tree_Create(sizeof(Simplex *), compareQuality);
sur->Vertices = Tree_Create(sizeof(Vertex *), compareVertex);
}
List_T *listVrtx = Tree2List(treeVrtx);
List_T *listSimp = Tree2List(treeSimp);
// ===================================
// TRANSFER TETGEN => GMSH
// ===================================
// Gets total number of vertices of Tetgen's mesh
int nbv = out->numberofpoints;
Msg(INFO, " NIKO : nbv=%i",nbv);
Tree_Delete(treeVrtx);
Tree_Delete(treeSimp);
for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
Surface *s;
List_Read(_vol->Surfaces, i, &s);
Tree_Delete(s->Simplexes);
}
Tree_Delete(_vol->Vertices);
Tree_Delete(THEM->Vertices);
_vol->Vertices = Tree_Create(sizeof(Simplex *), compareSimplex);
THEM->Vertices = Tree_Create(sizeof(Simplex *), compareSimplex);
THEM->MaxPointNum=0;
// Create new vertices
Vertex **vtable = (Vertex **)Malloc(nbv * sizeof(Vertex*));
for(int i = 0; i < nbv; i++) {
vtable[i] = Create_Vertex(++(THEM->MaxPointNum), out->pointlist[3*i + 0], out->pointlist[3*i + 1],out->pointlist[3*i + 2], 1., 0);
Tree_Add(_vol->Vertices, &vtable[i]);
Tree_Add(THEM->Vertices, &vtable[i]);
}
// Set input parameters
in.mesh_dim = 3;
in.firstnumber = 1;
Msg(INFO, "<NIKO>: nbTriFace=%i",out->numberoftrifaces);
if(!nbv) {
Msg(STATUS3, "<NIKO>: AUCUN ELEMENT");
return;}
in.numberofpoints = List_Nbr(listVrtx);
in.pointlist = new REAL[in.numberofpoints * 3];
in.pointmarkerlist = NULL;
in.numberoffacets = List_Nbr(listSimp);
in.facetlist = new tetgenio::facet[in.numberoffacets];
in.facetmarkerlist = new int[in.numberoffacets];
// Get total number of simplices of Tetgen's mesh
int nbe = out->numberoftetrahedra;
double lc_max = -1.0;
for(int i = 0; i < List_Nbr(listVrtx); i++){
Vertex *v;
List_Read(listVrtx, i, &v);
in.pointlist[i*3 + 0] = v->Pos.X;
in.pointlist[i*3 + 1] = v->Pos.Y;
in.pointlist[i*3 + 2] = v->Pos.Z;
if(v->lc > lc_max) lc_max = v->lc;
// Create new volume simplices
for(int i = 0; i < nbe; i++) {
Simplex *simp = Create_Simplex(vtable[out->tetrahedronlist[4*i + 0]-1], vtable[out->tetrahedronlist[4*i + 1]-1],
vtable[out->tetrahedronlist[4*i + 2]-1], vtable[out->tetrahedronlist[4*i + 3]-1]);
simp->iEnt = _vol->Num;
Tree_Add(_vol->Simplexes, &simp);
// also add a copy in the global simplex tree
Tree_Add(THEM->Simplexes, &simp);
// DELETE the vertices from global mesh
Tree_Suppress(THEM->Vertices, &v);
}
Free(vtable);
}
for(int i = 0; i < List_Nbr(listSimp); i++) {
Simplex *s;
List_Read(listSimp, i, &s);
tetgenio::facet *f = &in.facetlist[i];
tetgenio::init(f);
Tetgen::~Tetgen()
{
List_Delete(_surverts);
List_Delete(_volverts);
// Strange : this method is not seen by the compiler
// tetgenio::deinitialize();
}
f->numberofholes = 0;
f->numberofpolygons = 1;
tetgenio::polygon *p = f->polygonlist = new tetgenio::polygon[f->numberofpolygons];
tetgenio::init(p);
void Tetgen::MeshVolume()
{
tetrahedralize("pqV", in, out);
}
p->numberofvertices = 3;
p->vertexlist = new int[p->numberofvertices];
p->vertexlist[0] = List_ISearch(listVrtx, &s->V[0], compareVertex) +1;
p->vertexlist[1] = List_ISearch(listVrtx, &s->V[1], compareVertex) +1;
p->vertexlist[2] = List_ISearch(listVrtx, &s->V[2], compareVertex) +1;
void Tetgen::TransferVolumeMesh()
{
// Gets total number of vertices of Tetgen's mesh
int nbv = out->numberofpoints;
Msg(INFO, " NIKO : nbv=%i",nbv);
in.facetmarkerlist[i] = s->iEnt;
}
if(!nbv) {
Msg(STATUS3, " NIKO : AUCUN ELEMENT");
return;}
snprintf(opts, 128, "pqa%f%c", (float)CTX.mesh.quality,
(CTX.verbosity < 3)? 'Q': (CTX.verbosity > 6)? 'V': '\0');
Msg(STATUS3, "Meshing with volume constraint %f", (float)CTX.mesh.quality);
Vertex **vtable = (Vertex **)Malloc(nbv * sizeof(Vertex*));
tetrahedralize(opts, &in, &out);
// Get existing unmodified surface vertices
for(int i = 0; i < List_Nbr(_surverts); i++){
List_Read(_surverts, i, &vtable[i]);
Tree_Insert(_vol->Vertices, &vtable[i]);
Tree_Insert(THEM->Vertices, &vtable[i]);
}
// CHECK
// out.trifacemarkerlist != NULL
// out.numberofcorners == 4
// out.mesh_dim == 3
// Create new volume vertices
for(int i = List_Nbr(_surverts); i < nbv; i++) {
vtable[i] = Create_Vertex(++(THEM->MaxPointNum), out->pointlist[3*i + 0], out->pointlist[3*i + 1],out->pointlist[3*i + 2], 1., 0);
Tree_Add(_vol->Vertices, &vtable[i]);
Vertex **vtable = new Vertex*[out.numberofpoints];
for (int i = 0; i < out.numberofpoints; i++) {
vtable[i] = Create_Vertex(++(THEM->MaxPointNum),
out.pointlist[i * 3 + 0],
out.pointlist[i * 3 + 1],
out.pointlist[i * 3 + 2],
1., 0);
Tree_Add(vol->Vertices, &vtable[i]);
Tree_Add(THEM->Vertices, &vtable[i]);
}
// Get total number of simplices of Tetgen's mesh
int nbe = out->numberoftetrahedra;
// Create new volume simplices
for(int i = 0; i < nbe; i++) {
Msg(INFO, " NIKO out->tetrahedronlist: =%i",out->tetrahedronlist[4*i + 0]);
Simplex *simp = Create_Simplex(vtable[out->tetrahedronlist[4*i + 0]], vtable[out->tetrahedronlist[4*i + 1]],
vtable[out->tetrahedronlist[4*i + 2]], vtable[out->tetrahedronlist[4*i + 3]]);
simp->iEnt = _vol->Num;
Tree_Add(_vol->Simplexes, &simp);
// also add a copy in the global simplex tree
Tree_Add(THEM->Simplexes, &simp);
for(int j = 0; j < List_Nbr(vol->Surfaces); j++){
Surface *sur;
List_Read(vol->Surfaces, j, &sur);
for (int i = 0; i<out.numberoftrifaces; i++){
if (out.trifacemarkerlist[i] == sur->Num) {
Simplex *s = Create_Simplex(vtable[ out.trifacelist[i * 3 + 0]-1 ],
vtable[ out.trifacelist[i * 3 + 1]-1 ],
vtable[ out.trifacelist[i * 3 + 2]-1 ],
NULL);
s->iEnt = sur->Num;
Tree_Add(sur->Simplexes, &s);
}
Free(vtable);
}
int Mesh_Tetgen(Volume * v)
{
if(CTX.mesh.algo3d != DELAUNAY_TETGEN)
return 0;
Msg(GERROR, "Tetgen is not ready yet!");
return 0;
if(THEM->BGM.Typ == ONFILE){
Msg(GERROR, "Tetgen is not ready to be used with a background mesh");
return 0;
}
Msg(STATUS3, "Meshing volume %d", v->Num);
Tetgen tg(v);
//tg.MeshVolume();
//tg.TransferVolumeMesh();
for (int i = 0; i < out.numberoftetrahedra; i++) {
Simplex *s = Create_Simplex(vtable[ out.tetrahedronlist[i * 4 + 0] -1],
vtable[ out.tetrahedronlist[i * 4 + 1] -1],
vtable[ out.tetrahedronlist[i * 4 + 2] -1],
vtable[ out.tetrahedronlist[i * 4 + 3] -1]);
s->iEnt = vol->Num;
Tree_Add(vol->Simplexes, &s);
Tree_Add(THEM->Simplexes, &s);
}
List_Delete(listVrtx);
List_Delete(listSimp);
return 1;
}
......
$Id: CREDITS,v 1.27 2005-02-28 23:57:59 geuzaine Exp $
$Id: CREDITS,v 1.28 2005-08-22 00:29:11 geuzaine Exp $
Gmsh is copyright (C) 1997-2005
......@@ -20,8 +20,9 @@ ulg.ac.be> for transfinite mesh bug fixes; Laurent Stainier
<l.stainier at ulg.ac.be> for the eigenvalue solvers and for help with
the MacOS port and the tensor display code; Pierre Badel <badel at
freesurf.fr> for help with the GSL integration; Marc Ume <Marc.Ume at
digitalgraphics.be> for the original list code; Matt Gundry
<mjgundry at faa-engineers.com> for the Plot3d mesh format.
digitalgraphics.be> for the original list code; Matt Gundry <mjgundry
at faa-engineers.com> for the Plot3d mesh format; Jozef Vesely <vesely
at gjh.sk> for the Tetgen integration.
The AVL tree code (DataStr/avl.*) and the YUV image code
(Graphics/gl2yuv.*) are copyright (C) 1988-1993, 1995 The Regents of
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment