X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESHUtils%2FSMESH_MeshAlgos.cxx;h=1eae5f7818e23d205812d4330a7e8437c34349cb;hp=aa702977d6cea5bdc95d0eb264cd8f18d2254323;hb=d1bb1f5d44a2566316419a238a615bc4a69e7028;hpb=ab2525d9053d2d046d22e816b2b72229b9197cbe diff --git a/src/SMESHUtils/SMESH_MeshAlgos.cxx b/src/SMESHUtils/SMESH_MeshAlgos.cxx index aa702977d..1eae5f781 100644 --- a/src/SMESHUtils/SMESH_MeshAlgos.cxx +++ b/src/SMESHUtils/SMESH_MeshAlgos.cxx @@ -441,7 +441,8 @@ struct SMESH_ElementSearcherImpl: public SMESH_ElementSearcher { SMDS_Mesh* _mesh; SMDS_ElemIteratorPtr _meshPartIt; - ElementBndBoxTree* _ebbTree; + ElementBndBoxTree* _ebbTree [SMDSAbs_NbElementTypes]; + int _ebbTreeHeight[SMDSAbs_NbElementTypes]; SMESH_NodeSearcherImpl* _nodeSearcher; SMDSAbs_ElementType _elementType; double _tolerance; @@ -451,10 +452,21 @@ struct SMESH_ElementSearcherImpl: public SMESH_ElementSearcher SMESH_ElementSearcherImpl( SMDS_Mesh& mesh, double tol=-1, SMDS_ElemIteratorPtr elemIt=SMDS_ElemIteratorPtr()) - : _mesh(&mesh),_meshPartIt(elemIt),_ebbTree(0),_nodeSearcher(0),_tolerance(tol),_outerFacesFound(false) {} + : _mesh(&mesh),_meshPartIt(elemIt),_nodeSearcher(0),_tolerance(tol),_outerFacesFound(false) + { + for ( int i = 0; i < SMDSAbs_NbElementTypes; ++i ) + { + _ebbTree[i] = NULL; + _ebbTreeHeight[i] = -1; + } + _elementType = SMDSAbs_All; + } virtual ~SMESH_ElementSearcherImpl() { - if ( _ebbTree ) delete _ebbTree; _ebbTree = 0; + for ( int i = 0; i < SMDSAbs_NbElementTypes; ++i ) + { + delete _ebbTree[i]; _ebbTree[i] = NULL; + } if ( _nodeSearcher ) delete _nodeSearcher; _nodeSearcher = 0; } virtual int FindElementsByPoint(const gp_Pnt& point, @@ -479,6 +491,12 @@ struct SMESH_ElementSearcherImpl: public SMESH_ElementSearcher { return _outerFaces.empty() || _outerFaces.count(face); } + int getTreeHeight() + { + if ( _ebbTreeHeight[ _elementType ] < 0 ) + _ebbTreeHeight[ _elementType ] = _ebbTree[ _elementType ]->getHeight(); + return _ebbTreeHeight[ _elementType ]; + } struct TInters //!< data of intersection of the line and the mesh face (used in GetPointState()) { @@ -521,9 +539,9 @@ double SMESH_ElementSearcherImpl::getTolerance() double boxSize = _nodeSearcher->getTree()->maxSize(); _tolerance = 1e-8 * boxSize/* / meshInfo.NbNodes()*/; } - else if ( _ebbTree && meshInfo.NbElements() > 0 ) + else if ( _ebbTree[_elementType] && meshInfo.NbElements() > 0 ) { - double boxSize = _ebbTree->maxSize(); + double boxSize = _ebbTree[_elementType]->maxSize(); _tolerance = 1e-8 * boxSize/* / meshInfo.NbElements()*/; } if ( _tolerance == 0 ) @@ -544,10 +562,9 @@ double SMESH_ElementSearcherImpl::getTolerance() } else { - SMDS_ElemIteratorPtr elemIt = - _mesh->elementsIterator( SMDSAbs_ElementType( complexType )); + SMDS_ElemIteratorPtr elemIt = _mesh->elementsIterator( SMDSAbs_ElementType( complexType )); const SMDS_MeshElement* elem = elemIt->next(); - SMDS_ElemIteratorPtr nodeIt = elem->nodesIterator(); + SMDS_ElemIteratorPtr nodeIt = elem->nodesIterator(); SMESH_TNodeXYZ n1( nodeIt->next() ); elemSize = 0; while ( nodeIt->more() ) @@ -672,7 +689,7 @@ void SMESH_ElementSearcherImpl::findOuterBoundary(const SMDS_MeshElement* outerF outerFace2 = angle2Face.begin()->second; } } - // store the found outer face and add its links to continue seaching from + // store the found outer face and add its links to continue searching from if ( outerFace2 ) { _outerFaces.insert( outerFace2 ); @@ -716,6 +733,7 @@ FindElementsByPoint(const gp_Pnt& point, vector< const SMDS_MeshElement* >& foundElements) { foundElements.clear(); + _elementType = type; double tolerance = getTolerance(); @@ -749,13 +767,12 @@ FindElementsByPoint(const gp_Pnt& point, // ================================================================================= else // elements more complex than 0D { - if ( !_ebbTree || _elementType != type ) + if ( !_ebbTree[type] ) { - if ( _ebbTree ) delete _ebbTree; - _ebbTree = new ElementBndBoxTree( *_mesh, _elementType = type, _meshPartIt, tolerance ); + _ebbTree[_elementType] = new ElementBndBoxTree( *_mesh, type, _meshPartIt, tolerance ); } TIDSortedElemSet suspectElems; - _ebbTree->getElementsNearPoint( point, suspectElems ); + _ebbTree[ type ]->getElementsNearPoint( point, suspectElems ); TIDSortedElemSet::iterator elem = suspectElems.begin(); for ( ; elem != suspectElems.end(); ++elem ) if ( !SMESH_MeshAlgos::IsOut( *elem, point, tolerance )) @@ -777,29 +794,29 @@ SMESH_ElementSearcherImpl::FindClosestTo( const gp_Pnt& point, SMDSAbs_ElementType type ) { const SMDS_MeshElement* closestElem = 0; + _elementType = type; if ( type == SMDSAbs_Face || type == SMDSAbs_Volume ) { - if ( !_ebbTree || _elementType != type ) - { - if ( _ebbTree ) delete _ebbTree; - _ebbTree = new ElementBndBoxTree( *_mesh, _elementType = type, _meshPartIt ); - } + ElementBndBoxTree*& ebbTree = _ebbTree[ type ]; + if ( !ebbTree ) + ebbTree = new ElementBndBoxTree( *_mesh, type, _meshPartIt ); + TIDSortedElemSet suspectElems; - _ebbTree->getElementsNearPoint( point, suspectElems ); + ebbTree->getElementsNearPoint( point, suspectElems ); - if ( suspectElems.empty() && _ebbTree->maxSize() > 0 ) + if ( suspectElems.empty() && ebbTree->maxSize() > 0 ) { - gp_Pnt boxCenter = 0.5 * ( _ebbTree->getBox()->CornerMin() + - _ebbTree->getBox()->CornerMax() ); + gp_Pnt boxCenter = 0.5 * ( ebbTree->getBox()->CornerMin() + + ebbTree->getBox()->CornerMax() ); double radius = -1; - if ( _ebbTree->getBox()->IsOut( point.XYZ() )) - radius = point.Distance( boxCenter ) - 0.5 * _ebbTree->maxSize(); + if ( ebbTree->getBox()->IsOut( point.XYZ() )) + radius = point.Distance( boxCenter ) - 0.5 * ebbTree->maxSize(); if ( radius < 0 ) - radius = _ebbTree->maxSize() / pow( 2., _ebbTree->getHeight()) / 2; + radius = ebbTree->maxSize() / pow( 2., getTreeHeight()) / 2; while ( suspectElems.empty() ) { - _ebbTree->getElementsInSphere( point.XYZ(), radius, suspectElems ); + ebbTree->getElementsInSphere( point.XYZ(), radius, suspectElems ); radius *= 1.1; } } @@ -862,12 +879,14 @@ SMESH_ElementSearcherImpl::FindClosestTo( const gp_Pnt& point, TopAbs_State SMESH_ElementSearcherImpl::GetPointState(const gp_Pnt& point) { + _elementType = SMDSAbs_Face; + double tolerance = getTolerance(); - if ( !_ebbTree || _elementType != SMDSAbs_Face ) - { - if ( _ebbTree ) delete _ebbTree; - _ebbTree = new ElementBndBoxTree( *_mesh, _elementType = SMDSAbs_Face, _meshPartIt ); - } + + ElementBndBoxTree*& ebbTree = _ebbTree[ SMDSAbs_Face ]; + if ( !ebbTree ) + ebbTree = new ElementBndBoxTree( *_mesh, _elementType, _meshPartIt ); + // Algo: analyse transition of a line starting at the point through mesh boundary; // try three lines parallel to axis of the coordinate system and perform rough // analysis. If solution is not clear perform thorough analysis. @@ -883,7 +902,7 @@ TopAbs_State SMESH_ElementSearcherImpl::GetPointState(const gp_Pnt& point) gp_Lin line ( lineAxis ); TIDSortedElemSet suspectFaces; // faces possibly intersecting the line - _ebbTree->getElementsNearLine( lineAxis, suspectFaces ); + ebbTree->getElementsNearLine( lineAxis, suspectFaces ); // Intersect faces with the line @@ -987,16 +1006,19 @@ TopAbs_State SMESH_ElementSearcherImpl::GetPointState(const gp_Pnt& point) // skip tangent intersections int nbTgt = 0; - const SMDS_MeshElement* prevFace = u_int1->second._face; - while ( ok && u_int2->second._coincides ) + if ( u_int2 != u2inters.end() ) { - if ( SMESH_MeshAlgos::GetCommonNodes(prevFace , u_int2->second._face).empty() ) - ok = false; - else + const SMDS_MeshElement* prevFace = u_int1->second._face; + while ( ok && u_int2->second._coincides ) { - nbTgt++; - u_int2++; - ok = ( u_int2 != u2inters.end() ); + if ( SMESH_MeshAlgos::GetCommonNodes(prevFace , u_int2->second._face).empty() ) + ok = false; + else + { + nbTgt++; + u_int2++; + ok = ( u_int2 != u2inters.end() ); + } } } if ( !ok ) break; @@ -1088,13 +1110,13 @@ void SMESH_ElementSearcherImpl::GetElementsNearLine( const gp_Ax1& SMDSAbs_ElementType type, vector< const SMDS_MeshElement* >& foundElems) { - if ( !_ebbTree || _elementType != type ) - { - if ( _ebbTree ) delete _ebbTree; - _ebbTree = new ElementBndBoxTree( *_mesh, _elementType = type, _meshPartIt ); - } + _elementType = type; + ElementBndBoxTree*& ebbTree = _ebbTree[ type ]; + if ( !ebbTree ) + ebbTree = new ElementBndBoxTree( *_mesh, _elementType, _meshPartIt ); + TIDSortedElemSet suspectFaces; // elements possibly intersecting the line - _ebbTree->getElementsNearLine( line, suspectFaces ); + ebbTree->getElementsNearLine( line, suspectFaces ); foundElems.assign( suspectFaces.begin(), suspectFaces.end()); } @@ -1109,13 +1131,13 @@ void SMESH_ElementSearcherImpl::GetElementsInSphere( const gp_XYZ& SMDSAbs_ElementType type, vector< const SMDS_MeshElement* >& foundElems) { - if ( !_ebbTree || _elementType != type ) - { - if ( _ebbTree ) delete _ebbTree; - _ebbTree = new ElementBndBoxTree( *_mesh, _elementType = type, _meshPartIt ); - } + _elementType = type; + ElementBndBoxTree*& ebbTree = _ebbTree[ type ]; + if ( !ebbTree ) + ebbTree = new ElementBndBoxTree( *_mesh, _elementType, _meshPartIt ); + TIDSortedElemSet suspectFaces; // elements possibly intersecting the line - _ebbTree->getElementsInSphere( center, radius, suspectFaces ); + ebbTree->getElementsInSphere( center, radius, suspectFaces ); foundElems.assign( suspectFaces.begin(), suspectFaces.end() ); } @@ -1153,8 +1175,8 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin gp_Vec edge2( xyz[i+1], xyz[(i+2)%nbNodes] ); faceNorm += edge1 ^ edge2; } - double normSize = faceNorm.Magnitude(); - if ( normSize <= tol ) + double fNormSize = faceNorm.Magnitude(); + if ( fNormSize <= tol ) { // degenerated face: point is out if it is out of all face edges for ( i = 0; i < nbNodes; ++i ) @@ -1165,7 +1187,7 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin } return true; } - faceNorm /= normSize; + faceNorm /= fNormSize; // check if the point lays on face plane gp_Vec n2p( xyz[0], point ); @@ -1194,9 +1216,10 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin // to find intersections of the ray with the boundary. gp_Vec ray = n2p; gp_Vec plnNorm = ray ^ faceNorm; - normSize = plnNorm.Magnitude(); - if ( normSize <= tol ) return false; // point coincides with the first node - plnNorm /= normSize; + double n2pSize = plnNorm.Magnitude(); + if ( n2pSize <= tol ) return false; // point coincides with the first node + if ( n2pSize * n2pSize > fNormSize * 100 ) return true; // point is very far + plnNorm /= n2pSize; // for each node of the face, compute its signed distance to the cutting plane vector dist( nbNodes + 1); for ( i = 0; i < nbNodes; ++i ) @@ -1242,7 +1265,7 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin if ( rClosest > 0. && rClosest < 1. ) // not node intersection return out; - // ray pass through a face node; analyze transition through an adjacent edge + // the ray passes through a face node; analyze transition through an adjacent edge gp_Pnt p1 = xyz[ (rClosest == 0.) ? ((iClosest+nbNodes-1) % nbNodes) : (iClosest+1) ]; gp_Pnt p2 = xyz[ (rClosest == 0.) ? iClosest : ((iClosest+2) % nbNodes) ]; gp_Vec edgeAdjacent( p1, p2 ); @@ -1252,6 +1275,7 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin bool covexCorner = ( edgeNorm * edgeAdjacent * (rClosest==1. ? 1. : -1.)) < 0; return covexCorner ? (out || out2) : (out && out2); } + if ( element->GetType() == SMDSAbs_Edge ) // -------------------------------------------------- { // point is out of edge if it is NOT ON any straight part of edge @@ -1279,6 +1303,7 @@ bool SMESH_MeshAlgos::IsOut( const SMDS_MeshElement* element, const gp_Pnt& poin } return true; } + // Node or 0D element ------------------------------------------------------------------------- { gp_Vec n2p ( xyz[0], point );