X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FControls%2FSMESH_Controls.cxx;h=2fe0e08891904acfef1deea8223a59c85781cad5;hb=61bbd63ee9b5f8b85237b8a45c8f454b766c25ca;hp=42b5b6a1a366b3fea541a7ad325aa64194d89e5a;hpb=2de294b09ac8b9ace071a01db9cb4e235f1eadbb;p=modules%2Fsmesh.git diff --git a/src/Controls/SMESH_Controls.cxx b/src/Controls/SMESH_Controls.cxx index 42b5b6a1a..2fe0e0889 100644 --- a/src/Controls/SMESH_Controls.cxx +++ b/src/Controls/SMESH_Controls.cxx @@ -1,52 +1,57 @@ -// Copyright (C) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2012 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 +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS // -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. // -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. // -// You should have received a copy of the GNU Lesser General Public -// License along with this library; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA // -// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com // #include "SMESH_ControlsDef.hxx" -#include -#include +#include "SMDS_BallElement.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_Mesh.hxx" +#include "SMESH_OctreeNode.hxx" #include #include #include - -#include -#include -#include -#include -#include -#include -#include - #include #include #include - #include #include #include #include #include - +#include +#include +#include +#include +#include +#include +#include #include #include #include @@ -55,16 +60,11 @@ #include #include -#include "SMDS_Mesh.hxx" -#include "SMDS_Iterator.hxx" -#include "SMDS_MeshElement.hxx" -#include "SMDS_MeshNode.hxx" -#include "SMDS_VolumeTool.hxx" -#include "SMDS_QuadraticFaceOfNodes.hxx" -#include "SMDS_QuadraticEdge.hxx" +#include + +#include +#include -#include "SMESHDS_Mesh.hxx" -#include "SMESHDS_GroupBase.hxx" /* AUXILIARY METHODS @@ -284,24 +284,25 @@ long NumericalFunctor::GetPrecision() const void NumericalFunctor::SetPrecision( const long thePrecision ) { myPrecision = thePrecision; + myPrecisionValue = pow( 10., (double)( myPrecision ) ); } double NumericalFunctor::GetValue( long theId ) { + double aVal = 0; + myCurrElement = myMesh->FindElement( theId ); + TSequenceOfXYZ P; if ( GetPoints( theId, P )) - { - double aVal = GetValue( P ); - if ( myPrecision >= 0 ) - { - double prec = pow( 10., (double)( myPrecision ) ); - aVal = floor( aVal * prec + 0.5 ) / prec; - } - return aVal; - } + aVal = Round( GetValue( P )); - return 0.; + return aVal; +} + +double NumericalFunctor::Round( const double & aVal ) +{ + return ( myPrecision >= 0 ) ? floor( aVal * myPrecisionValue + 0.5 ) / myPrecisionValue : aVal; } //================================================================================ @@ -463,7 +464,7 @@ double MaxElementLength2D::GetValue( long theElementId ) aVal = Max(L1,Max(L2,L3)); break; } - else if( len == 8 ) { // quadratic quadrangles + else if( len == 8 || len == 9 ) { // 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 )); @@ -570,6 +571,12 @@ double MaxElementLength3D::GetValue( long theElementId ) aVal = Max(aVal,Max(Max(D1,D2),Max(D3,D4))); break; } + else if( len == 12 ) { // hexagonal prism + for ( int i1 = 1; i1 < 12; ++i1 ) + for ( int i2 = i1+1; i1 <= 12; ++i1 ) + aVal = Max( aVal, getDistance(P( i1 ),P( i2 ))); + break; + } else if( len == 10 ) { // quadratic tetras double L1 = getDistance(P( 1 ),P( 5 )) + getDistance(P( 5 ),P( 2 )); double L2 = getDistance(P( 2 ),P( 6 )) + getDistance(P( 6 ),P( 3 )); @@ -607,7 +614,7 @@ double MaxElementLength3D::GetValue( long theElementId ) aVal = Max(aVal,Max(Max(L7,L8),L9)); break; } - else if( len == 20 ) { // quadratic hexas + else if( len == 20 || len == 27 ) { // quadratic hexas 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 )); @@ -684,7 +691,7 @@ double MinimumAngle::GetValue( const TSequenceOfXYZ& P ) aMin = Min(aMin,A0); } - return aMin * 180.0 / PI; + return aMin * 180.0 / M_PI; } double MinimumAngle::GetBadRate( double Value, int nbNodes ) const @@ -704,6 +711,26 @@ SMDSAbs_ElementType MinimumAngle::GetType() const Class : AspectRatio Description : Functor for calculating aspect ratio */ +double AspectRatio::GetValue( long theId ) +{ + 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 )); + } + return aVal; +} + double AspectRatio::GetValue( const TSequenceOfXYZ& P ) { // According to "Mesh quality control" by Nadir Bouhamau referring to @@ -802,7 +829,7 @@ double AspectRatio::GetValue( const TSequenceOfXYZ& P ) return 0.; return alpha * L * C1 / C2; } - else if( nbNodes == 8 ){ // nbNodes==8 - quadratic quadrangle + else if( nbNodes == 8 || nbNodes == 9 ) { // nbNodes==8 - quadratic quadrangle // Compute lengths of the sides std::vector< double > aLen (4); aLen[0] = getDistance( P(1), P(3) ); @@ -852,9 +879,10 @@ double AspectRatio::GetValue( const TSequenceOfXYZ& P ) double AspectRatio::GetBadRate( double Value, int /*nbNodes*/ ) const { // the aspect ratio is in the range [1.0,infinity] + // < 1.0 = very bad, zero area // 1.0 = good // infinity = bad - return Value / 1000.; + return ( Value < 0.9 ) ? 1000 : Value / 1000.; } SMDSAbs_ElementType AspectRatio::GetType() const @@ -929,6 +957,28 @@ namespace{ } +double AspectRatio3D::GetValue( long theId ) +{ + double aVal = 0; + myCurrElement = myMesh->FindElement( theId ); + if ( myCurrElement && myCurrElement->GetVtkType() == VTK_TETRA ) + { + // 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() )) + aVal = Round( vtkMeshQuality::TetAspectRatio( avtkCell )); + } + else + { + TSequenceOfXYZ P; + if ( GetPoints( myCurrElement, P )) + aVal = Round( GetValue( P )); + } + return aVal; +} + double AspectRatio3D::GetValue( const TSequenceOfXYZ& P ) { double aQuality = 0.0; @@ -941,10 +991,11 @@ double AspectRatio3D::GetValue( const TSequenceOfXYZ& P ) 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 return aQuality; } - switch(nbNodes){ + switch(nbNodes) { case 4:{ double aLen[6] = { getDistance(P( 1 ),P( 2 )), // a @@ -1162,7 +1213,22 @@ double AspectRatio3D::GetValue( const TSequenceOfXYZ& P ) } break; } - } + 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); + } + { + gp_XYZ aXYZ[8] = {P( 2 ),P( 3 ),P( 5 ),P( 6 ),P( 8 ),P( 9 ),P( 11 ),P( 12 )}; + aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[8])),aQuality); + } + { + gp_XYZ aXYZ[8] = {P( 3 ),P( 4 ),P( 6 ),P( 1 ),P( 9 ),P( 10 ),P( 12 ),P( 7 )}; + aQuality = std::max(GetValue(TSequenceOfXYZ(&aXYZ[0],&aXYZ[8])),aQuality); + } + break; + } // switch(nbNodes) + if ( nbNodes > 4 ) { // avaluate aspect ratio of quadranle faces AspectRatio aspect2D; @@ -1230,12 +1296,12 @@ double Warping::ComputeA( const gp_XYZ& thePnt1, gp_XYZ N = GI.Crossed( GJ ); if ( N.Modulus() < gp::Resolution() ) - return PI / 2; + return M_PI / 2; N.Normalize(); double H = ( thePnt2 - theG ).Dot( N ); - return asin( fabs( H / L ) ) * 180. / PI; + return asin( fabs( H / L ) ) * 180. / M_PI; } double Warping::GetBadRate( double Value, int /*nbNodes*/ ) const @@ -1314,14 +1380,14 @@ double Skew::GetValue( const TSequenceOfXYZ& P ) return 0.; // Compute skew - static double PI2 = PI / 2.; + static double PI2 = M_PI / 2.; if ( P.size() == 3 ) { double A0 = fabs( PI2 - skewAngle( P( 3 ), P( 1 ), P( 2 ) ) ); double A1 = fabs( PI2 - skewAngle( P( 1 ), P( 2 ), P( 3 ) ) ); double A2 = fabs( PI2 - skewAngle( P( 2 ), P( 3 ), P( 1 ) ) ); - return Max( A0, Max( A1, A2 ) ) * 180. / PI; + return Max( A0, Max( A1, A2 ) ) * 180. / M_PI; } else { @@ -1338,7 +1404,7 @@ double Skew::GetValue( const TSequenceOfXYZ& P ) if ( A < Precision::Angular() ) return 0.; - return A * 180. / PI; + return A * 180. / M_PI; } } @@ -1392,7 +1458,7 @@ SMDSAbs_ElementType Area::GetType() const /* Class : Length - Description : Functor for calculating length off edge + Description : Functor for calculating length of edge */ double Length::GetValue( const TSequenceOfXYZ& P ) { @@ -1896,6 +1962,34 @@ void MultiConnection2D::GetValues(MValues& theValues){ } +/* + Class : BallDiameter + Description : Functor returning diameter of a ball element +*/ +double BallDiameter::GetValue( long theId ) +{ + double diameter = 0; + + if ( const SMDS_BallElement* ball = + dynamic_cast( myMesh->FindElement( theId ))) + { + diameter = ball->GetDiameter(); + } + return diameter; +} + +double BallDiameter::GetBadRate( double Value, int /*nbNodes*/ ) const +{ + // meaningless as it is not a quality control functor + return Value; +} + +SMDSAbs_ElementType BallDiameter::GetType() const +{ + return SMDSAbs_Ball; +} + + /* PREDICATES */ @@ -1956,11 +2050,13 @@ bool BareBorderVolume::IsSatisfy(long theElementId ) bool BareBorderFace::IsSatisfy(long theElementId ) { + bool ok = false; if ( const SMDS_MeshElement* face = myMesh->FindElement(theElementId)) + { if ( face->GetType() == SMDSAbs_Face ) { int nbN = face->NbCornerNodes(); - for ( int i = 0; i < nbN; ++i ) + for ( int i = 0; i < nbN && !ok; ++i ) { // check if a link is shared by another face const SMDS_MeshNode* n1 = face->GetNode( i ); @@ -1974,16 +2070,18 @@ bool BareBorderFace::IsSatisfy(long theElementId ) } if ( !isShared ) { - myLinkNodes.resize( 2 + face->IsQuadratic()); + const int iQuad = face->IsQuadratic(); + myLinkNodes.resize( 2 + iQuad); myLinkNodes[0] = n1; myLinkNodes[1] = n2; - if ( face->IsQuadratic() ) + if ( iQuad ) myLinkNodes[2] = face->GetNode( i+nbN ); - return !myMesh->FindElement( myLinkNodes, SMDSAbs_Edge, /*noMedium=*/false); + ok = !myMesh->FindElement( myLinkNodes, SMDSAbs_Edge, /*noMedium=*/false); } } } - return false; + } + return ok; } /* @@ -2039,6 +2137,106 @@ bool OverConstrainedFace::IsSatisfy(long theElementId ) return false; } +/* + Class : CoincidentNodes + Description : Predicate of Coincident nodes +*/ + +CoincidentNodes::CoincidentNodes() +{ + myToler = 1e-5; +} + +bool CoincidentNodes::IsSatisfy( long theElementId ) +{ + return myCoincidentIDs.Contains( theElementId ); +} + +SMDSAbs_ElementType CoincidentNodes::GetType() const +{ + return SMDSAbs_Node; +} + +void CoincidentNodes::SetMesh( const SMDS_Mesh* theMesh ) +{ + myMeshModifTracer.SetMesh( theMesh ); + if ( myMeshModifTracer.IsMeshModified() ) + { + TIDSortedNodeSet nodesToCheck; + SMDS_NodeIteratorPtr nIt = theMesh->nodesIterator(/*idInceasingOrder=*/true); + while ( nIt->more() ) + nodesToCheck.insert( nodesToCheck.end(), nIt->next() ); + + list< list< const SMDS_MeshNode*> > nodeGroups; + SMESH_OctreeNode::FindCoincidentNodes ( nodesToCheck, &nodeGroups, myToler ); + + myCoincidentIDs.Clear(); + list< list< const SMDS_MeshNode*> >::iterator groupIt = nodeGroups.begin(); + for ( ; groupIt != nodeGroups.end(); ++groupIt ) + { + list< const SMDS_MeshNode*>& coincNodes = *groupIt; + list< const SMDS_MeshNode*>::iterator n = coincNodes.begin(); + for ( ; n != coincNodes.end(); ++n ) + myCoincidentIDs.Add( (*n)->GetID() ); + } + } +} + +/* + Class : CoincidentElements + Description : Predicate of Coincident Elements + Note : This class is suitable only for visualization of Coincident Elements +*/ + +CoincidentElements::CoincidentElements() +{ + myMesh = 0; +} + +void CoincidentElements::SetMesh( const SMDS_Mesh* theMesh ) +{ + myMesh = theMesh; +} + +bool CoincidentElements::IsSatisfy( long theElementId ) +{ + if ( !myMesh ) return false; + + if ( const SMDS_MeshElement* e = myMesh->FindElement( theElementId )) + { + if ( e->GetType() != GetType() ) return false; + set< const SMDS_MeshNode* > elemNodes( e->begin_nodes(), e->end_nodes() ); + const int nbNodes = e->NbNodes(); + SMDS_ElemIteratorPtr invIt = (*elemNodes.begin())->GetInverseElementIterator( GetType() ); + while ( invIt->more() ) + { + const SMDS_MeshElement* e2 = invIt->next(); + if ( e2 == e || e2->NbNodes() != nbNodes ) continue; + + bool sameNodes = true; + for ( size_t i = 0; i < elemNodes.size() && sameNodes; ++i ) + sameNodes = ( elemNodes.count( e2->GetNode( i ))); + if ( sameNodes ) + return true; + } + } + return false; +} + +SMDSAbs_ElementType CoincidentElements1D::GetType() const +{ + return SMDSAbs_Edge; +} +SMDSAbs_ElementType CoincidentElements2D::GetType() const +{ + return SMDSAbs_Face; +} +SMDSAbs_ElementType CoincidentElements3D::GetType() const +{ + return SMDSAbs_Volume; +} + + /* Class : FreeBorders Description : Predicate for free borders @@ -2084,17 +2282,13 @@ bool FreeEdges::IsFreeEdge( const SMDS_MeshNode** theNodes, const int theFaceId TColStd_MapOfInteger aMap; for ( int i = 0; i < 2; i++ ) { - SMDS_ElemIteratorPtr anElemIter = theNodes[ i ]->GetInverseElementIterator(); + SMDS_ElemIteratorPtr anElemIter = theNodes[ i ]->GetInverseElementIterator(SMDSAbs_Face); while( anElemIter->more() ) { - const SMDS_MeshElement* anElem = anElemIter->next(); - if ( anElem != 0 && anElem->GetType() == SMDSAbs_Face ) + if ( const SMDS_MeshElement* anElem = anElemIter->next()) { - int anId = anElem->GetID(); - - if ( i == 0 ) - aMap.Add( anId ); - else if ( aMap.Contains( anId ) && anId != theFaceId ) + const int anId = anElem->GetID(); + if ( anId != theFaceId && !aMap.Add( anId )) return false; } } @@ -2119,7 +2313,7 @@ bool FreeEdges::IsSatisfy( long theId ) else { anIter = aFace->nodesIterator(); } - if ( anIter == 0 ) + if ( !anIter ) return false; int i = 0, nbNodes = aFace->NbNodes(); @@ -2197,14 +2391,11 @@ void FreeEdges::GetBoreders(TBorders& theBorders) long anId = aNode->GetID(); Border aBorder(anElemId,aNodeId[1],anId); aNodeId[1] = anId; - //std::cout<