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=eceb92e3657040b6f37d0ec777c370f534dce5bf;hb=refs%2Ftags%2FV9_7_0a1;hpb=6d32f944a0a115b6419184c50b57bf7c4eef5786 diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index eceb92e36..3b5483b13 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2019 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2021 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 @@ -151,8 +151,6 @@ bool SMESH_MesherHelper::IsQuadraticSubMesh(const TopoDS_Shape& aSh) // we can create quadratic elements only if all elements // created on sub-shapes of given shape are quadratic myCreateQuadratic = true; - mySeamShapeIds.clear(); - myDegenShapeIds.clear(); TopAbs_ShapeEnum subType( aSh.ShapeType()==TopAbs_FACE ? TopAbs_EDGE : TopAbs_FACE ); if ( aSh.ShapeType()==TopAbs_COMPOUND ) { @@ -741,7 +739,8 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face& F, if ( !C2d.IsNull() ) { double u = ( V == IthVertex( 0, edge )) ? f : l; uv = C2d->Value( u ); - uvOK = true; + gp_Pnt p = GetSurface( F )->Value( uv ); + uvOK = ( p.Distance( BRep_Tool::Pnt( V )) < getFaceMaxTol( F )); break; } } @@ -1673,30 +1672,8 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetMediumNode(const SMDS_MeshNode* n1, // get positions of the given nodes on shapes if ( pos.second == TopAbs_FACE ) { - F = TopoDS::Face(meshDS->IndexToShape( faceID = pos.first )); + F = TopoDS::Face( meshDS->IndexToShape( faceID = pos.first )); uv[0] = GetNodeUV(F,n1,n2, force3d ? 0 : &uvOK[0]); - if (( !force3d ) && - ( HasDegeneratedEdges() || GetSurface( F )->HasSingularities( 1e-7 ))) - { - // IPAL52850 (degen VERTEX not at singularity) - // project middle point to a surface - SMESH_TNodeXYZ p1( n1 ), p2( n2 ); - gp_Pnt pMid = 0.5 * ( p1 + p2 ); - Handle(ShapeAnalysis_Surface) projector = GetSurface( F ); - gp_Pnt2d uvMid; - if ( uvOK[0] ) - uvMid = projector->NextValueOfUV( uv[0], pMid, BRep_Tool::Tolerance( F )); - else - uvMid = projector->ValueOfUV( pMid, getFaceMaxTol( F )); - if ( projector->Gap() * projector->Gap() < ( p1 - p2 ).SquareModulus() / 4 ) - { - gp_Pnt pProj = projector->Value( uvMid ); - n12 = meshDS->AddNode( pProj.X(), pProj.Y(), pProj.Z() ); - meshDS->SetNodeOnFace( n12, faceID, uvMid.X(), uvMid.Y() ); - myTLinkNodeMap.insert( make_pair ( link, n12 )); - return n12; - } - } uv[1] = GetNodeUV(F,n2,n1, force3d ? 0 : &uvOK[1]); } else if ( pos.second == TopAbs_EDGE ) @@ -1729,26 +1706,43 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetMediumNode(const SMDS_MeshNode* n1, // nodes, else - medium between corresponding 3d points if( ! F.IsNull() ) { - //if ( uvOK[0] && uvOK[1] ) + if ( IsDegenShape( n1->getshapeId() )) { + if ( myParIndex & U_periodic ) uv[0].SetCoord( 1, uv[1].Coord( 1 )); + else uv[0].SetCoord( 2, uv[1].Coord( 2 )); + } + else if ( IsDegenShape( n2->getshapeId() )) { + if ( myParIndex & U_periodic ) uv[1].SetCoord( 1, uv[0].Coord( 1 )); + else uv[1].SetCoord( 2, uv[0].Coord( 2 )); + } + TopLoc_Location loc; + Handle(Geom_Surface) S = BRep_Tool::Surface(F,loc); + gp_XY UV = GetMiddleUV( S, uv[0], uv[1] ); + gp_Pnt P = S->Value( UV.X(), UV.Y() ).Transformed(loc); + + SMESH_TNodeXYZ p1( n1 ), p2( n2 ); + gp_Pnt pMid = 0.5 * ( p1 + p2 ); + double distMid = pMid.SquareDistance( P ); + double dist12 = ( p1 - p2 ).SquareModulus(); + Handle(ShapeAnalysis_Surface) surfInfo = GetSurface( F ); + if ( distMid > dist12 || + HasDegeneratedEdges() || + surfInfo->HasSingularities( 1e-7 ) ) { - if ( IsDegenShape( n1->getshapeId() )) { - if ( myParIndex & U_periodic ) uv[0].SetCoord( 1, uv[1].Coord( 1 )); - else uv[0].SetCoord( 2, uv[1].Coord( 2 )); - } - else if ( IsDegenShape( n2->getshapeId() )) { - if ( myParIndex & U_periodic ) uv[1].SetCoord( 1, uv[0].Coord( 1 )); - else uv[1].SetCoord( 2, uv[0].Coord( 2 )); - } - TopLoc_Location loc; - Handle(Geom_Surface) S = BRep_Tool::Surface(F,loc); - gp_XY UV = GetMiddleUV( S, uv[0], uv[1] ); - gp_Pnt P = S->Value( UV.X(), UV.Y() ).Transformed(loc); - n12 = meshDS->AddNode(P.X(), P.Y(), P.Z()); - // if ( mySetElemOnShape ) node is not elem! - meshDS->SetNodeOnFace(n12, faceID, UV.X(), UV.Y()); - myTLinkNodeMap.insert(make_pair(link,n12)); - return n12; + // IPAL52850 (degen VERTEX not at singularity) + // project middle point to a surface + gp_Pnt2d uvMid; + if ( uvOK[0] ) + uvMid = surfInfo->NextValueOfUV( uv[0], pMid, BRep_Tool::Tolerance( F )); + else + uvMid = surfInfo->ValueOfUV( pMid, getFaceMaxTol( F )); + if ( surfInfo->Gap() * surfInfo->Gap() < distMid ) + P = surfInfo->Value( uvMid ); } + n12 = meshDS->AddNode(P.X(), P.Y(), P.Z()); + // if ( mySetElemOnShape ) node is not elem! + meshDS->SetNodeOnFace(n12, faceID, UV.X(), UV.Y()); + myTLinkNodeMap.insert(make_pair(link,n12)); + return n12; } else if ( !E.IsNull() ) { @@ -1930,7 +1924,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(); @@ -1959,7 +1953,7 @@ 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(); @@ -1993,7 +1987,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(); @@ -2043,7 +2037,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(); @@ -2107,7 +2101,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(); @@ -2153,7 +2147,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(); @@ -2216,7 +2210,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(); @@ -2257,7 +2251,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; @@ -2307,7 +2301,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(); @@ -2426,8 +2420,8 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1, const SMDS_MeshNode* n10, const SMDS_MeshNode* n11, const SMDS_MeshNode* n12, - const int id, - bool force3d) + const smIdType id, + bool /*force3d*/) { SMESHDS_Mesh * meshDS = GetMeshDS(); SMDS_MeshVolume* elem = 0; @@ -2448,7 +2442,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(); @@ -2978,7 +2972,7 @@ bool SMESH_MesherHelper::IsReversedSubMesh (const TopoDS_Face& theFace) if ( !aSubMeshDSFace ) return isReversed; - // find an element on a bounday of theFace + // find an element on a boundary of theFace SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements(); const SMDS_MeshNode* nn[2]; while ( iteratorElem->more() ) // loop on elements on theFace @@ -3419,9 +3413,9 @@ 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 ) @@ -3797,7 +3791,7 @@ namespace { // Structures used by FixQuadraticElements() //================================================================================ /*! - * \brief Construct QFace from QLinks + * \brief Construct QFace from QLinks */ //================================================================================ @@ -3814,7 +3808,7 @@ namespace { // Structures used by FixQuadraticElements() gp_Vec v1( XYZ( l1->node2()), XYZ( l1->node1())); gp_Vec v2( XYZ( l2->node1()), XYZ( l2->node2())); if ( l1->node1() != l2->node1() && l1->node2() != l2->node2() ) - v1.Reverse(); + v1.Reverse(); _normal += v1 ^ v2; } double normSqSize = _normal.SquareMagnitude(); @@ -3825,6 +3819,8 @@ namespace { // Structures used by FixQuadraticElements() #ifdef _DEBUG_ _face = face; +#else + (void)face; // unused in release mode #endif } //================================================================================ @@ -4136,30 +4132,30 @@ namespace { // Structures used by FixQuadraticElements() */ //================================================================================ - bool QFace::IsSpoiled(const QLink* bentLink ) const - { - // code is valid for convex faces only - gp_XYZ gc(0,0,0); - for ( TIDSortedNodeSet::const_iterator n = begin(); n != end(); ++n ) - gc += XYZ( *n ) / double( size() ); - for ( size_t i = 0; i < _sides.size(); ++i ) - { - if ( _sides[i] == bentLink ) continue; - gp_Vec linkNorm = _normal ^ gp_Vec( XYZ(_sides[i]->node1()), XYZ(_sides[i]->node2())); - gp_Vec vecOut( gc, _sides[i]->MiddlePnt() ); - if ( linkNorm * vecOut < 0 ) - linkNorm.Reverse(); - double mag2 = linkNorm.SquareMagnitude(); - if ( mag2 > numeric_limits::min() ) - linkNorm /= sqrt( mag2 ); - gp_Vec vecBent ( _sides[i]->MiddlePnt(), bentLink->MediumPnt()); - gp_Vec vecStraight( _sides[i]->MiddlePnt(), bentLink->MiddlePnt()); - if ( vecBent * linkNorm > -0.1*vecStraight.Magnitude() ) - return true; - } - return false; + // bool QFace::IsSpoiled(const QLink* bentLink ) const + // { + // // code is valid for convex faces only + // gp_XYZ gc(0,0,0); + // for ( TIDSortedNodeSet::const_iterator n = begin(); n != end(); ++n ) + // gc += XYZ( *n ) / double( size() ); + // for ( size_t i = 0; i < _sides.size(); ++i ) + // { + // if ( _sides[i] == bentLink ) continue; + // gp_Vec linkNorm = _normal ^ gp_Vec( XYZ(_sides[i]->node1()), XYZ(_sides[i]->node2())); + // gp_Vec vecOut( gc, _sides[i]->MiddlePnt() ); + // if ( linkNorm * vecOut < 0 ) + // linkNorm.Reverse(); + // double mag2 = linkNorm.SquareMagnitude(); + // if ( mag2 > numeric_limits::min() ) + // linkNorm /= sqrt( mag2 ); + // gp_Vec vecBent ( _sides[i]->MiddlePnt(), bentLink->MediumPnt()); + // gp_Vec vecStraight( _sides[i]->MiddlePnt(), bentLink->MiddlePnt()); + // if ( vecBent * linkNorm > -0.1*vecStraight.Magnitude() ) + // return true; + // } + // return false; - } + // } //================================================================================ /*! @@ -4619,7 +4615,7 @@ namespace { // Structures used by FixQuadraticElements() if ( curvNorm * D2 > 0 ) continue; // convex edge } - catch ( Standard_Failure ) + catch ( Standard_Failure& ) { continue; } @@ -4651,6 +4647,7 @@ namespace { // Structures used by FixQuadraticElements() const SMDS_MeshElement* f = faceIt->next(); if ( !faceSM->Contains( f ) || f->NbNodes() < 6 || // check quadratic triangles only + f->NbNodes() > 7 || !checkedFaces.insert( f ).second ) continue; @@ -4733,7 +4730,7 @@ namespace { // Structures used by FixQuadraticElements() if ( concaveU || concaveV ) concaveFaces.push_back( face ); } - catch ( Standard_Failure ) + catch ( Standard_Failure& ) { concaveFaces.push_back( face ); } @@ -4862,7 +4859,12 @@ namespace { // Structures used by FixQuadraticElements() gp_Pnt pMedium = SMESH_TNodeXYZ( linkIt->second ); double hMedium = faceNorm * gp_Vec( pOnFace0, pMedium ).XYZ(); double hVol = faceNorm * gp_Vec( pOnFace0, pInSolid ).XYZ(); - isDistorted = ( Abs( hMedium ) > Abs( hVol * 0.75 )); + if ( Abs( hMedium ) > Abs( hVol * 0.75 )) + { + SMESH_TNodeXYZ pI( nOnFace[i]), pJ( nOnFace[j]); + double angle = gp_Vec( pI, pMedium ).Angle( gp_Vec( pI, pJ )); + isDistorted = ( angle > M_PI / 20 ); + } } } } @@ -5260,7 +5262,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, try { gp_Vec x = x01.Normalized() + x12.Normalized(); trsf.SetTransformation( gp_Ax3( gp::Origin(), link1->Normal(), x), gp_Ax3() ); - } catch ( Standard_Failure ) { + } catch ( Standard_Failure& ) { trsf.Invert(); } move.Transform(trsf); @@ -5505,10 +5507,12 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError, volExp.Set( *pentIt, /*ignoreCentralNodes=*/false ); } } -#ifdef _DEBUG_ - // avoid warning: defined but not used operator<<() - SMESH_Comment() << *links.begin() << *faces.begin(); -#endif + + if ( false ) + // avoid warning: defined but not used operator<<() + SMESH_Comment() << *links.begin() << *faces.begin(); + + return; } //================================================================================