X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESH%2FSMESH_MesherHelper.cxx;h=48244b0552b612ea6a29a226c6c6beae231da1c9;hp=6e5b08498b8addc06a1012b8bac388092f0612c3;hb=133cb381393ea728f6c6ff3a732738ecd1ad72d8;hpb=82d27750a524709b26cef3f8ca2fbbd2607973e7 diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index 6e5b08498..48244b055 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -259,6 +259,7 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh) { // look for a "seam" edge, a real seam or an edge on period boundary TopoDS_Edge edge = TopoDS::Edge( exp.Current() ); + const int edgeID = meshDS->ShapeToIndex( edge ); if ( myParIndex ) { BRep_Tool::UVPoints( edge, face, uv1, uv2 ); @@ -305,7 +306,6 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh) if ( isSeam ) { // store seam shape indices, negative if shape encounters twice - int edgeID = meshDS->ShapeToIndex( edge ); mySeamShapeIds.insert( IsSeamShape( edgeID ) ? -edgeID : edgeID ); for ( TopExp_Explorer v( edge, TopAbs_VERTEX ); v.More(); v.Next() ) { int vertexID = meshDS->ShapeToIndex( v.Current() ); @@ -315,10 +315,15 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh) } // look for a degenerated edge if ( SMESH_Algo::isDegenerated( edge )) { - myDegenShapeIds.insert( meshDS->ShapeToIndex( edge )); + myDegenShapeIds.insert( edgeID ); for ( TopExp_Explorer v( edge, TopAbs_VERTEX ); v.More(); v.Next() ) myDegenShapeIds.insert( meshDS->ShapeToIndex( v.Current() )); } + if ( !BRep_Tool::SameParameter( edge ) || + !BRep_Tool::SameRange( edge )) + { + setPosOnShapeValidity( edgeID, false ); + } } } } @@ -527,11 +532,11 @@ void SMESH_MesherHelper::ToFixNodeParameters(bool toFix) //======================================================================= -//function : GetUVOnSeam +//function : getUVOnSeam //purpose : Select UV on either of 2 pcurves of a seam edge, closest to the given UV //======================================================================= -gp_Pnt2d SMESH_MesherHelper::GetUVOnSeam( const gp_Pnt2d& uv1, const gp_Pnt2d& uv2 ) const +gp_Pnt2d SMESH_MesherHelper::getUVOnSeam( const gp_Pnt2d& uv1, const gp_Pnt2d& uv2 ) const { gp_Pnt2d result = uv1; for ( int i = U_periodic; i <= V_periodic ; ++i ) @@ -568,38 +573,34 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, const SMDS_PositionPtr Pos = n->GetPosition(); bool uvOK = false; - if(Pos->GetTypeOfPosition()==SMDS_TOP_FACE) + if ( Pos->GetTypeOfPosition() == SMDS_TOP_FACE ) { // node has position on face - const SMDS_FacePosition* fpos = - static_cast( Pos ); - uv.SetCoord(fpos->GetUParameter(),fpos->GetVParameter()); + const SMDS_FacePosition* fpos = static_cast( Pos ); + uv.SetCoord( fpos->GetUParameter(), fpos->GetVParameter() ); if ( check ) - uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*MaxTolerance( F )); + uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*getFaceMaxTol( F )); } - else if(Pos->GetTypeOfPosition()==SMDS_TOP_EDGE) + else if ( Pos->GetTypeOfPosition() == SMDS_TOP_EDGE ) { - // node has position on edge => it is needed to find - // corresponding edge from face, get pcurve for this - // edge and retrieve value from this pcurve - const SMDS_EdgePosition* epos = - static_cast( Pos ); - int edgeID = n->getshapeId(); - TopoDS_Edge E = TopoDS::Edge(GetMeshDS()->IndexToShape(edgeID)); + // node has position on EDGE => it is needed to find + // corresponding EDGE from FACE, get pcurve for this + // EDGE and retrieve value from this pcurve + const SMDS_EdgePosition* epos = static_cast( Pos ); + 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); - bool validU = ( f < u && u < l ); - if ( validU ) - uv = C2d->Value( u ); - else - uv.SetCoord( Precision::Infinite(),0.); + 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.); if ( check || !validU ) - uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*MaxTolerance( F ),/*force=*/ !validU ); + uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*getFaceMaxTol( F ),/*force=*/ !validU ); - // for a node on a seam edge select one of UVs on 2 pcurves - if ( n2 && IsSeamShape( edgeID ) ) + // for a node on a seam EDGE select one of UVs on 2 pcurves + if ( n2 && IsSeamShape( edgeID )) { - uv = GetUVOnSeam( uv, GetNodeUV( F, n2, 0, check )); + uv = getUVOnSeam( uv, GetNodeUV( F, n2, 0, check )); } else { // adjust uv to period @@ -611,23 +612,22 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, if ( isUPeriodic || isVPeriodic ) { Standard_Real UF,UL,VF,VL; S->Bounds(UF,UL,VF,VL); - if ( isUPeriodic ) - newUV.SetX( uv.X() + ShapeAnalysis::AdjustToPeriod(uv.X(),UF,UL)); - if ( isVPeriodic ) - newUV.SetY( uv.Y() + ShapeAnalysis::AdjustToPeriod(uv.Y(),VF,VL)); - } - if ( n2 ) - { - gp_Pnt2d uv2 = GetNodeUV( F, n2, 0, check ); - if ( isUPeriodic && Abs( uv.X()-uv2.X() ) < Abs( newUV.X()-uv2.X() )) - newUV.SetX( uv.X() ); - if ( isVPeriodic && Abs( uv.Y()-uv2.Y() ) < Abs( newUV.Y()-uv2.Y() )) - newUV.SetY( uv.Y() ); + if ( isUPeriodic ) newUV.SetX( uv.X() + ShapeAnalysis::AdjustToPeriod(uv.X(),UF,UL)); + if ( isVPeriodic ) newUV.SetY( uv.Y() + ShapeAnalysis::AdjustToPeriod(uv.Y(),VF,VL)); + + if ( n2 ) + { + gp_Pnt2d uv2 = GetNodeUV( F, n2, 0, check ); + if ( isUPeriodic && Abs( uv.X()-uv2.X() ) < Abs( newUV.X()-uv2.X() )) + newUV.SetX( uv.X() ); + if ( isVPeriodic && Abs( uv.Y()-uv2.Y() ) < Abs( newUV.Y()-uv2.Y() )) + newUV.SetY( uv.Y() ); + } } uv = newUV; } } - else if(Pos->GetTypeOfPosition()==SMDS_TOP_VERTEX) + else if ( Pos->GetTypeOfPosition() == SMDS_TOP_VERTEX ) { if ( int vertexID = n->getshapeId() ) { const TopoDS_Vertex& V = TopoDS::Vertex(GetMeshDS()->IndexToShape(vertexID)); @@ -646,7 +646,7 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, // get UV of a vertex closest to the node double dist = 1e100; gp_Pnt pn = XYZ( n ); - for ( TopExp_Explorer vert(F,TopAbs_VERTEX); !uvOK && vert.More(); vert.Next() ) { + for ( TopExp_Explorer vert( F,TopAbs_VERTEX ); !uvOK && vert.More(); vert.Next() ) { TopoDS_Vertex curV = TopoDS::Vertex( vert.Current() ); gp_Pnt p = BRep_Tool::Pnt( curV ); double curDist = p.SquareDistance( pn ); @@ -675,13 +675,23 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, } } } - if ( n2 && IsSeamShape( vertexID ) ) - uv = GetUVOnSeam( uv, GetNodeUV( F, n2, 0 )); + if ( n2 && IsSeamShape( vertexID )) + { + bool isSeam = ( myShape.IsSame( F )); + if ( !isSeam ) { + SMESH_MesherHelper h( *myMesh ); + h.SetSubShape( F ); + isSeam = IsSeamShape( vertexID ); + } + + if ( isSeam ) + uv = getUVOnSeam( uv, GetNodeUV( F, n2, 0 )); + } } } else { - uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*MaxTolerance( F )); + uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*getFaceMaxTol( F )); } if ( check ) @@ -750,7 +760,7 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, const_cast(n)->SetPosition ( SMDS_PositionPtr( new SMDS_FacePosition( U, V ))); } - else if ( uv.Modulus() > numeric_limits::min() ) + else if ( myShape.IsSame(F) && uv.Modulus() > numeric_limits::min() ) { setPosOnShapeValidity( shapeID, true ); } @@ -941,7 +951,7 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, int shapeID = n->getshapeId(); bool infinit = Precision::IsInfinite( u ); bool zero = ( u == 0. ); - if ( force || toCheckPosOnShape( shapeID ) || infinit || zero ) + if ( force || infinit || zero || toCheckPosOnShape( shapeID )) { TopLoc_Location loc; double f,l; Handle(Geom_Curve) curve = BRep_Tool::Curve( E,loc,f,l ); @@ -958,7 +968,7 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, gp_Pnt nodePnt = SMESH_TNodeXYZ( n ); if ( !loc.IsIdentity() ) nodePnt.Transform( loc.Transformation().Inverted() ); gp_Pnt curvPnt; - double dist = u; + double dist = 2*tol; if ( !infinit ) { curvPnt = curve->Value( u ); @@ -1037,9 +1047,12 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, //======================================================================= //function : GetMediumPos //purpose : Return index and type of the shape (EDGE or FACE only) to -// set a medium node on +// set a medium node on //param : useCurSubShape - if true, returns the shape set via SetSubShape() // if any +// calling GetMediumPos() with useCurSubShape=true is OK only for the +// case where the lower dim mesh is already constructed and converted to quadratic, +// else, nodes on EDGEs are assigned to FACE, for example. //======================================================================= std::pair @@ -1442,11 +1455,9 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetMediumNode(const SMDS_MeshNode* n1, TopoDS_Edge E; double u [2]; TopoDS_Face F; gp_XY uv[2]; bool uvOK[2] = { false, false }; + const bool useCurSubShape = ( !myShape.IsNull() && myShape.ShapeType() == TopAbs_EDGE ); - pair pos = GetMediumPos( n1, n2, mySetElemOnShape ); - // calling GetMediumPos() with useCurSubShape=mySetElemOnShape is OK only for the - // case where the lower dim mesh is already constructed, else, nodes on EDGEs are - // assigned to FACE, for example. + pair pos = GetMediumPos( n1, n2, useCurSubShape ); // get positions of the given nodes on shapes if ( pos.second == TopAbs_FACE ) @@ -2464,7 +2475,7 @@ namespace //======================================================================= //function : IsStructured -//purpose : Return true if 2D mesh on FACE is structured +//purpose : Return true if 2D mesh on FACE is a structured rectangle //======================================================================= bool SMESH_MesherHelper::IsStructured( SMESH_subMesh* faceSM ) @@ -2559,7 +2570,8 @@ bool SMESH_MesherHelper::IsStructured( SMESH_subMesh* faceSM ) //purpose : Return true if 2D mesh on FACE is ditorted //======================================================================= -bool SMESH_MesherHelper::IsDistorted2D( SMESH_subMesh* faceSM ) +bool SMESH_MesherHelper::IsDistorted2D( SMESH_subMesh* faceSM, + bool checkUV) { if ( !faceSM || faceSM->GetSubShape().ShapeType() != TopAbs_FACE ) return false; @@ -2574,9 +2586,10 @@ bool SMESH_MesherHelper::IsDistorted2D( SMESH_subMesh* faceSM ) if ( !smDS || smDS->NbElements() == 0 ) return false; SMDS_ElemIteratorPtr faceIt = smDS->GetElements(); - double prevArea2D = 0; + double prevArea = 0; vector< const SMDS_MeshNode* > nodes; vector< gp_XY > uv; + bool* toCheckUV = checkUV ? & checkUV : 0; while ( faceIt->more() && !haveBadFaces ) { const SMDS_MeshElement* face = faceIt->next(); @@ -2608,18 +2621,18 @@ bool SMESH_MesherHelper::IsDistorted2D( SMESH_subMesh* faceSM ) // get UVs uv.resize( nodes.size() ); for ( size_t i = 0; i < nodes.size(); ++i ) - uv[ i ] = helper.GetNodeUV( F, nodes[ i ], inFaceNode ); + uv[ i ] = helper.GetNodeUV( F, nodes[ i ], inFaceNode, toCheckUV ); // compare orientation of triangles + double faceArea = 0; for ( int iT = 0, nbT = nodes.size()-2; iT < nbT; ++iT ) { gp_XY v1 = uv[ iT+1 ] - uv[ 0 ]; gp_XY v2 = uv[ iT+2 ] - uv[ 0 ]; - double area2D = v2 ^ v1; - if (( haveBadFaces = ( area2D * prevArea2D < 0 ))) - break; - prevArea2D = area2D; + faceArea += v2 ^ v1; } + haveBadFaces = ( faceArea * prevArea < 0 ); + prevArea = faceArea; } return haveBadFaces; @@ -2804,6 +2817,28 @@ bool SMESH_MesherHelper::IsSubShape( const TopoDS_Shape& shape, SMESH_Mesh* aMes (shape.ShapeType() == TopAbs_COMPOUND && aMesh->GetMeshDS()->IsGroupOfSubShapes( shape )); } +//======================================================================= +//function : IsBlock +//purpose : +//======================================================================= + +bool SMESH_MesherHelper::IsBlock( const TopoDS_Shape& shape ) +{ + if ( shape.IsNull() ) + return false; + + TopoDS_Shell shell; + TopExp_Explorer exp( shape, TopAbs_SHELL ); + if ( !exp.More() ) return false; + shell = TopoDS::Shell( exp.Current() ); + if ( exp.Next(), exp.More() ) return false; + + TopoDS_Vertex v; + TopTools_IndexedMapOfOrientedShape map; + return SMESH_Block::FindBlockShapes( shell, v, v, map ); +} + + //================================================================================ /*! * \brief Return maximal tolerance of shape @@ -2824,6 +2859,24 @@ double SMESH_MesherHelper::MaxTolerance( const TopoDS_Shape& shape ) return tol; } +//================================================================================ +/*! + * \brief Return MaxTolerance( face ), probably cached + */ +//================================================================================ + +double SMESH_MesherHelper::getFaceMaxTol( const TopoDS_Shape& face ) const +{ + int faceID = GetMeshDS()->ShapeToIndex( face ); + + SMESH_MesherHelper* me = const_cast< SMESH_MesherHelper* >( this ); + double & tol = me->myFaceMaxTol.insert( make_pair( faceID, -1. )).first->second; + if ( tol < 0 ) + tol = MaxTolerance( face ); + + return tol; +} + //================================================================================ /*! * \brief Return an angle between two EDGEs sharing a common VERTEX with reference @@ -4582,6 +4635,8 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, // 3. Compute displacement of medium nodes // --------------------------------------- + SMESH_MesherHelper faceHlp(*myMesh); + // two loops on QFaces: the first is to treat boundary links, the second is for internal ones. TopLoc_Location loc; bool checkUV; @@ -4665,22 +4720,23 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, TopoDS_Face face; if ( !isInside ) { - // compute node displacement of end links of chain in parametric space of face + // compute node displacement of end links of chain in parametric space of FACE TChainLink& linkOnFace = *(++chain.begin()); const SMDS_MeshNode* nodeOnFace = linkOnFace->_mediumNode; TopoDS_Shape f = GetSubShapeByNode( nodeOnFace, GetMeshDS() ); if ( !f.IsNull() && f.ShapeType() == TopAbs_FACE ) { face = TopoDS::Face( f ); + faceHlp.SetSubShape( face ); Handle(Geom_Surface) surf = BRep_Tool::Surface(face,loc); bool isStraight[2]; for ( int is1 = 0; is1 < 2; ++is1 ) // move0 or move1 { TChainLink& link = is1 ? chain.back() : chain.front(); - gp_XY uvm = GetNodeUV( face, link->_mediumNode, nodeOnFace, &checkUV); - gp_XY uv1 = GetNodeUV( face, link->node1(), nodeOnFace, &checkUV); - gp_XY uv2 = GetNodeUV( face, link->node2(), nodeOnFace, &checkUV); - gp_XY uv12 = GetMiddleUV( surf, uv1, uv2); + gp_XY uvm = faceHlp.GetNodeUV( face, link->_mediumNode, nodeOnFace, &checkUV ); + gp_XY uv1 = faceHlp.GetNodeUV( face, link->node1(), nodeOnFace, &checkUV ); + gp_XY uv2 = faceHlp.GetNodeUV( face, link->node2(), nodeOnFace, &checkUV ); + gp_XY uv12 = faceHlp.GetMiddleUV( surf, uv1, uv2 ); // uvMove = uvm - uv12 gp_XY uvMove = applyIn2D(surf, uvm, uv12, gp_XY_Subtracted, /*inPeriod=*/false); ( is1 ? move1 : move0 ).SetCoord( uvMove.X(), uvMove.Y(), 0 ); @@ -4695,10 +4751,10 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, } // check if a chain is already fixed - gp_XY uvm = GetNodeUV( face, linkOnFace->_mediumNode, 0, &checkUV); - gp_XY uv1 = GetNodeUV( face, linkOnFace->node1(), nodeOnFace, &checkUV); - gp_XY uv2 = GetNodeUV( face, linkOnFace->node2(), nodeOnFace, &checkUV); - gp_XY uv12 = GetMiddleUV( surf, uv1, uv2); + gp_XY uvm = faceHlp.GetNodeUV( face, linkOnFace->_mediumNode, 0, &checkUV ); + gp_XY uv1 = faceHlp.GetNodeUV( face, linkOnFace->node1(), nodeOnFace, &checkUV ); + gp_XY uv2 = faceHlp.GetNodeUV( face, linkOnFace->node2(), nodeOnFace, &checkUV ); + gp_XY uv12 = faceHlp.GetMiddleUV( surf, uv1, uv2 ); if (( uvm - uv12 ).SquareModulus() > 1e-10 ) { MSG("Already fixed - ignore"); @@ -4743,8 +4799,8 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, else { // compute 3D displacement by 2D one Handle(Geom_Surface) s = BRep_Tool::Surface(face,loc); - gp_XY oldUV = GetNodeUV( face, (*link1)->_mediumNode, 0, &checkUV); - gp_XY newUV = applyIn2D( s, oldUV, gp_XY( move.X(),move.Y()), gp_XY_Added); + gp_XY oldUV = faceHlp.GetNodeUV( face, (*link1)->_mediumNode, 0, &checkUV ); + gp_XY newUV = applyIn2D( s, oldUV, gp_XY( move.X(),move.Y()), gp_XY_Added ); gp_Pnt newPnt = s->Value( newUV.X(), newUV.Y()); move = gp_Vec( XYZ((*link1)->_mediumNode), newPnt.Transformed(loc) ); if ( SMDS_FacePosition* nPos = @@ -4754,8 +4810,8 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, if ( (XYZ((*link1)->node1()) - XYZ((*link1)->node2())).SquareModulus() < move.SquareMagnitude()) { - gp_XY uv0 = GetNodeUV( face, (*link0)->_mediumNode, 0, &checkUV); - gp_XY uv2 = GetNodeUV( face, (*link2)->_mediumNode, 0, &checkUV); + 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: "<