diff --git a/Geo/GModel.cpp b/Geo/GModel.cpp index 308ecf0fbefc7cc0ec91e00ff43f6cdd09de17de..a6c8b05cb93d7fa789c5f33c8eba92616f255b88 100644 --- a/Geo/GModel.cpp +++ b/Geo/GModel.cpp @@ -27,14 +27,11 @@ #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" @@ -185,7 +182,7 @@ void GModel::setFactory(std::string name) if(_factory) delete _factory; _factory = 0; if(name == "Gmsh") { - Msg::Error("Gmsh factory not created yet"); + _factory = new GeoFactory(); } else if(name == "OpenCASCADE"){ #if defined(HAVE_OCC) @@ -1624,97 +1621,6 @@ 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) @@ -1722,6 +1628,12 @@ GFace* GModel::addPlanarFace (std::vector<std::vector<GEdge *> > edges){ return 0; } +GRegion* GModel::addVolume (std::vector<std::vector<GFace *> > faces){ + if(_factory) + return _factory->addVolume(this, faces); + return 0; +} + GEntity *GModel::revolve(GEntity *e, std::vector<double> p1, std::vector<double> p2, double angle) { @@ -2096,6 +2008,9 @@ void GModel::registerBindings(binding *b) cm->setDescription("Merge the file 'filename' in this model, the file can be " "in any format (guessed from the extension) known by gmsh."); cm->setArgNames("filename", NULL); + cm = cb->addMethod("setFactory", &GModel::setFactory); + cm->setDescription("Set the GModel factory: choose between 'Gmsh' or 'OpenCASCADE'."); + cm->setArgNames("name", NULL); cm = cb->addMethod("save", &GModel::save); cm->setDescription("Save this model in the file 'filename'. The content of the " "file depends on the format (guessed from the extension)."); @@ -2167,10 +2082,7 @@ void GModel::registerBindings(binding *b) 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 = cb->addMethod("addVolume", &GModel::addVolume); 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); diff --git a/Geo/GModel.h b/Geo/GModel.h index 41c7b0d30fe188ddf7c841eb8f20f5cecd4c4a74..f4445d0def335bfb0ced1e4d2ad14fb0b72530c9 100644 --- a/Geo/GModel.h +++ b/Geo/GModel.h @@ -392,8 +392,7 @@ 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); + GRegion* addVolume (std::vector<std::vector<GFace*> > faces); // create solid geometry primitives using the factory GEntity *addSphere(double cx, double cy, double cz, double radius); diff --git a/Geo/GModelFactory.cpp b/Geo/GModelFactory.cpp index 0b0928e9f9ad83c6b9c05bdc14097981cb609d8b..6bb549ec12dd5ba649e377581b8a512a3249a0ca 100644 --- a/Geo/GModelFactory.cpp +++ b/Geo/GModelFactory.cpp @@ -13,6 +13,13 @@ #include "OCCEdge.h" #include "OCCFace.h" #include "OCCRegion.h" +#include "ListUtils.h" +#include "Context.h" +#include "GVertex.h" +#include "gmshVertex.h" +#include "gmshEdge.h" +#include "gmshFace.h" +#include "gmshRegion.h" #include "BRepBuilderAPI_MakeVertex.hxx" #include "BRepOffsetAPI_MakePipe.hxx" #include "BRepBuilderAPI_MakeEdge.hxx" @@ -38,6 +45,150 @@ #include <TColStd_HArray1OfInteger.hxx> #include "MLine.h" + +GVertex *GeoFactory::addVertex(GModel *gm, double x, double y, double z, double lc) +{ + + int num = gm->getMaxElementaryNumber(0) + 1; + + x *= CTX::instance()->geom.scalingFactor; + y *= CTX::instance()->geom.scalingFactor; + z *= CTX::instance()->geom.scalingFactor; + lc *= CTX::instance()->geom.scalingFactor; + if(lc == 0.) lc = MAX_LC; // no mesh size given at the point + Vertex *p; + p = Create_Vertex(num, x, y, z, lc, 1.0); + Tree_Add(GModel::current()->getGEOInternals()->Points, &p); + p->Typ = MSH_POINT; + p->Num = num; + + GVertex *v = new gmshVertex(gm, p); + gm->add(v); + + return v; + +} + +GEdge *GeoFactory::addLine(GModel *gm, GVertex *start, GVertex *end) +{ + + int num = gm->getMaxElementaryNumber(1) + 1; + List_T *iList = List_Create(2, 2, sizeof(int)); + int tagBeg = start->tag(); + int tagEnd = end->tag(); + List_Add(iList, &tagBeg); + List_Add(iList, &tagEnd); + + Curve *c = Create_Curve(num, MSH_SEGM_LINE, 1, iList, NULL, + -1, -1, 0., 1.); + Tree_Add(GModel::current()->getGEOInternals()->Curves, &c); + CreateReversedCurve(c); + List_Delete(iList); + c->Typ = MSH_SEGM_LINE; + c->Num = num; + + GEdge *e = new gmshEdge(gm, c, start, end); + gm->add(e); + + return e; + +} + +GFace *GeoFactory::addPlanarFace(GModel *gm, 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 = gm->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 = gm->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(gm,s); + gm->add(gf); + + return gf; + +} + +GRegion* GeoFactory::addVolume (GModel *gm, 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 = gm->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 = gm->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(gm,v); + gm->add(gr); + + return gr; +} + + +//--------------------------------------------------------------------------------- GVertex *OCCFactory::addVertex(GModel *gm, double x, double y, double z, double lc) { if (!gm->_occ_internals) @@ -284,6 +435,12 @@ GEntity *OCCFactory::addSphere(GModel *gm, double xc, double yc, double zc, doub return getOCCRegionByNativePtr(gm, TopoDS::Solid(shape)); } +GRegion* OCCFactory::addVolume (GModel *gm, std::vector<std::vector<GFace *> > faces) +{ + Msg::Error("add Volume not implemented yet for OCCFactory"); + return 0; +} + GEntity *OCCFactory::addCylinder(GModel *gm, std::vector<double> p1, std::vector<double> p2, double radius) { diff --git a/Geo/GModelFactory.h b/Geo/GModelFactory.h index a3aa06fbb285e8f23ae50b64850ce654fae0b9e1..b4a25cf0f83ed59eaf436508a7e8bc42ff692c05 100644 --- a/Geo/GModelFactory.h +++ b/Geo/GModelFactory.h @@ -88,8 +88,57 @@ class GModelFactory { int createNewModel) = 0; virtual GModel *computeBooleanDifference(GModel *obj, GModel*tool, int createNewModel) = 0; + virtual GRegion* addVolume (GModel *gm, std::vector<std::vector<GFace *> > faces) = 0; }; +class GeoFactory : public GModelFactory { + public: + GeoFactory(){} + GVertex *addVertex(GModel *gm,double x, double y, double z, double lc); + virtual GEdge *addLine(GModel *gm,GVertex *v1, GVertex *v2); + GFace *addPlanarFace(GModel *gm, std::vector<std::vector<GEdge *> > edges); + GRegion* addVolume (GModel *gm, std::vector<std::vector<GFace *> > faces); + + //not implemented yet + GEdge *addCircleArc(GModel *gm,const arcCreationMethod &method, + GVertex *start, GVertex *end, + const SPoint3 &aPoint){}; + GEdge *addSpline(GModel *gm,const splineType &type, + GVertex *start, GVertex *end, + std::vector<std::vector<double> > controlPoints){}; + GEdge *addNURBS(GModel *gm, + GVertex *start, GVertex *end, + std::vector<std::vector<double> > controlPoints, + std::vector<double> knots, + std::vector<double> weights, + std::vector<int> multiplicity){}; + GEntity *revolve(GModel *gm, GEntity*,std::vector<double> p1, + std::vector<double> p2, double angle){}; + GEntity *extrude(GModel *gm, GEntity*,std::vector<double> p1, + std::vector<double> p2){}; + GEntity *addPipe(GModel *gm, GEntity *base, std::vector<GEdge *> wire){}; + GEntity *addSphere(GModel *gm,double cx, double cy, double cz, double radius){}; + GEntity *addCylinder(GModel *gm,std::vector<double> p1, std::vector<double> p2, + double radius){}; + std::vector<GFace *> addRuledFaces(GModel *gm, std::vector<std::vector<GEdge *> > edges){}; + GFace *addFace(GModel *gm, std::vector<GEdge *> edges, + std::vector< std::vector<double > > points){}; + GEntity *addTorus(GModel *gm,std::vector<double> p1, std::vector<double> p2, + double radius1, double radius2){}; + GEntity *addBlock(GModel *gm,std::vector<double> p1, std::vector<double> p2){}; + GEntity *addCone(GModel *gm,std::vector<double> p1, std::vector<double> p2, + double radius1, double radius2){}; + void translate(GModel *gm, std::vector<double> dx, int addToTheModel){}; + void rotate(GModel *gm, std::vector<double> p1,std::vector<double> p2, + double angle, int addToTheModel){}; + GModel *computeBooleanUnion(GModel *obj, GModel *tool, int createNewModel){}; + GModel *computeBooleanIntersection(GModel *obj, GModel *tool, int createNewModel){}; + GModel *computeBooleanDifference(GModel *obj, GModel *tool, int createNewModel){}; + void fillet(GModel *gm, std::vector<int> edges, double radius){}; + +}; + + #if defined(HAVE_OCC) class OCCFactory : public GModelFactory { @@ -133,6 +182,7 @@ class OCCFactory : public GModelFactory { GModel *computeBooleanIntersection(GModel *obj, GModel *tool, int createNewModel); GModel *computeBooleanDifference(GModel *obj, GModel *tool, int createNewModel); void fillet(GModel *gm, std::vector<int> edges, double radius); + GRegion* addVolume (GModel *gm, std::vector<std::vector<GFace *> > faces); }; #endif diff --git a/benchmarks/boolean/aorta_bound.lua b/benchmarks/boolean/aorta_bound.lua index d5eb1b92e6c932d8fe0057cf6f05c2e963f9a2d6..ca18233467acf99751e35ecac83e49e76ef9cbfe 100644 --- a/benchmarks/boolean/aorta_bound.lua +++ b/benchmarks/boolean/aorta_bound.lua @@ -6,12 +6,14 @@ 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) +options:numberSet('Mesh', 0, 'CharacteristicLengthExtendFromBoundary', 1) print'*** create GModel from stl ***' model = GModel() model:load ('aortaADAPT.stl') +model:setFactory('Gmsh') + print'*** create Topology ***' model:createTopology() @@ -24,14 +26,14 @@ 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}}) +f1 = model:addPlanarFace({{e1}}) +f2 = model:addPlanarFace({{e2}}) +f3 = model:addPlanarFace({{e3}}) +f4 = model:addPlanarFace({{e4}}) +f5 = model:addPlanarFace({{e5}}) print'*** add Volume ***' -r1 = model:addGeoVolume({{f0,f1,f2,f3,f4,f5}}); +r1 = model:addVolume({{f0,f1,f2,f3,f4,f5}}); print'*** meshing ***' model:mesh(3);