-// Copyright (C) 2007-2015 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2020 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
#include "SMESH_ControlsDef.hxx"
#include "SMDS_BallElement.hxx"
+#include "SMDS_FacePosition.hxx"
#include "SMDS_Iterator.hxx"
#include "SMDS_Mesh.hxx"
#include "SMDS_MeshElement.hxx"
#include "SMDS_MeshNode.hxx"
-#include "SMDS_QuadraticEdge.hxx"
-#include "SMDS_QuadraticFaceOfNodes.hxx"
#include "SMDS_VolumeTool.hxx"
#include "SMESHDS_GroupBase.hxx"
#include "SMESHDS_GroupOnFilter.hxx"
#include "SMESH_MeshAlgos.hxx"
#include "SMESH_OctreeNode.hxx"
+#include <GEOMUtils.hxx>
#include <Basics_Utils.hxx>
#include <BRepAdaptor_Surface.hxx>
+#include <BRepBndLib.hxx>
+#include <BRepBuilderAPI_Copy.hxx>
+#include <BRepClass3d_SolidClassifier.hxx>
#include <BRepClass_FaceClassifier.hxx>
#include <BRep_Tool.hxx>
+#include <GeomLib_IsPlanarSurface.hxx>
#include <Geom_CylindricalSurface.hxx>
#include <Geom_Plane.hxx>
#include <Geom_Surface.hxx>
#include <NCollection_Map.hxx>
#include <Precision.hxx>
+#include <ShapeAnalysis_Surface.hxx>
#include <TColStd_MapIteratorOfMapOfInteger.hxx>
#include <TColStd_MapOfInteger.hxx>
#include <TColStd_SequenceOfAsciiString.hxx>
v2.Magnitude() < gp::Resolution() ? 0 : v1.Angle( v2 );
}
+ inline double getCos2( const gp_XYZ& P1, const gp_XYZ& P2, const gp_XYZ& P3 )
+ {
+ gp_Vec v1( P1 - P2 ), v2( P3 - P2 );
+ double dot = v1 * v2, len1 = v1.SquareMagnitude(), len2 = v2.SquareMagnitude();
+
+ return ( dot < 0 || len1 < gp::Resolution() || len2 < gp::Resolution() ? -1 :
+ dot * dot / len1 / len2 );
+ }
+
inline double getArea( const gp_XYZ& P1, const gp_XYZ& P2, const gp_XYZ& P3 )
{
gp_Vec aVec1( P2 - P1 );
// Case 1 Case 2
// | | | | |
// | | | | |
- // +-----+------+ +-----+------+
+ // +-----+------+ +-----+------+
// | | | |
// | | | |
- // result sould be 2 in both cases
+ // result should be 2 in both cases
//
int aResult0 = 0, aResult1 = 0;
- // last node, it is a medium one in a quadratic edge
+ // last node, it is a medium one in a quadratic edge
const SMDS_MeshNode* aLastNode = anEdge->GetNode( anEdge->NbNodes() - 1 );
const SMDS_MeshNode* aNode0 = anEdge->GetNode( 0 );
const SMDS_MeshNode* aNode1 = anEdge->GetNode( 1 );
return false;
const SMDS_MeshElement* anElem = myMesh->FindElement( theId );
- if ( !anElem || anElem->GetType() != this->GetType() )
+ if ( !IsApplicable( anElem ))
return false;
return GetPoints( anElem, theRes );
theRes.setElement( anElem );
// Get nodes of the element
- SMDS_ElemIteratorPtr anIter;
-
- if ( anElem->IsQuadratic() ) {
- switch ( anElem->GetType() ) {
- case SMDSAbs_Edge:
- anIter = dynamic_cast<const SMDS_VtkEdge*>
- (anElem)->interlacedNodesElemIterator();
- break;
- case SMDSAbs_Face:
- anIter = dynamic_cast<const SMDS_VtkFace*>
- (anElem)->interlacedNodesElemIterator();
- break;
- default:
- anIter = anElem->nodesIterator();
- }
- }
- else {
- anIter = anElem->nodesIterator();
- }
-
+ SMDS_NodeIteratorPtr anIter= anElem->interlacedNodesIterator();
if ( anIter ) {
- double xyz[3];
+ SMESH_NodeXYZ p;
while( anIter->more() ) {
- if ( const SMDS_MeshNode* aNode = static_cast<const SMDS_MeshNode*>( anIter->next() ))
- {
- aNode->GetXYZ( xyz );
- theRes.push_back( gp_XYZ( xyz[0], xyz[1], xyz[2] ));
- }
+ if ( p.Set( anIter->next() ))
+ theRes.push_back( p );
}
}
return ( myPrecision >= 0 ) ? floor( aVal * myPrecisionValue + 0.5 ) / myPrecisionValue : aVal;
}
+//================================================================================
+/*!
+ * \brief Return true if a value can be computed for a given element.
+ * Some NumericalFunctor's are meaningful for elements of a certain
+ * geometry only.
+ */
+//================================================================================
+
+bool NumericalFunctor::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return element && element->GetType() == this->GetType();
+}
+
+bool NumericalFunctor::IsApplicable( long theElementId ) const
+{
+ return IsApplicable( myMesh->FindElement( theElementId ));
+}
+
//================================================================================
/*!
* \brief Return histogram of functor values
// }
// }
// { // polygons
-
+
// }
if( myPrecision >= 0 )
aVal = Max(aVal,Max(L7,L8));
break;
}
- case SMDSEntity_Quad_Penta: { // quadratic pentas
+ case SMDSEntity_Quad_Penta:
+ case SMDSEntity_BiQuad_Penta: { // quadratic pentas
double L1 = getDistance(P( 1 ),P( 7 )) + getDistance(P( 7 ),P( 2 ));
double L2 = getDistance(P( 2 ),P( 8 )) + getDistance(P( 8 ),P( 3 ));
double L3 = getDistance(P( 3 ),P( 9 )) + getDistance(P( 9 ),P( 1 ));
double MinimumAngle::GetValue( const TSequenceOfXYZ& P )
{
- double aMin;
-
- if (P.size() <3)
+ if ( P.size() < 3 )
return 0.;
- aMin = getAngle(P( P.size() ), P( 1 ), P( 2 ));
- aMin = Min(aMin,getAngle(P( P.size()-1 ), P( P.size() ), P( 1 )));
+ double aMaxCos2;
+
+ aMaxCos2 = getCos2( P( P.size() ), P( 1 ), P( 2 ));
+ aMaxCos2 = Max( aMaxCos2, getCos2( P( P.size()-1 ), P( P.size() ), P( 1 )));
for ( size_t i = 2; i < P.size(); i++ )
{
- double A0 = getAngle( P( i-1 ), P( i ), P( i+1 ) );
- aMin = Min(aMin,A0);
+ double A0 = getCos2( P( i-1 ), P( i ), P( i+1 ) );
+ aMaxCos2 = Max( aMaxCos2, A0 );
}
+ if ( aMaxCos2 < 0 )
+ return 0; // all nodes coincide
- return aMin * 180.0 / M_PI;
+ double cos = sqrt( aMaxCos2 );
+ if ( cos >= 1 ) return 0;
+ return acos( cos ) * 180.0 / M_PI;
}
double MinimumAngle::GetBadRate( double Value, int nbNodes ) const
{
double aVal = 0;
myCurrElement = myMesh->FindElement( theId );
- if ( myCurrElement && myCurrElement->GetVtkType() == VTK_QUAD )
- {
- // issue 21723
- vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myCurrElement->getMeshId()]->getGrid();
- if ( vtkCell* avtkCell = grid->GetCell( myCurrElement->getVtkId() ))
- aVal = Round( vtkMeshQuality::QuadAspectRatio( avtkCell ));
- }
- else
- {
- TSequenceOfXYZ P;
- if ( GetPoints( myCurrElement, P ))
- aVal = Round( GetValue( P ));
- }
+ TSequenceOfXYZ P;
+ if ( GetPoints( myCurrElement, P ))
+ aVal = Round( GetValue( P ));
return aVal;
}
if ( nbNodes == 3 ) {
// Compute lengths of the sides
- std::vector< double > aLen (nbNodes);
- for ( int i = 0; i < nbNodes - 1; i++ )
- aLen[ i ] = getDistance( P( i + 1 ), P( i + 2 ) );
- aLen[ nbNodes - 1 ] = getDistance( P( 1 ), P( nbNodes ) );
+ double aLen1 = getDistance( P( 1 ), P( 2 ));
+ double aLen2 = getDistance( P( 2 ), P( 3 ));
+ double aLen3 = getDistance( P( 3 ), P( 1 ));
// Q = alfa * h * p / S, where
//
// alfa = sqrt( 3 ) / 6
// h - length of the longest edge
// p - half perimeter
// S - triangle surface
- const double alfa = sqrt( 3. ) / 6.;
- double maxLen = Max( aLen[ 0 ], Max( aLen[ 1 ], aLen[ 2 ] ) );
- double half_perimeter = ( aLen[0] + aLen[1] + aLen[2] ) / 2.;
- double anArea = getArea( P( 1 ), P( 2 ), P( 3 ) );
+ const double alfa = sqrt( 3. ) / 6.;
+ double maxLen = Max( aLen1, Max( aLen2, aLen3 ));
+ double half_perimeter = ( aLen1 + aLen2 + aLen3 ) / 2.;
+ double anArea = getArea( P( 1 ), P( 2 ), P( 3 ));
if ( anArea <= theEps )
return theInf;
return alfa * maxLen * half_perimeter / anArea;
}
else if ( nbNodes == 6 ) { // quadratic triangles
// Compute lengths of the sides
- std::vector< double > aLen (3);
- aLen[0] = getDistance( P(1), P(3) );
- aLen[1] = getDistance( P(3), P(5) );
- aLen[2] = getDistance( P(5), P(1) );
- // Q = alfa * h * p / S, where
- //
- // alfa = sqrt( 3 ) / 6
- // h - length of the longest edge
- // p - half perimeter
- // S - triangle surface
- const double alfa = sqrt( 3. ) / 6.;
- double maxLen = Max( aLen[ 0 ], Max( aLen[ 1 ], aLen[ 2 ] ) );
- double half_perimeter = ( aLen[0] + aLen[1] + aLen[2] ) / 2.;
- double anArea = getArea( P(1), P(3), P(5) );
+ double aLen1 = getDistance( P( 1 ), P( 3 ));
+ double aLen2 = getDistance( P( 3 ), P( 5 ));
+ double aLen3 = getDistance( P( 5 ), P( 1 ));
+ // algo same as for the linear triangle
+ const double alfa = sqrt( 3. ) / 6.;
+ double maxLen = Max( aLen1, Max( aLen2, aLen3 ));
+ double half_perimeter = ( aLen1 + aLen2 + aLen3 ) / 2.;
+ double anArea = getArea( P( 1 ), P( 3 ), P( 5 ));
if ( anArea <= theEps )
return theInf;
return alfa * maxLen * half_perimeter / anArea;
}
else if( nbNodes == 4 ) { // quadrangle
// Compute lengths of the sides
- std::vector< double > aLen (4);
+ double aLen[4];
aLen[0] = getDistance( P(1), P(2) );
aLen[1] = getDistance( P(2), P(3) );
aLen[2] = getDistance( P(3), P(4) );
aLen[3] = getDistance( P(4), P(1) );
// Compute lengths of the diagonals
- std::vector< double > aDia (2);
+ double aDia[2];
aDia[0] = getDistance( P(1), P(3) );
aDia[1] = getDistance( P(2), P(4) );
// Compute areas of all triangles which can be built
// taking three nodes of the quadrangle
- std::vector< double > anArea (4);
+ double anArea[4];
anArea[0] = getArea( P(1), P(2), P(3) );
anArea[1] = getArea( P(1), P(2), P(4) );
anArea[2] = getArea( P(1), P(3), P(4) );
//
// alpha = sqrt( 1/32 )
// L = max( L1, L2, L3, L4, D1, D2 )
- // C1 = sqrt( ( L1^2 + L1^2 + L1^2 + L1^2 ) / 4 )
+ // C1 = sqrt( L1^2 + L1^2 + L1^2 + L1^2 )
// C2 = min( S1, S2, S3, S4 )
// Li - lengths of the edges
// Di - lengths of the diagonals
// Si - areas of the triangles
const double alpha = sqrt( 1 / 32. );
double L = Max( aLen[ 0 ],
- Max( aLen[ 1 ],
- Max( aLen[ 2 ],
- Max( aLen[ 3 ],
- Max( aDia[ 0 ], aDia[ 1 ] ) ) ) ) );
- double C1 = sqrt( ( aLen[0] * aLen[0] +
- aLen[1] * aLen[1] +
- aLen[2] * aLen[2] +
- aLen[3] * aLen[3] ) / 4. );
+ Max( aLen[ 1 ],
+ Max( aLen[ 2 ],
+ Max( aLen[ 3 ],
+ Max( aDia[ 0 ], aDia[ 1 ] ) ) ) ) );
+ double C1 = sqrt( aLen[0] * aLen[0] +
+ aLen[1] * aLen[1] +
+ aLen[2] * aLen[2] +
+ aLen[3] * aLen[3] );
double C2 = Min( anArea[ 0 ],
- Min( anArea[ 1 ],
- Min( anArea[ 2 ], anArea[ 3 ] ) ) );
+ Min( anArea[ 1 ],
+ Min( anArea[ 2 ], anArea[ 3 ] ) ) );
if ( C2 <= theEps )
return theInf;
return alpha * L * C1 / C2;
}
else if( nbNodes == 8 || nbNodes == 9 ) { // nbNodes==8 - quadratic quadrangle
// Compute lengths of the sides
- std::vector< double > aLen (4);
+ double aLen[4];
aLen[0] = getDistance( P(1), P(3) );
aLen[1] = getDistance( P(3), P(5) );
aLen[2] = getDistance( P(5), P(7) );
aLen[3] = getDistance( P(7), P(1) );
// Compute lengths of the diagonals
- std::vector< double > aDia (2);
+ double aDia[2];
aDia[0] = getDistance( P(1), P(5) );
aDia[1] = getDistance( P(3), P(7) );
// Compute areas of all triangles which can be built
// taking three nodes of the quadrangle
- std::vector< double > anArea (4);
+ double anArea[4];
anArea[0] = getArea( P(1), P(3), P(5) );
anArea[1] = getArea( P(1), P(3), P(7) );
anArea[2] = getArea( P(1), P(5), P(7) );
//
// alpha = sqrt( 1/32 )
// L = max( L1, L2, L3, L4, D1, D2 )
- // C1 = sqrt( ( L1^2 + L1^2 + L1^2 + L1^2 ) / 4 )
+ // C1 = sqrt( L1^2 + L1^2 + L1^2 + L1^2 )
// C2 = min( S1, S2, S3, S4 )
// Li - lengths of the edges
// Di - lengths of the diagonals
// Si - areas of the triangles
const double alpha = sqrt( 1 / 32. );
double L = Max( aLen[ 0 ],
- Max( aLen[ 1 ],
- Max( aLen[ 2 ],
- Max( aLen[ 3 ],
- Max( aDia[ 0 ], aDia[ 1 ] ) ) ) ) );
- double C1 = sqrt( ( aLen[0] * aLen[0] +
- aLen[1] * aLen[1] +
- aLen[2] * aLen[2] +
- aLen[3] * aLen[3] ) / 4. );
+ Max( aLen[ 1 ],
+ Max( aLen[ 2 ],
+ Max( aLen[ 3 ],
+ Max( aDia[ 0 ], aDia[ 1 ] ) ) ) ) );
+ double C1 = sqrt( aLen[0] * aLen[0] +
+ aLen[1] * aLen[1] +
+ aLen[2] * aLen[2] +
+ aLen[3] * aLen[3] );
double C2 = Min( anArea[ 0 ],
- Min( anArea[ 1 ],
- Min( anArea[ 2 ], anArea[ 3 ] ) ) );
+ Min( anArea[ 1 ],
+ Min( anArea[ 2 ], anArea[ 3 ] ) ) );
if ( C2 <= theEps )
return theInf;
return alpha * L * C1 / C2;
return 0;
}
+bool AspectRatio::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return ( NumericalFunctor::IsApplicable( element ) && !element->IsPoly() );
+}
+
double AspectRatio::GetBadRate( double Value, int /*nbNodes*/ ) const
{
// the aspect ratio is in the range [1.0,infinity]
// Action from CoTech | ACTION 31.3:
// EURIWARE BO: Homogenize the formulas used to calculate the Controls in SMESH to fit with
// those of ParaView. The library used by ParaView for those calculations can be reused in SMESH.
- vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myCurrElement->getMeshId()]->getGrid();
- if ( vtkCell* avtkCell = grid->GetCell( myCurrElement->getVtkId() ))
+ vtkUnstructuredGrid* grid = const_cast<SMDS_Mesh*>( myMesh )->GetGrid();
+ if ( vtkCell* avtkCell = grid->GetCell( myCurrElement->GetVtkID() ))
aVal = Round( vtkMeshQuality::TetAspectRatio( avtkCell ));
}
else
return aVal;
}
+bool AspectRatio3D::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return ( NumericalFunctor::IsApplicable( element ) && !element->IsPoly() );
+}
+
double AspectRatio3D::GetValue( const TSequenceOfXYZ& P )
{
double aQuality = 0.0;
int nbNodes = P.size();
- if(myCurrElement->IsQuadratic()) {
- if(nbNodes==10) nbNodes=4; // quadratic tetrahedron
+ if( myCurrElement->IsQuadratic() ) {
+ if (nbNodes==10) nbNodes=4; // quadratic tetrahedron
else if(nbNodes==13) nbNodes=5; // quadratic pyramid
else if(nbNodes==15) nbNodes=6; // quadratic pentahedron
else if(nbNodes==20) nbNodes=8; // quadratic hexahedron
- else if(nbNodes==27) nbNodes=8; // quadratic hexahedron
+ else if(nbNodes==27) nbNodes=8; // tri-quadratic hexahedron
else return aQuality;
}
case 5:{
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 2 ),P( 3 ),P( 5 )};
- aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4])),aQuality);
+ aQuality = GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4]));
}
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 3 ),P( 4 ),P( 5 )};
case 6:{
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 2 ),P( 4 ),P( 6 )};
- aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4])),aQuality);
+ aQuality = GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4]));
}
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 2 ),P( 4 ),P( 3 )};
case 8:{
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 2 ),P( 5 ),P( 3 )};
- aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4])),aQuality);
+ aQuality = GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[4]));
}
{
gp_XYZ aXYZ[4] = {P( 1 ),P( 2 ),P( 5 ),P( 4 )};
case 12:
{
gp_XYZ aXYZ[8] = {P( 1 ),P( 2 ),P( 4 ),P( 5 ),P( 7 ),P( 8 ),P( 10 ),P( 11 )};
- aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[8])),aQuality);
+ aQuality = GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[8]));
}
{
gp_XYZ aXYZ[8] = {P( 2 ),P( 3 ),P( 5 ),P( 6 ),P( 8 ),P( 9 ),P( 11 ),P( 12 )};
} // switch(nbNodes)
if ( nbNodes > 4 ) {
- // avaluate aspect ratio of quadranle faces
+ // evaluate aspect ratio of quadrangle faces
AspectRatio aspect2D;
SMDS_VolumeTool::VolumeType type = SMDS_VolumeTool::GetType( nbNodes );
int nbFaces = SMDS_VolumeTool::NbFaces( type );
if ( SMDS_VolumeTool::NbFaceNodes( type, i ) != 4 )
continue;
const int* pInd = SMDS_VolumeTool::GetFaceNodesIndices( type, i, true );
- for ( int p = 0; p < 4; ++p ) // loop on nodes of a quadranle face
+ for ( int p = 0; p < 4; ++p ) // loop on nodes of a quadrangle face
points( p + 1 ) = P( pInd[ p ] + 1 );
aQuality = std::max( aQuality, aspect2D.GetValue( points ));
}
*/
//================================================================================
+bool Warping::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return NumericalFunctor::IsApplicable( element ) && element->NbNodes() == 4;
+}
+
double Warping::GetValue( const TSequenceOfXYZ& P )
{
if ( P.size() != 4 )
*/
//================================================================================
+bool Taper::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return ( NumericalFunctor::IsApplicable( element ) && element->NbNodes() == 4 );
+}
+
double Taper::GetValue( const TSequenceOfXYZ& P )
{
if ( P.size() != 4 )
return v1.Magnitude() < gp::Resolution() || v2.Magnitude() < gp::Resolution() ? 0. : v1.Angle( v2 );
}
+bool Skew::IsApplicable( const SMDS_MeshElement* element ) const
+{
+ return ( NumericalFunctor::IsApplicable( element ) && element->NbNodes() <= 4 );
+}
+
double Skew::GetValue( const TSequenceOfXYZ& P )
{
if ( P.size() != 3 && P.size() != 4 )
//================================================================================
/*
- Class : Length2D
- Description : Functor for calculating minimal length of edge
+ Class : Length3D
+ Description : Functor for calculating minimal length of element edge
*/
//================================================================================
-double Length2D::GetValue( long theElementId )
-{
- TSequenceOfXYZ P;
-
- if ( GetPoints( theElementId, P ))
- {
- double aVal = 0;
- int len = P.size();
- SMDSAbs_EntityType aType = P.getElementEntity();
-
- switch (aType) {
- case SMDSEntity_Edge:
- if (len == 2)
- aVal = getDistance( P( 1 ), P( 2 ) );
- break;
- case SMDSEntity_Quad_Edge:
- if (len == 3) // quadratic edge
- aVal = getDistance(P( 1 ),P( 3 )) + getDistance(P( 3 ),P( 2 ));
- break;
- case SMDSEntity_Triangle:
- if (len == 3){ // triangles
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 1 ));
- aVal = Min(L1,Min(L2,L3));
- }
- break;
- case SMDSEntity_Quadrangle:
- if (len == 4){ // quadrangles
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 1 ));
- aVal = Min(Min(L1,L2),Min(L3,L4));
- }
- break;
- case SMDSEntity_Quad_Triangle:
- case SMDSEntity_BiQuad_Triangle:
- if (len >= 6){ // quadratic triangles
- double L1 = getDistance(P( 1 ),P( 2 )) + getDistance(P( 2 ),P( 3 ));
- double L2 = getDistance(P( 3 ),P( 4 )) + getDistance(P( 4 ),P( 5 ));
- double L3 = getDistance(P( 5 ),P( 6 )) + getDistance(P( 6 ),P( 1 ));
- aVal = Min(L1,Min(L2,L3));
- }
- break;
- case SMDSEntity_Quad_Quadrangle:
- case SMDSEntity_BiQuad_Quadrangle:
- if (len >= 8){ // quadratic quadrangles
- double L1 = getDistance(P( 1 ),P( 2 )) + getDistance(P( 2 ),P( 3 ));
- double L2 = getDistance(P( 3 ),P( 4 )) + getDistance(P( 4 ),P( 5 ));
- double L3 = getDistance(P( 5 ),P( 6 )) + getDistance(P( 6 ),P( 7 ));
- double L4 = getDistance(P( 7 ),P( 8 )) + getDistance(P( 8 ),P( 1 ));
- aVal = Min(Min(L1,L2),Min(L3,L4));
- }
- break;
- case SMDSEntity_Tetra:
- if (len == 4){ // tetrahedra
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 1 ));
- double L4 = getDistance(P( 1 ),P( 4 ));
- double L5 = getDistance(P( 2 ),P( 4 ));
- double L6 = getDistance(P( 3 ),P( 4 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- }
- break;
- case SMDSEntity_Pyramid:
- if (len == 5){ // piramids
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 1 ));
- double L5 = getDistance(P( 1 ),P( 5 ));
- double L6 = getDistance(P( 2 ),P( 5 ));
- double L7 = getDistance(P( 3 ),P( 5 ));
- double L8 = getDistance(P( 4 ),P( 5 ));
-
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(L7,L8));
- }
- break;
- case SMDSEntity_Penta:
- if (len == 6) { // pentaidres
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 1 ));
- double L4 = getDistance(P( 4 ),P( 5 ));
- double L5 = getDistance(P( 5 ),P( 6 ));
- double L6 = getDistance(P( 6 ),P( 4 ));
- double L7 = getDistance(P( 1 ),P( 4 ));
- double L8 = getDistance(P( 2 ),P( 5 ));
- double L9 = getDistance(P( 3 ),P( 6 ));
-
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(Min(L7,L8),L9));
- }
- break;
- case SMDSEntity_Hexa:
- if (len == 8){ // hexahedron
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 1 ));
- double L5 = getDistance(P( 5 ),P( 6 ));
- double L6 = getDistance(P( 6 ),P( 7 ));
- double L7 = getDistance(P( 7 ),P( 8 ));
- double L8 = getDistance(P( 8 ),P( 5 ));
- double L9 = getDistance(P( 1 ),P( 5 ));
- double L10= getDistance(P( 2 ),P( 6 ));
- double L11= getDistance(P( 3 ),P( 7 ));
- double L12= getDistance(P( 4 ),P( 8 ));
-
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(Min(L7,L8),Min(L9,L10)));
- aVal = Min(aVal,Min(L11,L12));
- }
- break;
- case SMDSEntity_Quad_Tetra:
- if (len == 10){ // quadratic tetraidrs
- double L1 = getDistance(P( 1 ),P( 5 )) + getDistance(P( 5 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 6 )) + getDistance(P( 6 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 7 )) + getDistance(P( 7 ),P( 1 ));
- double L4 = getDistance(P( 1 ),P( 8 )) + getDistance(P( 8 ),P( 4 ));
- double L5 = getDistance(P( 2 ),P( 9 )) + getDistance(P( 9 ),P( 4 ));
- double L6 = getDistance(P( 3 ),P( 10 )) + getDistance(P( 10 ),P( 4 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- }
- break;
- case SMDSEntity_Quad_Pyramid:
- if (len == 13){ // quadratic piramids
- double L1 = getDistance(P( 1 ),P( 6 )) + getDistance(P( 6 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 7 )) + getDistance(P( 7 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 8 )) + getDistance(P( 8 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 9 )) + getDistance(P( 9 ),P( 1 ));
- double L5 = getDistance(P( 1 ),P( 10 )) + getDistance(P( 10 ),P( 5 ));
- double L6 = getDistance(P( 2 ),P( 11 )) + getDistance(P( 11 ),P( 5 ));
- double L7 = getDistance(P( 3 ),P( 12 )) + getDistance(P( 12 ),P( 5 ));
- double L8 = getDistance(P( 4 ),P( 13 )) + getDistance(P( 13 ),P( 5 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(L7,L8));
- }
- break;
- case SMDSEntity_Quad_Penta:
- if (len == 15){ // quadratic pentaidres
- double L1 = getDistance(P( 1 ),P( 7 )) + getDistance(P( 7 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 8 )) + getDistance(P( 8 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 9 )) + getDistance(P( 9 ),P( 1 ));
- double L4 = getDistance(P( 4 ),P( 10 )) + getDistance(P( 10 ),P( 5 ));
- double L5 = getDistance(P( 5 ),P( 11 )) + getDistance(P( 11 ),P( 6 ));
- double L6 = getDistance(P( 6 ),P( 12 )) + getDistance(P( 12 ),P( 4 ));
- double L7 = getDistance(P( 1 ),P( 13 )) + getDistance(P( 13 ),P( 4 ));
- double L8 = getDistance(P( 2 ),P( 14 )) + getDistance(P( 14 ),P( 5 ));
- double L9 = getDistance(P( 3 ),P( 15 )) + getDistance(P( 15 ),P( 6 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(Min(L7,L8),L9));
- }
- break;
- case SMDSEntity_Quad_Hexa:
- case SMDSEntity_TriQuad_Hexa:
- if (len >= 20) { // quadratic hexaider
- double L1 = getDistance(P( 1 ),P( 9 )) + getDistance(P( 9 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 10 )) + getDistance(P( 10 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 11 )) + getDistance(P( 11 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 12 )) + getDistance(P( 12 ),P( 1 ));
- double L5 = getDistance(P( 5 ),P( 13 )) + getDistance(P( 13 ),P( 6 ));
- double L6 = getDistance(P( 6 ),P( 14 )) + getDistance(P( 14 ),P( 7 ));
- double L7 = getDistance(P( 7 ),P( 15 )) + getDistance(P( 15 ),P( 8 ));
- double L8 = getDistance(P( 8 ),P( 16 )) + getDistance(P( 16 ),P( 5 ));
- double L9 = getDistance(P( 1 ),P( 17 )) + getDistance(P( 17 ),P( 5 ));
- double L10= getDistance(P( 2 ),P( 18 )) + getDistance(P( 18 ),P( 6 ));
- double L11= getDistance(P( 3 ),P( 19 )) + getDistance(P( 19 ),P( 7 ));
- double L12= getDistance(P( 4 ),P( 20 )) + getDistance(P( 20 ),P( 8 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal,Min(Min(L7,L8),Min(L9,L10)));
- aVal = Min(aVal,Min(L11,L12));
- }
- break;
- case SMDSEntity_Polygon:
- if ( len > 1 ) {
- aVal = getDistance( P(1), P( P.size() ));
- for ( size_t i = 1; i < P.size(); ++i )
- aVal = Min( aVal, getDistance( P( i ), P( i+1 )));
- }
- break;
- case SMDSEntity_Quad_Polygon:
- if ( len > 2 ) {
- aVal = getDistance( P(1), P( P.size() )) + getDistance( P(P.size()), P( P.size()-1 ));
- for ( size_t i = 1; i < P.size()-1; i += 2 )
- aVal = Min( aVal, getDistance( P( i ), P( i+1 )) + getDistance( P( i+1 ), P( i+2 )));
- }
- break;
- case SMDSEntity_Hexagonal_Prism:
- if (len == 12) { // hexagonal prism
- double L1 = getDistance(P( 1 ),P( 2 ));
- double L2 = getDistance(P( 2 ),P( 3 ));
- double L3 = getDistance(P( 3 ),P( 4 ));
- double L4 = getDistance(P( 4 ),P( 5 ));
- double L5 = getDistance(P( 5 ),P( 6 ));
- double L6 = getDistance(P( 6 ),P( 1 ));
-
- double L7 = getDistance(P( 7 ), P( 8 ));
- double L8 = getDistance(P( 8 ), P( 9 ));
- double L9 = getDistance(P( 9 ), P( 10 ));
- double L10= getDistance(P( 10 ),P( 11 ));
- double L11= getDistance(P( 11 ),P( 12 ));
- double L12= getDistance(P( 12 ),P( 7 ));
-
- double L13 = getDistance(P( 1 ),P( 7 ));
- double L14 = getDistance(P( 2 ),P( 8 ));
- double L15 = getDistance(P( 3 ),P( 9 ));
- double L16 = getDistance(P( 4 ),P( 10 ));
- double L17 = getDistance(P( 5 ),P( 11 ));
- double L18 = getDistance(P( 6 ),P( 12 ));
- aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
- aVal = Min(aVal, Min(Min(Min(L7,L8),Min(L9,L10)),Min(L11,L12)));
- aVal = Min(aVal, Min(Min(Min(L13,L14),Min(L15,L16)),Min(L17,L18)));
- }
- break;
- case SMDSEntity_Polyhedra:
- {
- }
- break;
- default:
- return 0;
- }
-
- if (aVal < 0 ) {
- return 0.;
- }
-
- if ( myPrecision >= 0 )
- {
- double prec = pow( 10., (double)( myPrecision ) );
- aVal = floor( aVal * prec + 0.5 ) / prec;
- }
-
- return aVal;
-
- }
- return 0.;
-}
-
-double Length2D::GetBadRate( double Value, int /*nbNodes*/ ) const
-{
- // meaningless as it is not a quality control functor
- return Value;
-}
-
-SMDSAbs_ElementType Length2D::GetType() const
-{
- return SMDSAbs_Face;
-}
-
-Length2D::Value::Value(double theLength,long thePntId1, long thePntId2):
- myLength(theLength)
-{
- myPntId[0] = thePntId1; myPntId[1] = thePntId2;
- if(thePntId1 > thePntId2){
- myPntId[1] = thePntId1; myPntId[0] = thePntId2;
- }
-}
-
-bool Length2D::Value::operator<(const Length2D::Value& x) const
-{
- if(myPntId[0] < x.myPntId[0]) return true;
- if(myPntId[0] == x.myPntId[0])
- if(myPntId[1] < x.myPntId[1]) return true;
- return false;
-}
-
-void Length2D::GetValues(TValues& theValues)
+Length3D::Length3D():
+ Length2D ( SMDSAbs_Volume )
{
- TValues aValues;
- SMDS_FaceIteratorPtr anIter = myMesh->facesIterator();
- for(; anIter->more(); ){
- const SMDS_MeshFace* anElem = anIter->next();
-
- if(anElem->IsQuadratic()) {
- const SMDS_VtkFace* F =
- dynamic_cast<const SMDS_VtkFace*>(anElem);
- // use special nodes iterator
- SMDS_ElemIteratorPtr anIter = F->interlacedNodesElemIterator();
- long aNodeId[4];
- gp_Pnt P[4];
-
- double aLength;
- const SMDS_MeshElement* aNode;
- if(anIter->more()){
- aNode = anIter->next();
- const SMDS_MeshNode* aNodes = (SMDS_MeshNode*) aNode;
- P[0] = P[1] = gp_Pnt(aNodes->X(),aNodes->Y(),aNodes->Z());
- aNodeId[0] = aNodeId[1] = aNode->GetID();
- aLength = 0;
- }
- for(; anIter->more(); ){
- const SMDS_MeshNode* N1 = static_cast<const SMDS_MeshNode*> (anIter->next());
- P[2] = gp_Pnt(N1->X(),N1->Y(),N1->Z());
- aNodeId[2] = N1->GetID();
- aLength = P[1].Distance(P[2]);
- if(!anIter->more()) break;
- const SMDS_MeshNode* N2 = static_cast<const SMDS_MeshNode*> (anIter->next());
- P[3] = gp_Pnt(N2->X(),N2->Y(),N2->Z());
- aNodeId[3] = N2->GetID();
- aLength += P[2].Distance(P[3]);
- Value aValue1(aLength,aNodeId[1],aNodeId[2]);
- Value aValue2(aLength,aNodeId[2],aNodeId[3]);
- P[1] = P[3];
- aNodeId[1] = aNodeId[3];
- theValues.insert(aValue1);
- theValues.insert(aValue2);
- }
- aLength += P[2].Distance(P[0]);
- Value aValue1(aLength,aNodeId[1],aNodeId[2]);
- Value aValue2(aLength,aNodeId[2],aNodeId[0]);
- theValues.insert(aValue1);
- theValues.insert(aValue2);
- }
- else {
- SMDS_ElemIteratorPtr aNodesIter = anElem->nodesIterator();
- long aNodeId[2];
- gp_Pnt P[3];
-
- double aLength;
- const SMDS_MeshElement* aNode;
- if(aNodesIter->more()){
- aNode = aNodesIter->next();
- const SMDS_MeshNode* aNodes = (SMDS_MeshNode*) aNode;
- P[0] = P[1] = gp_Pnt(aNodes->X(),aNodes->Y(),aNodes->Z());
- aNodeId[0] = aNodeId[1] = aNode->GetID();
- aLength = 0;
- }
- for(; aNodesIter->more(); ){
- aNode = aNodesIter->next();
- const SMDS_MeshNode* aNodes = (SMDS_MeshNode*) aNode;
- long anId = aNode->GetID();
-
- P[2] = gp_Pnt(aNodes->X(),aNodes->Y(),aNodes->Z());
-
- aLength = P[1].Distance(P[2]);
-
- Value aValue(aLength,aNodeId[1],anId);
- aNodeId[1] = anId;
- P[1] = P[2];
- theValues.insert(aValue);
- }
-
- aLength = P[0].Distance(P[1]);
-
- Value aValue(aLength,aNodeId[0],aNodeId[1]);
- theValues.insert(aValue);
- }
- }
}
//================================================================================
/*
- Class : MultiConnection
- Description : Functor for calculating number of faces conneted to the edge
+ Class : Length2D
+ Description : Functor for calculating minimal length of element edge
*/
//================================================================================
-double MultiConnection::GetValue( const TSequenceOfXYZ& P )
-{
- return 0;
-}
-double MultiConnection::GetValue( long theId )
-{
- return getNbMultiConnection( myMesh, theId );
-}
-
-double MultiConnection::GetBadRate( double Value, int /*nbNodes*/ ) const
-{
- // meaningless as it is not quality control functor
- return Value;
-}
-
-SMDSAbs_ElementType MultiConnection::GetType() const
+Length2D::Length2D( SMDSAbs_ElementType type ):
+ myType ( type )
{
- return SMDSAbs_Edge;
}
-//================================================================================
-/*
- Class : MultiConnection2D
- Description : Functor for calculating number of faces conneted to the edge
-*/
-//================================================================================
-
-double MultiConnection2D::GetValue( const TSequenceOfXYZ& P )
+bool Length2D::IsApplicable( const SMDS_MeshElement* element ) const
{
- return 0;
+ return ( NumericalFunctor::IsApplicable( element ) &&
+ element->GetEntityType() != SMDSEntity_Polyhedra );
}
-double MultiConnection2D::GetValue( long theElementId )
+double Length2D::GetValue( const TSequenceOfXYZ& P )
{
- int aResult = 0;
-
- const SMDS_MeshElement* aFaceElem = myMesh->FindElement(theElementId);
- SMDSAbs_ElementType aType = aFaceElem->GetType();
+ double aVal = 0;
+ int len = P.size();
+ SMDSAbs_EntityType aType = P.getElementEntity();
switch (aType) {
- case SMDSAbs_Face:
- {
- int i = 0, len = aFaceElem->NbNodes();
- SMDS_ElemIteratorPtr anIter = aFaceElem->nodesIterator();
- if (!anIter) break;
-
- const SMDS_MeshNode *aNode, *aNode0;
- TColStd_MapOfInteger aMap, aMapPrev;
-
- for (i = 0; i <= len; i++) {
- aMapPrev = aMap;
- aMap.Clear();
-
- int aNb = 0;
- if (anIter->more()) {
- aNode = (SMDS_MeshNode*)anIter->next();
- } else {
- if (i == len)
- aNode = aNode0;
- else
- break;
- }
- if (!aNode) break;
- if (i == 0) aNode0 = aNode;
-
- SMDS_ElemIteratorPtr anElemIter = aNode->GetInverseElementIterator();
- while (anElemIter->more()) {
- const SMDS_MeshElement* anElem = anElemIter->next();
- if (anElem != 0 && anElem->GetType() == SMDSAbs_Face) {
- int anId = anElem->GetID();
-
- aMap.Add(anId);
- if (aMapPrev.Contains(anId)) {
- aNb++;
- }
- }
- }
- aResult = Max(aResult, aNb);
- }
- }
+ case SMDSEntity_Edge:
+ if (len == 2)
+ aVal = getDistance( P( 1 ), P( 2 ) );
break;
- default:
+ case SMDSEntity_Quad_Edge:
+ if (len == 3) // quadratic edge
+ aVal = getDistance(P( 1 ),P( 3 )) + getDistance(P( 3 ),P( 2 ));
+ break;
+ case SMDSEntity_Triangle:
+ if (len == 3){ // triangles
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 1 ));
+ aVal = Min(L1,Min(L2,L3));
+ }
+ break;
+ case SMDSEntity_Quadrangle:
+ if (len == 4){ // quadrangles
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 1 ));
+ aVal = Min(Min(L1,L2),Min(L3,L4));
+ }
+ break;
+ case SMDSEntity_Quad_Triangle:
+ case SMDSEntity_BiQuad_Triangle:
+ if (len >= 6){ // quadratic triangles
+ double L1 = getDistance(P( 1 ),P( 2 )) + getDistance(P( 2 ),P( 3 ));
+ double L2 = getDistance(P( 3 ),P( 4 )) + getDistance(P( 4 ),P( 5 ));
+ double L3 = getDistance(P( 5 ),P( 6 )) + getDistance(P( 6 ),P( 1 ));
+ aVal = Min(L1,Min(L2,L3));
+ }
+ break;
+ case SMDSEntity_Quad_Quadrangle:
+ case SMDSEntity_BiQuad_Quadrangle:
+ if (len >= 8){ // quadratic quadrangles
+ double L1 = getDistance(P( 1 ),P( 2 )) + getDistance(P( 2 ),P( 3 ));
+ double L2 = getDistance(P( 3 ),P( 4 )) + getDistance(P( 4 ),P( 5 ));
+ double L3 = getDistance(P( 5 ),P( 6 )) + getDistance(P( 6 ),P( 7 ));
+ double L4 = getDistance(P( 7 ),P( 8 )) + getDistance(P( 8 ),P( 1 ));
+ aVal = Min(Min(L1,L2),Min(L3,L4));
+ }
+ break;
+ case SMDSEntity_Tetra:
+ if (len == 4){ // tetrahedra
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 1 ));
+ double L4 = getDistance(P( 1 ),P( 4 ));
+ double L5 = getDistance(P( 2 ),P( 4 ));
+ double L6 = getDistance(P( 3 ),P( 4 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ }
+ break;
+ case SMDSEntity_Pyramid:
+ if (len == 5){ // pyramid
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 1 ));
+ double L5 = getDistance(P( 1 ),P( 5 ));
+ double L6 = getDistance(P( 2 ),P( 5 ));
+ double L7 = getDistance(P( 3 ),P( 5 ));
+ double L8 = getDistance(P( 4 ),P( 5 ));
+
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(L7,L8));
+ }
+ break;
+ case SMDSEntity_Penta:
+ if (len == 6) { // pentahedron
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 1 ));
+ double L4 = getDistance(P( 4 ),P( 5 ));
+ double L5 = getDistance(P( 5 ),P( 6 ));
+ double L6 = getDistance(P( 6 ),P( 4 ));
+ double L7 = getDistance(P( 1 ),P( 4 ));
+ double L8 = getDistance(P( 2 ),P( 5 ));
+ double L9 = getDistance(P( 3 ),P( 6 ));
+
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(Min(L7,L8),L9));
+ }
+ break;
+ case SMDSEntity_Hexa:
+ if (len == 8){ // hexahedron
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 1 ));
+ double L5 = getDistance(P( 5 ),P( 6 ));
+ double L6 = getDistance(P( 6 ),P( 7 ));
+ double L7 = getDistance(P( 7 ),P( 8 ));
+ double L8 = getDistance(P( 8 ),P( 5 ));
+ double L9 = getDistance(P( 1 ),P( 5 ));
+ double L10= getDistance(P( 2 ),P( 6 ));
+ double L11= getDistance(P( 3 ),P( 7 ));
+ double L12= getDistance(P( 4 ),P( 8 ));
+
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(Min(L7,L8),Min(L9,L10)));
+ aVal = Min(aVal,Min(L11,L12));
+ }
+ break;
+ case SMDSEntity_Quad_Tetra:
+ if (len == 10){ // quadratic tetrahedron
+ double L1 = getDistance(P( 1 ),P( 5 )) + getDistance(P( 5 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 6 )) + getDistance(P( 6 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 7 )) + getDistance(P( 7 ),P( 1 ));
+ double L4 = getDistance(P( 1 ),P( 8 )) + getDistance(P( 8 ),P( 4 ));
+ double L5 = getDistance(P( 2 ),P( 9 )) + getDistance(P( 9 ),P( 4 ));
+ double L6 = getDistance(P( 3 ),P( 10 )) + getDistance(P( 10 ),P( 4 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ }
+ break;
+ case SMDSEntity_Quad_Pyramid:
+ if (len == 13){ // quadratic pyramid
+ double L1 = getDistance(P( 1 ),P( 6 )) + getDistance(P( 6 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 7 )) + getDistance(P( 7 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 8 )) + getDistance(P( 8 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 9 )) + getDistance(P( 9 ),P( 1 ));
+ double L5 = getDistance(P( 1 ),P( 10 )) + getDistance(P( 10 ),P( 5 ));
+ double L6 = getDistance(P( 2 ),P( 11 )) + getDistance(P( 11 ),P( 5 ));
+ double L7 = getDistance(P( 3 ),P( 12 )) + getDistance(P( 12 ),P( 5 ));
+ double L8 = getDistance(P( 4 ),P( 13 )) + getDistance(P( 13 ),P( 5 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(L7,L8));
+ }
+ break;
+ case SMDSEntity_Quad_Penta:
+ case SMDSEntity_BiQuad_Penta:
+ if (len >= 15){ // quadratic pentahedron
+ double L1 = getDistance(P( 1 ),P( 7 )) + getDistance(P( 7 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 8 )) + getDistance(P( 8 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 9 )) + getDistance(P( 9 ),P( 1 ));
+ double L4 = getDistance(P( 4 ),P( 10 )) + getDistance(P( 10 ),P( 5 ));
+ double L5 = getDistance(P( 5 ),P( 11 )) + getDistance(P( 11 ),P( 6 ));
+ double L6 = getDistance(P( 6 ),P( 12 )) + getDistance(P( 12 ),P( 4 ));
+ double L7 = getDistance(P( 1 ),P( 13 )) + getDistance(P( 13 ),P( 4 ));
+ double L8 = getDistance(P( 2 ),P( 14 )) + getDistance(P( 14 ),P( 5 ));
+ double L9 = getDistance(P( 3 ),P( 15 )) + getDistance(P( 15 ),P( 6 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(Min(L7,L8),L9));
+ }
+ break;
+ case SMDSEntity_Quad_Hexa:
+ case SMDSEntity_TriQuad_Hexa:
+ if (len >= 20) { // quadratic hexahedron
+ double L1 = getDistance(P( 1 ),P( 9 )) + getDistance(P( 9 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 10 )) + getDistance(P( 10 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 11 )) + getDistance(P( 11 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 12 )) + getDistance(P( 12 ),P( 1 ));
+ double L5 = getDistance(P( 5 ),P( 13 )) + getDistance(P( 13 ),P( 6 ));
+ double L6 = getDistance(P( 6 ),P( 14 )) + getDistance(P( 14 ),P( 7 ));
+ double L7 = getDistance(P( 7 ),P( 15 )) + getDistance(P( 15 ),P( 8 ));
+ double L8 = getDistance(P( 8 ),P( 16 )) + getDistance(P( 16 ),P( 5 ));
+ double L9 = getDistance(P( 1 ),P( 17 )) + getDistance(P( 17 ),P( 5 ));
+ double L10= getDistance(P( 2 ),P( 18 )) + getDistance(P( 18 ),P( 6 ));
+ double L11= getDistance(P( 3 ),P( 19 )) + getDistance(P( 19 ),P( 7 ));
+ double L12= getDistance(P( 4 ),P( 20 )) + getDistance(P( 20 ),P( 8 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal,Min(Min(L7,L8),Min(L9,L10)));
+ aVal = Min(aVal,Min(L11,L12));
+ }
+ break;
+ case SMDSEntity_Polygon:
+ if ( len > 1 ) {
+ aVal = getDistance( P(1), P( P.size() ));
+ for ( size_t i = 1; i < P.size(); ++i )
+ aVal = Min( aVal, getDistance( P( i ), P( i+1 )));
+ }
+ break;
+ case SMDSEntity_Quad_Polygon:
+ if ( len > 2 ) {
+ aVal = getDistance( P(1), P( P.size() )) + getDistance( P(P.size()), P( P.size()-1 ));
+ for ( size_t i = 1; i < P.size()-1; i += 2 )
+ aVal = Min( aVal, getDistance( P( i ), P( i+1 )) + getDistance( P( i+1 ), P( i+2 )));
+ }
+ break;
+ case SMDSEntity_Hexagonal_Prism:
+ if (len == 12) { // hexagonal prism
+ double L1 = getDistance(P( 1 ),P( 2 ));
+ double L2 = getDistance(P( 2 ),P( 3 ));
+ double L3 = getDistance(P( 3 ),P( 4 ));
+ double L4 = getDistance(P( 4 ),P( 5 ));
+ double L5 = getDistance(P( 5 ),P( 6 ));
+ double L6 = getDistance(P( 6 ),P( 1 ));
+
+ double L7 = getDistance(P( 7 ), P( 8 ));
+ double L8 = getDistance(P( 8 ), P( 9 ));
+ double L9 = getDistance(P( 9 ), P( 10 ));
+ double L10= getDistance(P( 10 ),P( 11 ));
+ double L11= getDistance(P( 11 ),P( 12 ));
+ double L12= getDistance(P( 12 ),P( 7 ));
+
+ double L13 = getDistance(P( 1 ),P( 7 ));
+ double L14 = getDistance(P( 2 ),P( 8 ));
+ double L15 = getDistance(P( 3 ),P( 9 ));
+ double L16 = getDistance(P( 4 ),P( 10 ));
+ double L17 = getDistance(P( 5 ),P( 11 ));
+ double L18 = getDistance(P( 6 ),P( 12 ));
+ aVal = Min(Min(Min(L1,L2),Min(L3,L4)),Min(L5,L6));
+ aVal = Min(aVal, Min(Min(Min(L7,L8),Min(L9,L10)),Min(L11,L12)));
+ aVal = Min(aVal, Min(Min(Min(L13,L14),Min(L15,L16)),Min(L17,L18)));
+ }
+ break;
+ case SMDSEntity_Polyhedra:
+ {
+ }
+ break;
+ default:
+ return 0;
+ }
+
+ if (aVal < 0 ) {
+ return 0.;
+ }
+
+ if ( myPrecision >= 0 )
+ {
+ double prec = pow( 10., (double)( myPrecision ) );
+ aVal = floor( aVal * prec + 0.5 ) / prec;
+ }
+
+ return aVal;
+}
+
+double Length2D::GetBadRate( double Value, int /*nbNodes*/ ) const
+{
+ // meaningless as it is not a quality control functor
+ return Value;
+}
+
+SMDSAbs_ElementType Length2D::GetType() const
+{
+ return myType;
+}
+
+Length2D::Value::Value(double theLength,long thePntId1, long thePntId2):
+ myLength(theLength)
+{
+ myPntId[0] = thePntId1; myPntId[1] = thePntId2;
+ if(thePntId1 > thePntId2){
+ myPntId[1] = thePntId1; myPntId[0] = thePntId2;
+ }
+}
+
+bool Length2D::Value::operator<(const Length2D::Value& x) const
+{
+ if(myPntId[0] < x.myPntId[0]) return true;
+ if(myPntId[0] == x.myPntId[0])
+ if(myPntId[1] < x.myPntId[1]) return true;
+ return false;
+}
+
+void Length2D::GetValues(TValues& theValues)
+{
+ if ( myType == SMDSAbs_Face )
+ {
+ for ( SMDS_FaceIteratorPtr anIter = myMesh->facesIterator(); anIter->more(); )
+ {
+ const SMDS_MeshFace* anElem = anIter->next();
+ if ( anElem->IsQuadratic() )
+ {
+ // use special nodes iterator
+ SMDS_NodeIteratorPtr anIter = anElem->interlacedNodesIterator();
+ long aNodeId[4] = { 0,0,0,0 };
+ gp_Pnt P[4];
+
+ double aLength = 0;
+ if ( anIter->more() )
+ {
+ const SMDS_MeshNode* aNode = anIter->next();
+ P[0] = P[1] = SMESH_NodeXYZ( aNode );
+ aNodeId[0] = aNodeId[1] = aNode->GetID();
+ aLength = 0;
+ }
+ for ( ; anIter->more(); )
+ {
+ const SMDS_MeshNode* N1 = anIter->next();
+ P[2] = SMESH_NodeXYZ( N1 );
+ aNodeId[2] = N1->GetID();
+ aLength = P[1].Distance(P[2]);
+ if(!anIter->more()) break;
+ const SMDS_MeshNode* N2 = anIter->next();
+ P[3] = SMESH_NodeXYZ( N2 );
+ aNodeId[3] = N2->GetID();
+ aLength += P[2].Distance(P[3]);
+ Value aValue1(aLength,aNodeId[1],aNodeId[2]);
+ Value aValue2(aLength,aNodeId[2],aNodeId[3]);
+ P[1] = P[3];
+ aNodeId[1] = aNodeId[3];
+ theValues.insert(aValue1);
+ theValues.insert(aValue2);
+ }
+ aLength += P[2].Distance(P[0]);
+ Value aValue1(aLength,aNodeId[1],aNodeId[2]);
+ Value aValue2(aLength,aNodeId[2],aNodeId[0]);
+ theValues.insert(aValue1);
+ theValues.insert(aValue2);
+ }
+ else {
+ SMDS_NodeIteratorPtr aNodesIter = anElem->nodeIterator();
+ long aNodeId[2] = {0,0};
+ gp_Pnt P[3];
+
+ double aLength;
+ const SMDS_MeshElement* aNode;
+ if ( aNodesIter->more())
+ {
+ aNode = aNodesIter->next();
+ P[0] = P[1] = SMESH_NodeXYZ( aNode );
+ aNodeId[0] = aNodeId[1] = aNode->GetID();
+ aLength = 0;
+ }
+ for( ; aNodesIter->more(); )
+ {
+ aNode = aNodesIter->next();
+ long anId = aNode->GetID();
+
+ P[2] = SMESH_NodeXYZ( aNode );
+
+ aLength = P[1].Distance(P[2]);
+
+ Value aValue(aLength,aNodeId[1],anId);
+ aNodeId[1] = anId;
+ P[1] = P[2];
+ theValues.insert(aValue);
+ }
+
+ aLength = P[0].Distance(P[1]);
+
+ Value aValue(aLength,aNodeId[0],aNodeId[1]);
+ theValues.insert(aValue);
+ }
+ }
+ }
+ else
+ {
+ // not implemented
+ }
+}
+
+//================================================================================
+/*
+ Class : Deflection2D
+ Description : computes distance between a face center and an underlying surface
+*/
+//================================================================================
+
+double Deflection2D::GetValue( const TSequenceOfXYZ& P )
+{
+ if ( myMesh && P.getElement() )
+ {
+ // get underlying surface
+ if ( myShapeIndex != P.getElement()->getshapeId() )
+ {
+ mySurface.Nullify();
+ myShapeIndex = P.getElement()->getshapeId();
+ const TopoDS_Shape& S =
+ static_cast< const SMESHDS_Mesh* >( myMesh )->IndexToShape( myShapeIndex );
+ if ( !S.IsNull() && S.ShapeType() == TopAbs_FACE )
+ {
+ mySurface = new ShapeAnalysis_Surface( BRep_Tool::Surface( TopoDS::Face( S )));
+
+ GeomLib_IsPlanarSurface isPlaneCheck( mySurface->Surface() );
+ if ( isPlaneCheck.IsPlanar() )
+ myPlane.reset( new gp_Pln( isPlaneCheck.Plan() ));
+ else
+ myPlane.reset();
+ }
+ }
+ // project gravity center to the surface
+ if ( !mySurface.IsNull() )
+ {
+ gp_XYZ gc(0,0,0);
+ gp_XY uv(0,0);
+ int nbUV = 0;
+ for ( size_t i = 0; i < P.size(); ++i )
+ {
+ gc += P(i+1);
+
+ if ( SMDS_FacePositionPtr fPos = P.getElement()->GetNode( i )->GetPosition() )
+ {
+ uv.ChangeCoord(1) += fPos->GetUParameter();
+ uv.ChangeCoord(2) += fPos->GetVParameter();
+ ++nbUV;
+ }
+ }
+ gc /= P.size();
+ if ( nbUV ) uv /= nbUV;
+
+ double maxLen = MaxElementLength2D().GetValue( P );
+ double tol = 1e-3 * maxLen;
+ double dist;
+ if ( myPlane )
+ {
+ dist = myPlane->Distance( gc );
+ if ( dist < tol )
+ dist = 0;
+ }
+ else
+ {
+ if ( uv.X() != 0 && uv.Y() != 0 ) // faster way
+ mySurface->NextValueOfUV( uv, gc, tol, 0.5 * maxLen );
+ else
+ mySurface->ValueOfUV( gc, tol );
+ dist = mySurface->Gap();
+ }
+ return Round( dist );
+ }
+ }
+ return 0;
+}
+
+void Deflection2D::SetMesh( const SMDS_Mesh* theMesh )
+{
+ NumericalFunctor::SetMesh( dynamic_cast<const SMESHDS_Mesh* >( theMesh ));
+ myShapeIndex = -100;
+ myPlane.reset();
+}
+
+SMDSAbs_ElementType Deflection2D::GetType() const
+{
+ return SMDSAbs_Face;
+}
+
+double Deflection2D::GetBadRate( double Value, int /*nbNodes*/ ) const
+{
+ // meaningless as it is not quality control functor
+ return Value;
+}
+
+//================================================================================
+/*
+ Class : MultiConnection
+ Description : Functor for calculating number of faces conneted to the edge
+*/
+//================================================================================
+
+double MultiConnection::GetValue( const TSequenceOfXYZ& /*P*/ )
+{
+ return 0;
+}
+double MultiConnection::GetValue( long theId )
+{
+ return getNbMultiConnection( myMesh, theId );
+}
+
+double MultiConnection::GetBadRate( double Value, int /*nbNodes*/ ) const
+{
+ // meaningless as it is not quality control functor
+ return Value;
+}
+
+SMDSAbs_ElementType MultiConnection::GetType() const
+{
+ return SMDSAbs_Edge;
+}
+
+//================================================================================
+/*
+ Class : MultiConnection2D
+ Description : Functor for calculating number of faces conneted to the edge
+*/
+//================================================================================
+
+double MultiConnection2D::GetValue( const TSequenceOfXYZ& /*P*/ )
+{
+ return 0;
+}
+
+double MultiConnection2D::GetValue( long theElementId )
+{
+ int aResult = 0;
+
+ const SMDS_MeshElement* aFaceElem = myMesh->FindElement(theElementId);
+ SMDSAbs_ElementType aType = aFaceElem->GetType();
+
+ switch (aType) {
+ case SMDSAbs_Face:
+ {
+ int i = 0, len = aFaceElem->NbNodes();
+ SMDS_ElemIteratorPtr anIter = aFaceElem->nodesIterator();
+ if (!anIter) break;
+
+ const SMDS_MeshNode *aNode, *aNode0 = 0;
+ TColStd_MapOfInteger aMap, aMapPrev;
+
+ for (i = 0; i <= len; i++) {
+ aMapPrev = aMap;
+ aMap.Clear();
+
+ int aNb = 0;
+ if (anIter->more()) {
+ aNode = (SMDS_MeshNode*)anIter->next();
+ } else {
+ if (i == len)
+ aNode = aNode0;
+ else
+ break;
+ }
+ if (!aNode) break;
+ if (i == 0) aNode0 = aNode;
+
+ SMDS_ElemIteratorPtr anElemIter = aNode->GetInverseElementIterator();
+ while (anElemIter->more()) {
+ const SMDS_MeshElement* anElem = anElemIter->next();
+ if (anElem != 0 && anElem->GetType() == SMDSAbs_Face) {
+ int anId = anElem->GetID();
+
+ aMap.Add(anId);
+ if (aMapPrev.Contains(anId)) {
+ aNb++;
+ }
+ }
+ }
+ aResult = Max(aResult, aNb);
+ }
+ }
+ break;
+ default:
aResult = 0;
}
void MultiConnection2D::GetValues(MValues& theValues)
{
if ( !myMesh ) return;
- SMDS_FaceIteratorPtr anIter = myMesh->facesIterator();
- for(; anIter->more(); ){
- const SMDS_MeshFace* anElem = anIter->next();
- SMDS_ElemIteratorPtr aNodesIter;
- if ( anElem->IsQuadratic() )
- aNodesIter = dynamic_cast<const SMDS_VtkFace*>
- (anElem)->interlacedNodesElemIterator();
- else
- aNodesIter = anElem->nodesIterator();
- long aNodeId[3];
+ for ( SMDS_FaceIteratorPtr anIter = myMesh->facesIterator(); anIter->more(); )
+ {
+ const SMDS_MeshFace* anElem = anIter->next();
+ SMDS_NodeIteratorPtr aNodesIter = anElem->interlacedNodesIterator();
- //int aNbConnects=0;
- const SMDS_MeshNode* aNode0;
- const SMDS_MeshNode* aNode1;
+ const SMDS_MeshNode* aNode1 = anElem->GetNode( anElem->NbNodes() - 1 );
const SMDS_MeshNode* aNode2;
- if(aNodesIter->more()){
- aNode0 = (SMDS_MeshNode*) aNodesIter->next();
- aNode1 = aNode0;
- const SMDS_MeshNode* aNodes = (SMDS_MeshNode*) aNode1;
- aNodeId[0] = aNodeId[1] = aNodes->GetID();
- }
- for(; aNodesIter->more(); ) {
- aNode2 = (SMDS_MeshNode*) aNodesIter->next();
- long anId = aNode2->GetID();
- aNodeId[2] = anId;
-
- Value aValue(aNodeId[1],aNodeId[2]);
- MValues::iterator aItr = theValues.find(aValue);
- if (aItr != theValues.end()){
- aItr->second += 1;
- //aNbConnects = nb;
- }
- else {
- theValues[aValue] = 1;
- //aNbConnects = 1;
- }
- //cout << "NodeIds: "<<aNodeId[1]<<","<<aNodeId[2]<<" nbconn="<<aNbConnects<<endl;
- aNodeId[1] = aNodeId[2];
+ for ( ; aNodesIter->more(); )
+ {
+ aNode2 = aNodesIter->next();
+
+ Value aValue ( aNode1->GetID(), aNode2->GetID() );
+ MValues::iterator aItr = theValues.insert( std::make_pair( aValue, 0 )).first;
+ aItr->second++;
aNode1 = aNode2;
}
- Value aValue(aNodeId[0],aNodeId[2]);
- MValues::iterator aItr = theValues.find(aValue);
- if (aItr != theValues.end()) {
- aItr->second += 1;
- //aNbConnects = nb;
- }
- else {
- theValues[aValue] = 1;
- //aNbConnects = 1;
- }
- //cout << "NodeIds: "<<aNodeId[0]<<","<<aNodeId[2]<<" nbconn="<<aNbConnects<<endl;
}
-
+ return;
}
//================================================================================
double diameter = 0;
if ( const SMDS_BallElement* ball =
- dynamic_cast<const SMDS_BallElement*>( myMesh->FindElement( theId )))
+ myMesh->DownCast< SMDS_BallElement >( myMesh->FindElement( theId )))
{
diameter = ball->GetDiameter();
}
return SMDSAbs_Ball;
}
+//================================================================================
+/*
+ Class : NodeConnectivityNumber
+ Description : Functor returning number of elements connected to a node
+*/
+//================================================================================
+
+double NodeConnectivityNumber::GetValue( long theId )
+{
+ double nb = 0;
+
+ if ( const SMDS_MeshNode* node = myMesh->FindNode( theId ))
+ {
+ SMDSAbs_ElementType type;
+ if ( myMesh->NbVolumes() > 0 )
+ type = SMDSAbs_Volume;
+ else if ( myMesh->NbFaces() > 0 )
+ type = SMDSAbs_Face;
+ else if ( myMesh->NbEdges() > 0 )
+ type = SMDSAbs_Edge;
+ else
+ return 0;
+ nb = node->NbInverseElements( type );
+ }
+ return nb;
+}
+
+double NodeConnectivityNumber::GetBadRate( double Value, int /*nbNodes*/ ) const
+{
+ return Value;
+}
+
+SMDSAbs_ElementType NodeConnectivityNumber::GetType() const
+{
+ return SMDSAbs_Node;
+}
/*
PREDICATES
return false;
SMDS_VolumeTool vTool( myMesh->FindElement( theId ));
- return !vTool.IsForward();
+
+ bool isOk = true;
+ if ( vTool.IsPoly() )
+ {
+ isOk = true;
+ for ( int i = 0; i < vTool.NbFaces() && isOk; ++i )
+ isOk = vTool.IsFaceExternal( i );
+ }
+ else
+ {
+ isOk = vTool.IsForward();
+ }
+ return !isOk;
}
SMDSAbs_ElementType BadOrientedVolume::GetType() const
return SMDSAbs_Node;
}
+void CoincidentNodes::SetTolerance( const double theToler )
+{
+ if ( myToler != theToler )
+ {
+ SetMesh(0);
+ myToler = theToler;
+ }
+}
+
void CoincidentNodes::SetMesh( const SMDS_Mesh* theMesh )
{
myMeshModifTracer.SetMesh( theMesh );
if ( myMeshModifTracer.IsMeshModified() )
{
TIDSortedNodeSet nodesToCheck;
- SMDS_NodeIteratorPtr nIt = theMesh->nodesIterator(/*idInceasingOrder=*/true);
+ SMDS_NodeIteratorPtr nIt = theMesh->nodesIterator();
while ( nIt->more() )
nodesToCheck.insert( nodesToCheck.end(), nIt->next() );
bool FreeEdges::IsFreeEdge( const SMDS_MeshNode** theNodes, const int theFaceId )
{
- TColStd_MapOfInteger aMap;
- for ( int i = 0; i < 2; i++ )
+ SMDS_ElemIteratorPtr anElemIter = theNodes[ 0 ]->GetInverseElementIterator(SMDSAbs_Face);
+ while( anElemIter->more() )
{
- SMDS_ElemIteratorPtr anElemIter = theNodes[ i ]->GetInverseElementIterator(SMDSAbs_Face);
- while( anElemIter->more() )
+ if ( const SMDS_MeshElement* anElem = anElemIter->next())
{
- if ( const SMDS_MeshElement* anElem = anElemIter->next())
- {
- const int anId = anElem->GetID();
- if ( anId != theFaceId && !aMap.Add( anId ))
- return false;
- }
+ const int anId = anElem->GetID();
+ if ( anId != theFaceId && anElem->GetNodeIndex( theNodes[1] ) >= 0 )
+ return false;
}
}
return true;
void FreeEdges::GetBoreders(TBorders& theBorders)
{
TBorders aRegistry;
- SMDS_FaceIteratorPtr anIter = myMesh->facesIterator();
- for(; anIter->more(); ){
+ for ( SMDS_FaceIteratorPtr anIter = myMesh->facesIterator(); anIter->more(); )
+ {
const SMDS_MeshFace* anElem = anIter->next();
long anElemId = anElem->GetID();
- SMDS_ElemIteratorPtr aNodesIter;
- if ( anElem->IsQuadratic() )
- aNodesIter = static_cast<const SMDS_VtkFace*>(anElem)->
- interlacedNodesElemIterator();
- else
- aNodesIter = anElem->nodesIterator();
- long aNodeId[2];
- const SMDS_MeshElement* aNode;
- if(aNodesIter->more()){
- aNode = aNodesIter->next();
- aNodeId[0] = aNodeId[1] = aNode->GetID();
- }
- for(; aNodesIter->more(); ){
- aNode = aNodesIter->next();
- long anId = aNode->GetID();
- Border aBorder(anElemId,aNodeId[1],anId);
- aNodeId[1] = anId;
- UpdateBorders(aBorder,aRegistry,theBorders);
+ SMDS_NodeIteratorPtr aNodesIter = anElem->interlacedNodesIterator();
+ if ( !aNodesIter->more() ) continue;
+ long aNodeId[2] = {0,0};
+ aNodeId[0] = anElem->GetNode( anElem->NbNodes()-1 )->GetID();
+ for ( ; aNodesIter->more(); )
+ {
+ aNodeId[1] = aNodesIter->next()->GetID();
+ Border aBorder( anElemId, aNodeId[0], aNodeId[1] );
+ UpdateBorders( aBorder, aRegistry, theBorders );
+ aNodeId[0] = aNodeId[1];
}
- Border aBorder(anElemId,aNodeId[0],aNodeId[1]);
- UpdateBorders(aBorder,aRegistry,theBorders);
}
}
SMDS_MeshElement* aVol = (SMDS_MeshElement*)volItr->next();
TItrMapOfVolume itr = mapOfVol.insert( std::make_pair( aVol, 0 )).first;
(*itr).second++;
- }
+ }
}
int nbVol = 0;
TItrMapOfVolume volItr = mapOfVol.begin();
TItrMapOfVolume volEnd = mapOfVol.end();
for ( ; volItr != volEnd; ++volItr )
if ( (*volItr).second >= nbNode )
- nbVol++;
- // face is not free if number of volumes constructed on thier nodes more than one
+ nbVol++;
+ // face is not free if number of volumes constructed on their nodes more than one
return (nbVol < 2);
}
//================================================================================
/*
Class : GroupColor
- Description : Functor for check color of group to whic mesh element belongs to
+ Description : Functor for check color of group to which mesh element belongs to
*/
//================================================================================
{
// keep elements of myType
const SMDS_MeshElement* element = eIt->next();
- if ( element->GetType() == myType )
+ if ( myType == SMDSAbs_All || element->GetType() == myType )
myOkIDs.insert( myOkIDs.end(), element->GetID() );
// enqueue nodes of the element
double dot = v1 * v2; // cos * |v1| * |v2|
double l1 = v1.SquareMagnitude();
double l2 = v2.SquareMagnitude();
- return (( dot * cos >= 0 ) &&
+ return (( dot * cos >= 0 ) &&
( dot * dot ) / l1 / l2 >= ( cos * cos ));
}
}
{
char c = aStr.Value( i );
if ( !isdigit( c ) && c != ',' && c != '-' )
- aStr.SetValue( i, ' ');
+ aStr.SetValue( i, ',');
}
aStr.RemoveAll( ' ' );
myPredicate = thePredicate;
}
-void Filter::GetElementsId( const SMDS_Mesh* theMesh,
- PredicatePtr thePredicate,
- TIdSequence& theSequence )
+void Filter::GetElementsId( const SMDS_Mesh* theMesh,
+ PredicatePtr thePredicate,
+ TIdSequence& theSequence,
+ SMDS_ElemIteratorPtr theElements )
{
theSequence.clear();
thePredicate->SetMesh( theMesh );
- SMDS_ElemIteratorPtr elemIt = theMesh->elementsIterator( thePredicate->GetType() );
- if ( elemIt ) {
- while ( elemIt->more() ) {
- const SMDS_MeshElement* anElem = elemIt->next();
- long anId = anElem->GetID();
- if ( thePredicate->IsSatisfy( anId ) )
- theSequence.push_back( anId );
+ if ( !theElements )
+ theElements = theMesh->elementsIterator( thePredicate->GetType() );
+
+ if ( theElements ) {
+ while ( theElements->more() ) {
+ const SMDS_MeshElement* anElem = theElements->next();
+ if ( thePredicate->GetType() == SMDSAbs_All ||
+ thePredicate->GetType() == anElem->GetType() )
+ {
+ long anId = anElem->GetID();
+ if ( thePredicate->IsSatisfy( anId ) )
+ theSequence.push_back( anId );
+ }
}
}
}
void Filter::GetElementsId( const SMDS_Mesh* theMesh,
- Filter::TIdSequence& theSequence )
+ Filter::TIdSequence& theSequence,
+ SMDS_ElemIteratorPtr theElements )
{
- GetElementsId(theMesh,myPredicate,theSequence);
+ GetElementsId(theMesh,myPredicate,theSequence,theElements);
}
/*
void ManifoldPart::getFacesByLink( const ManifoldPart::Link& theLink,
ManifoldPart::TVectorOfFacePtr& theFaces ) const
{
- std::set<SMDS_MeshCell *> aSetOfFaces;
+
// take all faces that shared first node
- SMDS_ElemIteratorPtr anItr = theLink.myNode1->facesIterator();
- for ( ; anItr->more(); )
- {
- SMDS_MeshFace* aFace = (SMDS_MeshFace*)anItr->next();
- if ( !aFace )
- continue;
- aSetOfFaces.insert( aFace );
- }
+ SMDS_ElemIteratorPtr anItr = theLink.myNode1->GetInverseElementIterator( SMDSAbs_Face );
+ SMDS_StdIterator< const SMDS_MeshElement*, SMDS_ElemIteratorPtr > faces( anItr ), facesEnd;
+ std::set<const SMDS_MeshElement *> aSetOfFaces( faces, facesEnd );
+
// take all faces that shared second node
- anItr = theLink.myNode2->facesIterator();
+ anItr = theLink.myNode2->GetInverseElementIterator( SMDSAbs_Face );
// find the common part of two sets
for ( ; anItr->more(); )
{
- SMDS_MeshFace* aFace = (SMDS_MeshFace*)anItr->next();
- if ( aSetOfFaces.count( aFace ) )
- theFaces.push_back( aFace );
+ const SMDS_MeshElement* aFace = anItr->next();
+ if ( aSetOfFaces.count( aFace ))
+ theFaces.push_back( (SMDS_MeshFace*) aFace );
}
}
return myGroup ? myGroup->GetType() : SMDSAbs_All;
}
-/*
- ElementsOnSurface
-*/
+//================================================================================
+// ElementsOnSurface
+//================================================================================
ElementsOnSurface::ElementsOnSurface()
{
void ElementsOnSurface::SetTolerance( const double theToler )
{
if ( myToler != theToler )
- myIds.Clear();
- myToler = theToler;
+ {
+ myToler = theToler;
+ process();
+ }
}
double ElementsOnSurface::GetTolerance() const
}
-/*
- ElementsOnShape
-*/
+//================================================================================
+// ElementsOnShape
+//================================================================================
-ElementsOnShape::ElementsOnShape()
- : //myMesh(0),
- myType(SMDSAbs_All),
- myToler(Precision::Confusion()),
- myAllNodesFlag(false)
+namespace {
+ const int theIsCheckedFlag = 0x0000100;
+}
+
+struct ElementsOnShape::Classifier
+{
+ Classifier() { mySolidClfr = 0; myFlags = 0; }
+ ~Classifier();
+ void Init(const TopoDS_Shape& s, double tol, const Bnd_B3d* box = 0 );
+ bool IsOut(const gp_Pnt& p) { return SetChecked( true ), (this->*myIsOutFun)( p ); }
+ TopAbs_ShapeEnum ShapeType() const { return myShape.ShapeType(); }
+ const TopoDS_Shape& Shape() const { return myShape; }
+ const Bnd_B3d* GetBndBox() const { return & myBox; }
+ double Tolerance() const { return myTol; }
+ bool IsChecked() { return myFlags & theIsCheckedFlag; }
+ bool IsSetFlag( int flag ) const { return myFlags & flag; }
+ void SetChecked( bool is ) { is ? SetFlag( theIsCheckedFlag ) : UnsetFlag( theIsCheckedFlag ); }
+ void SetFlag ( int flag ) { myFlags |= flag; }
+ void UnsetFlag( int flag ) { myFlags &= ~flag; }
+
+private:
+ bool isOutOfSolid (const gp_Pnt& p);
+ bool isOutOfBox (const gp_Pnt& p);
+ bool isOutOfFace (const gp_Pnt& p);
+ bool isOutOfEdge (const gp_Pnt& p);
+ bool isOutOfVertex(const gp_Pnt& p);
+ bool isOutOfNone (const gp_Pnt& /*p*/) { return true; }
+ bool isBox (const TopoDS_Shape& s);
+
+ TopoDS_Shape prepareSolid( const TopoDS_Shape& theSolid );
+
+ bool (Classifier::* myIsOutFun)(const gp_Pnt& p);
+ BRepClass3d_SolidClassifier* mySolidClfr; // ptr because of a run-time forbidden copy-constructor
+ Bnd_B3d myBox;
+ GeomAPI_ProjectPointOnSurf myProjFace;
+ GeomAPI_ProjectPointOnCurve myProjEdge;
+ gp_Pnt myVertexXYZ;
+ TopoDS_Shape myShape;
+ double myTol;
+ int myFlags;
+};
+
+struct ElementsOnShape::OctreeClassifier : public SMESH_Octree
+{
+ OctreeClassifier( const std::vector< ElementsOnShape::Classifier* >& classifiers );
+ OctreeClassifier( const OctreeClassifier* otherTree,
+ const std::vector< ElementsOnShape::Classifier >& clsOther,
+ std::vector< ElementsOnShape::Classifier >& cls );
+ void GetClassifiersAtPoint( const gp_XYZ& p,
+ std::vector< ElementsOnShape::Classifier* >& classifiers );
+ size_t GetSize();
+
+protected:
+ OctreeClassifier() {}
+ SMESH_Octree* newChild() const { return new OctreeClassifier; }
+ void buildChildrenData();
+ Bnd_B3d* buildRootBox();
+
+ std::vector< ElementsOnShape::Classifier* > myClassifiers;
+};
+
+
+ElementsOnShape::ElementsOnShape():
+ myOctree(0),
+ myType(SMDSAbs_All),
+ myToler(Precision::Confusion()),
+ myAllNodesFlag(false)
{
}
clearClassifiers();
}
+Predicate* ElementsOnShape::clone() const
+{
+ size_t size = sizeof( *this );
+ if ( myOctree )
+ size += myOctree->GetSize();
+ if ( !myClassifiers.empty() )
+ size += sizeof( myClassifiers[0] ) * myClassifiers.size();
+ if ( !myWorkClassifiers.empty() )
+ size += sizeof( myWorkClassifiers[0] ) * myWorkClassifiers.size();
+ if ( size > 1e+9 ) // 1G
+ {
+#ifdef _DEBUG_
+ std::cout << "Avoid ElementsOnShape::clone(), too large: " << size << " bytes " << std::endl;
+#endif
+ return 0;
+ }
+
+ ElementsOnShape* cln = new ElementsOnShape();
+ cln->SetAllNodes ( myAllNodesFlag );
+ cln->SetTolerance( myToler );
+ cln->SetMesh ( myMeshModifTracer.GetMesh() );
+ cln->myShape = myShape; // avoid creation of myClassifiers
+ cln->SetShape ( myShape, myType );
+ cln->myClassifiers.resize( myClassifiers.size() );
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ cln->myClassifiers[ i ].Init( BRepBuilderAPI_Copy( myClassifiers[ i ].Shape()),
+ myToler, myClassifiers[ i ].GetBndBox() );
+ if ( myOctree ) // copy myOctree
+ {
+ cln->myOctree = new OctreeClassifier( myOctree, myClassifiers, cln->myClassifiers );
+ }
+ return cln;
+}
+
SMDSAbs_ElementType ElementsOnShape::GetType() const
{
return myType;
void ElementsOnShape::SetShape (const TopoDS_Shape& theShape,
const SMDSAbs_ElementType theType)
{
+ bool shapeChanges = ( myShape != theShape );
myType = theType;
myShape = theShape;
if ( myShape.IsNull() ) return;
- TopTools_IndexedMapOfShape shapesMap;
- TopAbs_ShapeEnum shapeTypes[4] = { TopAbs_SOLID, TopAbs_FACE, TopAbs_EDGE, TopAbs_VERTEX };
- TopExp_Explorer sub;
- for ( int i = 0; i < 4; ++i )
+ if ( shapeChanges )
{
- if ( shapesMap.IsEmpty() )
- for ( sub.Init( myShape, shapeTypes[i] ); sub.More(); sub.Next() )
- shapesMap.Add( sub.Current() );
- if ( i > 0 )
- for ( sub.Init( myShape, shapeTypes[i], shapeTypes[i-1] ); sub.More(); sub.Next() )
- shapesMap.Add( sub.Current() );
- }
+ // find most complex shapes
+ TopTools_IndexedMapOfShape shapesMap;
+ TopAbs_ShapeEnum shapeTypes[4] = { TopAbs_SOLID, TopAbs_FACE, TopAbs_EDGE, TopAbs_VERTEX };
+ TopExp_Explorer sub;
+ for ( int i = 0; i < 4; ++i )
+ {
+ if ( shapesMap.IsEmpty() )
+ for ( sub.Init( myShape, shapeTypes[i] ); sub.More(); sub.Next() )
+ shapesMap.Add( sub.Current() );
+ if ( i > 0 )
+ for ( sub.Init( myShape, shapeTypes[i], shapeTypes[i-1] ); sub.More(); sub.Next() )
+ shapesMap.Add( sub.Current() );
+ }
- clearClassifiers();
- myClassifiers.resize( shapesMap.Extent() );
- for ( int i = 0; i < shapesMap.Extent(); ++i )
- myClassifiers[ i ] = new TClassifier( shapesMap( i+1 ), myToler );
+ clearClassifiers();
+ myClassifiers.resize( shapesMap.Extent() );
+ for ( int i = 0; i < shapesMap.Extent(); ++i )
+ myClassifiers[ i ].Init( shapesMap( i+1 ), myToler );
+ }
if ( theType == SMDSAbs_Node )
{
void ElementsOnShape::clearClassifiers()
{
- for ( size_t i = 0; i < myClassifiers.size(); ++i )
- delete myClassifiers[ i ];
+ // for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ // delete myClassifiers[ i ];
myClassifiers.clear();
+
+ delete myOctree;
+ myOctree = 0;
}
-bool ElementsOnShape::IsSatisfy (long elemId)
+bool ElementsOnShape::IsSatisfy( long elemId )
{
- const SMDS_Mesh* mesh = myMeshModifTracer.GetMesh();
- const SMDS_MeshElement* elem =
- ( myType == SMDSAbs_Node ? mesh->FindNode( elemId ) : mesh->FindElement( elemId ));
- if ( !elem || myClassifiers.empty() )
+ if ( myClassifiers.empty() )
+ return false;
+
+ const SMDS_Mesh* mesh = myMeshModifTracer.GetMesh();
+ if ( myType == SMDSAbs_Node )
+ return IsSatisfy( mesh->FindNode( elemId ));
+ return IsSatisfy( mesh->FindElement( elemId ));
+}
+
+bool ElementsOnShape::IsSatisfy (const SMDS_MeshElement* elem)
+{
+ if ( !elem )
return false;
bool isSatisfy = myAllNodesFlag, isNodeOut;
gp_XYZ centerXYZ (0, 0, 0);
- SMDS_ElemIteratorPtr aNodeItr = elem->nodesIterator();
- while (aNodeItr->more() && (isSatisfy == myAllNodesFlag))
+ if ( !myOctree && myClassifiers.size() > 5 )
+ {
+ myWorkClassifiers.resize( myClassifiers.size() );
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ myWorkClassifiers[ i ] = & myClassifiers[ i ];
+ myOctree = new OctreeClassifier( myWorkClassifiers );
+
+ SMESHUtils::FreeVector( myWorkClassifiers );
+ }
+
+ for ( int i = 0, nb = elem->NbNodes(); i < nb && (isSatisfy == myAllNodesFlag); ++i )
{
- SMESH_TNodeXYZ aPnt( aNodeItr->next() );
+ SMESH_TNodeXYZ aPnt( elem->GetNode( i ));
centerXYZ += aPnt;
isNodeOut = true;
if ( !getNodeIsOut( aPnt._node, isNodeOut ))
{
- for ( size_t i = 0; i < myClassifiers.size() && isNodeOut; ++i )
- isNodeOut = myClassifiers[i]->IsOut( aPnt );
+ if ( myOctree )
+ {
+ myWorkClassifiers.clear();
+ myOctree->GetClassifiersAtPoint( aPnt, myWorkClassifiers );
+ for ( size_t i = 0; i < myWorkClassifiers.size(); ++i )
+ myWorkClassifiers[i]->SetChecked( false );
+
+ for ( size_t i = 0; i < myWorkClassifiers.size() && isNodeOut; ++i )
+ if ( !myWorkClassifiers[i]->IsChecked() )
+ isNodeOut = myWorkClassifiers[i]->IsOut( aPnt );
+ }
+ else
+ {
+ for ( size_t i = 0; i < myClassifiers.size() && isNodeOut; ++i )
+ isNodeOut = myClassifiers[i].IsOut( aPnt );
+ }
setNodeIsOut( aPnt._node, isNodeOut );
}
isSatisfy = !isNodeOut;
}
// Check the center point for volumes MantisBug 0020168
- if (isSatisfy &&
- myAllNodesFlag &&
- myClassifiers[0]->ShapeType() == TopAbs_SOLID)
+ if ( isSatisfy &&
+ myAllNodesFlag &&
+ myClassifiers[0].ShapeType() == TopAbs_SOLID )
{
centerXYZ /= elem->NbNodes();
isSatisfy = false;
- for ( size_t i = 0; i < myClassifiers.size() && !isSatisfy; ++i )
- isSatisfy = ! myClassifiers[i]->IsOut( centerXYZ );
+ if ( myOctree )
+ {
+ myWorkClassifiers.clear();
+ myOctree->GetClassifiersAtPoint( centerXYZ, myWorkClassifiers );
+ for ( size_t i = 0; i < myWorkClassifiers.size() && !isSatisfy; ++i )
+ isSatisfy = ! myWorkClassifiers[i]->IsOut( centerXYZ );
+ }
+ else
+ {
+ for ( size_t i = 0; i < myClassifiers.size() && !isSatisfy; ++i )
+ isSatisfy = ! myClassifiers[i].IsOut( centerXYZ );
+ }
}
return isSatisfy;
}
-TopAbs_ShapeEnum ElementsOnShape::TClassifier::ShapeType() const
-{
- return myShape.ShapeType();
-}
+//================================================================================
+/*!
+ * \brief Check and optionally return a satisfying shape
+ */
+//================================================================================
-bool ElementsOnShape::TClassifier::IsOut(const gp_Pnt& p)
+bool ElementsOnShape::IsSatisfy (const SMDS_MeshNode* node,
+ TopoDS_Shape* okShape)
{
- return (this->*myIsOutFun)( p );
+ if ( !node )
+ return false;
+
+ if ( !myOctree && myClassifiers.size() > 5 )
+ {
+ myWorkClassifiers.resize( myClassifiers.size() );
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ myWorkClassifiers[ i ] = & myClassifiers[ i ];
+ myOctree = new OctreeClassifier( myWorkClassifiers );
+ }
+
+ bool isNodeOut = true;
+
+ if ( okShape || !getNodeIsOut( node, isNodeOut ))
+ {
+ SMESH_NodeXYZ aPnt = node;
+ if ( myOctree )
+ {
+ myWorkClassifiers.clear();
+ myOctree->GetClassifiersAtPoint( aPnt, myWorkClassifiers );
+
+ for ( size_t i = 0; i < myWorkClassifiers.size(); ++i )
+ myWorkClassifiers[i]->SetChecked( false );
+
+ for ( size_t i = 0; i < myWorkClassifiers.size(); ++i )
+ if ( !myWorkClassifiers[i]->IsChecked() &&
+ !myWorkClassifiers[i]->IsOut( aPnt ))
+ {
+ isNodeOut = false;
+ if ( okShape )
+ *okShape = myWorkClassifiers[i]->Shape();
+ break;
+ }
+ }
+ else
+ {
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ if ( !myClassifiers[i].IsOut( aPnt ))
+ {
+ isNodeOut = false;
+ if ( okShape )
+ *okShape = myClassifiers[i].Shape();
+ break;
+ }
+ }
+ setNodeIsOut( node, isNodeOut );
+ }
+
+ return !isNodeOut;
}
-void ElementsOnShape::TClassifier::Init (const TopoDS_Shape& theShape, double theTol)
+void ElementsOnShape::Classifier::Init( const TopoDS_Shape& theShape,
+ double theTol,
+ const Bnd_B3d* theBox )
{
myShape = theShape;
myTol = theTol;
+ myFlags = 0;
+
+ bool isShapeBox = false;
switch ( myShape.ShapeType() )
{
- case TopAbs_SOLID: {
- if ( isBox( theShape ))
+ case TopAbs_SOLID:
+ {
+ if (( isShapeBox = isBox( theShape )))
{
- myIsOutFun = & ElementsOnShape::TClassifier::isOutOfBox;
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfBox;
}
else
{
- mySolidClfr.Load(theShape);
- myIsOutFun = & ElementsOnShape::TClassifier::isOutOfSolid;
+ mySolidClfr = new BRepClass3d_SolidClassifier( prepareSolid( theShape ));
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfSolid;
}
break;
}
- case TopAbs_FACE: {
+ case TopAbs_FACE:
+ {
Standard_Real u1,u2,v1,v2;
Handle(Geom_Surface) surf = BRep_Tool::Surface( TopoDS::Face( theShape ));
- surf->Bounds( u1,u2,v1,v2 );
- myProjFace.Init(surf, u1,u2, v1,v2, myTol );
- myIsOutFun = & ElementsOnShape::TClassifier::isOutOfFace;
+ if ( surf.IsNull() )
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfNone;
+ else
+ {
+ surf->Bounds( u1,u2,v1,v2 );
+ myProjFace.Init(surf, u1,u2, v1,v2, myTol );
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfFace;
+ }
break;
}
- case TopAbs_EDGE: {
+ case TopAbs_EDGE:
+ {
Standard_Real u1, u2;
- Handle(Geom_Curve) curve = BRep_Tool::Curve( TopoDS::Edge(theShape), u1, u2);
- myProjEdge.Init(curve, u1, u2);
- myIsOutFun = & ElementsOnShape::TClassifier::isOutOfEdge;
+ Handle(Geom_Curve) curve = BRep_Tool::Curve( TopoDS::Edge( theShape ), u1, u2);
+ if ( curve.IsNull() )
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfNone;
+ else
+ {
+ myProjEdge.Init(curve, u1, u2);
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfEdge;
+ }
break;
}
- case TopAbs_VERTEX:{
+ case TopAbs_VERTEX:
+ {
myVertexXYZ = BRep_Tool::Pnt( TopoDS::Vertex( theShape ) );
- myIsOutFun = & ElementsOnShape::TClassifier::isOutOfVertex;
+ myIsOutFun = & ElementsOnShape::Classifier::isOutOfVertex;
break;
}
default:
- throw SALOME_Exception("Programmer error in usage of ElementsOnShape::TClassifier");
+ throw SALOME_Exception("Programmer error in usage of ElementsOnShape::Classifier");
+ }
+
+ if ( !isShapeBox )
+ {
+ if ( theBox )
+ {
+ myBox = *theBox;
+ }
+ else
+ {
+ Bnd_Box box;
+ if ( myShape.ShapeType() == TopAbs_FACE )
+ {
+ BRepAdaptor_Surface SA( TopoDS::Face( myShape ), /*useBoundaries=*/false );
+ if ( SA.GetType() == GeomAbs_BSplineSurface )
+ BRepBndLib::AddOptimal( myShape, box,
+ /*useTriangulation=*/true, /*useShapeTolerance=*/true );
+ }
+ if ( box.IsVoid() )
+ BRepBndLib::Add( myShape, box );
+ myBox.Clear();
+ myBox.Add( box.CornerMin() );
+ myBox.Add( box.CornerMax() );
+ gp_XYZ halfSize = 0.5 * ( box.CornerMax().XYZ() - box.CornerMin().XYZ() );
+ for ( int iDim = 1; iDim <= 3; ++iDim )
+ {
+ double x = halfSize.Coord( iDim );
+ halfSize.SetCoord( iDim, x + Max( myTol, 1e-2 * x ));
+ }
+ myBox.SetHSize( halfSize );
+ }
}
}
-bool ElementsOnShape::TClassifier::isOutOfSolid (const gp_Pnt& p)
+ElementsOnShape::Classifier::~Classifier()
+{
+ delete mySolidClfr; mySolidClfr = 0;
+}
+
+TopoDS_Shape ElementsOnShape::Classifier::prepareSolid( const TopoDS_Shape& theSolid )
{
- mySolidClfr.Perform( p, myTol );
- return ( mySolidClfr.State() != TopAbs_IN && mySolidClfr.State() != TopAbs_ON );
+ // try to limit tolerance of theSolid down to myTol (issue #19026)
+
+ // check if tolerance of theSolid is more than myTol
+ bool tolIsOk = true; // max tolerance is at VERTEXes
+ for ( TopExp_Explorer exp( theSolid, TopAbs_VERTEX ); exp.More() && tolIsOk; exp.Next() )
+ tolIsOk = ( myTol >= BRep_Tool::Tolerance( TopoDS::Vertex( exp.Current() )));
+ if ( tolIsOk )
+ return theSolid;
+
+ // make a copy to prevent the original shape from changes
+ TopoDS_Shape resultShape = BRepBuilderAPI_Copy( theSolid );
+
+ if ( !GEOMUtils::FixShapeTolerance( resultShape, TopAbs_SHAPE, myTol ))
+ return theSolid;
+ return resultShape;
+}
+
+bool ElementsOnShape::Classifier::isOutOfSolid( const gp_Pnt& p )
+{
+ if ( isOutOfBox( p )) return true;
+ mySolidClfr->Perform( p, myTol );
+ return ( mySolidClfr->State() != TopAbs_IN && mySolidClfr->State() != TopAbs_ON );
}
-bool ElementsOnShape::TClassifier::isOutOfBox (const gp_Pnt& p)
+bool ElementsOnShape::Classifier::isOutOfBox( const gp_Pnt& p )
{
return myBox.IsOut( p.XYZ() );
}
-bool ElementsOnShape::TClassifier::isOutOfFace (const gp_Pnt& p)
+bool ElementsOnShape::Classifier::isOutOfFace( const gp_Pnt& p )
{
+ if ( isOutOfBox( p )) return true;
myProjFace.Perform( p );
if ( myProjFace.IsDone() && myProjFace.LowerDistance() <= myTol )
{
// check relatively to the face
- Quantity_Parameter u, v;
+ Standard_Real u, v;
myProjFace.LowerDistanceParameters(u, v);
gp_Pnt2d aProjPnt (u, v);
BRepClass_FaceClassifier aClsf ( TopoDS::Face( myShape ), aProjPnt, myTol );
return true;
}
-bool ElementsOnShape::TClassifier::isOutOfEdge (const gp_Pnt& p)
+bool ElementsOnShape::Classifier::isOutOfEdge( const gp_Pnt& p )
{
+ if ( isOutOfBox( p )) return true;
myProjEdge.Perform( p );
return ! ( myProjEdge.NbPoints() > 0 && myProjEdge.LowerDistance() <= myTol );
}
-bool ElementsOnShape::TClassifier::isOutOfVertex(const gp_Pnt& p)
+bool ElementsOnShape::Classifier::isOutOfVertex( const gp_Pnt& p )
{
return ( myVertexXYZ.Distance( p ) > myTol );
}
-bool ElementsOnShape::TClassifier::isBox (const TopoDS_Shape& theShape)
+bool ElementsOnShape::Classifier::isBox(const TopoDS_Shape& theShape )
{
TopTools_IndexedMapOfShape vMap;
TopExp::MapShapes( theShape, TopAbs_VERTEX, vMap );
return true;
}
+ElementsOnShape::
+OctreeClassifier::OctreeClassifier( const std::vector< ElementsOnShape::Classifier* >& classifiers )
+ :SMESH_Octree( new SMESH_TreeLimit )
+{
+ myClassifiers = classifiers;
+ compute();
+}
+
+ElementsOnShape::
+OctreeClassifier::OctreeClassifier( const OctreeClassifier* otherTree,
+ const std::vector< ElementsOnShape::Classifier >& clsOther,
+ std::vector< ElementsOnShape::Classifier >& cls )
+ :SMESH_Octree( new SMESH_TreeLimit )
+{
+ myBox = new Bnd_B3d( *otherTree->getBox() );
+
+ if (( myIsLeaf = otherTree->isLeaf() ))
+ {
+ myClassifiers.resize( otherTree->myClassifiers.size() );
+ for ( size_t i = 0; i < otherTree->myClassifiers.size(); ++i )
+ {
+ int ind = otherTree->myClassifiers[i] - & clsOther[0];
+ myClassifiers[ i ] = & cls[ ind ];
+ }
+ }
+ else if ( otherTree->myChildren )
+ {
+ myChildren = new SMESH_Tree< Bnd_B3d, 8 > * [ 8 ];
+ for ( int i = 0; i < nbChildren(); i++ )
+ myChildren[i] =
+ new OctreeClassifier( static_cast<const OctreeClassifier*>( otherTree->myChildren[i]),
+ clsOther, cls );
+ }
+}
+
+void ElementsOnShape::
+OctreeClassifier::GetClassifiersAtPoint( const gp_XYZ& point,
+ std::vector< ElementsOnShape::Classifier* >& result )
+{
+ if ( getBox()->IsOut( point ))
+ return;
+
+ if ( isLeaf() )
+ {
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ if ( !myClassifiers[i]->GetBndBox()->IsOut( point ))
+ result.push_back( myClassifiers[i] );
+ }
+ else
+ {
+ for (int i = 0; i < nbChildren(); i++)
+ ((OctreeClassifier*) myChildren[i])->GetClassifiersAtPoint( point, result );
+ }
+}
+
+size_t ElementsOnShape::OctreeClassifier::GetSize()
+{
+ size_t res = sizeof( *this );
+ if ( !myClassifiers.empty() )
+ res += sizeof( myClassifiers[0] ) * myClassifiers.size();
+
+ if ( !isLeaf() )
+ for (int i = 0; i < nbChildren(); i++)
+ res += ((OctreeClassifier*) myChildren[i])->GetSize();
+
+ return res;
+}
+
+void ElementsOnShape::OctreeClassifier::buildChildrenData()
+{
+ // distribute myClassifiers among myChildren
+
+ const int childFlag[8] = { 0x0000001,
+ 0x0000002,
+ 0x0000004,
+ 0x0000008,
+ 0x0000010,
+ 0x0000020,
+ 0x0000040,
+ 0x0000080 };
+ int nbInChild[8] = { 0, 0, 0, 0, 0, 0, 0, 0 };
+
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ {
+ for ( int j = 0; j < nbChildren(); j++ )
+ {
+ if ( !myClassifiers[i]->GetBndBox()->IsOut( *myChildren[j]->getBox() ))
+ {
+ myClassifiers[i]->SetFlag( childFlag[ j ]);
+ ++nbInChild[ j ];
+ }
+ }
+ }
+
+ for ( int j = 0; j < nbChildren(); j++ )
+ {
+ OctreeClassifier* child = static_cast<OctreeClassifier*>( myChildren[ j ]);
+ child->myClassifiers.resize( nbInChild[ j ]);
+ for ( size_t i = 0; nbInChild[ j ] && i < myClassifiers.size(); ++i )
+ {
+ if ( myClassifiers[ i ]->IsSetFlag( childFlag[ j ]))
+ {
+ --nbInChild[ j ];
+ child->myClassifiers[ nbInChild[ j ]] = myClassifiers[ i ];
+ myClassifiers[ i ]->UnsetFlag( childFlag[ j ]);
+ }
+ }
+ }
+ SMESHUtils::FreeVector( myClassifiers );
+
+ // define if a child isLeaf()
+ for ( int i = 0; i < nbChildren(); i++ )
+ {
+ OctreeClassifier* child = static_cast<OctreeClassifier*>( myChildren[ i ]);
+ child->myIsLeaf = ( child->myClassifiers.size() <= 5 ||
+ child->maxSize() < child->myClassifiers[0]->Tolerance() );
+ }
+}
+
+Bnd_B3d* ElementsOnShape::OctreeClassifier::buildRootBox()
+{
+ Bnd_B3d* box = new Bnd_B3d;
+ for ( size_t i = 0; i < myClassifiers.size(); ++i )
+ box->Add( *myClassifiers[i]->GetBndBox() );
+ return box;
+}
/*
Class : BelongToGeom
BelongToGeom::BelongToGeom()
: myMeshDS(NULL),
- myType(SMDSAbs_All),
+ myType(SMDSAbs_NbElementTypes),
myIsSubshape(false),
myTolerance(Precision::Confusion())
{}
+Predicate* BelongToGeom::clone() const
+{
+ BelongToGeom* cln = 0;
+ if ( myElementsOnShapePtr )
+ if ( ElementsOnShape* eos = static_cast<ElementsOnShape*>( myElementsOnShapePtr->clone() ))
+ {
+ cln = new BelongToGeom( *this );
+ cln->myElementsOnShapePtr.reset( eos );
+ }
+ return cln;
+}
+
void BelongToGeom::SetMesh( const SMDS_Mesh* theMesh )
{
- myMeshDS = dynamic_cast<const SMESHDS_Mesh*>(theMesh);
- init();
+ if ( myMeshDS != theMesh )
+ {
+ myMeshDS = dynamic_cast<const SMESHDS_Mesh*>(theMesh);
+ init();
+ }
+ if ( myElementsOnShapePtr )
+ myElementsOnShapePtr->SetMesh( myMeshDS );
}
void BelongToGeom::SetGeom( const TopoDS_Shape& theShape )
{
- myShape = theShape;
- init();
+ if ( myShape != theShape )
+ {
+ myShape = theShape;
+ init();
+ }
}
static bool IsSubShape (const TopTools_IndexedMapOfShape& theMap,
- const TopoDS_Shape& theShape)
+ const TopoDS_Shape& theShape)
{
if (theMap.Contains(theShape)) return true;
void BelongToGeom::init()
{
- if (!myMeshDS || myShape.IsNull()) return;
+ if ( !myMeshDS || myShape.IsNull() ) return;
// is sub-shape of main shape?
TopoDS_Shape aMainShape = myMeshDS->ShapeToMesh();
}
else {
TopTools_IndexedMapOfShape aMap;
- TopExp::MapShapes(aMainShape, aMap);
- myIsSubshape = IsSubShape(aMap, myShape);
+ TopExp::MapShapes( aMainShape, aMap );
+ myIsSubshape = IsSubShape( aMap, myShape );
+ if ( myIsSubshape )
+ {
+ aMap.Clear();
+ TopExp::MapShapes( myShape, aMap );
+ mySubShapesIDs.Clear();
+ for ( int i = 1; i <= aMap.Extent(); ++i )
+ {
+ int subID = myMeshDS->ShapeToIndex( aMap( i ));
+ if ( subID > 0 )
+ mySubShapesIDs.Add( subID );
+ }
+ }
}
//if (!myIsSubshape) // to be always ready to check an element not bound to geometry
{
- myElementsOnShapePtr.reset(new ElementsOnShape());
- myElementsOnShapePtr->SetTolerance(myTolerance);
- myElementsOnShapePtr->SetAllNodes(true); // "belong", while false means "lays on"
- myElementsOnShapePtr->SetMesh(myMeshDS);
- myElementsOnShapePtr->SetShape(myShape, myType);
- }
-}
-
-static bool IsContains( const SMESHDS_Mesh* theMeshDS,
- const TopoDS_Shape& theShape,
- const SMDS_MeshElement* theElem,
- TopAbs_ShapeEnum theFindShapeEnum,
- TopAbs_ShapeEnum theAvoidShapeEnum = TopAbs_SHAPE )
-{
- TopExp_Explorer anExp( theShape,theFindShapeEnum,theAvoidShapeEnum );
-
- while( anExp.More() )
- {
- const TopoDS_Shape& aShape = anExp.Current();
- if( SMESHDS_SubMesh* aSubMesh = theMeshDS->MeshElements( aShape ) ){
- if( aSubMesh->Contains( theElem ) )
- return true;
- }
- anExp.Next();
+ if ( !myElementsOnShapePtr )
+ myElementsOnShapePtr.reset( new ElementsOnShape() );
+ myElementsOnShapePtr->SetTolerance( myTolerance );
+ myElementsOnShapePtr->SetAllNodes( true ); // "belong", while false means "lays on"
+ myElementsOnShapePtr->SetMesh( myMeshDS );
+ myElementsOnShapePtr->SetShape( myShape, myType );
}
- return false;
}
bool BelongToGeom::IsSatisfy (long theId)
return myElementsOnShapePtr->IsSatisfy(theId);
}
- // Case of submesh
+ // Case of sub-mesh
+
if (myType == SMDSAbs_Node)
{
- if( const SMDS_MeshNode* aNode = myMeshDS->FindNode( theId ) )
+ if ( const SMDS_MeshNode* aNode = myMeshDS->FindNode( theId ))
{
if ( aNode->getshapeId() < 1 )
return myElementsOnShapePtr->IsSatisfy(theId);
-
- const SMDS_PositionPtr& aPosition = aNode->GetPosition();
- SMDS_TypeOfPosition aTypeOfPosition = aPosition->GetTypeOfPosition();
- switch( aTypeOfPosition )
- {
- case SMDS_TOP_VERTEX : return ( IsContains( myMeshDS,myShape,aNode,TopAbs_VERTEX ));
- case SMDS_TOP_EDGE : return ( IsContains( myMeshDS,myShape,aNode,TopAbs_EDGE ));
- case SMDS_TOP_FACE : return ( IsContains( myMeshDS,myShape,aNode,TopAbs_FACE ));
- case SMDS_TOP_3DSPACE: return ( IsContains( myMeshDS,myShape,aNode,TopAbs_SOLID ) ||
- IsContains( myMeshDS,myShape,aNode,TopAbs_SHELL ));
- default:;
- }
+ else
+ return mySubShapesIDs.Contains( aNode->getshapeId() );
}
}
else
{
if ( const SMDS_MeshElement* anElem = myMeshDS->FindElement( theId ))
{
- if ( anElem->getshapeId() < 1 )
- return myElementsOnShapePtr->IsSatisfy(theId);
-
- if( myType == SMDSAbs_All )
- {
- return ( IsContains( myMeshDS,myShape,anElem,TopAbs_EDGE ) ||
- IsContains( myMeshDS,myShape,anElem,TopAbs_FACE ) ||
- IsContains( myMeshDS,myShape,anElem,TopAbs_SOLID )||
- IsContains( myMeshDS,myShape,anElem,TopAbs_SHELL ));
- }
- else if( myType == anElem->GetType() )
+ if ( myType == SMDSAbs_All || anElem->GetType() == myType )
{
- switch( myType )
- {
- case SMDSAbs_Edge : return ( IsContains( myMeshDS,myShape,anElem,TopAbs_EDGE ));
- case SMDSAbs_Face : return ( IsContains( myMeshDS,myShape,anElem,TopAbs_FACE ));
- case SMDSAbs_Volume: return ( IsContains( myMeshDS,myShape,anElem,TopAbs_SOLID )||
- IsContains( myMeshDS,myShape,anElem,TopAbs_SHELL ));
- default:;
- }
+ if ( anElem->getshapeId() < 1 )
+ return myElementsOnShapePtr->IsSatisfy(theId);
+ else
+ return mySubShapesIDs.Contains( anElem->getshapeId() );
}
}
}
void BelongToGeom::SetType (SMDSAbs_ElementType theType)
{
- myType = theType;
- init();
+ if ( myType != theType )
+ {
+ myType = theType;
+ init();
+ }
}
SMDSAbs_ElementType BelongToGeom::GetType() const
void BelongToGeom::SetTolerance (double theTolerance)
{
myTolerance = theTolerance;
- if (!myIsSubshape)
- init();
+ init();
}
double BelongToGeom::GetTolerance()
/*
Class : LyingOnGeom
Description : Predicate for verifying whether entiy lying or partially lying on
- specified geometrical support
+ specified geometrical support
*/
LyingOnGeom::LyingOnGeom()
: myMeshDS(NULL),
- myType(SMDSAbs_All),
+ myType(SMDSAbs_NbElementTypes),
myIsSubshape(false),
myTolerance(Precision::Confusion())
{}
+Predicate* LyingOnGeom::clone() const
+{
+ LyingOnGeom* cln = 0;
+ if ( myElementsOnShapePtr )
+ if ( ElementsOnShape* eos = static_cast<ElementsOnShape*>( myElementsOnShapePtr->clone() ))
+ {
+ cln = new LyingOnGeom( *this );
+ cln->myElementsOnShapePtr.reset( eos );
+ }
+ return cln;
+}
+
void LyingOnGeom::SetMesh( const SMDS_Mesh* theMesh )
{
- myMeshDS = dynamic_cast<const SMESHDS_Mesh*>(theMesh);
- init();
+ if ( myMeshDS != theMesh )
+ {
+ myMeshDS = dynamic_cast<const SMESHDS_Mesh*>(theMesh);
+ init();
+ }
+ if ( myElementsOnShapePtr )
+ myElementsOnShapePtr->SetMesh( myMeshDS );
}
void LyingOnGeom::SetGeom( const TopoDS_Shape& theShape )
{
- myShape = theShape;
- init();
+ if ( myShape != theShape )
+ {
+ myShape = theShape;
+ init();
+ }
}
void LyingOnGeom::init()
mySubShapesIDs.Add( subID );
}
}
- else
+ // else // to be always ready to check an element not bound to geometry
{
- myElementsOnShapePtr.reset(new ElementsOnShape());
- myElementsOnShapePtr->SetTolerance(myTolerance);
- myElementsOnShapePtr->SetAllNodes(false); // lays on, while true means "belong"
- myElementsOnShapePtr->SetMesh(myMeshDS);
- myElementsOnShapePtr->SetShape(myShape, myType);
+ if ( !myElementsOnShapePtr )
+ myElementsOnShapePtr.reset( new ElementsOnShape() );
+ myElementsOnShapePtr->SetTolerance( myTolerance );
+ myElementsOnShapePtr->SetAllNodes( false ); // lays on, while true means "belong"
+ myElementsOnShapePtr->SetMesh( myMeshDS );
+ myElementsOnShapePtr->SetShape( myShape, myType );
}
}
if ( mySubShapesIDs.Contains( elem->getshapeId() ))
return true;
- if ( elem->GetType() != SMDSAbs_Node )
+ if (( elem->GetType() != SMDSAbs_Node ) &&
+ ( myType == SMDSAbs_All || elem->GetType() == myType ))
{
SMDS_ElemIteratorPtr nodeItr = elem->nodesIterator();
while ( nodeItr->more() )
void LyingOnGeom::SetType( SMDSAbs_ElementType theType )
{
- myType = theType;
- init();
+ if ( myType != theType )
+ {
+ myType = theType;
+ init();
+ }
}
SMDSAbs_ElementType LyingOnGeom::GetType() const
void LyingOnGeom::SetTolerance (double theTolerance)
{
myTolerance = theTolerance;
- if (!myIsSubshape)
- init();
+ init();
}
double LyingOnGeom::GetTolerance()
return myTolerance;
}
-bool LyingOnGeom::Contains( const SMESHDS_Mesh* theMeshDS,
- const TopoDS_Shape& theShape,
- const SMDS_MeshElement* theElem,
- TopAbs_ShapeEnum theFindShapeEnum,
- TopAbs_ShapeEnum theAvoidShapeEnum )
-{
- // if (IsContains(theMeshDS, theShape, theElem, theFindShapeEnum, theAvoidShapeEnum))
- // return true;
-
- // TopTools_MapOfShape aSubShapes;
- // TopExp_Explorer exp( theShape, theFindShapeEnum, theAvoidShapeEnum );
- // for ( ; exp.More(); exp.Next() )
- // {
- // const TopoDS_Shape& aShape = exp.Current();
- // if ( !aSubShapes.Add( aShape )) continue;
-
- // if ( SMESHDS_SubMesh* aSubMesh = theMeshDS->MeshElements( aShape ))
- // {
- // if ( aSubMesh->Contains( theElem ))
- // return true;
-
- // SMDS_ElemIteratorPtr nodeItr = theElem->nodesIterator();
- // while ( nodeItr->more() )
- // {
- // const SMDS_MeshElement* aNode = nodeItr->next();
- // if ( aSubMesh->Contains( aNode ))
- // return true;
- // }
- // }
- // }
- return false;
-}
-
TSequenceOfXYZ::TSequenceOfXYZ(): myElem(0)
{}