refactor the netgen code in a dedicated class

Mesh/3D_Mesh.cpp

-// $Id: 3D_Mesh.cpp,v 1.61 2004-06-26 17:58:14 geuzaine Exp $
+// $Id: 3D_Mesh.cpp,v 1.62 2004-06-26 21:34:50 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -897,9 +897,7 @@ void Maillage_Volume(void *data, void *dum)
   }
   else if(MeshTransfiniteVolume(v)) {
   }
-  else if(v->Typ != 99999){
-    if (Mesh_Netgen(v)) {
-    }
+  else if((v->Typ != 99999) && Mesh_Netgen(v)) {
   }
   else if((v->Typ == 99999) && (CTX.mesh.algo3d != FRONTAL_NETGEN)) {
     

Mesh/3D_Mesh_Netgen.cpp

-// $Id: 3D_Mesh_Netgen.cpp,v 1.1 2004-06-26 17:58:14 geuzaine Exp $
+// $Id: 3D_Mesh_Netgen.cpp,v 1.2 2004-06-26 21:34:50 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -30,12 +30,18 @@
 
 #if !defined(HAVE_NETGEN)
 
-int Mesh_Netgen(Surface * s)
+int Mesh_Netgen(Volume * v)
 {
+  if(CTX.mesh.algo3d == FRONTAL_NETGEN)
     Msg(GERROR, "Netgen is not compiled in this version of Gmsh");
   return 0;
 }
 
+void Optimize_Netgen(Volume * v)
+{
+  Msg(GERROR, "Netgen is not compiled in this version of Gmsh");
+}
+
 #else
 
 extern "C"
@@ -46,58 +52,65 @@ extern "C"
 extern Context_T CTX;
 extern Mesh *THEM;
 
-int Mesh_Netgen(Volume * v)
+class Netgen{
+ private:
+  List_T *_vertices;
+  Volume *_vol;
+  Ng_Mesh *_ngmesh;
+ public:
+  Netgen(Volume *vol, int importVolumeMesh = 0);
+  Netgen(Surface *s, int importSurfaceMesh = 0);
+  ~Netgen();
+  void MeshVolume();
+  void OptimizeVolume();
+};
+
+Netgen::Netgen(Volume *vol, int importVolumeMesh)
+  : _vol(vol)
 {
-  if(CTX.mesh.algo3d != FRONTAL_NETGEN)
-    return 0;
-
-  if(THEM->BGM.Typ == ONFILE){
-    Msg(GERROR, "Netgen is not ready to be used with a background mesh");
-    return 0;
-  }
-
-  List_T *lSurfaceVertexList = List_Create(100, 100, sizeof(Vertex*));
-  
-  // ===================================
-  //         GMSH => NETGEN 
-  // ===================================
-  
   // creates Netgen mesh structure
   Ng_Init();
-  Ng_Mesh *lNetgenMesh = Ng_NewMesh();
+  _ngmesh = Ng_NewMesh();
 
-  // Transfert all surface vertices to Netgen
-  for(int i = 0; i < List_Nbr(v->Surfaces); i++) {
+  if(importVolumeMesh){
+    _vertices = Tree2List(_vol->Vertices);
+  }
+  else{
+    // Transfert all surface vertices we must *not* add 2 times the
+    // same vertex (the same vertex can belong to several surfaces)
+    _vertices = List_Create(100, 100, sizeof(Vertex*));
+    for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
       Surface *s;
-    List_Read(v->Surfaces, i, &s);
-    List_T *lVertexList = Tree2List(s->Vertices);
-    for(int j = 0; j < List_Nbr(lVertexList); j++) {
-      Vertex *lVertex;
-      List_Read(lVertexList, j, &lVertex);
-      // We must *not* add 2 times the same vertex (the same vertex
-      // can belong to several surfaces)
-      if(List_ISearchSeq(lSurfaceVertexList, &lVertex, compareVertex) < 0) {
-	List_Add(lSurfaceVertexList, &lVertex);
-	double tmp[3];
-	tmp[0] = lVertex->Pos.X;
-	tmp[1] = lVertex->Pos.Y;
-	tmp[2] = lVertex->Pos.Z;
-      	Ng_AddPoint(lNetgenMesh, tmp);
+      List_Read(_vol->Surfaces, i, &s);
+      List_T *vlist = Tree2List(s->Vertices);
+      for(int j = 0; j < List_Nbr(vlist); j++) 
+	List_Insert(_vertices, List_Pointer(vlist, j), compareVertex);
+      List_Delete(vlist);
     }
   }
-    List_Delete(lVertexList) ;
+  List_Sort(_vertices, compareVertex);
+
+  // Transfer the vertices
+  for(int i = 0; i < List_Nbr(_vertices); i++){
+    Vertex *v;
+    List_Read(_vertices, i, &v);
+    double tmp[3];
+    tmp[0] = v->Pos.X;
+    tmp[1] = v->Pos.Y;
+    tmp[2] = v->Pos.Z;
+    Ng_AddPoint(_ngmesh, tmp);
   }
 
-  // Transfert all surface simplices to Netgen
-  for(int i = 0; i < List_Nbr(v->Surfaces); i++) {
+  // Transfert all surface simplices
+  for(int i = 0; i < List_Nbr(_vol->Surfaces); i++) {
     Surface *s;
+    List_Read(_vol->Surfaces, i, &s);
     int sign;
-    List_Read(v->Surfaces, i, &s);
-    List_Read(v->SurfacesOrientations, i, &sign);
-    List_T *lSimplexeList = Tree2List(s->Simplexes);
-    for(int j = 0; j < List_Nbr(lSimplexeList); j++) {
-      Simplex *lSimplex;
-      List_Read(lSimplexeList, j, &lSimplex);
+    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;
@@ -109,71 +122,116 @@ int Mesh_Netgen(Volume * v)
 	index[1] = 2;
 	index[2] = 1;
       }
-      tmp[0] = 1 + List_ISearchSeq(lSurfaceVertexList, &lSimplex->V[index[0]],
-				   compareVertex);
-      tmp[1] = 1 + List_ISearchSeq(lSurfaceVertexList, &lSimplex->V[index[1]],
-				   compareVertex);
-      tmp[2] = 1 + List_ISearchSeq(lSurfaceVertexList, &lSimplex->V[index[2]],
-				   compareVertex);
-      Ng_AddSurfaceElement(lNetgenMesh, NG_TRIG, tmp);
+      tmp[0] = 1 + List_ISearch(_vertices, &simp->V[index[0]], compareVertex);
+      tmp[1] = 1 + List_ISearch(_vertices, &simp->V[index[1]], compareVertex);
+      tmp[2] = 1 + List_ISearch(_vertices, &simp->V[index[2]], compareVertex);
+      Ng_AddSurfaceElement(_ngmesh, NG_TRIG, tmp);
     }
-    List_Delete(lSimplexeList);
+    List_Delete(simplist);
+  }
+  
+  // Transfer the volume elements
+  if(importVolumeMesh){
+    List_T *simplist = Tree2List(_vol->Simplexes);
+    for(int i = 0; i < List_Nbr(simplist); i++) {
+      Simplex *simp;
+      List_Read(simplist, i, &simp);
+      int tmp[4];
+      tmp[0] = 1 + List_ISearch(_vertices, &simp->V[0], compareVertex);
+      tmp[1] = 1 + List_ISearch(_vertices, &simp->V[1], compareVertex);
+      tmp[2] = 1 + List_ISearch(_vertices, &simp->V[2], compareVertex);
+      tmp[3] = 1 + List_ISearch(_vertices, &simp->V[3], compareVertex);
+      Ng_AddVolumeElement(_ngmesh, NG_TET, tmp);
+    }
+    List_Delete(simplist);    
+  }
+}
+
+Netgen::Netgen(Surface *sur, int importSurfaceMesh)
+{
+  // todo
 }
 
-  // ===================================
-  //            MESHING
-  // ===================================
+Netgen::~Netgen()
+{
+  List_Delete(_vertices);
+  Ng_DeleteMesh(_ngmesh);
+  Ng_Exit();
+}
 
-  // define Netgen meshing option
+void Netgen::MeshVolume()
+{
   Ng_Meshing_Parameters mp;
   mp.maxh = 1;
   mp.fineness = 1;
   mp.secondorder = 0;
-
-  // generate Netgen volume mesh
-  Msg(STATUS3, "Meshing volume %d", v->Num);
-  Ng_GenerateVolumeMesh(lNetgenMesh, &mp);
-  
-  // ===================================
-  //            NETGEN => GMSH
-  // ===================================
+  Ng_GenerateVolumeMesh(_ngmesh, &mp);
 
   // Gets total number of vertices of Netgen's mesh
-  int lNbOfNGPoints = Ng_GetNP(lNetgenMesh);
-  Vertex **vtable = (Vertex **)Malloc(lNbOfNGPoints * sizeof(Vertex*));
+  int nbv = Ng_GetNP(_ngmesh);
+  Vertex **vtable = (Vertex **)Malloc(nbv * sizeof(Vertex*));
   
-  // Writes existing surface vertices
-  for(int i = 0; i < List_Nbr(lSurfaceVertexList); i++)
-    List_Read(lSurfaceVertexList, i, &vtable[i]);
+  // Get existing vertices
+  for(int i = 0; i < List_Nbr(_vertices); i++)
+    List_Read(_vertices, i, &vtable[i]);
 
-  // Writes and transfers new volume vertices
-  for(int i = List_Nbr(lSurfaceVertexList); i < lNbOfNGPoints; i++) {
+  // Create new volume vertices
+  for(int i = List_Nbr(_vertices); i < nbv; i++) {
     double tmp[3];
-    Ng_GetPoint(lNetgenMesh, i+1, tmp);
+    Ng_GetPoint(_ngmesh, i+1, tmp);
     vtable[i] = Create_Vertex(++(THEM->MaxPointNum), tmp[0], tmp[1], tmp[2], 1., 0);
-    Tree_Add(v->Vertices, &vtable[i]);
+    Tree_Add(_vol->Vertices, &vtable[i]);
     Tree_Add(THEM->Vertices, &vtable[i]);
   }
 
   // Get total number of simplices of Netgen's mesh
-  int lNbOfNGElements = Ng_GetNE(lNetgenMesh);
-  // Transfers new volume simplices
-  for(int i = 1; i <= lNbOfNGElements; i++) {
+  int nbe = Ng_GetNE(_ngmesh);
+
+  // Create new volume simplices
+  for(int i = 0; i < nbe; i++) {
     int tmp[4];
-    Ng_GetVolumeElement(lNetgenMesh, i, tmp);
-    Simplex *lSimplex = Create_Simplex(vtable[tmp[0]-1], vtable[tmp[1]-1],
+    Ng_GetVolumeElement(_ngmesh, i+1, tmp);
+    Simplex *simp = Create_Simplex(vtable[tmp[0]-1], vtable[tmp[1]-1],
 				   vtable[tmp[2]-1], vtable[tmp[3]-1]);
-    lSimplex->iEnt = v->Num;
-    Tree_Add(v->Simplexes, &lSimplex);
+    simp->iEnt = _vol->Num;
+    Tree_Add(_vol->Simplexes, &simp);
   }
   
-  List_Delete(lSurfaceVertexList);
   Free(vtable);
+}
 
-  Ng_DeleteMesh(lNetgenMesh);
-  Ng_Exit();
+void Netgen::OptimizeVolume()
+{
+  // remove the pure volume vertices from THEM->Vertices and v->Vertices
+  // remove the tets from THEM->Simplexes
+  // reset v->Simplexes
+  //RemoveIllegalElements(*_ngmesh);
+  //OptimizeVolume(mparam, *_ngmesh, NULL);
+  // transfer vertices and tets back into v and THEM
+}
+
+int Mesh_Netgen(Volume * v)
+{
+  if(CTX.mesh.algo3d != FRONTAL_NETGEN)
+    return 0;
+
+  if(THEM->BGM.Typ == ONFILE){
+    Msg(GERROR, "Netgen is not ready to be used with a background mesh");
+    return 0;
+  }
+
+  Msg(STATUS3, "Meshing volume %d", v->Num);
+  Netgen ng(v);
+  ng.MeshVolume();
   
   return 1;
 }
 
+void Optimize_Netgen(Volume * v)
+{
+  Msg(STATUS3, "Optimizing volume %d", v->Num);
+  Netgen ng(v, 1);
+  ng.OptimizeVolume();
+}
+
 #endif // !HAVE_NETGEN

Mesh/Mesh.h

@@ -469,6 +469,7 @@ int AlgorithmeMaillage2DAnisotropeModeJF(Surface * s);
 void Maillage_Automatique_VieuxCode(Surface * pS, Mesh * m, int ori);
 int Mesh_Shewchuk(Surface *s);
 int Mesh_Netgen(Volume * v);
+void Optimize_Netgen(Volume * v);
 
 int Calcule_Contours(Surface * s);
 void Link_Simplexes(List_T * Sim, Tree_T * Tim);

tutorial/t5.geo

@@ -116,8 +116,7 @@ Function CheeseHole
 
   theloops[t] = newreg ; 
 
-  Surface Loop(theloops[t]) = {l8+1, l5+1, l1+1, l2+1, -(l3+1), -(l7+1),
-			       l6+1, l4+1};
+  Surface Loop(theloops[t]) = {l8+1, l5+1, l1+1, l2+1, l3+1, l7+1, l6+1, l4+1};
 
   thehole = newreg ; 
   Volume(thehole) = theloops[t] ;