X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_RadialPrism_3D.cxx;h=f2d2970afd48e095143ec341a08a9167a21d68d6;hb=382f2cc4abb2ee8553a911aeb27348e96c39d197;hp=1814680af4c91b05a457b0d91a2ffb6c3fe2f8dc;hpb=c2a18a423efe58a6b74ae9c9cfd811c084bb3952;p=modules%2Fsmesh.git diff --git a/src/StdMeshers/StdMeshers_RadialPrism_3D.cxx b/src/StdMeshers/StdMeshers_RadialPrism_3D.cxx index 1814680af..f2d2970af 100644 --- a/src/StdMeshers/StdMeshers_RadialPrism_3D.cxx +++ b/src/StdMeshers/StdMeshers_RadialPrism_3D.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2013 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 @@ -6,7 +6,7 @@ // 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. +// version 2.1 of the License, or (at your option) any later version. // // This library is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of @@ -20,7 +20,7 @@ // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com // -// SMESH SMESH : implementaion of SMESH idl descriptions +// SMESH SMESH : implementation of SMESH idl descriptions // File : StdMeshers_RadialPrism_3D.cxx // Module : SMESH // Created : Fri Oct 20 11:37:07 2006 @@ -29,8 +29,6 @@ #include "StdMeshers_RadialPrism_3D.hxx" -#include - #include "StdMeshers_ProjectionUtils.hxx" #include "StdMeshers_NumberOfLayers.hxx" #include "StdMeshers_LayerDistribution.hxx" @@ -57,11 +55,7 @@ #include #include #include -#if OCC_VERSION_LARGE > 0x06050400 #include -#else -#include -#endif using namespace std; @@ -75,8 +69,8 @@ namespace ProjectionUtils = StdMeshers_ProjectionUtils; //purpose : //======================================================================= -StdMeshers_RadialPrism_3D::StdMeshers_RadialPrism_3D(int hypId, int studyId, SMESH_Gen* gen) - :SMESH_3D_Algo(hypId, studyId, gen) +StdMeshers_RadialPrism_3D::StdMeshers_RadialPrism_3D(int hypId, SMESH_Gen* gen) + :SMESH_3D_Algo(hypId, gen) { _name = "RadialPrism_3D"; _shapeType = (1 << TopAbs_SOLID); // 1 bit per shape type @@ -165,15 +159,11 @@ bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& a myHelper = new SMESH_MesherHelper( aMesh ); myHelper->IsQuadraticSubMesh( aShape ); // to delete helper at exit from Compute() - std::auto_ptr helperDeleter( myHelper ); + std::unique_ptr helperDeleter( myHelper ); // get 2 shells TopoDS_Solid solid = TopoDS::Solid( aShape ); -#if OCC_VERSION_LARGE > 0x06050400 TopoDS_Shell outerShell = BRepClass3d::OuterShell( solid ); -#else - TopoDS_Shell outerShell = BRepTools::OuterShell( solid ); -#endif TopoDS_Shape innerShell; int nbShells = 0; for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells ) @@ -187,13 +177,13 @@ bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& a // ---------------------------------- ProjectionUtils::TShapeShapeMap shape2ShapeMaps[2]; - if ( !ProjectionUtils::FindSubShapeAssociation( innerShell, &aMesh, - outerShell, &aMesh, - shape2ShapeMaps[0]) - && - !ProjectionUtils::FindSubShapeAssociation( innerShell.Reversed(), &aMesh, - outerShell, &aMesh, - shape2ShapeMaps[1])) + bool mapOk1 = ProjectionUtils::FindSubShapeAssociation( innerShell, &aMesh, + outerShell, &aMesh, + shape2ShapeMaps[0]); + bool mapOk2 = ProjectionUtils::FindSubShapeAssociation( innerShell.Reversed(), &aMesh, + outerShell, &aMesh, + shape2ShapeMaps[1]); + if ( !mapOk1 && !mapOk2 ) return error(COMPERR_BAD_SHAPE,"Topology of inner and outer shells seems different" ); int iMap; @@ -239,7 +229,7 @@ bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& a } // Find matching nodes of in and out faces - TNodeNodeMap nodeIn2OutMap; + ProjectionUtils::TNodeNodeMap nodeIn2OutMap; if ( ! ProjectionUtils::FindMatchingNodesOnFaces( inFace, &aMesh, outFace, &aMesh, shape2ShapeMap, nodeIn2OutMap )) return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Mesh on faces #") @@ -287,7 +277,7 @@ bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& a /*! * \brief Create a column of nodes from outNode to inNode * \param n2ColMap - map of node columns to add a created column - * \param outNode - botton node of a column + * \param outNode - bottom node of a column * \param inNode - top node of a column * \retval const TNodeColumn* - a new column pointer */ @@ -346,7 +336,7 @@ public: const int myID = -1000; TNodeDistributor* myHyp = dynamic_cast( aMesh.GetHypothesis( myID )); if ( !myHyp ) - myHyp = new TNodeDistributor( myID, 0, aMesh.GetGen() ); + myHyp = new TNodeDistributor( myID, aMesh.GetGen() ); return myHyp; } // ----------------------------------------------------------------------------- @@ -384,8 +374,8 @@ public: } protected: // ----------------------------------------------------------------------------- - TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen) - : StdMeshers_Regular_1D( hypId, studyId, gen) + TNodeDistributor( int hypId, SMESH_Gen* gen) + : StdMeshers_Regular_1D( hypId, gen) { } // ----------------------------------------------------------------------------- @@ -439,11 +429,7 @@ bool StdMeshers_RadialPrism_3D::Evaluate(SMESH_Mesh& aMesh, { // get 2 shells TopoDS_Solid solid = TopoDS::Solid( aShape ); -#if OCC_VERSION_LARGE > 0x06050400 TopoDS_Shell outerShell = BRepClass3d::OuterShell( solid ); -#else - TopoDS_Shell outerShell = BRepTools::OuterShell( solid ); -#endif TopoDS_Shape innerShell; int nbShells = 0; for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells ) @@ -598,3 +584,56 @@ bool StdMeshers_RadialPrism_3D::Evaluate(SMESH_Mesh& aMesh, return true; } + +//================================================================================ +/*! + * \brief Return true if the algorithm can mesh this shape + * \param [in] aShape - shape to check + * \param [in] toCheckAll - if true, this check returns OK if all shapes are OK, + * else, returns OK if at least one shape is OK + */ +//================================================================================ + +bool StdMeshers_RadialPrism_3D::IsApplicable( const TopoDS_Shape & aShape, bool toCheckAll ) +{ + int nbFoundSolids = 0; + for (TopExp_Explorer exp( aShape, TopAbs_SOLID ); exp.More(); exp.Next(), ++nbFoundSolids ) + { + TopoDS_Shape shell[2]; + int nbShells = 0; + for ( TopoDS_Iterator It (exp.Current()); It.More(); It.Next() ) + { + nbShells++; + if ( nbShells > 2 ) { + if ( toCheckAll ) return false; + break; + } + shell[ nbShells-1 ] = It.Value(); + } + if ( nbShells != 2 ) { + if ( toCheckAll ) return false; + continue; + } + + int nbFaces1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_FACE, 0 ); + int nbFaces2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_FACE, 0 ); + if ( nbFaces1 != nbFaces2 ){ + if( toCheckAll ) return false; + continue; + } + int nbEdges1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_EDGE, 0 ); + int nbEdges2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_EDGE, 0 ); + if ( nbEdges1 != nbEdges2 ){ + if( toCheckAll ) return false; + continue; + } + int nbVertices1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_VERTEX, 0 ); + int nbVertices2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_VERTEX, 0 ); + if ( nbVertices1 != nbVertices2 ){ + if( toCheckAll ) return false; + continue; + } + if ( !toCheckAll ) return true; + } + return ( toCheckAll && nbFoundSolids != 0); +}