X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESH%2FSMESH_MesherHelper.cxx;h=3b5483b13e49c1a72d5677277d0a49067b6a063d;hp=ed97f5cc202afc547aa077a25315b56ab80eda44;hb=HEAD;hpb=624e57e97d4aaea94f4318200b82f7756f27453a diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index ed97f5cc2..49a115281 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2020 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2024 CEA, EDF, OPEN CASCADE // // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS @@ -28,7 +28,7 @@ #include "SMDS_EdgePosition.hxx" #include "SMDS_FaceOfNodes.hxx" -#include "SMDS_FacePosition.hxx" +#include "SMDS_FacePosition.hxx" #include "SMDS_IteratorOnIterators.hxx" #include "SMDS_VolumeTool.hxx" #include "SMESHDS_Mesh.hxx" @@ -51,6 +51,7 @@ #include #include #include +#include #include #include #include @@ -99,7 +100,7 @@ SMESH_MesherHelper::SMESH_MesherHelper(SMESH_Mesh& theMesh) //======================================================================= //function : ~SMESH_MesherHelper -//purpose : +//purpose : //======================================================================= SMESH_MesherHelper::~SMESH_MesherHelper() @@ -424,7 +425,7 @@ bool SMESH_MesherHelper::GetNodeUVneedInFaceNode(const TopoDS_Face& F) const //======================================================================= //function : IsMedium -//purpose : +//purpose : //======================================================================= bool SMESH_MesherHelper::IsMedium(const SMDS_MeshNode* node, @@ -582,7 +583,7 @@ bool SMESH_MesherHelper::toCheckPosOnShape(int shapeID ) const void SMESH_MesherHelper::setPosOnShapeValidity(int shapeID, bool ok ) const { - std::map< int,bool >::iterator sh_ok = + std::map< int,bool >::iterator sh_ok = ((SMESH_MesherHelper*)this)->myNodePosShapesValidity.insert( make_pair( shapeID, ok)).first; if ( !ok ) sh_ok->second = ok; @@ -590,7 +591,7 @@ void SMESH_MesherHelper::setPosOnShapeValidity(int shapeID, bool ok ) const //======================================================================= //function : ToFixNodeParameters -//purpose : Enables fixing node parameters on EDGEs and FACEs in +//purpose : Enables fixing node parameters on EDGEs and FACEs in // GetNodeU(...,check=true), GetNodeUV(...,check=true), CheckNodeUV() and // CheckNodeU() in case if a node lies on a shape set via SetSubShape(). // Default is False @@ -658,10 +659,10 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, // corresponding EDGE from FACE, get pcurve for this // EDGE and retrieve value from this pcurve SMDS_EdgePositionPtr epos = pos; - const int edgeID = n->getshapeId(); - const TopoDS_Edge& E = TopoDS::Edge( GetMeshDS()->IndexToShape( edgeID )); + const int edgeID = n->getshapeId(); + const TopoDS_Edge& E = TopoDS::Edge( GetMeshDS()->IndexToShape( edgeID )); double f, l, u = epos->GetUParameter(); - Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface( E, F, f, l ); + Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface( E, F, f, l ); bool validU = ( !C2d.IsNull() && ( f < u ) && ( u < l )); if ( validU ) uv = C2d->Value( u ); else uv.SetCoord( Precision::Infinite(),0.); @@ -769,6 +770,17 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, if ( isSeam ) uv = getUVOnSeam( uv, GetNodeUV( F, n2, 0 )); } + else if ( myParIndex && n2 ) + { + gp_Pnt2d oldUV = uv; + gp_Pnt2d uv2 = GetNodeUV( F, n2, 0 ); + if ( myParIndex & 1 ) + uv.SetX( uv.X() + ShapeAnalysis::AdjustToPeriod( uv.X(), myPar1[0], myPar2[0])); + if ( myParIndex & 2 ) + uv.SetY( uv.Y() + ShapeAnalysis::AdjustToPeriod( uv.Y(), myPar1[1], myPar2[1])); + if ( uv.SquareDistance( uv2 ) > oldUV.SquareDistance( uv2 )) + uv = oldUV; + } } } else @@ -864,6 +876,51 @@ GeomAPI_ProjectPointOnSurf& SMESH_MesherHelper::GetProjector(const TopoDS_Face& return *( i_proj->second ); } +//======================================================================= +//function : GetProjector +//purpose : Return projector initialized by given face, which is returned +//======================================================================= + +GeomAPI_ProjectPointOnSurf& SMESH_MesherHelper::GetProjector(const TopoDS_Face& F, + double tol ) const +{ + Handle(Geom_Surface) surface = BRep_Tool::Surface( F ); + int faceID = GetMeshDS()->ShapeToIndex( F ); + TID2ProjectorOnSurf& i2proj = const_cast< TID2ProjectorOnSurf&>( myFace2Projector ); + TID2ProjectorOnSurf::iterator i_proj = i2proj.find( faceID ); + if ( i_proj == i2proj.end() ) + { + if ( tol == 0 ) tol = BRep_Tool::Tolerance( F ); + double U1, U2, V1, V2; + surface->Bounds(U1, U2, V1, V2); + GeomAPI_ProjectPointOnSurf* proj = new GeomAPI_ProjectPointOnSurf(); + proj->Init( surface, U1, U2, V1, V2, tol ); + i_proj = i2proj.insert( make_pair( faceID, proj )).first; + } + return *( i_proj->second ); +} + +//======================================================================= +//function : GetPCProjector +//purpose : Return projector initialized by given EDGE +//======================================================================= + +GeomAPI_ProjectPointOnCurve& SMESH_MesherHelper::GetPCProjector(const TopoDS_Edge& E ) const +{ + int edgeID = GetMeshDS()->ShapeToIndex( E ); + TID2ProjectorOnCurve& i2proj = const_cast< TID2ProjectorOnCurve&>( myEdge2Projector ); + TID2ProjectorOnCurve::iterator i_proj = i2proj.insert( make_pair( edgeID, nullptr )).first; + if ( !i_proj->second ) + { + double f,l; + Handle(Geom_Curve) curve = BRep_Tool::Curve( E,f,l ); + i_proj->second = new GeomAPI_ProjectPointOnCurve(); + i_proj->second->Init( curve, f, l ); + } + GeomAPI_ProjectPointOnCurve* projector = i_proj->second; + return *projector; +} + //======================================================================= //function : GetSurface //purpose : Return a cached ShapeAnalysis_Surface of a FACE @@ -886,7 +943,7 @@ namespace { gp_XY AverageUV(const gp_XY& uv1, const gp_XY& uv2) { return ( uv1 + uv2 ) / 2.; } gp_XY_FunPtr(Added); // define gp_XY_Added pointer to function calling gp_XY::Added(gp_XY) - gp_XY_FunPtr(Subtracted); + gp_XY_FunPtr(Subtracted); } //======================================================================= @@ -910,9 +967,9 @@ gp_XY SMESH_MesherHelper::ApplyIn2D(Handle(Geom_Surface) surface, return fun(uv1,uv2); // move uv2 not far than half-period from uv1 - double u2 = + double u2 = uv2.X()+(isUPeriodic ? ShapeAnalysis::AdjustByPeriod(uv2.X(),uv1.X(),surface->UPeriod()) :0); - double v2 = + double v2 = uv2.Y()+(isVPeriodic ? ShapeAnalysis::AdjustByPeriod(uv2.Y(),uv1.Y(),surface->VPeriod()) :0); // execute operation @@ -981,8 +1038,8 @@ gp_XY SMESH_MesherHelper::GetMiddleUV(const Handle(Geom_Surface)& surface, //======================================================================= gp_XY SMESH_MesherHelper::GetCenterUV(const gp_XY& uv1, - const gp_XY& uv2, - const gp_XY& uv3, + const gp_XY& uv2, + const gp_XY& uv3, const gp_XY& uv12, const gp_XY& uv23, const gp_XY& uv31, @@ -1110,27 +1167,17 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, { setPosOnShapeValidity( shapeID, false ); // u incorrect, project the node to the curve - int edgeID = GetMeshDS()->ShapeToIndex( E ); - TID2ProjectorOnCurve& i2proj = const_cast< TID2ProjectorOnCurve&>( myEdge2Projector ); - TID2ProjectorOnCurve::iterator i_proj = - i2proj.insert( make_pair( edgeID, (GeomAPI_ProjectPointOnCurve*) 0 )).first; - if ( !i_proj->second ) - { - i_proj->second = new GeomAPI_ProjectPointOnCurve(); - i_proj->second->Init( curve, f, l ); - } - GeomAPI_ProjectPointOnCurve* projector = i_proj->second; - projector->Perform( nodePnt ); - if ( projector->NbPoints() < 1 ) - { - MESSAGE( "SMESH_MesherHelper::CheckNodeU() failed to project" ); - return false; - } - Standard_Real U = projector->LowerDistanceParameter(); - u = double( U ); - curvPnt = curve->Value( u ); - dist = nodePnt.Distance( curvPnt ); + //GeomAPI_ProjectPointOnCurve& projector = GetPCProjector( E ); -- bug in OCCT-7.5.3p1 + GeomAdaptor_Curve curveAd( curve, f, l ); + ShapeAnalysis_Curve projector; + dist = projector.Project( curveAd, nodePnt, tol, curvPnt, u, false ); + // if ( projector.NbPoints() < 1 ) + // { + // MESSAGE( "SMESH_MesherHelper::CheckNodeU() failed to project" ); + // return false; + // } if ( distXYZ ) { + curvPnt = curve->Value( u ); curvPnt.Transform( loc ); distXYZ[0] = dist; distXYZ[1] = curvPnt.X(); distXYZ[2] = curvPnt.Y(); distXYZ[3]=curvPnt.Z(); @@ -1143,7 +1190,7 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, // store the fixed U on the edge if ( myShape.IsSame(E) && shapeID == myShapeID && myFixNodeParameters ) const_cast(n)->SetPosition - ( SMDS_PositionPtr( new SMDS_EdgePosition( U ))); + ( SMDS_PositionPtr( new SMDS_EdgePosition( u ))); } else if ( fabs( u ) > numeric_limits::min() ) { @@ -1323,7 +1370,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1, TBiQuad keyOfMap(n1,n2,n3,n4); std::map::iterator itMapCentralNode; itMapCentralNode = myMapWithCentralNode.find( keyOfMap ); - if ( itMapCentralNode != myMapWithCentralNode.end() ) + if ( itMapCentralNode != myMapWithCentralNode.end() ) { return (*itMapCentralNode).second; } @@ -1338,9 +1385,9 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1, std::map< int, int > faceId2nbNodes; std::map< int, int > ::iterator itMapWithIdFace; - + SMESHDS_Mesh* meshDS = GetMeshDS(); - + // check if a face lies on a FACE, i.e. its all corner nodes lie either on the FACE or // on sub-shapes of the FACE if ( GetMesh()->HasShapeToMesh() ) @@ -1498,7 +1545,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1, TBiQuad keyOfMap(n1,n2,n3); std::map::iterator itMapCentralNode; itMapCentralNode = myMapWithCentralNode.find( keyOfMap ); - if ( itMapCentralNode != myMapWithCentralNode.end() ) + if ( itMapCentralNode != myMapWithCentralNode.end() ) { return (*itMapCentralNode).second; } @@ -1513,9 +1560,9 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1, std::map< int, int > faceId2nbNodes; std::map< int, int > ::iterator itMapWithIdFace; - + SMESHDS_Mesh* meshDS = GetMeshDS(); - + // check if a face lies on a FACE, i.e. its all corner nodes lie either on the FACE or // on sub-shapes of the FACE if ( GetMesh()->HasShapeToMesh() ) @@ -1924,7 +1971,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::getMediumNodeOnComposedWire(const SMDS_ //purpose : Creates a node //======================================================================= -SMDS_MeshNode* SMESH_MesherHelper::AddNode(double x, double y, double z, int ID, +SMDS_MeshNode* SMESH_MesherHelper::AddNode(double x, double y, double z, smIdType ID, double u, double v) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -1953,11 +2000,11 @@ SMDS_MeshNode* SMESH_MesherHelper::AddNode(double x, double y, double z, int ID, SMDS_MeshEdge* SMESH_MesherHelper::AddEdge(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); - + SMDS_MeshEdge* edge = 0; if (myCreateQuadratic) { const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,force3d); @@ -1987,7 +2034,7 @@ SMDS_MeshEdge* SMESH_MesherHelper::AddEdge(const SMDS_MeshNode* n1, SMDS_MeshFace* SMESH_MesherHelper::AddFace(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2, const SMDS_MeshNode* n3, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2037,7 +2084,7 @@ SMDS_MeshFace* SMESH_MesherHelper::AddFace(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2, const SMDS_MeshNode* n3, const SMDS_MeshNode* n4, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2101,7 +2148,7 @@ SMDS_MeshFace* SMESH_MesherHelper::AddFace(const SMDS_MeshNode* n1, //======================================================================= SMDS_MeshFace* SMESH_MesherHelper::AddPolygonalFace (const vector& nodes, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2147,7 +2194,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n4, const SMDS_MeshNode* n5, const SMDS_MeshNode* n6, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2210,7 +2257,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2, const SMDS_MeshNode* n3, const SMDS_MeshNode* n4, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2251,7 +2298,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n3, const SMDS_MeshNode* n4, const SMDS_MeshNode* n5, - const int id, + const smIdType id, const bool force3d) { SMDS_MeshVolume* elem = 0; @@ -2301,7 +2348,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n6, const SMDS_MeshNode* n7, const SMDS_MeshNode* n8, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2420,7 +2467,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n10, const SMDS_MeshNode* n11, const SMDS_MeshNode* n12, - const int id, + const smIdType id, bool /*force3d*/) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2442,7 +2489,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, SMDS_MeshVolume* SMESH_MesherHelper::AddPolyhedralVolume (const std::vector& nodes, const std::vector& quantities, - const int id, + const smIdType id, const bool force3d) { SMESHDS_Mesh * meshDS = GetMeshDS(); @@ -2724,33 +2771,30 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap & theParam2 theParam2ColumnMap.begin()->second.size() == prevNbRows + expectNbRows ); } -namespace +//================================================================================ +/*! + * \brief Return true if a node is at a corner of a 2D structured mesh of FACE + */ +//================================================================================ + +bool SMESH_MesherHelper::IsCornerOfStructure( const SMDS_MeshNode* n, + const SMESHDS_SubMesh* faceSM, + SMESH_MesherHelper& faceAnalyser ) { - //================================================================================ - /*! - * \brief Return true if a node is at a corner of a 2D structured mesh of FACE - */ - //================================================================================ + int nbFacesInSM = 0; + if ( n ) { + SMDS_ElemIteratorPtr fIt = n->GetInverseElementIterator( SMDSAbs_Face ); + while ( fIt->more() ) + nbFacesInSM += faceSM->Contains( fIt->next() ); + } + if ( nbFacesInSM == 1 ) + return true; - bool isCornerOfStructure( const SMDS_MeshNode* n, - const SMESHDS_SubMesh* faceSM, - SMESH_MesherHelper& faceAnalyser ) + if ( nbFacesInSM == 2 && n->GetPosition()->GetTypeOfPosition() == SMDS_TOP_VERTEX ) { - int nbFacesInSM = 0; - if ( n ) { - SMDS_ElemIteratorPtr fIt = n->GetInverseElementIterator( SMDSAbs_Face ); - while ( fIt->more() ) - nbFacesInSM += faceSM->Contains( fIt->next() ); - } - if ( nbFacesInSM == 1 ) - return true; - - if ( nbFacesInSM == 2 && n->GetPosition()->GetTypeOfPosition() == SMDS_TOP_VERTEX ) - { - return faceAnalyser.IsRealSeam( n->getshapeId() ); - } - return false; + return faceAnalyser.IsRealSeam( n->getshapeId() ); } + return false; } //======================================================================= @@ -2785,7 +2829,7 @@ bool SMESH_MesherHelper::IsStructured( SMESH_subMesh* faceSM ) int nbRemainEdges = nbEdgesInWires.front(); do { TopoDS_Vertex V = IthVertex( 0, edges.front() ); - isCorner = isCornerOfStructure( SMESH_Algo::VertexNode( V, meshDS ), + isCorner = IsCornerOfStructure( SMESH_Algo::VertexNode( V, meshDS ), fSM, faceAnalyser); if ( !isCorner ) { edges.splice( edges.end(), edges, edges.begin() ); @@ -2826,7 +2870,7 @@ bool SMESH_MesherHelper::IsStructured( SMESH_subMesh* faceSM ) for ( ; n != nodes.end(); ++n ) { ++nbEdges; - if ( isCornerOfStructure( *n, fSM, faceAnalyser )) { + if ( IsCornerOfStructure( *n, fSM, faceAnalyser )) { nbEdgesInSide.push_back( nbEdges ); nbEdges = 0; } @@ -3162,7 +3206,7 @@ bool SMESH_MesherHelper::IsSubShape( const TopoDS_Shape& shape, //======================================================================= //function : IsSubShape -//purpose : +//purpose : //======================================================================= bool SMESH_MesherHelper::IsSubShape( const TopoDS_Shape& shape, SMESH_Mesh* aMesh ) @@ -3177,7 +3221,7 @@ bool SMESH_MesherHelper::IsSubShape( const TopoDS_Shape& shape, SMESH_Mesh* aMes //======================================================================= //function : IsBlock -//purpose : +//purpose : //======================================================================= bool SMESH_MesherHelper::IsBlock( const TopoDS_Shape& shape ) @@ -3235,6 +3279,16 @@ double SMESH_MesherHelper::getFaceMaxTol( const TopoDS_Shape& face ) const return tol; } +bool CheckAlmostZero(gp_Vec & vec1,gp_Vec & vec2, gp_Vec & vecref) +{ + auto v1 = gp_Dir(vec1); + auto v2 = gp_Dir(vec2); + auto vref = gp_Dir(vecref); + auto XYZ = v1.Crossed (v2); + double cond = XYZ.X()*vref.X()+XYZ.Y()*vref.Y()+XYZ.Z()*vref.Z(); + return (Abs(cond) <= 1e-12); +} + //================================================================================ /*! * \brief Return an angle between two EDGEs sharing a common VERTEX with reference @@ -3279,9 +3333,7 @@ double SMESH_MesherHelper::GetAngle( const TopoDS_Edge & theE1, vecRef = du ^ dv; if ( ++nbLoops > 10 ) { -#ifdef _DEBUG_ - cout << "SMESH_MesherHelper::GetAngle(): Captured in a sigularity" << endl; -#endif + MESSAGE("SMESH_MesherHelper::GetAngle(): Captured in a singularity"); return angle; } } @@ -3299,7 +3351,9 @@ double SMESH_MesherHelper::GetAngle( const TopoDS_Edge & theE1, if ( theE2.Orientation() /*GetSubShapeOri( F, theE2 )*/ == TopAbs_REVERSED ) vec2.Reverse(); angle = vec1.AngleWithRef( vec2, vecRef ); - + if ( angle < 0. && CheckAlmostZero(vec1,vec2,vecRef)) + angle*=-1; + if ( Abs ( angle ) >= 0.99 * M_PI ) { BRep_Tool::Range( theE1, f, l ); @@ -3360,7 +3414,7 @@ TopoDS_Vertex SMESH_MesherHelper::IthVertex( const bool is2nd, //================================================================================ /*! - * \brief Return type of shape contained in a group + * \brief Return type of shape contained in a group * \param group - a shape of type TopAbs_COMPOUND * \param avoidCompound - not to return TopAbs_COMPOUND */ @@ -3413,20 +3467,20 @@ TopoDS_Shape SMESH_MesherHelper::GetShapeOfHypothesis( const SMESHDS_Hypothesis SMESH_MesherHelper:: MType SMESH_MesherHelper::IsQuadraticMesh() { - int NbAllEdgsAndFaces=0; - int NbQuadFacesAndEdgs=0; - int NbFacesAndEdges=0; + smIdType NbAllEdgsAndFaces=0; + smIdType NbQuadFacesAndEdgs=0; + smIdType NbFacesAndEdges=0; //All faces and edges NbAllEdgsAndFaces = myMesh->NbEdges() + myMesh->NbFaces(); if ( NbAllEdgsAndFaces == 0 ) return SMESH_MesherHelper::LINEAR; - + //Quadratic faces and edges NbQuadFacesAndEdgs = myMesh->NbEdges(ORDER_QUADRATIC) + myMesh->NbFaces(ORDER_QUADRATIC); //Linear faces and edges NbFacesAndEdges = myMesh->NbEdges(ORDER_LINEAR) + myMesh->NbFaces(ORDER_LINEAR); - + if (NbAllEdgsAndFaces == NbQuadFacesAndEdgs) { //Quadratic mesh return SMESH_MesherHelper::QUADRATIC; @@ -3701,9 +3755,7 @@ namespace { // Structures used by FixQuadraticElements() mutable vector< const QLink* > _sides; mutable bool _sideIsAdded[4]; // added in chain of links gp_Vec _normal; -#ifdef _DEBUG_ mutable const SMDS_MeshElement* _face; -#endif QFace( const vector< const QLink*>& links, const SMDS_MeshElement* face=0 ); @@ -3817,11 +3869,8 @@ namespace { // Structures used by FixQuadraticElements() else _normal.SetCoord(1e-33,0,0); -#ifdef _DEBUG_ - _face = face; -#else - (void)face; // unused in release mode -#endif + if (SALOME::VerbosityActivated()) + _face = face; } //================================================================================ /*! @@ -4076,7 +4125,7 @@ namespace { // Structures used by FixQuadraticElements() } } else if ( _sides.size() < 4 ) - return thePrevLen; + return thePrevLen; // propagate to adjacent faces till limit step or boundary double len1 = thePrevLen + (theLink->MiddlePnt() - _sides[iL1]->MiddlePnt()).Modulus(); @@ -4166,7 +4215,7 @@ namespace { // Structures used by FixQuadraticElements() void QLink::SetContinuesFaces() const { // x0 x - QLink, [-|] - QFace, v - volume - // v0 | v1 + // v0 | v1 // | Between _faces of link x2 two vertical faces are continues // x1----x2-----x3 and two horizontal faces are continues. We set vertical faces // | to _faces[0] and _faces[1] and horizontal faces to @@ -4273,7 +4322,7 @@ namespace { // Structures used by FixQuadraticElements() } return isStraight; } - + //================================================================================ /*! * \brief Move medium nodes of vertical links of pentahedrons adjacent by side faces @@ -4402,13 +4451,13 @@ namespace { // Structures used by FixQuadraticElements() while ( startLink != linksEnd) // loop on columns { // We suppose we have a rectangular structure like shown here. We have found a - // corner of the rectangle (startCorner) and a boundary link sharing - // |/ |/ | the startCorner (startLink). We are going to loop on rows of the - // --o---o---o structure making several chains at once. One chain (columnChain) - // |\ | /| starts at startLink and continues upward (we look at the structure - // \ | \ | / | from such point that startLink is on the bottom of the structure). - // \| \|/ | While going upward we also fill horizontal chains (rowChains) we - // --o---o---o encounter. + // corner of the rectangle (startCorner) and a boundary link sharing + // |/ |/ | the startCorner (startLink). We are going to loop on rows of the + // --o---o---o structure making several chains at once. One chain (columnChain) + // |\ | /| starts at startLink and continues upward (we look at the structure + // \ | \ | / | from such point that startLink is on the bottom of the structure). + // \| \|/ | While going upward we also fill horizontal chains (rowChains) we + // --o---o---o encounter. // /|\ |\ | // / | \ | \ | startCorner // | \| \|,' @@ -4671,7 +4720,7 @@ namespace { // Structures used by FixQuadraticElements() continue; gp_XYZ edgeDir = SMESH_TNodeXYZ( nOnEdge[0] ) - SMESH_TNodeXYZ( nOnEdge[1] ); gp_XYZ edgeNorm = faceNorm ^ edgeDir; - n = theHelper.GetMediumNode( nOnEdge[0], nOnEdge[1], true ); // find n, not create + n = theHelper.GetMediumNode( nOnEdge[0], nOnEdge[1], true ); // find n, not create gp_XYZ pN0 = SMESH_TNodeXYZ( nOnEdge[0] ); gp_XYZ pMedium = SMESH_TNodeXYZ( n ); // on-edge node location gp_XYZ pFaceN = SMESH_TNodeXYZ( nOnFace ); // on-face node location @@ -4919,12 +4968,11 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, if ( !myMesh->HasShapeToMesh() ) return; SetSubShape( myMesh->GetShapeToMesh() ); -#ifdef _DEBUG_ int nbSolids = 0; TopTools_IndexedMapOfShape solids; TopExp::MapShapes(myShape,TopAbs_SOLID,solids); nbSolids = solids.Extent(); -#endif + TopTools_MapOfShape faces; // faces not in solid or in not meshed solid for ( TopExp_Explorer f(myShape,TopAbs_FACE,TopAbs_SOLID); f.More(); f.Next() ) { faces.Add( f.Current() ); // not in solid @@ -4935,27 +4983,39 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, faces.Add( f.Current() ); // in not meshed solid } else { // fix nodes in the solid and its faces -#ifdef _DEBUG_ MSG("FIX SOLID " << nbSolids-- << " #" << GetMeshDS()->ShapeToIndex(s.Current())); -#endif + SMESH_MesherHelper h(*myMesh); h.SetSubShape( s.Current() ); h.ToFixNodeParameters(true); - h.FixQuadraticElements( compError, false ); + try { + OCC_CATCH_SIGNALS; + h.FixQuadraticElements( compError, false ); + } + catch(...) { + if ( compError && compError->myComment.empty() ) + compError->myComment = "SMESH_MesherHelper::FixQuadraticElements() failed"; + } } } // fix nodes on geom faces -#ifdef _DEBUG_ - int nbfaces = nbSolids; - nbfaces = faces.Extent(); /*avoid "unused varianbles": */ nbfaces++, nbfaces--; -#endif + + int nbfaces = faces.Extent(); + for ( TopTools_MapIteratorOfMapOfShape fIt( faces ); fIt.More(); fIt.Next() ) { MESSAGE("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key())); MSG("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key())); SMESH_MesherHelper h(*myMesh); h.SetSubShape( fIt.Key() ); h.ToFixNodeParameters(true); - h.FixQuadraticElements( compError, true); + try { + OCC_CATCH_SIGNALS; + h.FixQuadraticElements( compError, true); + } + catch(...) { + if ( compError && compError->myComment.empty() ) + compError->myComment = "SMESH_MesherHelper::FixQuadraticElements() failed"; + } } //perf_print_all_meters(1); if ( compError && compError->myName == EDITERR_NO_MEDIUM_ON_GEOM ) @@ -5044,13 +5104,15 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, // hasRectFaces = hasRectFaces || // ( volTool.GetVolumeType() == SMDS_VolumeTool::QUAD_HEXA || // volTool.GetVolumeType() == SMDS_VolumeTool::QUAD_PENTA ); -#ifdef _DEBUG_ - if ( nbN == 6 ) - pFace->_face = GetMeshDS()->FindFace(faceNodes[0],faceNodes[2],faceNodes[4]); - else - pFace->_face = GetMeshDS()->FindFace(faceNodes[0],faceNodes[2], - faceNodes[4],faceNodes[6] ); -#endif + + if (SALOME::VerbosityActivated()) + { + if ( nbN == 6 ) + pFace->_face = GetMeshDS()->FindFace(faceNodes[0],faceNodes[2],faceNodes[4]); + else + pFace->_face = GetMeshDS()->FindFace(faceNodes[0],faceNodes[2], + faceNodes[4],faceNodes[6] ); + } } // collect pyramid apexes for further correction if ( vol->NbCornerNodes() == 5 ) @@ -5163,8 +5225,10 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, while ( len < numeric_limits::min() ) { // remove degenerated link if ( savedChain.empty() ) savedChain = chain; link1 = chain.erase( link1 ); - if ( link1 == chain.end() ) + if ( link1 == chain.end() ) { + link1 = --chain.end(); break; + } len = ((*link0)->MiddlePnt() - (*link1)->MiddlePnt()).Modulus(); } chainLen += len; @@ -5181,6 +5245,9 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, linkPos.push_back( chainLen ); } } + if ( chain.begin() == --chain.end() ) // chain.size() == 1 + continue; + gp_Vec move0 = chain.front()->_nodeMove; gp_Vec move1 = chain.back ()->_nodeMove; @@ -5277,20 +5344,22 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, move = gp_Vec( XYZ((*link1)->_mediumNode), newPnt.Transformed(loc) ); if ( SMDS_FacePositionPtr nPos = (*link1)->_mediumNode->GetPosition()) nPos->SetParameters( newUV.X(), newUV.Y() ); -#ifdef _DEBUG_ - if ( (XYZ((*link1)->node1()) - XYZ((*link1)->node2())).SquareModulus() < - move.SquareMagnitude()) + + if (SALOME::VerbosityActivated()) { - gp_XY uv0 = faceHlp.GetNodeUV( face, (*link0)->_mediumNode, 0, &checkUV ); - gp_XY uv2 = faceHlp.GetNodeUV( face, (*link2)->_mediumNode, 0, &checkUV ); - MSG( "TOO LONG MOVE \t" << - "uv0: "<node1()) - XYZ((*link1)->node2())).SquareModulus() < + move.SquareMagnitude()) + { + gp_XY uv0 = faceHlp.GetNodeUV( face, (*link0)->_mediumNode, 0, &checkUV ); + gp_XY uv2 = faceHlp.GetNodeUV( face, (*link2)->_mediumNode, 0, &checkUV ); + MSG( "TOO LONG MOVE \t" << + "uv0: "<Move( move, /*sum=*/false, /*is2dFixed=*/true ); } MSG( "Move " << (*link1)->_mediumNode->GetID() << " following " @@ -5377,7 +5446,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, { uv[ i ] = GetNodeUV( F, nodes[i], nodes[8], &checkUV ); // as this method is used after mesh generation, UV of nodes is not - // updated according to bending links, so we update + // updated according to bending links, so we update if ( i > 3 && nodes[i]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE ) CheckNodeUV( F, nodes[i], uv[ i ], 2*tol, /*force=*/true ); } @@ -5412,7 +5481,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, { uv[ i ] = GetNodeUV( F, nodes[i], nodes[(i+1)%3], &uvOK ); // as this method is used after mesh generation, UV of nodes is not - // updated according to bending links, so we update + // updated according to bending links, so we update if ( nodes[i]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE ) CheckNodeUV( F, nodes[i], uv[ i ], 2*tol, /*force=*/true ); } @@ -5479,16 +5548,16 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, pointsOnShapes[ SMESH_Block::ID_Ex11 ] = SMESH_TNodeXYZ( hexNodes[ 13 ] ); pointsOnShapes[ SMESH_Block::ID_E0y1 ] = SMESH_TNodeXYZ( hexNodes[ 12 ] ); pointsOnShapes[ SMESH_Block::ID_E1y1 ] = SMESH_TNodeXYZ( hexNodes[ 14 ] ); - pointsOnShapes[ SMESH_Block::ID_E00z ] = SMESH_TNodeXYZ( hexNodes[ 16 ] ); - pointsOnShapes[ SMESH_Block::ID_E10z ] = SMESH_TNodeXYZ( hexNodes[ 19 ] ); - pointsOnShapes[ SMESH_Block::ID_E01z ] = SMESH_TNodeXYZ( hexNodes[ 17 ] ); + pointsOnShapes[ SMESH_Block::ID_E00z ] = SMESH_TNodeXYZ( hexNodes[ 16 ] ); + pointsOnShapes[ SMESH_Block::ID_E10z ] = SMESH_TNodeXYZ( hexNodes[ 19 ] ); + pointsOnShapes[ SMESH_Block::ID_E01z ] = SMESH_TNodeXYZ( hexNodes[ 17 ] ); pointsOnShapes[ SMESH_Block::ID_E11z ] = SMESH_TNodeXYZ( hexNodes[ 18 ] ); pointsOnShapes[ SMESH_Block::ID_Fxy0 ] = SMESH_TNodeXYZ( hexNodes[ 20 ] ); pointsOnShapes[ SMESH_Block::ID_Fxy1 ] = SMESH_TNodeXYZ( hexNodes[ 25 ] ); - pointsOnShapes[ SMESH_Block::ID_Fx0z ] = SMESH_TNodeXYZ( hexNodes[ 21 ] ); - pointsOnShapes[ SMESH_Block::ID_Fx1z ] = SMESH_TNodeXYZ( hexNodes[ 23 ] ); - pointsOnShapes[ SMESH_Block::ID_F0yz ] = SMESH_TNodeXYZ( hexNodes[ 24 ] ); + pointsOnShapes[ SMESH_Block::ID_Fx0z ] = SMESH_TNodeXYZ( hexNodes[ 21 ] ); + pointsOnShapes[ SMESH_Block::ID_Fx1z ] = SMESH_TNodeXYZ( hexNodes[ 23 ] ); + pointsOnShapes[ SMESH_Block::ID_F0yz ] = SMESH_TNodeXYZ( hexNodes[ 24 ] ); pointsOnShapes[ SMESH_Block::ID_F1yz ] = SMESH_TNodeXYZ( hexNodes[ 22 ] ); gp_XYZ nCenterParams(0.5, 0.5, 0.5), nCenterCoords; @@ -5525,8 +5594,6 @@ void SMESH_MesherHelper::WriteShape(const TopoDS_Shape& s) { const char* name = "/tmp/shape.brep"; BRepTools::Write( s, name ); -#ifdef _DEBUG_ - std::cout << name << std::endl; -#endif + MESSAGE(name); }