-// Copyright (C) 2007-2021 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
double Warping::GetValue( const TSequenceOfXYZ& P )
{
- if ( P.size() != 4 )
- return 0;
-
- gp_XYZ G = ( P( 1 ) + P( 2 ) + P( 3 ) + P( 4 ) ) / 4.;
-
- double A1 = ComputeA( P( 1 ), P( 2 ), P( 3 ), G );
- double A2 = ComputeA( P( 2 ), P( 3 ), P( 4 ), G );
- double A3 = ComputeA( P( 3 ), P( 4 ), P( 1 ), G );
- double A4 = ComputeA( P( 4 ), P( 1 ), P( 2 ), G );
-
- double val = Max( Max( A1, A2 ), Max( A3, A4 ) );
-
- const double eps = 0.1; // val is in degrees
-
- return val < eps ? 0. : val;
+ return ComputeValue(P);
}
double Warping::ComputeA( const gp_XYZ& thePnt1,
return asin( fabs( H / L ) ) * 180. / M_PI;
}
+double Warping::ComputeValue(const TSequenceOfXYZ& thePoints) const
+{
+ if (thePoints.size() != 4)
+ return 0;
+
+ gp_XYZ G = (thePoints(1) + thePoints(2) + thePoints(3) + thePoints(4)) / 4.;
+
+ double A1 = ComputeA(thePoints(1), thePoints(2), thePoints(3), G);
+ double A2 = ComputeA(thePoints(2), thePoints(3), thePoints(4), G);
+ double A3 = ComputeA(thePoints(3), thePoints(4), thePoints(1), G);
+ double A4 = ComputeA(thePoints(4), thePoints(1), thePoints(2), G);
+
+ double val = Max(Max(A1, A2), Max(A3, A4));
+
+ const double eps = 0.1; // val is in degrees
+
+ return val < eps ? 0. : val;
+}
+
double Warping::GetBadRate( double Value, int /*nbNodes*/ ) const
{
// the warp is in the range [0.0,PI/2]
}
+//================================================================================
+/*
+ Class : Warping3D
+ Description : Functor for calculating warping
+*/
+//================================================================================
+
+bool Warping3D::IsApplicable(const SMDS_MeshElement* element) const
+{
+ return NumericalFunctor::IsApplicable(element);//&& element->NbNodes() == 4;
+}
+
+double Warping3D::GetValue(long theId)
+{
+ double aVal = 0;
+ myCurrElement = myMesh->FindElement(theId);
+ if (myCurrElement)
+ {
+ WValues aValues;
+ ProcessVolumeELement(aValues);
+ for (const auto& aValue: aValues)
+ {
+ aVal = Max(aVal, aValue.myWarp);
+ }
+ }
+ return aVal;
+}
+
+double Warping3D::GetValue(const TSequenceOfXYZ& P)
+{
+ return ComputeValue(P);
+}
+
+SMDSAbs_ElementType Warping3D::GetType() const
+{
+ return SMDSAbs_Volume;
+}
+
+bool Warping3D::Value::operator<(const Warping3D::Value& x) const
+{
+ if (myPntIds.size() != x.myPntIds.size())
+ return myPntIds.size() < x.myPntIds.size();
+
+ for (int anInd = 0; anInd < myPntIds.size(); ++anInd)
+ if (myPntIds[anInd] != x.myPntIds[anInd])
+ return myPntIds[anInd] != x.myPntIds[anInd];
+
+ return false;
+}
+
+// Compute value on each face of volume
+void Warping3D::ProcessVolumeELement(WValues& theValues)
+{
+ SMDS_VolumeTool aVTool(myCurrElement);
+ double aCoord[3];
+ for (int aFaceID = 0; aFaceID < aVTool.NbFaces(); ++aFaceID)
+ {
+ TSequenceOfXYZ aPoints;
+ std::set<const SMDS_MeshNode*> aNodes;
+ std::vector<long> aNodeIds;
+ const SMDS_MeshNode** aNodesPtr = aVTool.GetFaceNodes(aFaceID);
+
+ if (aNodesPtr)
+ {
+ for (int i = 0; i < aVTool.NbFaceNodes(aFaceID); ++i)
+ {
+ aNodesPtr[i]->GetXYZ(aCoord);
+ aPoints.push_back(gp_XYZ{ aCoord[0], aCoord[1], aCoord[2] });
+ aNodeIds.push_back(aNodesPtr[i]->GetID());
+ }
+ double aWarp = GetValue(aPoints);
+ Value aVal{ aWarp, aNodeIds };
+
+ theValues.push_back(aVal);
+ }
+ }
+}
+
+void Warping3D::GetValues(WValues& theValues)
+{
+ for (SMDS_VolumeIteratorPtr anIter = myMesh->volumesIterator(); anIter->more(); )
+ {
+ myCurrElement = anIter->next();
+ ProcessVolumeELement(theValues);
+ }
+}
+
//================================================================================
/*
Class : Taper
return SMDSAbs_Node;
}
+//================================================================================
+/*
+ Class : ScaledJacobian
+ Description : Functor returning the ScaledJacobian for volumetric elements
+*/
+//================================================================================
+
+double ScaledJacobian::GetValue( long theElementId )
+{
+ if ( theElementId && myMesh ) {
+ SMDS_VolumeTool aVolumeTool;
+ if ( aVolumeTool.Set( myMesh->FindElement( theElementId )))
+ return aVolumeTool.GetScaledJacobian();
+ }
+ return 0;
+
+ /*
+ //VTK version not used because lack of implementation for HEXAGONAL_PRISM.
+ //Several mesh quality measures implemented in vtkMeshQuality can be accessed left here as reference
+ double aVal = 0;
+ myCurrElement = myMesh->FindElement( theElementId );
+ if ( myCurrElement )
+ {
+ VTKCellType cellType = myCurrElement->GetVtkType();
+ vtkUnstructuredGrid* grid = const_cast<SMDS_Mesh*>( myMesh )->GetGrid();
+ vtkCell* avtkCell = grid->GetCell( myCurrElement->GetVtkID() );
+ switch ( cellType )
+ {
+ case VTK_QUADRATIC_TETRA:
+ case VTK_TETRA:
+ aVal = Round( vtkMeshQuality::TetScaledJacobian( avtkCell ));
+ break;
+ case VTK_QUADRATIC_HEXAHEDRON:
+ case VTK_HEXAHEDRON:
+ aVal = Round( vtkMeshQuality::HexScaledJacobian( avtkCell ));
+ break;
+ case VTK_QUADRATIC_WEDGE:
+ case VTK_WEDGE: //Pentahedron
+ aVal = Round( vtkMeshQuality::WedgeScaledJacobian( avtkCell ));
+ break;
+ case VTK_QUADRATIC_PYRAMID:
+ case VTK_PYRAMID:
+ aVal = Round( vtkMeshQuality::PyramidScaledJacobian( avtkCell ));
+ break;
+ case VTK_HEXAGONAL_PRISM:
+ case VTK_POLYHEDRON:
+ default:
+ break;
+ }
+ }
+ return aVal;
+ */
+}
+
+double ScaledJacobian::GetBadRate( double Value, int /*nbNodes*/ ) const
+{
+ return Value;
+}
+
+SMDSAbs_ElementType ScaledJacobian::GetType() const
+{
+ return SMDSAbs_Volume;
+}
+
/*
PREDICATES
*/
struct ElementsOnShape::Classifier
{
- Classifier() { mySolidClfr = 0; myFlags = 0; }
+ Classifier(): mySolidClfr(0), myProjFace(0), myProjEdge(0), myFlags(0) { myU = myV = 1e100; }
~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 ); }
void SetChecked( bool is ) { is ? SetFlag( theIsCheckedFlag ) : UnsetFlag( theIsCheckedFlag ); }
void SetFlag ( int flag ) { myFlags |= flag; }
void UnsetFlag( int flag ) { myFlags &= ~flag; }
+ void GetParams( double & u, double & v ) const { u = myU; v = myV; }
private:
bool isOutOfSolid (const gp_Pnt& p);
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
+ BRepClass3d_SolidClassifier* mySolidClfr;
Bnd_B3d myBox;
- GeomAPI_ProjectPointOnSurf myProjFace;
- GeomAPI_ProjectPointOnCurve myProjEdge;
+ GeomAPI_ProjectPointOnSurf* myProjFace;
+ GeomAPI_ProjectPointOnCurve* myProjEdge;
gp_Pnt myVertexXYZ;
TopoDS_Shape myShape;
double myTol;
+ double myU, myV; // result of isOutOfFace() and isOutOfEdge()
int myFlags;
};
size += sizeof( myWorkClassifiers[0] ) * myWorkClassifiers.size();
if ( size > 1e+9 ) // 1G
{
-#ifdef _DEBUG_
+
+ if (SALOME::VerbosityActivated())
std::cout << "Avoid ElementsOnShape::clone(), too large: " << size << " bytes " << std::endl;
-#endif
+
return 0;
}
isNodeOut = false;
if ( okShape )
*okShape = myWorkClassifiers[i]->Shape();
+ myWorkClassifiers[i]->GetParams( myU, myV );
break;
}
}
isNodeOut = false;
if ( okShape )
*okShape = myClassifiers[i].Shape();
+ myClassifiers[i].GetParams( myU, myV );
break;
}
}
else
{
surf->Bounds( u1,u2,v1,v2 );
- myProjFace.Init(surf, u1,u2, v1,v2, myTol );
+ myProjFace = new GeomAPI_ProjectPointOnSurf;
+ myProjFace->Init( surf, u1,u2, v1,v2, myTol );
myIsOutFun = & ElementsOnShape::Classifier::isOutOfFace;
}
break;
myIsOutFun = & ElementsOnShape::Classifier::isOutOfNone;
else
{
- myProjEdge.Init(curve, u1, u2);
+ myProjEdge = new GeomAPI_ProjectPointOnCurve;
+ myProjEdge->Init( curve, u1, u2 );
myIsOutFun = & ElementsOnShape::Classifier::isOutOfEdge;
}
break;
ElementsOnShape::Classifier::~Classifier()
{
delete mySolidClfr; mySolidClfr = 0;
+ delete myProjFace; myProjFace = 0;
+ delete myProjEdge; myProjEdge = 0;
}
TopoDS_Shape ElementsOnShape::Classifier::prepareSolid( const TopoDS_Shape& theSolid )
bool ElementsOnShape::Classifier::isOutOfFace( const gp_Pnt& p )
{
if ( isOutOfBox( p )) return true;
- myProjFace.Perform( p );
- if ( myProjFace.IsDone() && myProjFace.LowerDistance() <= myTol )
+ myProjFace->Perform( p );
+ if ( myProjFace->IsDone() && myProjFace->LowerDistance() <= myTol )
{
// check relatively to the face
- Standard_Real u, v;
- myProjFace.LowerDistanceParameters(u, v);
- gp_Pnt2d aProjPnt (u, v);
+ myProjFace->LowerDistanceParameters( myU, myV );
+ gp_Pnt2d aProjPnt( myU, myV );
BRepClass_FaceClassifier aClsf ( TopoDS::Face( myShape ), aProjPnt, myTol );
if ( aClsf.State() == TopAbs_IN || aClsf.State() == TopAbs_ON )
return false;
bool ElementsOnShape::Classifier::isOutOfEdge( const gp_Pnt& p )
{
if ( isOutOfBox( p )) return true;
- myProjEdge.Perform( p );
- return ! ( myProjEdge.NbPoints() > 0 && myProjEdge.LowerDistance() <= myTol );
+ myProjEdge->Perform( p );
+ bool isOn = ( myProjEdge->NbPoints() > 0 && myProjEdge->LowerDistance() <= myTol );
+ if ( isOn )
+ myU = myProjEdge->LowerDistanceParameter();
+ return !isOn;
}
bool ElementsOnShape::Classifier::isOutOfVertex( const gp_Pnt& p )
