X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_FaceSide.cxx;h=0ebea37408a71c5d3e62c1036016a2f96b5a8b5e;hp=0c7c334f20cd7e4b1fbf33411756b737a3a01612;hb=382f2cc4abb2ee8553a911aeb27348e96c39d197;hpb=b6e8c17211f3a65d729b5ffe7269c789a042ae24 diff --git a/src/StdMeshers/StdMeshers_FaceSide.cxx b/src/StdMeshers/StdMeshers_FaceSide.cxx index 0c7c334f2..0ebea3740 100644 --- a/src/StdMeshers/StdMeshers_FaceSide.cxx +++ b/src/StdMeshers/StdMeshers_FaceSide.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2020 CEA/DEN, EDF R&D, OPEN CASCADE // // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS @@ -32,10 +32,11 @@ #include "SMESHDS_Mesh.hxx" #include "SMESHDS_SubMesh.hxx" #include "SMESH_Algo.hxx" +#include "SMESH_Block.hxx" +#include "SMESH_ComputeError.hxx" #include "SMESH_Mesh.hxx" +#include "SMESH_MeshEditor.hxx" #include "SMESH_MesherHelper.hxx" -#include "SMESH_ComputeError.hxx" -#include "SMESH_Block.hxx" #include #include @@ -63,7 +64,7 @@ using namespace std; * \brief Constructor of a side of one edge * \param theFace - the face * \param theEdge - the edge - */ + */ //================================================================================ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, @@ -71,11 +72,15 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, SMESH_Mesh* theMesh, const bool theIsForward, const bool theIgnoreMediumNodes, + SMESH_MesherHelper* theFaceHelper, SMESH_ProxyMesh::Ptr theProxyMesh) { std::list edges(1,theEdge); - *this = StdMeshers_FaceSide( theFace, edges, theMesh, theIsForward, - theIgnoreMediumNodes, theProxyMesh ); + StdMeshers_FaceSide tmp( theFace, edges, theMesh, theIsForward, + theIgnoreMediumNodes, theFaceHelper, theProxyMesh ); + *this = tmp; + + tmp.myHelper = NULL; } //================================================================================ @@ -84,12 +89,13 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, */ //================================================================================ -StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, - std::list& theEdges, - SMESH_Mesh* theMesh, - const bool theIsForward, - const bool theIgnoreMediumNodes, - SMESH_ProxyMesh::Ptr theProxyMesh) +StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, + const std::list& theEdges, + SMESH_Mesh* theMesh, + const bool theIsForward, + const bool theIgnoreMediumNodes, + SMESH_MesherHelper* theFaceHelper, + SMESH_ProxyMesh::Ptr theProxyMesh) { int nbEdges = theEdges.size(); myEdge.resize ( nbEdges ); @@ -108,13 +114,19 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, myMissingVertexNodes = false; myIgnoreMediumNodes = theIgnoreMediumNodes; myDefaultPnt2d = gp_Pnt2d( 1e+100, 1e+100 ); + myHelper = NULL; if ( !myProxyMesh ) myProxyMesh.reset( new SMESH_ProxyMesh( *theMesh )); + if ( theFaceHelper && theFaceHelper->GetSubShape() == myFace ) + { + myHelper = new SMESH_MesherHelper( * myProxyMesh->GetMesh() ); + myHelper->CopySubShapeInfo( *theFaceHelper ); + } if ( nbEdges == 0 ) return; SMESHDS_Mesh* meshDS = myProxyMesh->GetMeshDS(); int nbDegen = 0; - std::list::iterator edge = theEdges.begin(); + std::list::const_iterator edge = theEdges.begin(); for ( int index = 0; edge != theEdges.end(); ++index, ++edge ) { int i = theIsForward ? index : nbEdges-index-1; @@ -155,12 +167,55 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, myC3dAdaptor[i].Load( C3d, 0, 0.5 * BRep_Tool::Tolerance( V )); } } + else if ( myEdgeLength[i] > DBL_MIN ) + { + Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],myFirst[i], myLast[i] ); + myC3dAdaptor[i].Load( C3d, myFirst[i], myLast[i] ); + if ( myEdge[i].Orientation() == TopAbs_REVERSED ) + std::swap( myFirst[i], myLast[i] ); + } + // reverse a proxy sub-mesh if ( !theIsForward ) reverseProxySubmesh( myEdge[i] ); } // loop on edges + // orient seam edges (#19982) + const double tol = Precision::Confusion(); + if ( NbEdges() > 1 && !myC2d[0].IsNull() ) + for ( int i = 0; i < NbEdges(); ++i ) + { + int iPrev = SMESH_MesherHelper::WrapIndex( i - 1, NbEdges() ); + if ( !BRep_Tool::IsClosed( myEdge[i], myFace ) || !myC2d[iPrev] ) + continue; + gp_Pnt2d pLastPrev = myC2d[iPrev]->Value( myLast[iPrev] ); + gp_Pnt2d pFirst = myC2d[i]->Value( myFirst[i] ); + if ( pLastPrev.IsEqual( pFirst, tol )) + continue; // OK + pFirst = myC2d[i]->Value( myLast[i] ); + if ( pLastPrev.IsEqual( pFirst, tol )) + { + std::swap( myFirst[i], myLast[i] ); + continue; + } + TopoDS_Edge E = myEdge[i]; + E.Reverse(); + Handle(Geom2d_Curve) c2dRev = BRep_Tool::CurveOnSurface( E, myFace, myFirst[i], myLast[i] ); + pFirst = c2dRev->Value( myFirst[i] ); + if ( pLastPrev.IsEqual( pFirst, tol )) + { + myC2d[i] = c2dRev; + continue; + } + pFirst = c2dRev->Value( myLast[i] ); + if ( pLastPrev.IsEqual( pFirst, tol )) + { + myC2d[i] = c2dRev; + std::swap( myFirst[i], myLast[i] ); + } + } + // count nodes and segments NbPoints( /*update=*/true ); @@ -203,6 +258,7 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const StdMeshers_FaceSide* theSide, myLast.push_back ( theULast ); myNormPar.push_back ( 1. ); myIsUniform.push_back( true ); + myHelper = NULL; myLength = 0; myProxyMesh = theSide->myProxyMesh; myDefaultPnt2d = *thePnt2d1; @@ -260,6 +316,7 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(UVPtStructVec& theSideNodes, } myFace = theFace; + myHelper = NULL; myPoints = theSideNodes; myNbPonits = myPoints.size(); myNbSegments = myNbPonits + 1; @@ -312,6 +369,17 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(UVPtStructVec& theSideNodes, myEdgeLength.resize( 1, myLength ); } +//================================================================================ +/*! + * \brief Destructor + */ +//================================================================================ + +StdMeshers_FaceSide::~StdMeshers_FaceSide() +{ + delete myHelper; myHelper = NULL; +} + //================================================================================ /* * Return info on nodes on the side @@ -326,12 +394,12 @@ const std::vector& StdMeshers_FaceSide::GetUVPtStruct(bool isXCons if ( NbEdges() == 0 ) return myPoints; StdMeshers_FaceSide* me = const_cast< StdMeshers_FaceSide* >( this ); - SMESH_MesherHelper eHelper( *myProxyMesh->GetMesh() ); - SMESH_MesherHelper fHelper( *myProxyMesh->GetMesh() ); - fHelper.SetSubShape( myFace ); - bool paramOK; + bool paramOK = true; double eps = 1e-100; + SMESH_MesherHelper eHelper( *myProxyMesh->GetMesh() ); + SMESH_MesherHelper& fHelper = *FaceHelper(); + // sort nodes of all edges by putting them into a map map< double, const SMDS_MeshNode*> u2node; @@ -417,6 +485,7 @@ const std::vector& StdMeshers_FaceSide::GetUVPtStruct(bool isXCons for ( size_t j = 0; j < u2nodeVec.size(); ++j ) u2node.insert( u2node.end(), u2nodeVec[j] ); } + continue; } // loop on myEdge's // Add 2nd VERTEX node for a last EDGE @@ -517,8 +586,7 @@ const std::vector& StdMeshers_FaceSide::GetUVPtStruct(bool isXCons uvPt.normParam = u_node->first; uvPt.x = uvPt.y = uvPt.normParam; // -- U ---------------------------------------------- - const SMDS_EdgePosition* epos = - dynamic_cast(uvPt.node->GetPosition()); + SMDS_EdgePositionPtr epos = uvPt.node->GetPosition(); if ( epos && uvPt.node->getshapeId() == myEdgeID[iE] ) { uvPt.param = epos->GetUParameter(); } @@ -618,9 +686,8 @@ std::vector StdMeshers_FaceSide::GetOrderedNodes(int theEd if ( NbEdges() == 0 ) return resultNodes; //SMESHDS_Mesh* meshDS = myProxyMesh->GetMeshDS(); - SMESH_MesherHelper eHelper( *myProxyMesh->GetMesh() ); - SMESH_MesherHelper fHelper( *myProxyMesh->GetMesh() ); - fHelper.SetSubShape( myFace ); + SMESH_MesherHelper eHelper( *myProxyMesh->GetMesh() ); + SMESH_MesherHelper& fHelper = * FaceHelper(); bool paramOK = true; // Sort nodes of all edges putting them into a map @@ -634,7 +701,7 @@ std::vector StdMeshers_FaceSide::GetOrderedNodes(int theEd iE = theEdgeInd % NbEdges(); iEnd = iE + 1; } - for ( iE = 0; iE < iEnd; ++iE ) + for ( ; iE < iEnd; ++iE ) { double prevNormPar = ( iE == 0 ? 0 : myNormPar[ iE-1 ]); // normalized param @@ -1014,8 +1081,7 @@ int StdMeshers_FaceSide::NbPoints(const bool update) const } } - SMESH_MesherHelper helper( *myProxyMesh->GetMesh() ); - helper.SetSubShape( myFace ); + SMESH_MesherHelper* helper = FaceHelper(); std::set< const SMDS_MeshNode* > vNodes; const int nbV = NbEdges() + !IsClosed(); @@ -1023,8 +1089,8 @@ int StdMeshers_FaceSide::NbPoints(const bool update) const if ( const SMDS_MeshNode* n = VertexNode( i )) { if ( !vNodes.insert( n ).second && - ( helper.IsRealSeam ( n->getshapeId() ) || - helper.IsDegenShape( n->getshapeId() ))) + ( helper->IsRealSeam ( n->getshapeId() ) || + helper->IsDegenShape( n->getshapeId() ))) me->myNbPonits++; } else @@ -1110,6 +1176,8 @@ void StdMeshers_FaceSide::dump(const char* msg) const MESSAGE_ADD ( "\tF: "<myLast[i] ? -1. : 1.); + double aLen3dU = r * myEdgeLength[i] * ( myFirst[i] > myLast[i] ? -1. : 1. ); GCPnts_AbscissaPoint AbPnt ( const_cast( myC3dAdaptor[i]), aLen3dU, myFirst[i] ); if( AbPnt.IsDone() ) { @@ -1241,9 +1309,14 @@ TSideVector StdMeshers_FaceSide::GetFaceWires(const TopoDS_Face& theFace, SMESH_Mesh & theMesh, const bool theIgnoreMediumNodes, TError & theError, + SMESH_MesherHelper* theFaceHelper, SMESH_ProxyMesh::Ptr theProxyMesh, const bool theCheckVertexNodes) { + SMESH_MesherHelper helper( theMesh ); + if ( theFaceHelper && theFaceHelper->GetSubShape() == theFace ) + helper.CopySubShapeInfo( *theFaceHelper ); + list< TopoDS_Edge > edges, internalEdges; list< int > nbEdgesInWires; int nbWires = SMESH_Block::GetOrderedEdges (theFace, edges, nbEdgesInWires); @@ -1287,7 +1360,7 @@ TSideVector StdMeshers_FaceSide::GetFaceWires(const TopoDS_Face& theFace, StdMeshers_FaceSide* wire = new StdMeshers_FaceSide( theFace, wireEdges, &theMesh, /*isForward=*/true, theIgnoreMediumNodes, - theProxyMesh ); + &helper, theProxyMesh ); wires[ iW ] = StdMeshers_FaceSidePtr( wire ); from = to; } @@ -1295,7 +1368,7 @@ TSideVector StdMeshers_FaceSide::GetFaceWires(const TopoDS_Face& theFace, { StdMeshers_FaceSide* wire = new StdMeshers_FaceSide( theFace, internalEdges.back(), &theMesh, /*isForward=*/true, theIgnoreMediumNodes, - theProxyMesh ); + &helper, theProxyMesh ); wires.push_back( StdMeshers_FaceSidePtr( wire )); internalEdges.pop_back(); } @@ -1351,3 +1424,20 @@ bool StdMeshers_FaceSide::IsClosed() const { return myEdge.empty() ? false : FirstVertex().IsSame( LastVertex() ); } + +//================================================================================ +/*! + * \brief Return a helper initialized with the FACE + */ +//================================================================================ + +SMESH_MesherHelper* StdMeshers_FaceSide::FaceHelper() const +{ + StdMeshers_FaceSide* me = const_cast< StdMeshers_FaceSide* >( this ); + if ( !myHelper && myProxyMesh ) + { + me->myHelper = new SMESH_MesherHelper( *myProxyMesh->GetMesh() ); + me->myHelper->SetSubShape( myFace ); + } + return me->myHelper; +}