Added simple possibility to create geo model from lua using old Geo and not...

Added simple possibility to create geo model from lua using old Geo and not OCC, added example in bench/boolean. 

Geo/GModel.cpp

@@ -27,10 +27,14 @@
 #include "boundaryLayerEdge.h"
 #include "boundaryLayerVertex.h"
 #include "gmshSurface.h"
+#include "ListUtils.h"
+#include "Geo.h"
 #include "SmoothData.h"
 #include "Context.h"
 #include "OS.h"
 #include "GEdgeLoop.h"
+#include "gmshFace.h"
+#include "gmshRegion.h"
 #include "MVertexPositionSet.h"
 #include "OpenFile.h"
 #include "CreateFile.h"
@@ -1620,6 +1624,98 @@ GFace* GModel::addFace (std::vector<GEdge *> edges, std::vector< std::vector<dou
   return 0;
 }
 
+GFace* GModel::addGeoPlanarFace (std::vector<std::vector<GEdge *> > edges){
+  
+  //create line loops
+   std::vector<EdgeLoop *> vecLoops;
+   int nLoops = edges.size();
+   for (int i=0; i< nLoops; i++){
+     int numl = getMaxElementaryNumber(1) + i;
+     while (FindEdgeLoop(numl)){
+       numl++;
+       if (!FindEdgeLoop(numl)) break;
+     }
+     int nl=(int)edges[i].size();
+     List_T *iListl = List_Create(nl, nl, sizeof(int));
+     for(int j = 0; j < nl; j++){
+       int numEdge = edges[i][j]->tag();
+       List_Add(iListl, &numEdge);
+     }
+     sortEdgesInLoop(numl, iListl);
+     EdgeLoop *l = Create_EdgeLoop(numl, iListl);
+     vecLoops.push_back(l);
+     Tree_Add(GModel::current()->getGEOInternals()->EdgeLoops, &l);
+     l->Num = numl;
+     List_Delete(iListl);
+   }
+ 
+  //create plane surfaces
+  int numf = getMaxElementaryNumber(2) + 1;
+  Surface *s = Create_Surface(numf, MSH_SURF_PLAN);
+  List_T *iList = List_Create(nLoops, nLoops, sizeof(int));
+  for (int i=0; i< vecLoops.size(); i++){
+    int numl = vecLoops[i]->Num;
+    List_Add(iList, &numl);
+  }
+  setSurfaceGeneratrices(s, iList);
+  End_Surface(s);
+  Tree_Add(GModel::current()->getGEOInternals()->Surfaces, &s);
+  s->Typ= MSH_SURF_PLAN;
+  s->Num = numf;
+  List_Delete(iList);
+ 
+  //gmsh surface
+  GFace *gf = new gmshFace(this,s);
+  add(gf);
+
+  return gf;
+
+}
+
+GRegion* GModel::addGeoVolume (std::vector<std::vector<GFace *> > faces){
+  
+  //surface loop
+   std::vector<SurfaceLoop *> vecLoops;
+   int nLoops = faces.size();
+   for (int i=0; i< nLoops; i++){
+     int numfl = getMaxElementaryNumber(2) + 1;
+     while (FindSurfaceLoop(numfl)){
+       numfl++;
+       if (!FindSurfaceLoop(numfl)) break;
+     }
+     int nl=(int)faces[i].size();
+     List_T *iListl = List_Create(nl, nl, sizeof(int));
+     for(int j = 0; j < nl; j++){
+       int numFace = faces[i][j]->tag();
+       List_Add(iListl, &numFace);
+     }
+     SurfaceLoop *l = Create_SurfaceLoop(numfl, iListl);
+     vecLoops.push_back(l);
+     Tree_Add(GModel::current()->getGEOInternals()->SurfaceLoops, &l);
+     List_Delete(iListl);
+   }
+
+   //volume
+   int numv = getMaxElementaryNumber(3) + 1;
+   Volume *v = Create_Volume(numv, MSH_VOLUME);
+   List_T *iList = List_Create(nLoops, nLoops, sizeof(int));
+   for (int i=0; i< vecLoops.size(); i++){
+     int numl = vecLoops[i]->Num;
+     List_Add(iList, &numl);
+   }
+   setVolumeSurfaces(v, iList);
+   List_Delete(iList);
+   Tree_Add(GModel::current()->getGEOInternals()->Volumes, &v);
+   v->Typ = MSH_VOLUME;
+   v->Num = numv;
+   
+   //gmsh volume
+  GRegion *gr = new gmshRegion(this,v);
+  add(gr);
+
+  return gr;
+}
+
 GFace* GModel::addPlanarFace (std::vector<std::vector<GEdge *> > edges){
   if(_factory)
     return _factory->addPlanarFace(this, edges);
@@ -2066,11 +2162,17 @@ void GModel::registerBindings(binding *b)
   cm->setDescription("creates a face that is constraint by edges and points");
   cm->setArgNames("{list of edges}","{{x,y,z},{x,y,z},...}",NULL);
   cm = cb->addMethod("addRuledFaces", &GModel::addRuledFaces);
-  cm->setDescription("creates ruled faces that contains a list of wires");
+  cm->setDescription("create ruled faces that contains a list of wires");
   cm->setArgNames("{{list of edges},{list of edges},...}",NULL);
   cm = cb->addMethod("addPlanarFace", &GModel::addPlanarFace);
   cm->setDescription("creates a planar face that contains a list of wires");
   cm->setArgNames("{{list of edges},{list of edges},...}",NULL);
+  cm = cb->addMethod("addGeoPlanarFace", &GModel::addGeoPlanarFace);
+  cm->setDescription("creates a planar face that contains a list with the edges");
+  cm->setArgNames("{{list of edges},{list of edges},...}",NULL);
+  cm = cb->addMethod("addGeoVolume", &GModel::addGeoVolume);
+  cm->setDescription("creates a Volume bounded by a list of faces");
+  cm->setArgNames("{{list of faces},{list of faces},...}",NULL);
   cm = cb->addMethod("addPipe", &GModel::addPipe);
   cm->setDescription("creates a pipe with a base and a wire for the direction");
   cm->setArgNames("ge","{list of edges}",NULL);
@@ -2137,6 +2239,10 @@ void GModel::registerBindings(binding *b)
                      "is chosen. Returns the number.");
   cm->setArgNames("physicalName","dim","number",NULL);
 
