// Copyright (C) 2012-2015 ALNEOS
-// Copyright (C) 2016 EDF R&D
+// Copyright (C) 2016-2021 EDF R&D
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
// See http://www.alneos.com/ or email : contact@alneos.fr
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
#include "GMSHPlugin_Mesher.hxx"
#include "GMSHPlugin_Hypothesis_2D.hxx"
#include <TopTools_MapOfShape.hxx>
#include <TopoDS.hxx>
+#if GMSH_MAJOR_VERSION >=4
+#include <GmshGlobal.h>
+#include <Context.h>
+#endif
+
using namespace std;
+namespace
+{
+ struct ShapeBounds
+ {
+ SBoundingBox3d _bounds;
+ TopoDS_Shape _shape;
+ };
+
+ //================================================================================
+ /*!
+ * \brief Retrieve ShapeBounds from a compound GEdge
+ */
+ //================================================================================
+
+ bool getBoundsOfShapes( GEdge* gEdge,
+ std::vector< ShapeBounds > & topoEdges )
+ {
+ topoEdges.clear();
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ for ( size_t i = 0; i < gEdge->compound.size(); ++i )
+ {
+ GEdge* gE = static_cast< GEdge* >( gEdge->compound[ i ]);
+ topoEdges.push_back( ShapeBounds{ gE->bounds(), *((TopoDS_Edge*)gE->getNativePtr()) });
+ }
+#endif
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION <3
+ for ( size_t i = 0; i < gEdge->_compound.size(); ++i )
+ {
+ GEdge* gE = static_cast< GEdge* >( gEdge->_compound[ i ]);
+ topoEdges.push_back( ShapeBounds{ gE->bounds(), *((TopoDS_Edge*)gE->getNativePtr()) });
+ }
+#endif
+#if GMSH_MAJOR_VERSION <4
+ if ( gEdge->geomType() == GEntity::CompoundCurve )
+ {
+ std::vector<GEdge*> gEdges = ((GEdgeCompound*)gEdge)->getCompounds();
+ for ( size_t i = 0; i < gEdges.size(); ++i )
+ {
+ GEdge* gE = gEdges[ i ];
+ topoEdges.push_back( ShapeBounds{ gE->bounds(), *((TopoDS_Edge*)gE->getNativePtr()) });
+ }
+ }
+#endif
+ return topoEdges.size();
+ }
+
+ //================================================================================
+ /*!
+ * \brief Retrieve ShapeBounds from a compound GFace
+ */
+ //================================================================================
+
+ bool getBoundsOfShapes( GFace* gFace,
+ std::vector< ShapeBounds > & topoFaces )
+ {
+ topoFaces.clear();
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ for ( size_t i = 0; i < gFace->compound.size(); ++i )
+ {
+ GFace* gF = static_cast< GFace* >( gFace->compound[ i ]);
+ topoFaces.push_back( ShapeBounds{ gF->bounds(), *((TopoDS_Face*)gF->getNativePtr()) });
+ }
+#endif
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION <3
+ for ( size_t i = 0; i < gFace->_compound.size(); ++i )
+ {
+ GFace* gF = static_cast< GFace* >( gFace->_compound[ i ]);
+ topoFaces.push_back( ShapeBounds{ gF->bounds(), *((TopoDS_Face*)gF->getNativePtr()) });
+ }
+#endif
+#if GMSH_MAJOR_VERSION <4
+ if ( gFace->geomType() == GEntity::CompoundSurface )
+ {
+ std::list<GFace*> gFaces = ((GFaceCompound*)gFace)->getCompounds();
+ for ( std::list<GFace*>::const_iterator itl = gFaces.begin();itl != gFaces.end(); ++itl )
+ {
+ GFace* gF = *itl;
+ topoFaces.push_back( ShapeBounds{ gF->bounds(), *((TopoDS_Face*)gF->getNativePtr()) });
+ }
+ }
+#endif
+ return topoFaces.size();
+ }
+ //================================================================================
+ /*!
+ * \brief Find a shape whose bounding box includes a given point
+ */
+ //================================================================================
+
+ TopoDS_Shape getShapeAtPoint( const SPoint3& point, const std::vector< ShapeBounds > & shapes )
+ {
+ TopoDS_Shape shape;
+ float distmin = std::numeric_limits<float>::max();
+ for ( size_t i = 0; i < shapes.size(); ++i )
+ {
+ float dist = GMSHPlugin_Mesher::DistBoundingBox( shapes[i]._bounds, point );
+ if (dist < distmin)
+ {
+ shape = shapes[i]._shape;
+ distmin = dist;
+ if ( distmin == 0. )
+ break;
+ }
+ }
+ return shape;
+ }
+}
+
//=============================================================================
/*!
*
//=============================================================================
GMSHPlugin_Mesher::GMSHPlugin_Mesher (SMESH_Mesh* mesh,
- const TopoDS_Shape& aShape)
+ const TopoDS_Shape& aShape)
: _mesh (mesh),
_shape (aShape)
-{
+{
// il faudra peut être mettre un truc par defaut si l'utilisateur ne rentre rien en para
//defaultParameters();
}
_secondOrder = false;
_useIncomplElem = true;
_is2d = false;
+ _compounds.clear();
}
}
printf("We chose _useIncomplElem %d \n", (_useIncomplElem)?1:0);
printf("We are in dimension %d \n", (_is2d)?2:3);
//*/
-
+
std::map <int,double> mapAlgo2d;
- mapAlgo2d[0]=2; mapAlgo2d[1]=1; mapAlgo2d[2]=5; mapAlgo2d[3]=6; mapAlgo2d[4]=8; mapAlgo2d[5]=9;
+ mapAlgo2d[0]=2; // Automatic
+ mapAlgo2d[1]=1; // MeshAdapt
+ mapAlgo2d[2]=5; // Delaunay
+ mapAlgo2d[3]=6; // Frontal-Delaunay
+ mapAlgo2d[4]=8; // DelQuad (Frontal-Delaunay for Quads)
+ mapAlgo2d[5]=9; // Packing of parallelograms
+
std::map <int,double> mapAlgo3d;
- mapAlgo3d[0]=1; mapAlgo3d[1]=4; mapAlgo3d[2]=5; mapAlgo3d[3]=6; mapAlgo3d[4]=7; mapAlgo3d[5]=9;
+ mapAlgo3d[0]=1; // Delaunay
+ mapAlgo3d[1]=4; // Frontal
+ mapAlgo3d[2]=7; // MMG3D
+ mapAlgo3d[3]=9; // R-tree
int ok;
ok = GmshSetOption("Mesh", "Algorithm" , mapAlgo2d[_algo2d]) ;
ASSERT(ok);
if ( !_is2d)
- {
- ok = GmshSetOption("Mesh", "Algorithm3D" , mapAlgo2d[_algo3d]) ;
+ {
+ ok = GmshSetOption("Mesh", "Algorithm3D" , mapAlgo3d[_algo3d]) ;
ASSERT(ok);
- }
+ }
ok = GmshSetOption("Mesh", "RecombinationAlgorithm" , (double)_recomb2DAlgo) ;
ASSERT(ok);
ok = GmshSetOption("Mesh", "RecombineAll" , (_recombineAll)?1.:0.) ;
ok = GmshSetOption("Mesh", "SubdivisionAlgorithm" , (double)_subdivAlgo) ;
ASSERT(ok);
ok = GmshSetOption("Mesh", "RemeshAlgorithm" , (double)_remeshAlgo) ;
- ASSERT(ok);
+ //ASSERT(ok);
ok = GmshSetOption("Mesh", "RemeshParametrization" , (double)_remeshPara) ;
- ASSERT(ok);
+ //ASSERT(ok);
ok = GmshSetOption("Mesh", "Smoothing" , (double)_smouthSteps) ;
- ASSERT(ok);
+ //ASSERT(ok);
ok = GmshSetOption("Mesh", "CharacteristicLengthFactor", _sizeFactor) ;
- ASSERT(ok);
+ //ASSERT(ok);
ok = GmshSetOption("Mesh", "CharacteristicLengthMin" , _minSize) ;
ASSERT(ok);
ok = GmshSetOption("Mesh", "CharacteristicLengthMax" , _maxSize) ;
ok = GmshSetOption("Mesh", "ElementOrder" , (_secondOrder)?2.:1.) ;
ASSERT(ok);
if (_secondOrder)
- {
+ {
ok = GmshSetOption("Mesh", "SecondOrderIncomplete" ,(_useIncomplElem)?1.:0.);
ASSERT(ok);
- }
+ }
}
//================================================================================
void GMSHPlugin_Mesher::CreateGmshCompounds()
{
MESSAGE("GMSHPlugin_Mesher::CreateGmshCompounds");
-
+
SMESH_Gen_i* smeshGen_i = SMESH_Gen_i::GetSMESHGen();
-
+
OCC_Internals* occgeo = _gModel->getOCCInternals();
-
+ bool toSynchronize = false;
+
for(std::set<std::string>::const_iterator its = _compounds.begin();its != _compounds.end(); ++its )
{
GEOM::GEOM_Object_var aGeomObj;
TopoDS_Shape geomShape = TopoDS_Shape();
- SALOMEDS::SObject_var aSObj = SMESH_Gen_i::getStudyServant()->FindObjectID( (*its).c_str() );
+ SALOMEDS::SObject_var aSObj = SMESH_Gen_i::GetSMESHGen()->getStudyServant()->FindObjectID( (*its).c_str() );
SALOMEDS::GenericAttribute_var anAttr;
if (!aSObj->_is_nil() && aSObj->FindAttribute(anAttr, "AttributeIOR"))
{
SALOMEDS::AttributeIOR_var anIOR = SALOMEDS::AttributeIOR::_narrow(anAttr);
CORBA::String_var aVal = anIOR->Value();
- CORBA::Object_var obj = SMESH_Gen_i::getStudyServant()->ConvertIORToObject(aVal);
+ CORBA::Object_var obj = SMESH_Gen_i::GetSMESHGen()->getStudyServant()->ConvertIORToObject(aVal);
aGeomObj = GEOM::GEOM_Object::_narrow(obj);
}
- if ( !aGeomObj->_is_nil() )
- geomShape = smeshGen_i->GeomObjectToShape( aGeomObj.in() );
-
+ geomShape = smeshGen_i->GeomObjectToShape( aGeomObj.in() );
+ if ( geomShape.IsNull() )
+ continue;
+
TopAbs_ShapeEnum geomType = geomShape.ShapeType();
if ( geomType == TopAbs_COMPOUND)// voir s'il ne faut pas mettre une erreur dans le cas contraire
{
MESSAGE("shapeType == TopAbs_COMPOUND");
TopoDS_Iterator it(geomShape);
+ if ( !it.More() )
+ continue;
TopAbs_ShapeEnum shapeType = it.Value().ShapeType();
+#if GMSH_MAJOR_VERSION >=4
+ std::vector< std::pair< int, int > > dimTags;
+ for ( ; it.More(); it.Next())
+ {
+ const TopoDS_Shape& topoShape = it.Value();
+ ASSERT(topoShape.ShapeType() == shapeType);
+ if ( _mesh->GetMeshDS()->ShapeToIndex( topoShape ) > 0 )
+ occgeo->importShapes( &topoShape, false, dimTags );
+ else
+ {
+ TopoDS_Shape face = TopExp_Explorer( _shape, shapeType ).Current();
+ SMESH_subMesh* sm = _mesh->GetSubMesh( face );
+ sm->GetComputeError() =
+ SMESH_ComputeError::New
+ ( COMPERR_WARNING, "Compound shape does not belong to the main geometry. Ignored");
+ }
+ }
+ std::vector<int> tags;
+ int dim = ( shapeType == TopAbs_EDGE ) ? 1 : 2;
+ for ( size_t i = 0; i < dimTags.size(); ++i )
+ {
+ if ( dimTags[i].first == dim )
+ tags.push_back( dimTags[i].second );
+ }
+ if ( !tags.empty() )
+ {
+ _gModel->getGEOInternals()->setCompoundMesh( dim, tags );
+ toSynchronize = true;
+ }
+#else
// compound of edges
if (shapeType == TopAbs_EDGE)
{
ASSERT(topoShape.ShapeType() == shapeType);
curve->compound.push_back(occgeo->addEdgeToModel(_gModel, (TopoDS_Edge&)topoShape)->tag());
}
+ toSynchronize = true;
Tree_Add(_gModel->getGEOInternals()->Curves, &curve);
//_gModel->importGEOInternals();
}
ASSERT(topoShape.ShapeType() == shapeType);
surface->compound.push_back(occgeo->addFaceToModel(_gModel, (TopoDS_Face&)topoShape)->tag());
}
+ toSynchronize = true;
Tree_Add(_gModel->getGEOInternals()->Surfaces, &surface);
- _gModel->getGEOInternals()->synchronize(_gModel);
}
+#endif
+ if ( toSynchronize )
+ _gModel->getGEOInternals()->synchronize(_gModel);
}
}
}
+//================================================================================
+/*!
+ * \brief For a compound mesh set the mesh components to be transmitted to SMESH
+ */
+//================================================================================
+
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+void GMSHPlugin_Mesher::SetCompoundMeshVisibility()
+{
+
+ // Loop over all faces, if the face belongs to a compound entry then
+ // for all (boundary) edges whithin the face visibility is set to 0,
+ // if the face doesn't belong to a compound entry then visibility is
+ // set to 1 for all its (boundary) edges. Later, in FillSMesh() func
+ // getVisibility() (returns either 1 or 0) is used to decide weather
+ // the mesh of edge should be transmitted to SMESH or not.
+
+ for ( GModel::fiter itF = _gModel->firstFace(); itF != _gModel->lastFace(); ++itF )
+ {
+ std::vector< GEdge *> faceEdges = (*itF)->edges();
+
+ for ( auto itE = faceEdges.begin(); itE != faceEdges.end(); ++itE )
+ {
+ if ( ((*itF)->compound.size()) )
+ (*itE)->setVisibility(0);
+ else
+ (*itE)->setVisibility(1);
+ }
+ }
+
+
+ // Loop over all edges, if the edge belongs to a compound entry then
+ // for all (boundary) vertices whithin the edge visibility is set to
+ // 0, if the edge doesn't belong to a compound entry then visibility
+ // is set to 1 for all its (boundary) vertices. Later, in FillSMesh()
+ // func getVisibility() (returns either 1 or 0) is used to decide we-
+ // ather the mesh of vertices should be transmitted to SMESH or not.
+
+ for ( GModel::eiter itE = _gModel->firstEdge(); itE != _gModel->lastEdge(); ++itE )
+ {
+ std::vector<GVertex *> bndVerticies = (*itE)->vertices();
+
+ for( auto itV = bndVerticies.begin(); itV != bndVerticies.end(); ++itV )
+ {
+ if(((*itE)->compound.size()))
+ (*itV)->setVisibility(0);
+ else
+ (*itV)->setVisibility(1);
+ }
+ }
+
+}
+#endif
+
//================================================================================
/*!
* \brief Write mesh from GModel instance to SMESH instance
void GMSHPlugin_Mesher::FillSMesh()
{
- MESSAGE("GMSHPlugin_Mesher::FillSMesh");
- // GET MESH FROM SMESH
SMESHDS_Mesh* meshDS = _mesh->GetMeshDS();
-
+
// ADD 0D ELEMENTS
- for(GModel::viter it = _gModel->firstVertex(); it != _gModel->lastVertex(); ++it){
+ for ( GModel::viter it = _gModel->firstVertex(); it != _gModel->lastVertex(); ++it)
+ {
GVertex *gVertex = *it;
-
+
// GET topoVertex CORRESPONDING TO gVertex
TopoDS_Vertex topoVertex = *((TopoDS_Vertex*)gVertex->getNativePtr());
-
+
if (gVertex->getVisibility() == 0) // belongs to a compound
{
SMESH_subMesh* sm = _mesh->GetSubMesh(topoVertex);
sm->SetIsAlwaysComputed(true); // prevent from displaying errors
continue;
}
-
+
// FILL SMESH FOR topoVertex
//nodes
for(unsigned int i = 0; i < gVertex->mesh_vertices.size(); i++)
meshDS->SetNodeOnVertex( node, topoVertex );
}
}
+ // WE DON'T ADD 0D ELEMENTS because it does not follow the salome meshers philosophy
//elements
- for(unsigned int i = 0; i < gVertex->getNumMeshElements(); i++)
- {
- MElement *e = gVertex->getMeshElement(i);
- std::vector<MVertex*> verts;
- e->getVertices(verts);
- ASSERT(verts.size()==1);
- //SMDS_Mesh0DElement* zeroDElement = 0;
- // WE DON'T ADD 0D ELEMENTS because it does not follow the salome meshers philosophy
- //zeroDElement = meshDS->Add0DElementWithID(verts[0]->getNum(),e->getNum());
- //meshDS->SetMeshElementOnShape(zeroDElement, topoVertex);
- }
+ // for(unsigned int i = 0; i < gVertex->getNumMeshElements(); i++)
+ // {
+ // MElement *e = gVertex->getMeshElement(i);
+ // std::vector<MVertex*> verts;
+ // e->getVertices(verts);
+ // ASSERT(verts.size()==1);
+ // SMDS_Mesh0DElement* zeroDElement = 0;
+ // zeroDElement = meshDS->Add0DElementWithID(verts[0]->getNum(),e->getNum());
+ // meshDS->SetMeshElementOnShape(zeroDElement, topoVertex);
+ // }
}
-
+
// ADD 1D ELEMENTS
- for(GModel::eiter it = _gModel->firstEdge(); it != _gModel->lastEdge(); ++it){
+ for(GModel::eiter it = _gModel->firstEdge(); it != _gModel->lastEdge(); ++it)
+ {
GEdge *gEdge = *it;
-
- // GET topoEdge CORRESPONDING TO gEdge (map if compound)
+
+ // GET topoEdge CORRESPONDING TO gEdge
TopoDS_Edge topoEdge;
- std::map<GEdge*,TopoDS_Edge> topoEdges;
-
- if(gEdge->geomType() != GEntity::CompoundCurve)
+ std::vector< ShapeBounds > topoEdges;
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ if(gEdge->haveParametrization())
+#else
+ if ( gEdge->geomType() != GEntity::CompoundCurve )
+#endif
{
topoEdge = *((TopoDS_Edge*)gEdge->getNativePtr());
if (gEdge->getVisibility() == 0) // belongs to a compound
continue;
}
}
- else
- {
- // compound case, map construction GEdge/TopoDS_Edge
- std::vector<GEdge*> gEdges = ((GEdgeCompound*)gEdge)->getCompounds();
- for(std::vector<GEdge*>::const_iterator itv = gEdges.begin();itv != gEdges.end(); ++itv)
- {
- topoEdges.insert( pair<GEdge*,TopoDS_Edge>(*itv,*((TopoDS_Edge*)(*itv)->getNativePtr())) );
- }
- }
-
+ bool isCompound = getBoundsOfShapes( gEdge, topoEdges );
+
// FILL SMESH FOR topoEdge
//nodes
- for(unsigned int i = 0; i < gEdge->mesh_vertices.size(); i++)
+ for ( size_t i = 0; i < gEdge->mesh_vertices.size(); i++ )
{
MVertex *v = gEdge->mesh_vertices[i];
- if(v->getIndex() >= 0)
+ if ( v->getIndex() >= 0 )
{
SMDS_MeshNode *node = meshDS->AddNodeWithID(v->x(),v->y(),v->z(),v->getNum());
-
- if (topoEdges.size() != 0)
- {
- // get back corresponding topoEdge in compound case
- topoEdge.Nullify();
- SPoint3 point = v->point();
- float distmin = std::numeric_limits<float>::max();
- for(map<GEdge*,TopoDS_Edge>::const_iterator itm=topoEdges.begin();itm != topoEdges.end(); itm++)
- {
- SBoundingBox3d bounds = (*itm).first->bounds();
- float dist = DistBoundingBox(bounds,point);
- if (dist < distmin)
- {
- topoEdge = (*itm).second;
- distmin = dist;
- if (distmin == 0.)
- break;
- }
- }
- }
-
+
+ if ( isCompound )
+ topoEdge = TopoDS::Edge( getShapeAtPoint( v->point(), topoEdges ));
+
meshDS->SetNodeOnEdge( node, topoEdge );
}
}
+ }
+
+ for ( GModel::eiter it = _gModel->firstEdge(); it != _gModel->lastEdge(); ++it )
+ {
+ GEdge *gEdge = *it;
+ if ( gEdge->getVisibility() == 0) // belongs to a compound
+ continue;
+
+ TopoDS_Edge topoEdge;
+ std::vector< ShapeBounds > topoEdges;
+ bool isCompound = getBoundsOfShapes( gEdge, topoEdges );
+ if ( !isCompound )
+ topoEdge = *((TopoDS_Edge*)gEdge->getNativePtr());
+
//elements
- for(unsigned int i = 0; i < gEdge->getNumMeshElements(); i++)
+ std::vector<MVertex*> verts(3);
+ for ( size_t i = 0; i < gEdge->getNumMeshElements(); i++ )
{
MElement *e = gEdge->getMeshElement(i);
- std::vector<MVertex*> verts;
+ verts.clear();
e->getVertices(verts);
- SMDS_MeshEdge* edge = 0;
-
- if (topoEdges.size() !=0)
- {
- // get back corresponding topoEdge in compound case
- topoEdge.Nullify();
- SPoint3 point = e->barycenter();
- float distmin = std::numeric_limits<float>::max();
- for(map<GEdge*,TopoDS_Edge>::const_iterator itm=topoEdges.begin();itm != topoEdges.end(); itm++)
- {
- SBoundingBox3d bounds = (*itm).first->bounds();
- float dist = DistBoundingBox(bounds,point);
- if (dist < distmin)
- {
- topoEdge = (*itm).second;
- distmin = dist;
- if (distmin == 0.)
- break;
- }
- }
- }
-
+
// if a node wasn't set, it is assigned here
- for (unsigned j = 0; j < verts.size(); j++)
+ for ( size_t j = 0; j < verts.size(); j++ )
{
- if(verts[j]->onWhat()->getVisibility() == 0)
+ if ( verts[j]->onWhat()->getVisibility() == 0 )
{
SMDS_MeshNode *node = meshDS->AddNodeWithID(verts[i]->x(),verts[j]->y(),verts[j]->z(),verts[j]->getNum());
meshDS->SetNodeOnEdge( node, topoEdge );
verts[j]->setEntity(gEdge);
}
}
-
+
+ SMDS_MeshEdge* edge = 0;
switch (verts.size())
{
case 2:
ASSERT(false);
continue;
}
- meshDS->SetMeshElementOnShape(edge, topoEdge);
+ if ( isCompound )
+ topoEdge = TopoDS::Edge( getShapeAtPoint( e->barycenter(), topoEdges ));
+
+ meshDS->SetMeshElementOnShape( edge, topoEdge );
}
}
// ADD 2D ELEMENTS
- for(GModel::fiter it = _gModel->firstFace(); it != _gModel->lastFace(); ++it){
+ for ( GModel::fiter it = _gModel->firstFace(); it != _gModel->lastFace(); ++it)
+ {
GFace *gFace = *it;
-
- // GET topoFace CORRESPONDING TO gFace (map if compound)
+
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ // Gmsh since version 4.3 is now producing extra surface and mesh when
+ // compounds are involved. Since in Gmsh meshing procedure needs acess
+ // to each of the original topology and the meshed topology. Hence we
+ // bypass the additional mesh in case of compounds. Note, similar cri-
+ // teria also occus in the following 'for' loop.
+ if ( _compounds.size() && gFace->geomType() == GEntity::DiscreteSurface )
+ continue;
+#endif
+
+ // GET topoFace CORRESPONDING TO gFace
TopoDS_Face topoFace;
- std::map<GFace*,TopoDS_Face> topoFaces;
-
- if(gFace->geomType() != GEntity::CompoundSurface)
+ std::vector< ShapeBounds > topoFaces;
+
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ if(gFace->haveParametrization())
+#else
+ if ( gFace->geomType() != GEntity::CompoundSurface )
+#endif
{
topoFace = *((TopoDS_Face*)gFace->getNativePtr());
if (gFace->getVisibility() == 0) // belongs to a compound
continue;
}
}
- else
- {
- // compound case, map construction GFace/TopoDS_Face
- std::list<GFace*> gFaces = ((GFaceCompound*)gFace)->getCompounds();
- for(std::list<GFace*>::const_iterator itl = gFaces.begin();itl != gFaces.end(); ++itl)
- {
- topoFaces.insert( pair<GFace*,TopoDS_Face>(*itl,*((TopoDS_Face*)(*itl)->getNativePtr())) );
- }
- }
-
+ bool isCompound = getBoundsOfShapes( gFace, topoFaces );
+
// FILL SMESH FOR topoFace
//nodes
- for(unsigned int i = 0; i < gFace->mesh_vertices.size(); i++)
+ for ( size_t i = 0; i < gFace->mesh_vertices.size(); i++ )
{
MVertex *v = gFace->mesh_vertices[i];
- if(v->getIndex() >= 0)
+ if ( v->getIndex() >= 0 )
{
- SMDS_MeshNode *node = meshDS->AddNodeWithID(v->x(),v->y(),v->z(),v->getNum());
-
- if (topoFaces.size() != 0)
- {
- // get back corresponding topoFace in compound case
- topoFace.Nullify();
- SPoint3 point = v->point();
- float distmin = std::numeric_limits<float>::max();
- for(map<GFace*,TopoDS_Face>::const_iterator itm=topoFaces.begin();itm != topoFaces.end(); itm++)
- {
- SBoundingBox3d bounds = (*itm).first->bounds();
- float dist = DistBoundingBox(bounds,point);
- if (dist < distmin)
- {
- topoFace = (*itm).second;
- distmin = dist;
- if (distmin == 0.)
- break;
- }
- }
- }
-
+ SMDS_MeshNode *node = meshDS->AddNodeWithID(v->x(),v->y(),v->z(),v->getNum());
+
+ if ( isCompound )
+ topoFace = TopoDS::Face( getShapeAtPoint( v->point(), topoFaces ));
+
meshDS->SetNodeOnFace( node, topoFace );
}
}
+ }
+
+ for ( GModel::fiter it = _gModel->firstFace(); it != _gModel->lastFace(); ++it)
+ {
+ GFace *gFace = *it;
+
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ if ( _compounds.size() && gFace->geomType() == GEntity::DiscreteSurface )
+ continue;
+
+ bool isCompound = (!gFace->haveParametrization());
+#else
+ bool isCompound = ( gFace->geomType() == GEntity::CompoundSurface );
+#endif
+
+ if ( !isCompound && gFace->getVisibility() == 0 )
+ continue; // belongs to a compound
+
+ TopoDS_Face topoFace;
+ std::vector< ShapeBounds > topoFaces;
+ if ( isCompound )
+ getBoundsOfShapes( gFace, topoFaces );
+ else
+ topoFace = *((TopoDS_Face*)gFace->getNativePtr());
+
//elements
- for(unsigned int i = 0; i < gFace->getNumMeshElements(); i++)
+ std::vector<MVertex*> verts;
+ for ( size_t i = 0; i < gFace->getNumMeshElements(); i++ )
{
MElement *e = gFace->getMeshElement(i);
- std::vector<MVertex*> verts;
+ verts.clear();
e->getVertices(verts);
SMDS_MeshFace* face = 0;
-
- if (topoFaces.size() !=0)
- {
- // get back corresponding topoFace in compound case
- topoFace.Nullify();
- SPoint3 point = e->barycenter();
- float distmin = std::numeric_limits<float>::max();
- for(map<GFace*,TopoDS_Face>::const_iterator itm=topoFaces.begin();itm != topoFaces.end(); itm++)
- {
- SBoundingBox3d bounds = (*itm).first->bounds();
- float dist = DistBoundingBox(bounds,point);
- if (dist < distmin)
- {
- topoFace = (*itm).second;
- distmin = dist;
- if (distmin == 0.)
- break;
- }
- }
- }
-
+
// if a node wasn't set, it is assigned here
- for (unsigned j = 0; j < verts.size(); j++)
+ for ( size_t j = 0; j < verts.size(); j++)
{
if(verts[j]->onWhat()->getVisibility() == 0)
{
ASSERT(false);
continue;
}
+
+ if ( isCompound )
+ topoFace = TopoDS::Face( getShapeAtPoint( e->barycenter(), topoFaces ));
+
meshDS->SetMeshElementOnShape(face, topoFace);
}
}
-
+
// ADD 3D ELEMENTS
- for(GModel::riter it = _gModel->firstRegion(); it != _gModel->lastRegion(); ++it){
+ for ( GModel::riter it = _gModel->firstRegion(); it != _gModel->lastRegion(); ++it)
+ {
GRegion *gRegion = *it;
if (gRegion->getVisibility() == 0)
continue;
-
+
// GET topoSolid CORRESPONDING TO gRegion
TopoDS_Solid topoSolid = *((TopoDS_Solid*)gRegion->getNativePtr());
-
- //printf("volume %d (%d) contains %d elements\n", meshDS->ShapeToIndex(topoSolid), gRegion->tag(), gRegion->getNumMeshElements());
+
// FILL SMESH FOR topoSolid
-
+
//nodes
for(unsigned int i = 0; i < gRegion->mesh_vertices.size(); i++)
{
meshDS->SetNodeInVolume( node, topoSolid );
}
}
-
+
//elements
+ std::vector<MVertex*> verts;
for(unsigned int i = 0; i < gRegion->getNumMeshElements(); i++)
{
MElement *e = gRegion->getMeshElement(i);
- std::vector<MVertex*> verts;
+ verts.clear();
e->getVertices(verts);
SMDS_MeshVolume* volume = 0;
switch (verts.size()){
meshDS->SetMeshElementOnShape(volume, topoSolid);
}
}
-
+
//return 0;
}
*/
//================================================================================
-float GMSHPlugin_Mesher::DistBoundingBox(SBoundingBox3d& bounds,SPoint3& point)
+float GMSHPlugin_Mesher::DistBoundingBox(const SBoundingBox3d& bounds, const SPoint3& point)
{
SPoint3 min = bounds.min();
SPoint3 max = bounds.max();
-
+
float x,y,z;
-
+
if (point.x() < min.x())
x = min.x()-point.x();
else if (point.x() > max.x())
x = point.x()-max.x();
else
x = 0.;
-
+
if (point.y() < min.y())
y = min.y()-point.y();
else if (point.y() > max.y())
y = point.y()-max.y();
else
y = 0.;
-
+
if (point.z() < min.z())
z = min.z()-point.z();
else if (point.z() > max.z())
z = point.z()-max.z();
else
z = 0.;
-
- return sqrt(x*x+y*y+z*z);
+
+ return x*x+y*y+z*z;
}
//================================================================================
/*!
{
//MESSAGE("level="<< level.c_str() << ", msg=" << msg.c_str()<< "\n");
printf("level=%s msg=%s\n", level.c_str(), msg.c_str());
-
+
if(level == "Fatal" || level == "Error")
{
std::ostringstream oss;
//comment << SMESH_Comment(oss.str);
//std::string str = oss.str();
//smError.reset( new SMESH_ComputeError( str ));
-
+
// plutot que de faire de la merde ici, on pourait simplement
// remplir une liste pour dire sur quelles entités gmsh se plante
// (puis faire le fillsmesh)
bool GMSHPlugin_Mesher::Compute()
{
MESSAGE("GMSHPlugin_Mesher::Compute");
-
+
int err = 0;
-
+
GmshInitialize();
SetGmshOptions();
_gModel = new GModel();
err = 1;
MESSAGE("Unrecoverable error during Generation of Gmsh Mesh");
}
-
+
if (!err)
{
+#if GMSH_MAJOR_VERSION < 4
if (_compounds.size() > 0) _gModel->setCompoundVisibility();
+#endif
+#if GMSH_MAJOR_VERSION >=4 && GMSH_MINOR_VERSION >=3
+ if (_compounds.size() > 0)
+ SetCompoundMeshVisibility();
+#endif
FillSMesh();
}
delete _gModel;
