X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_ViscousLayers.cxx;h=27afe5d0211b61fa20c56c15f6c0702b5cbd8ae7;hp=0c40b47438c6ccf930cbf9e6aed0dfbfdd3ae30b;hb=efbf393dc4a81161ddfeab879b2db6e17aeceb8b;hpb=bc95c31002da130977f9cdbb0a5482741c529104 diff --git a/src/StdMeshers/StdMeshers_ViscousLayers.cxx b/src/StdMeshers/StdMeshers_ViscousLayers.cxx index 0c40b4743..27afe5d02 100644 --- a/src/StdMeshers/StdMeshers_ViscousLayers.cxx +++ b/src/StdMeshers/StdMeshers_ViscousLayers.cxx @@ -30,6 +30,7 @@ #include "SMDS_SetIterator.hxx" #include "SMESHDS_Group.hxx" #include "SMESHDS_Hypothesis.hxx" +#include "SMESHDS_Mesh.hxx" #include "SMESH_Algo.hxx" #include "SMESH_ComputeError.hxx" #include "SMESH_ControlsDef.hxx" @@ -49,6 +50,7 @@ #include #include #include +//#include #include #include #include @@ -93,7 +95,7 @@ #include #ifdef _DEBUG_ -#define __myDEBUG +//#define __myDEBUG //#define __NOT_INVALIDATE_BAD_SMOOTH #endif @@ -457,7 +459,7 @@ namespace VISCOUS_3D void SmoothPos( const vector< double >& segLen, const double tol ); int Smooth(const int step, const bool isConcaveFace, bool findBest); int Smooth(const int step, bool findBest, vector< _LayerEdge* >& toSmooth ); - int CheckNeiborsOnBoundary(vector< _LayerEdge* >* badNeibors = 0 ); + int CheckNeiborsOnBoundary(vector< _LayerEdge* >* badNeibors = 0, bool * needSmooth = 0 ); void SmoothWoCheck(); bool SmoothOnEdge(Handle(ShapeAnalysis_Surface)& surface, const TopoDS_Face& F, @@ -1000,6 +1002,7 @@ namespace VISCOUS_3D { gp_XYZ _xyz; // coord of a point inflated from EDGE w/o smooth double _len; // length reached at previous inflation step + double _param; // on EDGE _2NearEdges _2edges; // 2 neighbor _LayerEdge's double Distance( const OffPnt& p ) const { return ( _xyz - p._xyz ).Modulus(); } }; @@ -1009,6 +1012,7 @@ namespace VISCOUS_3D _LayerEdge _leOnV[2]; // _LayerEdge's holding normal to the EDGE at VERTEXes size_t _iSeg[2]; // index of segment where extreme tgt node is projected _EdgesOnShape& _eos; + double _curveLen; // length of the EDGE static Handle(Geom_Curve) CurveForSmooth( const TopoDS_Edge& E, _EdgesOnShape& eos, @@ -2887,7 +2891,7 @@ bool _ViscousBuilder::findShapesToSmooth( _SolidData& data ) { TGeomID iV = getMeshDS()->ShapeToIndex( vIt.Value() ); vector<_LayerEdge*>& eV = edgesByGeom[ iV ]._edges; - if ( eV.empty() ) continue; + if ( eV.empty() || eV[0]->Is( _LayerEdge::MULTI_NORMAL )) continue; gp_Vec eDir = getEdgeDir( E, TopoDS::Vertex( vIt.Value() )); double angle = eDir.Angle( eV[0]->_normal ); double cosin = Cos( angle ); @@ -3880,11 +3884,11 @@ gp_XYZ _OffsetPlane::GetCommonPoint(bool& isFound) const if ( NbLines() == 2 ) { gp_Vec lPerp0 = _lines[0].Direction().XYZ() ^ _plane.Axis().Direction().XYZ(); - gp_Vec l0l1 = _lines[1].Location().XYZ() - _lines[0].Location().XYZ(); double dot01 = lPerp0 * _lines[1].Direction().XYZ(); if ( Abs( dot01 ) > std::numeric_limits::min() ) { - double u1 = - ( lPerp0 * l0l1 ) / dot01; + gp_Vec l0l1 = _lines[1].Location().XYZ() - _lines[0].Location().XYZ(); + double u1 = - ( lPerp0 * l0l1 ) / dot01; p = ( _lines[1].Location().XYZ() + _lines[1].Direction().XYZ() * u1 ); isFound = true; } @@ -4267,7 +4271,10 @@ bool _ViscousBuilder::inflate(_SolidData& data) if ( data._edgesOnShape[i].ShapeType() == TopAbs_VERTEX && data._edgesOnShape[i]._edges.size() > 0 && data._edgesOnShape[i]._edges[0]->Is( _LayerEdge::MULTI_NORMAL )) + { + data._edgesOnShape[i]._edges[0]->Unset( _LayerEdge::BLOCKED ); data._edgesOnShape[i]._edges[0]->Block( data ); + } const double safeFactor = ( 2*data._maxThickness < data._geomSize ) ? 1 : theThickToIntersection; @@ -4414,7 +4421,7 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, bool moved, improved; double vol; - vector< _LayerEdge* > movedEdges, badSmooEdges; + vector< _LayerEdge* > movedEdges, badEdges; vector< _EdgesOnShape* > eosC1; // C1 continues shapes vector< bool > isConcaveFace; @@ -4517,13 +4524,13 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, makeOffsetSurface( *eosC1[ iEOS ], helper ); } - int step = 0, stepLimit = 5, badNb = 0; + int step = 0, stepLimit = 5, nbBad = 0; while (( ++step <= stepLimit ) || improved ) { dumpFunction(SMESH_Comment("smooth")<Unset( _LayerEdge::SMOOTHED ); if ( movedEdges[i]->Smooth( step, findBest, movedEdges ) > 0 ) - badSmooEdges.push_back( movedEdges[i] ); + badEdges.push_back( movedEdges[i] ); } #else bool findBest = ( step == stepLimit || isConcaveFace[ iEOS ]); @@ -4542,18 +4549,18 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, { edges[i]->Unset( _LayerEdge::SMOOTHED ); if ( edges[i]->Smooth( step, findBest, false ) > 0 ) - badSmooEdges.push_back( eos._edges[i] ); + badEdges.push_back( eos._edges[i] ); } } #endif - badNb = badSmooEdges.size(); + nbBad = badEdges.size(); - if ( badNb > 0 ) - debugMsg(SMESH_Comment("badNb = ") << badNb ); + if ( nbBad > 0 ) + debugMsg(SMESH_Comment("nbBad = ") << nbBad ); - if ( !badSmooEdges.empty() && step >= stepLimit / 2 ) + if ( !badEdges.empty() && step >= stepLimit / 2 ) { - if ( badSmooEdges[0]->Is( _LayerEdge::ON_CONCAVE_FACE )) + if ( badEdges[0]->Is( _LayerEdge::ON_CONCAVE_FACE )) stepLimit = 9; // resolve hard smoothing situation around concave VERTEXes @@ -4562,26 +4569,26 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, vector< _EdgesOnShape* > & eosCoVe = eosC1[ iEOS ]->_eosConcaVer; for ( size_t i = 0; i < eosCoVe.size(); ++i ) eosCoVe[i]->_edges[0]->MoveNearConcaVer( eosCoVe[i], eosC1[ iEOS ], - step, badSmooEdges ); + step, badEdges ); } - // look for the best smooth of _LayerEdge's neighboring badSmooEdges - badNb = 0; - for ( size_t i = 0; i < badSmooEdges.size(); ++i ) + // look for the best smooth of _LayerEdge's neighboring badEdges + nbBad = 0; + for ( size_t i = 0; i < badEdges.size(); ++i ) { - _LayerEdge* ledge = badSmooEdges[i]; + _LayerEdge* ledge = badEdges[i]; for ( size_t iN = 0; iN < ledge->_neibors.size(); ++iN ) { ledge->_neibors[iN]->Unset( _LayerEdge::SMOOTHED ); - badNb += ledge->_neibors[iN]->Smooth( step, true, /*findBest=*/true ); + nbBad += ledge->_neibors[iN]->Smooth( step, true, /*findBest=*/true ); } ledge->Unset( _LayerEdge::SMOOTHED ); - badNb += ledge->Smooth( step, true, /*findBest=*/true ); + nbBad += ledge->Smooth( step, true, /*findBest=*/true ); } - debugMsg(SMESH_Comment("badNb = ") << badNb ); + debugMsg(SMESH_Comment("nbBad = ") << nbBad ); } - if ( badNb == oldBadNb && - badNb > 0 && + if ( nbBad == oldBadNb && + nbBad > 0 && step < stepLimit ) // smooth w/o chech of validity { dumpFunctionEnd(); @@ -4595,11 +4602,11 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, stepLimit++; } - improved = ( badNb < oldBadNb ); + improved = ( nbBad < oldBadNb ); dumpFunctionEnd(); - if (( step % 3 == 1 ) || ( badNb > 0 && step >= stepLimit / 2 )) + if (( step % 3 == 1 ) || ( nbBad > 0 && step >= stepLimit / 2 )) for ( size_t iEOS = 0; iEOS < eosC1.size(); ++iEOS ) { putOnOffsetSurface( *eosC1[ iEOS ], infStep, step, /*moveAll=*/step == 1 ); @@ -4610,13 +4617,13 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, // project -- to prevent intersections or fix bad simplices for ( size_t iEOS = 0; iEOS < eosC1.size(); ++iEOS ) { - if ( ! eosC1[ iEOS ]->_eosConcaVer.empty() || badNb > 0 ) + if ( ! eosC1[ iEOS ]->_eosConcaVer.empty() || nbBad > 0 ) putOnOffsetSurface( *eosC1[ iEOS ], infStep ); } - if ( !badSmooEdges.empty() ) + if ( !badEdges.empty() ) { - badSmooEdges.clear(); + badEdges.clear(); for ( size_t iEOS = 0; iEOS < eosC1.size(); ++iEOS ) { for ( size_t i = 0; i < eosC1[ iEOS ]->_edges.size(); ++i ) @@ -4624,8 +4631,8 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, if ( !eosC1[ iEOS ]->_sWOL.IsNull() ) continue; _LayerEdge* edge = eosC1[ iEOS ]->_edges[i]; - edge->CheckNeiborsOnBoundary( & badSmooEdges ); - if ( badNb > 0 ) + edge->CheckNeiborsOnBoundary( & badEdges ); + if ( nbBad > 0 ) { SMESH_TNodeXYZ tgtXYZ = edge->_nodes.back(); const gp_XYZ& prevXYZ = edge->PrevCheckPos(); @@ -4636,7 +4643,7 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, << " "<< tgtXYZ._node->GetID() << " "<< edge->_simplices[j]._nPrev->GetID() << " "<< edge->_simplices[j]._nNext->GetID() << " )" ); - badSmooEdges.push_back( edge ); + badEdges.push_back( edge ); break; } } @@ -4644,9 +4651,9 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, } // try to fix bad simplices by removing the last inflation step of some _LayerEdge's - badNb = invalidateBadSmooth( data, helper, badSmooEdges, eosC1, infStep ); + nbBad = invalidateBadSmooth( data, helper, badEdges, eosC1, infStep ); - if ( badNb > 0 ) + if ( nbBad > 0 ) return false; } @@ -4664,14 +4671,14 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, !eos._sWOL.IsNull() ) continue; - badSmooEdges.clear(); + badEdges.clear(); for ( size_t i = 0; i < eos._edges.size(); ++i ) { _LayerEdge* edge = eos._edges[i]; if ( edge->_nodes.size() < 2 ) continue; SMESH_TNodeXYZ tgtXYZ = edge->_nodes.back(); - //const gp_XYZ& prevXYZ = edge->PrevCheckPos(); - const gp_XYZ& prevXYZ = edge->PrevPos(); + const gp_XYZ& prevXYZ = edge->PrevCheckPos(); + //const gp_XYZ& prevXYZ = edge->PrevPos(); for ( size_t j = 0; j < edge->_simplices.size(); ++j ) if ( !edge->_simplices[j].IsForward( &prevXYZ, &tgtXYZ, vol )) { @@ -4679,15 +4686,15 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, << " "<< tgtXYZ._node->GetID() << " "<< edge->_simplices[j]._nPrev->GetID() << " "<< edge->_simplices[j]._nNext->GetID() << " )" ); - badSmooEdges.push_back( edge ); + badEdges.push_back( edge ); break; } } // try to fix bad simplices by removing the last inflation step of some _LayerEdge's eosC1[0] = &eos; - int badNb = invalidateBadSmooth( data, helper, badSmooEdges, eosC1, infStep ); - if ( badNb > 0 ) + int nbBad = invalidateBadSmooth( data, helper, badEdges, eosC1, infStep ); + if ( nbBad > 0 ) return false; } @@ -4720,12 +4727,49 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, eos._edges[i]->Is( _LayerEdge::MULTI_NORMAL )) continue; if ( eos._edges[i]->FindIntersection( *searcher, dist, data._epsilon, eos, &intFace )) + { return false; + // commented due to "Illegal hash-positionPosition" error in NETGEN + // on Debian60 on viscous_layers_01/B2 case + // Collision; try to deflate _LayerEdge's causing it + // badEdges.clear(); + // badEdges.push_back( eos._edges[i] ); + // eosC1[0] = & eos; + // int nbBad = invalidateBadSmooth( data, helper, badEdges, eosC1, infStep ); + // if ( nbBad > 0 ) + // return false; + + // badEdges.clear(); + // if ( _EdgesOnShape* eof = data.GetShapeEdges( intFace->getshapeId() )) + // { + // if ( const _TmpMeshFace* f = dynamic_cast< const _TmpMeshFace*>( intFace )) + // { + // const SMDS_MeshElement* srcFace = + // eof->_subMesh->GetSubMeshDS()->GetElement( f->getIdInShape() ); + // SMDS_ElemIteratorPtr nIt = srcFace->nodesIterator(); + // while ( nIt->more() ) + // { + // const SMDS_MeshNode* srcNode = static_cast( nIt->next() ); + // TNode2Edge::iterator n2e = data._n2eMap.find( srcNode ); + // if ( n2e != data._n2eMap.end() ) + // badEdges.push_back( n2e->second ); + // } + // eosC1[0] = eof; + // nbBad = invalidateBadSmooth( data, helper, badEdges, eosC1, infStep ); + // if ( nbBad > 0 ) + // return false; + // } + // } + // if ( eos._edges[i]->FindIntersection( *searcher, dist, data._epsilon, eos, &intFace )) + // return false; + // else + // continue; + } if ( !intFace ) { SMESH_Comment msg("Invalid? normal at node "); msg << eos._edges[i]->_nodes[0]->GetID(); debugMsg( msg ); - continue; + continue; } const bool isShorterDist = ( distToIntersection > dist ); @@ -4846,18 +4890,34 @@ int _ViscousBuilder::invalidateBadSmooth( _SolidData& data, dumpFunction(SMESH_Comment("invalidateBadSmooth")<<"_S"<_shapeID<<"_InfStep"<NbSteps() < 2 || edge->Is( _LayerEdge::MARKED )) + if ( edge->NbSteps() < 2 /*|| edge->Is( _LayerEdge::MARKED )*/) continue; + _EdgesOnShape* eos = data.GetShapeEdges( edge ); edge->InvalidateStep( edge->NbSteps(), *eos, /*restoreLength=*/true ); edge->Block( data ); - edge->Set( _LayerEdge::MARKED ); + //edge->Set( _LayerEdge::MARKED ); + + // look for _LayerEdge's of bad _simplices + SMESH_TNodeXYZ tgtXYZ = edge->_nodes.back(); + const gp_XYZ& prevXYZ1 = edge->PrevCheckPos(); + const gp_XYZ& prevXYZ2 = edge->PrevPos(); + for ( size_t j = 0; j < edge->_simplices.size(); ++j ) + { + if (( edge->_simplices[j].IsForward( &prevXYZ1, &tgtXYZ, vol )) && + ( &prevXYZ1 == &prevXYZ2 || edge->_simplices[j].IsForward( &prevXYZ2, &tgtXYZ, vol ))) + continue; + for ( size_t iN = 0; iN < edge->_neibors.size(); ++iN ) + if ( edge->_simplices[j].Includes( edge->_neibors[iN]->_nodes.back() )) + badSmooEdges.push_back( edge->_neibors[iN] ); + } if ( eos->ShapeType() == TopAbs_VERTEX ) { @@ -4874,22 +4934,14 @@ int _ViscousBuilder::invalidateBadSmooth( _SolidData& data, } eoe->_edgeSmoother->Perform( data, surface, F, helper ); } - } - // look for _LayerEdge's of bad _simplices - SMESH_TNodeXYZ tgtXYZ = edge->_nodes.back(); - const gp_XYZ& prevXYZ = edge->PrevCheckPos(); - for ( size_t j = 0; j < edge->_simplices.size(); ++j ) - { - if ( edge->_simplices[j].IsForward( &prevXYZ, &tgtXYZ, vol )) - continue; - for ( size_t iN = 0; iN < edge->_neibors.size(); ++iN ) - if ( edge->_simplices[j].Includes( edge->_neibors[iN]->_nodes.back() )) - badSmooEdges.push_back( edge->_neibors[iN] ); } - } + } // loop on badSmooEdges + - int badNb = 0; + // check result of invalidation + + int nbBad = 0; for ( size_t iEOS = 0; iEOS < eosC1.size(); ++iEOS ) { for ( size_t i = 0; i < eosC1[ iEOS ]->_edges.size(); ++i ) @@ -4901,7 +4953,7 @@ int _ViscousBuilder::invalidateBadSmooth( _SolidData& data, for ( size_t j = 0; j < edge->_simplices.size(); ++j ) if ( !edge->_simplices[j].IsForward( &prevXYZ, &tgtXYZ, vol )) { - ++badNb; + ++nbBad; debugMsg("Bad simplex remains ( " << edge->_nodes[0]->GetID() << " "<< tgtXYZ._node->GetID() << " "<< edge->_simplices[j]._nPrev->GetID() @@ -4911,7 +4963,7 @@ int _ViscousBuilder::invalidateBadSmooth( _SolidData& data, } dumpFunctionEnd(); - return badNb; + return nbBad; } //================================================================================ @@ -5153,8 +5205,6 @@ bool _Smoother1D::smoothAnalyticEdge( _SolidData& data, gp_XYZ newPos; for ( size_t i = iFrom; i < iTo; ++i ) { - if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue; - _LayerEdge* edge = _eos._edges[i]; SMDS_MeshNode* tgtNode = const_cast( edge->_nodes.back() ); newPos = p0 * ( 1. - _leParams[i] ) + p1 * _leParams[i]; @@ -5167,12 +5217,20 @@ bool _Smoother1D::smoothAnalyticEdge( _SolidData& data, lineDir * ( curPos - pSrc0 )); newPos = curPos + lineDir * shift / lineDir.SquareModulus(); } + if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) + { + SMESH_TNodeXYZ pSrc( edge->_nodes[0] ); + double curThick = pSrc.SquareDistance( tgtNode ); + double newThink = ( pSrc - newPos ).SquareModulus(); + if ( newThink > curThick ) + continue; + } edge->_pos.back() = newPos; tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() ); dumpMove( tgtNode ); } } - else + else // 2D { _LayerEdge* e0 = getLEdgeOnV( 0 ); _LayerEdge* e1 = getLEdgeOnV( 1 ); @@ -5312,7 +5370,7 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, if ( _offPoints.empty() ) return false; - // move _offPoints to a new position + // move _offPoints to positions along normals of _LayerEdge's _LayerEdge* e[2] = { getLEdgeOnV(0), getLEdgeOnV(1) }; if ( e[0]->Is( _LayerEdge::NORMAL_UPDATED )) setNormalOnV( 0, helper ); @@ -5333,7 +5391,7 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, _offPoints[i]._len = avgLen; } - double fTol; + double fTol = 0; if ( !surface.IsNull() ) // project _offPoints to the FACE { fTol = 100 * BRep_Tool::Tolerance( F ); @@ -5360,7 +5418,8 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, int i = _iSeg[ is2nd ]; int di = is2nd ? -1 : +1; bool projected = false; - double uOnSeg, uOnSegDiff, uOnSegBestDiff = Precision::Infinite(); + double uOnSeg, uOnSegDiff, uOnSegBestDiff = Precision::Infinite(), uOnSegPrevDiff = 0; + int nbWorse = 0; do { gp_Vec v0p( _offPoints[i]._xyz, pExtreme[ is2nd ] ); gp_Vec v01( _offPoints[i]._xyz, _offPoints[i+1]._xyz ); @@ -5373,6 +5432,12 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, pProj[ is2nd ] = _offPoints[i]._xyz + ( v01 * uOnSeg ).XYZ(); uOnSegBestDiff = uOnSegDiff; } + else if ( uOnSegDiff > uOnSegPrevDiff ) + { + if ( ++nbWorse > 3 ) // avoid projection to the middle of a closed EDGE + break; + } + uOnSegPrevDiff = uOnSegDiff; i += di; } while ( !projected && @@ -5462,10 +5527,14 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, void _Smoother1D::prepare(_SolidData& data) { - // sort _LayerEdge's by position on the EDGE const TopoDS_Edge& E = TopoDS::Edge( _eos._shape ); + _curveLen = SMESH_Algo::EdgeLength( E ); + + // sort _LayerEdge's by position on the EDGE data.SortOnEdge( E, _eos._edges ); + SMESH_MesherHelper& helper = data.GetHelper(); + // compute normalized param of _eos._edges on EDGE _leParams.resize( _eos._edges.size() + 1 ); { @@ -5482,10 +5551,46 @@ void _Smoother1D::prepare(_SolidData& data) pPrev = p; } double fullLen = _leParams.back() + pPrev.Distance( SMESH_TNodeXYZ( getLEdgeOnV(1)->_nodes[0])); - for ( size_t i = 0; i < _leParams.size(); ++i ) + for ( size_t i = 0; i < _leParams.size()-1; ++i ) _leParams[i] = _leParams[i+1] / fullLen; } + // find intersection of neighbor _LayerEdge's to limit _maxLen + // according to EDGE curvature (IPAL52648) + _LayerEdge* e0 = _eos._edges[0]; + for ( size_t i = 1; i < _eos._edges.size(); ++i ) + { + _LayerEdge* ei = _eos._edges[i]; + gp_XYZ plnNorm = e0->_normal ^ ei->_normal; + gp_XYZ perp0 = e0->_normal ^ plnNorm; + double dot0i = perp0 * ei->_normal; + if ( Abs( dot0i ) > std::numeric_limits::min() ) + { + SMESH_TNodeXYZ srci( ei->_nodes[0] ), src0( e0->_nodes[0] ); + double ui = ( perp0 * ( src0 - srci )) / dot0i; + if ( ui > 0 ) + { + ei->_maxLen = Min( ei->_maxLen, 0.75 * ui / ei->_lenFactor ); + if ( ei->_maxLen < ei->_len ) + { + ei->InvalidateStep( ei->NbSteps(), _eos, /*restoreLength=*/true ); + ei->SetNewLength( ei->_maxLen, _eos, helper ); + ei->Block( data ); + } + gp_Pnt pi = srci + ei->_normal * ui; + double u0 = pi.Distance( src0 ); + e0->_maxLen = Min( e0->_maxLen, 0.75 * u0 / e0->_lenFactor ); + if ( e0->_maxLen < e0->_len ) + { + e0->InvalidateStep( e0->NbSteps(), _eos, /*restoreLength=*/true ); + e0->SetNewLength( e0->_maxLen, _eos, helper ); + e0->Block( data ); + } + } + } + e0 = ei; + } + if ( isAnalytic() ) return; @@ -7459,7 +7564,7 @@ void _LayerEdge::SmoothWoCheck() */ //================================================================================ -int _LayerEdge::CheckNeiborsOnBoundary( vector< _LayerEdge* >* badNeibors ) +int _LayerEdge::CheckNeiborsOnBoundary( vector< _LayerEdge* >* badNeibors, bool * needSmooth ) { if ( ! Is( NEAR_BOUNDARY )) return 0; @@ -7471,6 +7576,9 @@ int _LayerEdge::CheckNeiborsOnBoundary( vector< _LayerEdge* >* badNeibors ) _LayerEdge* eN = _neibors[iN]; if ( eN->_nodes[0]->getshapeId() == _nodes[0]->getshapeId() ) continue; + if ( needSmooth ) + *needSmooth |= ( eN->Is( _LayerEdge::BLOCKED ) || eN->Is( _LayerEdge::NORMAL_UPDATED )); + SMESH_TNodeXYZ curPosN ( eN->_nodes.back() ); SMESH_TNodeXYZ prevPosN( eN->_nodes[0] ); for ( size_t i = 0; i < eN->_simplices.size(); ++i ) @@ -7506,7 +7614,7 @@ int _LayerEdge::CheckNeiborsOnBoundary( vector< _LayerEdge* >* badNeibors ) int _LayerEdge::Smooth(const int step, bool findBest, vector< _LayerEdge* >& toSmooth ) { - if ( !Is( MOVED ) || Is( SMOOTHED )) + if ( !Is( MOVED ) || Is( SMOOTHED ) || Is( BLOCKED )) return 0; // shape of simplices not changed if ( _simplices.size() < 2 ) return 0; // _LayerEdge inflated along EDGE or FACE @@ -7527,13 +7635,14 @@ int _LayerEdge::Smooth(const int step, bool findBest, vector< _LayerEdge* >& toS } int nbBad = _simplices.size() - nbOkBefore; + bool bndNeedSmooth = false; if ( nbBad == 0 ) - nbBad = CheckNeiborsOnBoundary(); + nbBad = CheckNeiborsOnBoundary( 0, & bndNeedSmooth ); if ( nbBad > 0 ) Set( DISTORTED ); // evaluate min angle - if ( nbBad == 0 && !findBest ) + if ( nbBad == 0 && !findBest && !bndNeedSmooth ) { size_t nbGoodAngles = _simplices.size(); double angle; @@ -8581,13 +8690,19 @@ void _LayerEdge::Block( _SolidData& data ) minDist = Min( pSrc.SquareDistance( pTgtN ), minDist ); minDist = Min( pTgt.SquareDistance( pSrcN ), minDist ); double newMaxLen = edge->_maxLen + 0.5 * Sqrt( minDist ); + if ( edge->_nodes[0]->getshapeId() == neibor->_nodes[0]->getshapeId() ) + { + newMaxLen *= edge->_lenFactor / neibor->_lenFactor; + } if ( neibor->_maxLen > newMaxLen ) { neibor->_maxLen = newMaxLen; if ( neibor->_maxLen < neibor->_len ) { _EdgesOnShape* eos = data.GetShapeEdges( neibor ); - neibor->InvalidateStep( neibor->NbSteps(), *eos, /*restoreLength=*/true ); + while ( neibor->_len > neibor->_maxLen && + neibor->NbSteps() > 1 ) + neibor->InvalidateStep( neibor->NbSteps(), *eos, /*restoreLength=*/true ); neibor->SetNewLength( neibor->_maxLen, *eos, data.GetHelper() ); } queue.push( neibor ); @@ -8635,7 +8750,9 @@ void _LayerEdge::InvalidateStep( size_t curStep, const _EdgesOnShape& eos, bool dumpMove( n ); if ( restoreLength ) + { _len -= ( nXYZ.XYZ() - curXYZ ).Modulus() / _lenFactor; + } } } @@ -8648,7 +8765,7 @@ void _LayerEdge::InvalidateStep( size_t curStep, const _EdgesOnShape& eos, bool void _LayerEdge::SmoothPos( const vector< double >& segLen, const double tol ) { //return; - if ( Is( NORMAL_UPDATED ) || _pos.size() <= 2 ) + if ( /*Is( NORMAL_UPDATED ) ||*/ _pos.size() <= 2 ) return; // find the 1st smoothed _pos @@ -8665,30 +8782,54 @@ void _LayerEdge::SmoothPos( const vector< double >& segLen, const double tol ) { // if ( segLen[ iSmoothed ] / segLen.back() < 0.5 ) // return; - for ( size_t i = Max( 1, iSmoothed-1 ); i < _pos.size()-1; ++i ) - { - gp_XYZ midPos = 0.5 * ( _pos[i-1] + _pos[i+1] ); - gp_XYZ newPos = 0.5 * ( midPos + _pos[i] ); - _pos[i] = newPos; - double midLen = 0.5 * ( segLen[i-1] + segLen[i+1] ); - double newLen = 0.5 * ( midLen + segLen[i] ); - const_cast< double& >( segLen[i] ) = newLen; - } - } - else - { - for ( size_t i = 1; i < _pos.size()-1; ++i ) - { - if ((int) i < iSmoothed && ( segLen[i] / segLen.back() < 0.5 )) - continue; - - double wgt = segLen[i] / segLen.back(); - gp_XYZ normPos = _pos[0] + _normal * wgt * _len; - gp_XYZ tgtPos = ( 1 - wgt ) * _pos[0] + wgt * _pos.back(); - gp_XYZ newPos = ( 1 - wgt ) * normPos + wgt * tgtPos; - _pos[i] = newPos; + gp_XYZ normal = _normal; + if ( Is( NORMAL_UPDATED )) + for ( size_t i = 1; i < _pos.size(); ++i ) + { + normal = _pos[i] - _pos[0]; + double size = normal.Modulus(); + if ( size > RealSmall() ) + { + normal /= size; + break; + } + } + const double r = 0.2; + for ( int iter = 0; iter < 3; ++iter ) + { + double minDot = 1; + for ( size_t i = Max( 1, iSmoothed-1-iter ); i < _pos.size()-1; ++i ) + { + gp_XYZ midPos = 0.5 * ( _pos[i-1] + _pos[i+1] ); + gp_XYZ newPos = ( 1-r ) * midPos + r * _pos[i]; + _pos[i] = newPos; + double midLen = 0.5 * ( segLen[i-1] + segLen[i+1] ); + double newLen = ( 1-r ) * midLen + r * segLen[i]; + const_cast< double& >( segLen[i] ) = newLen; + // check angle between normal and (_pos[i+1], _pos[i] ) + gp_XYZ posDir = _pos[i+1] - _pos[i]; + double size = posDir.Modulus(); + if ( size > RealSmall() ) + minDot = Min( minDot, ( normal * posDir ) / size ); + } + if ( minDot > 0.5 ) + break; } } + // else + // { + // for ( size_t i = 1; i < _pos.size()-1; ++i ) + // { + // if ((int) i < iSmoothed && ( segLen[i] / segLen.back() < 0.5 )) + // continue; + + // double wgt = segLen[i] / segLen.back(); + // gp_XYZ normPos = _pos[0] + _normal * wgt * _len; + // gp_XYZ tgtPos = ( 1 - wgt ) * _pos[0] + wgt * _pos.back(); + // gp_XYZ newPos = ( 1 - wgt ) * normPos + wgt * tgtPos; + // _pos[i] = newPos; + // } + // } } //================================================================================ @@ -8707,7 +8848,7 @@ bool _ViscousBuilder::refine(_SolidData& data) TopoDS_Edge geomEdge; TopoDS_Face geomFace; TopLoc_Location loc; - double f,l, u; + double f,l, u = 0; gp_XY uv; vector< gp_XYZ > pos3D; bool isOnEdge; @@ -8830,7 +8971,12 @@ bool _ViscousBuilder::refine(_SolidData& data) { swapped = false; for ( size_t j = 1; j < edge._pos.size(); ++j ) - if ( segLen[j] < segLen[j-1] ) + if ( segLen[j] > segLen.back() ) + { + segLen.erase( segLen.begin() + j ); + edge._pos.erase( edge._pos.begin() + j ); + } + else if ( segLen[j] < segLen[j-1] ) { std::swap( segLen[j], segLen[j-1] ); std::swap( edge._pos[j], edge._pos[j-1] ); @@ -8839,8 +8985,8 @@ bool _ViscousBuilder::refine(_SolidData& data) } } // smooth a path formed by edge._pos - if (( smoothed ) && - ( eos.ShapeType() == TopAbs_FACE || edge.Is( _LayerEdge::SMOOTHED_C1 ))) + if (( smoothed ) /*&& + ( eos.ShapeType() == TopAbs_FACE || edge.Is( _LayerEdge::SMOOTHED_C1 ))*/) edge.SmoothPos( segLen, preci ); } else if ( eos._isRegularSWOL ) // usual SWOL @@ -9392,7 +9538,7 @@ bool _ViscousBuilder::shrink() // ================== bool shrinked = true; - int badNb, shriStep=0, smooStep=0; + int nbBad, shriStep=0, smooStep=0; _SmoothNode::SmoothType smoothType = isConcaveFace ? _SmoothNode::ANGULAR : _SmoothNode::LAPLACIAN; SMESH_Comment errMsg; @@ -9423,22 +9569,22 @@ bool _ViscousBuilder::shrink() // ----------------- int nbNoImpSteps = 0; bool moved = true; - badNb = 1; - while (( nbNoImpSteps < 5 && badNb > 0) && moved) + nbBad = 1; + while (( nbNoImpSteps < 5 && nbBad > 0) && moved) { dumpFunction(SMESH_Comment("shrinkFace")<first<<"_st"<<++smooStep); // debug - int oldBadNb = badNb; - badNb = 0; + int oldBadNb = nbBad; + nbBad = 0; moved = false; // '% 5' minimizes NB FUNCTIONS on viscous_layers_00/B2 case _SmoothNode::SmoothType smooTy = ( smooStep % 5 ) ? smoothType : _SmoothNode::LAPLACIAN; for ( size_t i = 0; i < nodesToSmooth.size(); ++i ) { - moved |= nodesToSmooth[i].Smooth( badNb, surface, helper, refSign, + moved |= nodesToSmooth[i].Smooth( nbBad, surface, helper, refSign, smooTy, /*set3D=*/isConcaveFace); } - if ( badNb < oldBadNb ) + if ( nbBad < oldBadNb ) nbNoImpSteps = 0; else nbNoImpSteps++; @@ -9447,7 +9593,7 @@ bool _ViscousBuilder::shrink() } errMsg.clear(); - if ( badNb > 0 ) + if ( nbBad > 0 ) errMsg << "Can't shrink 2D mesh on face " << f2sd->first; if ( shriStep > 200 ) errMsg << "Infinite loop at shrinking 2D mesh on face " << f2sd->first; @@ -9494,7 +9640,7 @@ bool _ViscousBuilder::shrink() // dumpFunction(SMESH_Comment("shrinkFace")<first<<"_st"<<++smooStep); // debug // for ( size_t i = 0; i < nodesToSmooth.size(); ++i ) // { - // nodesToSmooth[i].Smooth( badNb,surface,helper,refSign, + // nodesToSmooth[i].Smooth( nbBad,surface,helper,refSign, // _SmoothNode::LAPLACIAN,/*set3D=*/false); // } // } @@ -9681,7 +9827,7 @@ bool _ViscousBuilder::shrink() dumpFunction(SMESH_Comment("shrinkFace")<first<<"_st"<<++smooStep); // debug for ( size_t i = 0; i < nodesToSmooth.size(); ++i ) { - nodesToSmooth[i].Smooth( badNb,surface,helper,refSign, + nodesToSmooth[i].Smooth( nbBad,surface,helper,refSign, smoothType,/*set3D=*/st==1 ); } dumpFunctionEnd(); @@ -10094,7 +10240,7 @@ bool _LayerEdge::SetNewLength2d( Handle(Geom_Surface)& surface, */ //================================================================================ -bool _SmoothNode::Smooth(int& badNb, +bool _SmoothNode::Smooth(int& nbBad, Handle(Geom_Surface)& surface, SMESH_MesherHelper& helper, const double refSign, @@ -10175,7 +10321,7 @@ bool _SmoothNode::Smooth(int& badNb, if ( nbOkAfter < nbOkBefore ) { - badNb += _simplices.size() - nbOkBefore; + nbBad += _simplices.size() - nbOkBefore; return false; } @@ -10193,7 +10339,7 @@ bool _SmoothNode::Smooth(int& badNb, dumpMove( _node ); } - badNb += _simplices.size() - nbOkAfter; + nbBad += _simplices.size() - nbOkAfter; return ( (tgtUV-newPos).SquareModulus() > 1e-10 ); }