-// Copyright (C) 2007-2019 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2024 CEA, EDF, OPEN CASCADE
//
// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
void getElementsInBox ( const Bnd_B3d& box, TElemSeq& foundElems );
void getElementsInSphere ( const gp_XYZ& center, const double radius, TElemSeq& foundElems );
ElementBndBoxTree* getLeafAtPoint( const gp_XYZ& point );
+ int getNbElements();
protected:
ElementBndBoxTree() {}
void init(const SMDS_MeshElement* elem, double tolerance);
};
std::vector< ElementBox* > _elements;
+ //std::string _name;
typedef ObjectPool< ElementBox > TElementBoxPool;
double tolerance)
:SMESH_Octree( new LimitAndPool() )
{
- int nbElems = mesh.GetMeshInfo().NbElements( elemType );
+ smIdType nbElems = mesh.GetMeshInfo().NbElements( elemType );
_elements.reserve( nbElems );
TElementBoxPool& elBoPool = getLimitAndPool()->_elBoPool;
-#ifdef _DEBUG_
- if ( theElemIt && !theElemIt->more() )
+ if (SALOME::VerbosityActivated() && theElemIt && !theElemIt->more() )
std::cout << "WARNING: ElementBndBoxTree constructed on empty iterator!" << std::endl;
-#endif
SMDS_ElemIteratorPtr elemIt = theElemIt ? theElemIt : mesh.elementsIterator( elemType );
while ( elemIt->more() )
//================================================================================
/*!
- * \brief Redistrubute element boxes among children
+ * \brief Redistribute element boxes among children
*/
//================================================================================
if ( child->isLeaf() && child->_elements.capacity() > child->_elements.size() )
SMESHUtils::CompactVector( child->_elements );
+
+ // child->_name = _name + char('0' + j);
+ // if ( child->isLeaf() && child->_elements.size() )
+ // {
+ // cout << child->_name << " ";
+ // for ( size_t i = 0; i < child->_elements.size(); ++i )
+ // cout << child->_elements[i]->_element->GetID() << " ";
+ // cout << endl;
+ // }
}
}
return 0;
}
+ //================================================================================
+ /*!
+ * \brief Return number of elements
+ */
+ //================================================================================
+
+ int ElementBndBoxTree::getNbElements()
+ {
+ int nb = 0;
+ if ( isLeaf() )
+ {
+ nb = _elements.size();
+ }
+ else
+ {
+ for (int i = 0; i < 8; i++)
+ nb += ((ElementBndBoxTree*) myChildren[i])->getNbElements();
+ }
+ return nb;
+ }
+
//================================================================================
/*!
* \brief Construct the element box
{
delete _ebbTree[i]; _ebbTree[i] = NULL;
}
- if ( _nodeSearcher ) delete _nodeSearcher; _nodeSearcher = 0;
+ if ( _nodeSearcher ) delete _nodeSearcher;
+ _nodeSearcher = 0;
}
virtual int FindElementsByPoint(const gp_Pnt& point,
SMDSAbs_ElementType type,
GC_MakeSegment edge( SMESH_TNodeXYZ( face->GetNode( i )),
SMESH_TNodeXYZ( face->GetNode( (i+1)%nbNodes) ));
anExtCC.Init( lineCurve, edge.Value() );
+ if ( !anExtCC.Extrema().IsDone() ||
+ anExtCC.Extrema().IsParallel() )
+ continue;
if ( anExtCC.NbExtrema() > 0 && anExtCC.LowerDistance() <= tol)
{
Standard_Real pl, pe;
/*!
* \brief Find an element of given type most close to the given point
*
- * WARNING: Only face search is implemeneted so far
+ * WARNING: Only edge, face and volume search is implemented so far
*/
//=======================================================================
const SMDS_MeshElement* prevFace = u_int1->second._face;
while ( ok && u_int2->second._coincides )
{
- if ( SMESH_MeshAlgos::GetCommonNodes(prevFace , u_int2->second._face).empty() )
+ if ( SMESH_MeshAlgos::NbCommonNodes(prevFace , u_int2->second._face) == 0 )
ok = false;
else
{
minDist = d;
}
}
+ if ( minDist > radius )
+ {
+ ElementBndBoxTree::TElemSeq elems2;
+ ebbTree->getElementsInSphere( p, minDist, elems2 );
+ for ( e = elems2.begin(); e != elems2.end(); ++e )
+ {
+ if ( elems.count( *e ))
+ continue;
+ double d = SMESH_MeshAlgos::GetDistance( *e, point, &proj );
+ if ( d < minDist )
+ {
+ bestProj = proj;
+ elem = *e;
+ minDist = d;
+ }
+ }
+ }
if ( closestElem ) *closestElem = elem;
return bestProj;
try {
tgtCS = gp_Ax3( xyz[0], OZ, OX );
}
- catch ( Standard_Failure ) {
+ catch ( Standard_Failure& ) {
return badDistance;
}
trsf.SetTransformation( tgtCS );
typedef std::pair< bool, const SMDS_MeshNode* > TIsSharpAndMedium;
typedef NCollection_DataMap< SMESH_TLink, TIsSharpAndMedium, SMESH_TLink > TLinkSharpMap;
- TLinkSharpMap linkIsSharp( theMesh->NbFaces() );
+ TLinkSharpMap linkIsSharp;
+ Standard_Integer nbBuckets = FromSmIdType<Standard_Integer>( theMesh->NbFaces() );
+ if ( nbBuckets > 0 )
+ linkIsSharp.ReSize( nbBuckets );
+
TIsSharpAndMedium sharpMedium( true, 0 );
bool & isSharp = sharpMedium.first;
const SMDS_MeshNode* & nMedium = sharpMedium.second;
typedef std::vector< const SMDS_MeshElement* > TFaceVec;
typedef NCollection_DataMap< SMESH_TLink, TFaceVec, SMESH_TLink > TFacesByLinks;
- TFacesByLinks facesByLink( theMesh->NbFaces() );
+ TFacesByLinks facesByLink;
+ Standard_Integer nbBuckets = FromSmIdType<Standard_Integer>( theMesh->NbFaces() );
+ if ( nbBuckets > 0 )
+ facesByLink.ReSize( nbBuckets );
std::vector< const SMDS_MeshNode* > faceNodes;
for ( SMDS_FaceIteratorPtr faceIt = theMesh->facesIterator(); faceIt->more(); )
return ok;
}
-//=======================================================================
-//function : GetCommonNodes
-//purpose : Return nodes common to two elements
-//=======================================================================
+//================================================================================
+/*!
+ * \brief Return nodes common to two elements
+ */
+//================================================================================
+
+int SMESH_MeshAlgos::NbCommonNodes(const SMDS_MeshElement* e1,
+ const SMDS_MeshElement* e2)
+{
+ int nb = 0;
+ for ( int i = 0 ; i < e1->NbNodes(); ++i )
+ nb += ( e2->GetNodeIndex( e1->GetNode( i )) >= 0 );
+ return nb;
+}
+
+//================================================================================
+/*!
+ * \brief Return nodes common to two elements
+ */
+//================================================================================
std::vector< const SMDS_MeshNode*> SMESH_MeshAlgos::GetCommonNodes(const SMDS_MeshElement* e1,
const SMDS_MeshElement* e2)
common.push_back( e1->GetNode( i ));
return common;
}
+
+//================================================================================
+/*!
+ * \brief Return true if a node is on a boundary of 2D mesh.
+ * Optionally returns two neighboring boundary nodes (or more in non-manifold mesh)
+ */
+//================================================================================
+
+bool SMESH_MeshAlgos::IsOn2DBoundary( const SMDS_MeshNode* theNode,
+ std::vector< const SMDS_MeshNode*> * theNeibors )
+{
+ typedef NCollection_DataMap< SMESH_TLink, int, SMESH_TLink > TLinkCountMap;
+ TLinkCountMap linkCountMap( 10 );
+
+ int nbFreeLinks = 0;
+ for ( SMDS_ElemIteratorPtr fIt = theNode->GetInverseElementIterator(SMDSAbs_Face); fIt->more(); )
+ {
+ const SMDS_MeshElement* face = fIt->next();
+ const int nbCorners = face->NbCornerNodes();
+
+ int iN = face->GetNodeIndex( theNode );
+ int iPrev = ( iN - 1 + nbCorners ) % nbCorners;
+ int iNext = ( iN + 1 ) % nbCorners;
+
+ for ( int i : { iPrev, iNext } )
+ {
+ SMESH_TLink link( theNode, face->GetNode( i ));
+ int* count = linkCountMap.ChangeSeek( link );
+ if ( count ) ++( *count );
+ else linkCountMap.Bind( link, 1 );
+
+ if ( !count ) ++nbFreeLinks;
+ else --nbFreeLinks;
+ }
+ }
+
+ if ( theNeibors )
+ {
+ theNeibors->clear();
+ theNeibors->reserve( nbFreeLinks );
+ for ( TLinkCountMap::Iterator linkIt( linkCountMap ); linkIt.More(); linkIt.Next() )
+ if ( linkIt.Value() == 1 )
+ {
+ theNeibors->push_back( linkIt.Key().node1() );
+ if ( theNeibors->back() == theNode )
+ theNeibors->back() = linkIt.Key().node2();
+ }
+ }
+ return nbFreeLinks > 0;
+}
+
//================================================================================
/*!
* \brief Return true if node1 encounters first in the face and node2, after
{
return new SMESH_ElementSearcherImpl( mesh, tolerance, elemIt );
}
+
+
+//================================================================================
+/*!
+ * \brief Intersect a ray with a convex volume
+ * \param [in] ray - the ray
+ * \param [in] rayLen - ray length
+ * \param [in] vol - the volume
+ * \param [out] tMin - return a ray parameter where the ray enters the volume
+ * \param [out] tMax - return a ray parameter where the ray exit the volume
+ * \param [out] iFacetMin - facet index where the ray enters the volume
+ * \param [out] iFacetMax - facet index where the ray exit the volume
+ * \return bool - true if the ray intersects the volume
+ */
+//================================================================================
+
+bool SMESH_MeshAlgos::IntersectRayVolume( const gp_Ax1& ray,
+ const double rayLen,
+ const SMDS_MeshElement* vol,
+ double & tMin,
+ double & tMax,
+ int & iFacetMin,
+ int & iFacetMax)
+{
+ /* Ray-Convex Polyhedron Intersection Test by Eric Haines, erich@eye.com
+ *
+ * This test checks the ray against each face of a polyhedron, checking whether
+ * the set of intersection points found for each ray-plane intersection
+ * overlaps the previous intersection results. If there is no overlap (i.e.
+ * no line segment along the ray that is inside the polyhedron), then the
+ * ray misses and returns false; else true.
+ */
+ SMDS_VolumeTool vTool;
+ if ( !vTool.Set( vol ))
+ return false;
+
+ tMin = -Precision::Infinite() ;
+ tMax = rayLen ;
+
+ /* Test each plane in polyhedron */
+ for ( int iF = 0; iF < vTool.NbFaces(); ++iF )
+ {
+ const SMDS_MeshNode** fNodes = vTool.GetFaceNodes( iF );
+ gp_XYZ normal;
+ vTool.GetFaceNormal( iF,
+ normal.ChangeCoord(1),
+ normal.ChangeCoord(2),
+ normal.ChangeCoord(3) );
+ double D = - ( normal * SMESH_NodeXYZ( fNodes[0] ));
+
+ /* Compute intersection point T and sidedness */
+ double vd = ray.Direction().XYZ() * normal;
+ double vn = ray.Location().XYZ() * normal + D;
+ if ( vd == 0.0 ) {
+ /* ray is parallel to plane - check if ray origin is inside plane's
+ half-space */
+ if ( vn > 0.0 )
+ /* ray origin is outside half-space */
+ return false;
+ }
+ else
+ {
+ /* ray not parallel - get distance to plane */
+ double t = -vn / vd ;
+ if ( vd < 0.0 )
+ {
+ /* front face - T is a near point */
+ if ( t > tMax ) return false;
+ if ( t > tMin ) {
+ /* hit near face */
+ tMin = t ;
+ iFacetMin = iF;
+ }
+ }
+ else
+ {
+ /* back face - T is a far point */
+ if ( t < tMin ) return false;
+ if ( t < tMax ) {
+ /* hit far face */
+ tMax = t ;
+ iFacetMax = iF;
+ }
+ }
+ }
+ }
+
+ /* survived all tests */
+ /* Note: if ray originates on polyhedron, may want to change 0.0 to some
+ * epsilon to avoid intersecting the originating face.
+ */
+ if ( tMin >= 0.0 ) {
+ /* outside, hitting front face */
+ return true;
+ }
+ else
+ {
+ if ( tMax < rayLen ) {
+ /* inside, hitting back face */
+ return true;
+ }
+ else
+ {
+ /* inside, but back face beyond tmax */
+ return false;
+ }
+ }
+}