+  cm = cb->addMethod("createTopology", &GModel::createTopologyFromMesh);
+  cm->setDescription("build the topology for a given mesh.");
+  cm->setArgNames(NULL);
+
   cm = cb->addMethod("createBoundaryLayer", &GModel::createBoundaryLayer);
   cm->setDescription("create a boundary layer using a given field for the "
                      "extrusion height.");

Geo/GModel.h

@@ -392,6 +392,8 @@ class GModel
   void addRuledFaces (std::vector<std::vector<GEdge *> > edges);
   GFace* addFace (std::vector<GEdge *> edges, std::vector< std::vector<double > > points);
   GFace* addPlanarFace (std::vector<std::vector<GEdge *> > edges);
+  GFace* addGeoPlanarFace (std::vector<std::vector<GEdge*> > edges);
+  GRegion* addGeoVolume (std::vector<std::vector<GFace*> > faces);
 
   // create solid geometry primitives using the factory
   GEntity *addSphere(double cx, double cy, double cz, double radius);

Geo/GModelFactory.cpp

@@ -36,6 +36,7 @@
 #include <TColgp_HArray1OfPnt.hxx>
 #include <TColStd_HArray1OfReal.hxx>
 #include <TColStd_HArray1OfInteger.hxx>
+#include "MLine.h"
 
 GVertex *OCCFactory::addVertex(GModel *gm, double x, double y, double z, double lc)
 {
@@ -559,7 +560,6 @@ std::vector<GFace *> OCCFactory::addRuledFaces(GModel *gm,
     BRepBuilderAPI_MakeWire wire_maker;
     for (unsigned j=0;j<wires[i].size();j++) {
       GEdge *ge = wires[i][j];
-      printf("edge %d\n",ge->tag());
       OCCEdge *occe = dynamic_cast<OCCEdge*>(ge);
       if (occe){
 	wire_maker.Add(occe->getTopoDS_Edge());
@@ -610,6 +610,7 @@ extern void computeMeanPlane(const std::vector<SPoint3> &points, mean_plane &mea
 
 GFace *OCCFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> > wires)
 {
+
   std::set<GVertex*> verts;
   for (unsigned i = 0; i < wires.size(); i++) {
     for (unsigned j = 0; j < wires[i].size(); j++) {
@@ -623,6 +624,7 @@ GFace *OCCFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> >
   for ( ; it != verts.end(); ++it){
     points.push_back(SPoint3((*it)->x(), (*it)->y(), (*it)->z()));
   }
+ 
   mean_plane meanPlane;
   computeMeanPlane(points, meanPlane);
 
@@ -644,9 +646,10 @@ GFace *OCCFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> >
   }
 
   aGenerator.Build();
-  
   TopoDS_Shape aResult = aGenerator.Shape();
+
   return gm->_occ_internals->addFaceToModel(gm, TopoDS::Face(aResult));
+
 }
 
 GEntity *OCCFactory::addPipe(GModel *gm, GEntity *base, std::vector<GEdge *> wire)

Geo/GModelIO_Geo.cpp

@@ -44,6 +44,7 @@ int GModel::readGEO(const std::string &name)
 
 int GModel::exportDiscreteGEOInternals()
 {
+
   if(_geo_internals) delete _geo_internals;
   _geo_internals = new GEO_Internals;
 

Geo/GModelIO_OCC.cpp

@@ -55,7 +55,6 @@ void OCC_Internals::buildShapeFromLists(TopoDS_Shape _shape)
 void OCC_Internals::addShapeToLists(TopoDS_Shape _shape)
 {
   TopExp_Explorer exp0, exp1, exp2, exp3, exp4, exp5;
-  
   for(exp0.Init(_shape, TopAbs_SOLID); exp0.More(); exp0.Next()){
     TopoDS_Solid solid = TopoDS::Solid(exp0.Current());
     if(somap.FindIndex(TopoDS::Solid(exp0.Current())) < 1){
@@ -505,28 +504,32 @@ void OCC_Internals::buildGModel(GModel *model)
   // building geom vertices
   int nvertices = vmap.Extent();
   for(int i = 1; i <= nvertices; i++){
+    int num = model->getMaxElementaryNumber(0) + 1;
     if (!getOCCVertexByNativePtr(model,TopoDS::Vertex(vmap(i))))
-      model->add(new OCCVertex(model, i, TopoDS::Vertex(vmap(i))));
+      model->add(new OCCVertex(model,num , TopoDS::Vertex(vmap(i))));
   }
   // building geom edges
   int nedges = emap.Extent();
   for(int i = 1; i <= nedges; i++){
     int i1 = vmap.FindIndex(TopExp::FirstVertex(TopoDS::Edge(emap(i)))); 
     int i2 = vmap.FindIndex(TopExp::LastVertex(TopoDS::Edge(emap(i))));
+    int num = model->getMaxElementaryNumber(1) + 1;
     if (!getOCCEdgeByNativePtr(model,TopoDS::Edge(emap(i))))
-      model->add(new OCCEdge(model, TopoDS::Edge(emap(i)), i, model->getVertexByTag(i1), model->getVertexByTag(i2)));
+      model->add(new OCCEdge(model, TopoDS::Edge(emap(i)), num, model->getVertexByTag(i1), model->getVertexByTag(i2)));
   }
   // building geom faces
   int nfaces = fmap.Extent();
   for(int i = 1; i <= nfaces; i++){
+    int num = model->getMaxElementaryNumber(2) + 1;
     if (!getOCCFaceByNativePtr(model,TopoDS::Face(fmap(i))))
-      model->add(new OCCFace(model, TopoDS::Face(fmap(i)), i));
+      model->add(new OCCFace(model, TopoDS::Face(fmap(i)), num));
   }
   // building geom regions
   int nvolumes = somap.Extent();
   for(int i = 1; i <= nvolumes; i++){
+    int num = model->getMaxElementaryNumber(3) + 1;
     if (!getOCCRegionByNativePtr(model,TopoDS::Solid(somap(i))))
-      model->add(new OCCRegion(model, TopoDS::Solid(somap(i)), i));
+      model->add(new OCCRegion(model, TopoDS::Solid(somap(i)), num));
   }
 }
 

benchmarks/boolean/aorta_bound.lua

+
+
+options = gmshOptions()
+options:initOptions()
+
+options:numberSet('Mesh', 0, 'Algorithm3D', 4.0)
+options:numberSet('Mesh', 0, 'Optimize', 1.0)
+--options:numberSet('Mesh', 0, 'CharacteristicLengthFactor', 0.2)
+options:numberSet('Mesh', 0, 'Mesh.CharacteristicLengthExtendFromBoundary',  1)
+
+print'*** create GModel from stl ***'
+model = GModel()
+model:load ('aortaADAPT.stl')
+
+print'*** create Topology ***'
+model:createTopology()
+
+e1 = model:getEdgeByTag(1)
+e2 = model:getEdgeByTag(2)
+e3 = model:getEdgeByTag(3)
+e4 = model:getEdgeByTag(4)
+e5 = model:getEdgeByTag(5)
+
+print'*** add Planar face ***'
+f0 = model:getFaceByTag(1) 
+
+f1 = model:addGeoPlanarFace({{e1}})
+f2 = model:addGeoPlanarFace({{e2}}) 
+f3 = model:addGeoPlanarFace({{e3}}) 
+f4 = model:addGeoPlanarFace({{e4}}) 
+f5 = model:addGeoPlanarFace({{e5}}) 
+
+print'*** add Volume ***'
+r1 = model:addGeoVolume({{f0,f1,f2,f3,f4,f5}});
+
+print'*** meshing ***'
+model:mesh(3);
+model:save("aorta_bound_lua.msh")
+
+
+

benchmarks/boolean/nurbs.lua

@@ -20,11 +20,11 @@ v7 = g:addVertex(1.5,-1,0,.1);
 --	   {0,0,0,0,0.25,0.5,0.75,1,1,1,1},{1,1,1,1,1,1,1},{4,1,1,1,1,1,4}); 
 
 --g:addBezier (v1,v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}}); 
---g:addNURBS (v1,v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}},
---            {0,0.25,0.5,0.75,1},{1,1,1,1,1,1,1},{4,1,1,1,4}); 							   
+
+g:addNURBS (v1, v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}},          {0,0.25,0.5,0.75,1},{1,1,1,1,1,1,1},{4,1,1,1,4}); 							   
 --g:addNURBS (v1,v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}},	   	   {0,0.5,1},{1,1,1,1,1,1,1},{4,3,4}); 							   
 
-e1 = g:addNURBS (v1,v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}},	   
-	   {0,1/3.,2./3.,1},{1,1,1,1,1,1,1},{5,1,1,5}); 							   
+--e1 = g:addNURBS (v1,v2, {{v3:x(),v3:y(),0},{v4:x(),v4:y(),0},{v5:x(),v5:y(),0},{v6:x(),v6:y(),0},{v7:x(),v7:y(),0}},	   
+--	   {0,1/3.,2./3.,1},{1,1,1,1,1,1,1},{5,1,1,5}); 							   
 
-f1 = g:extrude(e1, {0,0,0}, {0,0,1})
+--f1 = g:extrude(e1, {0,0,0}, {0,0,1})