X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FSMESH%2FSMESH_MesherHelper.cxx;h=05a27b16ac3077b3f1951e99f87ef5bebb47656a;hb=27cd471609eed565450afd8c9cd3c1125be68abd;hp=bf029c15025e98696808fabc41850847e1c6a09e;hpb=52be03b5646945ae8dd24dda4e037df30a925bd6;p=modules%2Fsmesh.git diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index bf029c150..05a27b16a 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -479,7 +479,7 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, if ( !loc.IsIdentity() ) nodePnt.Transform( loc.Transformation().Inverted() ); if ( Precision::IsInfinite( uv.X() ) || Precision::IsInfinite( uv.Y() ) || - nodePnt.Distance( surface->Value( uv.X(), uv.Y() )) > tol ) + nodePnt.Distance( surface->Value( uv.X(), uv.Y() )) > 2*tol ) { // uv incorrect, project the node to surface GeomAPI_ProjectPointOnSurf& projector = GetProjector( F, loc, tol ); @@ -491,12 +491,12 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, } Quantity_Parameter U,V; projector.LowerDistanceParameters(U,V); - if ( nodePnt.Distance( surface->Value( U, V )) > tol ) + uv.SetCoord( U,V ); + if ( nodePnt.Distance( surface->Value( U, V )) > 2*tol ) { MESSAGE( "SMESH_MesherHelper::CheckNodeUV(), invalid projection" ); return false; } - uv.SetCoord( U,V ); } else if ( uv.Modulus() > numeric_limits::min() ) { @@ -676,12 +676,13 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, return false; } Quantity_Parameter U = projector.LowerDistanceParameter(); + u = double( U ); if ( nodePnt.Distance( curve->Value( U )) > tol ) { MESSAGE( "SMESH_MesherHelper::CheckNodeU(), invalid projection" ); return false; } - u = double( U ); + //u = double( U ); } else if ( fabs( u ) > numeric_limits::min() ) { @@ -1388,6 +1389,8 @@ double SMESH_MesherHelper::GetOtherParam(const double param) const return fabs(param-myPar1[i]) < fabs(param-myPar2[i]) ? myPar2[i] : myPar1[i]; } +//#include + //======================================================================= namespace { // Structures used by FixQuadraticElements() //======================================================================= @@ -1397,7 +1400,12 @@ namespace { // Structures used by FixQuadraticElements() #define MSG(txt) __DMP__(txt< < 1/15 * + return middleNodeMove2 < 1/15./15. * linkLen2; + } struct QFace; // --------------------------------------- @@ -1434,8 +1442,10 @@ namespace { // Structures used by FixQuadraticElements() { _nodeMove += move; _nbMoves += sum ? (_nbMoves==0) : 1; } gp_XYZ Move() const { return _nodeMove.XYZ() / _nbMoves; } bool IsMoved() const { return (_nbMoves > 0 && !IsStraight()); } - bool IsStraight() const { return _nodeMove.SquareMagnitude() <= straightTol2; } - + bool IsStraight() const + { return isStraightLink( (XYZ(node1())-XYZ(node2())).SquareModulus(), + _nodeMove.SquareMagnitude()); + } bool operator<(const QLink& other) const { return (node1()->GetID() == other.node1()->GetID() ? node2()->GetID() < other.node2()->GetID() : @@ -1457,7 +1467,7 @@ namespace { // Structures used by FixQuadraticElements() TChainLink(const QLink* qlink=0):_qlink(qlink) { _qfaces[0] = _qfaces[1] = 0; } - void SetFace(const QFace* face) { int iF = _qfaces[0] ? 1 : 0; _qfaces[iF]=face; } + void SetFace(const QFace* face) const { int iF = _qfaces[0] ? 1 : 0; _qfaces[iF]=face; } bool IsBoundary() const { return !_qfaces[1]; } @@ -1635,6 +1645,7 @@ namespace { // Structures used by FixQuadraticElements() if ( _sides.size() != 4 ) { // triangle - visit all my continous faces MSGBEG( *this ); + TLinkSet links; list< const QFace* > faces( 1, this ); while ( !faces.empty() ) { const QFace* face = faces.front(); @@ -1642,12 +1653,13 @@ namespace { // Structures used by FixQuadraticElements() if ( !face->_sideIsAdded[i] && face->_sides[i] ) { face->_sideIsAdded[i] = true; // find a face side in the chain - TChain::iterator chLink = chain.begin(); - for ( ; chLink != chain.end(); ++chLink ) - if ( chLink->_qlink == face->_sides[i] ) - break; - if ( chLink == chain.end() ) - chLink = chain.insert( chain.begin(), TChainLink(face->_sides[i])); + TLinkInSet chLink = links.insert( TChainLink(face->_sides[i])).first; +// TChain::iterator chLink = chain.begin(); +// for ( ; chLink != chain.end(); ++chLink ) +// if ( chLink->_qlink == face->_sides[i] ) +// break; +// if ( chLink == chain.end() ) +// chLink = chain.insert( chain.begin(), TChainLink(face->_sides[i])); // add a face to a chained link and put a continues face in the queue chLink->SetFace( face ); if ( face->_sides[i]->MediumPos() >= pos ) @@ -1659,6 +1671,7 @@ namespace { // Structures used by FixQuadraticElements() } if ( error < ERR_TRI ) error = ERR_TRI; + chain.insert( chain.end(), links.begin(),links.end() ); return false; } _sideIsAdded[iSide] = true; // not to add this link to chain again @@ -1724,7 +1737,7 @@ namespace { // Structures used by FixQuadraticElements() TLinkInSet link = links.find( _sides[iL] ); if ( link == linksEnd ) continue; if ( (*link)->MediumPos() > SMDS_TOP_FACE ) - continue; // We work on faces here, don't go into a volume + continue; // We work on faces here, don't go inside a solid // check link if ( link->IsBoundary() ) { @@ -1867,7 +1880,8 @@ namespace { // Structures used by FixQuadraticElements() // propagate to adjacent faces till limit step or boundary double len1 = thePrevLen + (theLink->MiddlePnt() - _sides[iL1]->MiddlePnt()).Modulus(); double len2 = thePrevLen + (theLink->MiddlePnt() - _sides[iL2]->MiddlePnt()).Modulus(); - gp_Vec linkDir1, linkDir2; + gp_Vec linkDir1(0,0,0); // initialize to avoid valgrind error ("Conditional jump...") + gp_Vec linkDir2(0,0,0); try { OCC_CATCH_SIGNALS; if ( f1 ) @@ -2103,7 +2117,7 @@ namespace { // Structures used by FixQuadraticElements() enum TSplitTriaResult { _OK, _NO_CORNERS, _FEW_ROWS, _MANY_ROWS, _NO_SIDELINK, _BAD_MIDQUAD, _NOT_RECT, - _NO_MIDQUAD, _NO_UPTRIA, _BAD_SET_SIZE, _BAD_CORNER, _BAD_START, _NO_BOTLINK }; + _NO_MIDQUAD, _NO_UPTRIA, _BAD_SET_SIZE, _BAD_CORNER, _BAD_START, _NO_BOTLINK, _TWISTED_CHAIN }; TSplitTriaResult splitTrianglesIntoChains( TChain & allLinks, vector< TChain> & resultChains, @@ -2169,6 +2183,8 @@ namespace { // Structures used by FixQuadraticElements() const QFace* botTria = botLink->_qfaces[0]; // bottom triangle bound by botLink if ( !botTria ) { // the column ends + if ( botLink == startLink ) + return _TWISTED_CHAIN; // issue 0020951 linkSet.erase( botLink ); if ( iRow != rowChains.size() ) return _FEW_ROWS; // different nb of rows in columns @@ -2228,8 +2244,10 @@ namespace { // Structures used by FixQuadraticElements() // next bottom link ends at the new corner linkSet.erase( botLink ); botLink = upTria->GetLinkByNode( linkSet, (isCase2 ? *sideLink : *midQuadLink), corner ); - if ( botLink == linksEnd || botLink == (isCase2 ? midQuadLink : sideLink)) + if ( botLink == linksEnd || botLink == midQuadLink || botLink == sideLink) return _NO_BOTLINK; + if ( midQuadLink == startLink || sideLink == startLink ) + return _TWISTED_CHAIN; // issue 0020951 linkSet.erase( midQuadLink ); linkSet.erase( sideLink ); @@ -2286,37 +2304,49 @@ namespace { // Structures used by FixQuadraticElements() void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) { - // apply algorithm to solids or geom faces + // 0. Apply algorithm to solids or geom faces // ---------------------------------------------- if ( myShape.IsNull() ) { 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() ); + faces.Add( f.Current() ); // not in solid } for ( TopExp_Explorer s(myShape,TopAbs_SOLID); s.More(); s.Next() ) { if ( myMesh->GetSubMesh( s.Current() )->IsEmpty() ) { // get faces of solid for ( TopExp_Explorer f( s.Current(), TopAbs_FACE); f.More(); f.Next() ) - faces.Add( f.Current() ); + faces.Add( f.Current() ); // in not meshed solid } else { // fix nodes in the solid and its faces + MSG("FIX SOLID " << nbSolids-- << " #" << GetMeshDS()->ShapeToIndex(s.Current())); SMESH_MesherHelper h(*myMesh); h.SetSubShape( s.Current() ); h.FixQuadraticElements(false); } } // fix nodes on geom faces +#ifdef _DEBUG_ + int nbfaces = faces.Extent(); +#endif for ( TopTools_MapIteratorOfMapOfShape fIt( faces ); fIt.More(); fIt.Next() ) { + MSG("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key())); SMESH_MesherHelper h(*myMesh); h.SetSubShape( fIt.Key() ); h.FixQuadraticElements(true); } + //perf_print_all_meters(1); return; } - // Find out type of elements and get iterator on them + // 1. Find out type of elements and get iterator on them // --------------------------------------------------- SMDS_ElemIteratorPtr elemIt; @@ -2335,7 +2365,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) if ( !elemIt || !elemIt->more() || elemType < SMDSAbs_Face ) return; - // Fill in auxiliary data structures + // 2. Fill in auxiliary data structures // ---------------------------------- set< QLink > links; @@ -2344,7 +2374,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) set< QFace >::iterator pFace; bool isCurved = false; - bool hasRectFaces = false; + //bool hasRectFaces = false; set nbElemNodeSet; if ( elemType == SMDSAbs_Volume ) @@ -2377,9 +2407,9 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) if ( pFace->NbVolumes() == 0 ) pFace->AddSelfToLinks(); pFace->SetVolume( vol ); - hasRectFaces = hasRectFaces || - ( volTool.GetVolumeType() == SMDS_VolumeTool::QUAD_HEXA || - volTool.GetVolumeType() == SMDS_VolumeTool::QUAD_PENTA ); +// 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]); @@ -2415,13 +2445,13 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) // store QFace pFace = faces.insert( QFace( faceLinks )).first; pFace->AddSelfToLinks(); - hasRectFaces = ( hasRectFaces || nbN == 4 ); + //hasRectFaces = ( hasRectFaces || nbN == 4 ); } } if ( !isCurved ) return; // no curved edges of faces - // Compute displacement of medium nodes + // 3. Compute displacement of medium nodes // ------------------------------------- // two loops on faces: the first is to treat boundary links, the second is for internal ones @@ -2512,6 +2542,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) { face = TopoDS::Face( f ); 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(); @@ -2524,9 +2555,11 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) ( is1 ? move1 : move0 ).SetCoord( uvMove.X(), uvMove.Y(), 0 ); if ( !is1 ) // correct nodeOnFace for move1 (issue 0020919) nodeOnFace = (*(++chain.rbegin()))->_mediumNode; + isStraight[is1] = isStraightLink( (uv2-uv1).SquareModulus(),uvMove.SquareModulus()); } - if ( move0.SquareMagnitude() < straightTol2 && - move1.SquareMagnitude() < straightTol2 ) { +// if ( move0.SquareMagnitude() < straightTol2 && +// move1.SquareMagnitude() < straightTol2 ) { + if ( isStraight[0] && isStraight[1] ) { MSG("2D straight - ignore"); continue; // straight - no need to move nodes of internal links } @@ -2598,7 +2631,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly) } // loop on faces } - // Move nodes + // 4. Move nodes // ----------- for ( pLink = links.begin(); pLink != links.end(); ++pLink ) {