1 // Copyright (C) 2007-2020 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 // SMESH SMESH : implementation of SMESH idl descriptions
24 // File : StdMeshers_RadialPrism_3D.cxx
26 // Created : Fri Oct 20 11:37:07 2006
27 // Author : Edward AGAPOV (eap)
30 #include "StdMeshers_RadialPrism_3D.hxx"
32 #include "StdMeshers_ProjectionUtils.hxx"
33 #include "StdMeshers_NumberOfLayers.hxx"
34 #include "StdMeshers_LayerDistribution.hxx"
35 #include "StdMeshers_Prism_3D.hxx"
36 #include "StdMeshers_Regular_1D.hxx"
38 #include "SMDS_MeshNode.hxx"
39 #include "SMESHDS_SubMesh.hxx"
40 #include "SMESH_Gen.hxx"
41 #include "SMESH_Mesh.hxx"
42 #include "SMESH_MesherHelper.hxx"
43 #include "SMESH_subMesh.hxx"
45 #include "utilities.h"
47 #include <BRepAdaptor_Curve.hxx>
48 #include <BRepBuilderAPI_MakeEdge.hxx>
49 #include <BRep_Tool.hxx>
50 #include <TopExp_Explorer.hxx>
51 #include <TopTools_DataMapIteratorOfDataMapOfShapeShape.hxx>
52 #include <TopTools_MapOfShape.hxx>
54 #include <TopoDS_Shell.hxx>
55 #include <TopoDS_Solid.hxx>
58 #include <BRepClass3d.hxx>
62 #define RETURN_BAD_RESULT(msg) { MESSAGE(")-: Error: " << msg); return false; }
63 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
65 namespace ProjectionUtils = StdMeshers_ProjectionUtils;
67 //=======================================================================
68 //function : StdMeshers_RadialPrism_3D
70 //=======================================================================
72 StdMeshers_RadialPrism_3D::StdMeshers_RadialPrism_3D(int hypId, SMESH_Gen* gen)
73 :SMESH_3D_Algo(hypId, gen)
75 _name = "RadialPrism_3D";
76 _shapeType = (1 << TopAbs_SOLID); // 1 bit per shape type
78 _compatibleHypothesis.push_back("LayerDistribution");
79 _compatibleHypothesis.push_back("NumberOfLayers");
81 myDistributionHypo = 0;
84 //================================================================================
88 //================================================================================
90 StdMeshers_RadialPrism_3D::~StdMeshers_RadialPrism_3D()
93 //=======================================================================
94 //function : CheckHypothesis
96 //=======================================================================
98 bool StdMeshers_RadialPrism_3D::CheckHypothesis(SMESH_Mesh& aMesh,
99 const TopoDS_Shape& aShape,
100 SMESH_Hypothesis::Hypothesis_Status& aStatus)
102 // check aShape that must have 2 shells
104 if ( ProjectionUtils::Count( aShape, TopAbs_SOLID, 0 ) != 1 ||
105 ProjectionUtils::Count( aShape, TopAbs_SHELL, 0 ) != 2 )
107 aStatus = HYP_BAD_GEOMETRY;
112 myDistributionHypo = 0;
114 list <const SMESHDS_Hypothesis * >::const_iterator itl;
116 const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
117 if ( hyps.size() == 0 )
119 aStatus = SMESH_Hypothesis::HYP_MISSING;
120 return false; // can't work with no hypothesis
123 if ( hyps.size() > 1 )
125 aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
129 const SMESHDS_Hypothesis *theHyp = hyps.front();
131 string hypName = theHyp->GetName();
133 if (hypName == "NumberOfLayers")
135 myNbLayerHypo = static_cast<const StdMeshers_NumberOfLayers *>(theHyp);
136 aStatus = SMESH_Hypothesis::HYP_OK;
139 if (hypName == "LayerDistribution")
141 myDistributionHypo = static_cast<const StdMeshers_LayerDistribution *>(theHyp);
142 aStatus = SMESH_Hypothesis::HYP_OK;
145 aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
149 //=======================================================================
152 //=======================================================================
154 bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape)
157 SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
159 myHelper = new SMESH_MesherHelper( aMesh );
160 myHelper->IsQuadraticSubMesh( aShape );
161 // to delete helper at exit from Compute()
162 std::unique_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
165 TopoDS_Solid solid = TopoDS::Solid( aShape );
166 TopoDS_Shell outerShell = BRepClass3d::OuterShell( solid );
167 TopoDS_Shape innerShell;
169 for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
170 if ( !outerShell.IsSame( It.Value() ))
171 innerShell = It.Value();
173 return error(COMPERR_BAD_SHAPE, SMESH_Comment("Must be 2 shells but not ")<<nbShells);
175 // ----------------------------------
176 // Associate sub-shapes of the shells
177 // ----------------------------------
179 ProjectionUtils::TShapeShapeMap shape2ShapeMaps[2];
180 bool mapOk1 = ProjectionUtils::FindSubShapeAssociation( innerShell, &aMesh,
183 bool mapOk2 = ProjectionUtils::FindSubShapeAssociation( innerShell.Reversed(), &aMesh,
186 if ( !mapOk1 && !mapOk2 )
187 return error(COMPERR_BAD_SHAPE,"Topology of inner and outer shells seems different" );
190 if ( shape2ShapeMaps[0].Extent() == shape2ShapeMaps[1].Extent() )
192 // choose an assiciation by shortest distance between VERTEXes
193 double dist1 = 0, dist2 = 0;
194 TopTools_DataMapIteratorOfDataMapOfShapeShape ssIt( shape2ShapeMaps[0]._map1to2 );
195 for (; ssIt.More(); ssIt.Next() )
197 if ( ssIt.Key().ShapeType() != TopAbs_VERTEX ) continue;
198 gp_Pnt pIn = BRep_Tool::Pnt( TopoDS::Vertex( ssIt.Key() ));
199 gp_Pnt pOut1 = BRep_Tool::Pnt( TopoDS::Vertex( ssIt.Value() ));
200 gp_Pnt pOut2 = BRep_Tool::Pnt( TopoDS::Vertex( shape2ShapeMaps[1]( ssIt.Key() )));
201 dist1 += pIn.SquareDistance( pOut1 );
202 dist2 += pIn.SquareDistance( pOut2 );
204 iMap = ( dist1 < dist2 ) ? 0 : 1;
208 iMap = ( shape2ShapeMaps[0].Extent() > shape2ShapeMaps[1].Extent() ) ? 0 : 1;
210 ProjectionUtils::TShapeShapeMap& shape2ShapeMap = shape2ShapeMaps[iMap];
212 // ------------------
214 // ------------------
216 TNode2ColumnMap node2columnMap;
217 myLayerPositions.clear();
219 for ( exp.Init( outerShell, TopAbs_FACE ); exp.More(); exp.Next() )
221 // Corresponding sub-shapes
222 TopoDS_Face outFace = TopoDS::Face( exp.Current() );
224 if ( !shape2ShapeMap.IsBound( outFace, /*isOut=*/true )) {
225 return error(SMESH_Comment("Corresponding inner face not found for face #" )
226 << meshDS->ShapeToIndex( outFace ));
228 inFace = TopoDS::Face( shape2ShapeMap( outFace, /*isOut=*/true ));
231 // Find matching nodes of in and out faces
232 ProjectionUtils::TNodeNodeMap nodeIn2OutMap;
233 if ( ! ProjectionUtils::FindMatchingNodesOnFaces( inFace, &aMesh, outFace, &aMesh,
234 shape2ShapeMap, nodeIn2OutMap ))
235 return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Mesh on faces #")
236 << meshDS->ShapeToIndex( outFace ) << " and "
237 << meshDS->ShapeToIndex( inFace ) << " seems different" );
241 SMDS_ElemIteratorPtr faceIt = meshDS->MeshElements( inFace )->GetElements();
242 while ( faceIt->more() ) // loop on faces on inFace
244 const SMDS_MeshElement* face = faceIt->next();
245 if ( !face || face->GetType() != SMDSAbs_Face )
247 int nbNodes = face->NbNodes();
248 if ( face->IsQuadratic() )
251 // find node columns for each node
252 vector< const TNodeColumn* > columns( nbNodes );
253 for ( int i = 0; i < nbNodes; ++i )
255 const SMDS_MeshNode* nIn = face->GetNode( i );
256 TNode2ColumnMap::iterator n_col = node2columnMap.find( nIn );
257 if ( n_col != node2columnMap.end() ) {
258 columns[ i ] = & n_col->second;
261 TNodeNodeMap::iterator nInOut = nodeIn2OutMap.find( nIn );
262 if ( nInOut == nodeIn2OutMap.end() )
263 RETURN_BAD_RESULT("No matching node for "<< nIn->GetID() <<
264 " in face "<< face->GetID());
265 columns[ i ] = makeNodeColumn( node2columnMap, nIn, nInOut->second );
269 StdMeshers_Prism_3D::AddPrisms( columns, myHelper );
271 } // loop on faces of out shell
276 //================================================================================
278 * \brief Create a column of nodes from outNode to inNode
279 * \param n2ColMap - map of node columns to add a created column
280 * \param outNode - bottom node of a column
281 * \param inNode - top node of a column
282 * \retval const TNodeColumn* - a new column pointer
284 //================================================================================
286 TNodeColumn* StdMeshers_RadialPrism_3D::makeNodeColumn( TNode2ColumnMap& n2ColMap,
287 const SMDS_MeshNode* outNode,
288 const SMDS_MeshNode* inNode)
290 SMESHDS_Mesh * meshDS = myHelper->GetMeshDS();
291 int shapeID = myHelper->GetSubShapeID();
293 if ( myLayerPositions.empty() ) {
294 gp_Pnt pIn = gpXYZ( inNode ), pOut = gpXYZ( outNode );
295 computeLayerPositions( pIn, pOut );
297 int nbSegments = myLayerPositions.size() + 1;
299 TNode2ColumnMap::iterator n_col =
300 n2ColMap.insert( make_pair( outNode, TNodeColumn() )).first;
301 TNodeColumn & column = n_col->second;
302 column.resize( nbSegments + 1 );
303 column.front() = outNode;
304 column.back() = inNode;
306 gp_XYZ p1 = gpXYZ( outNode );
307 gp_XYZ p2 = gpXYZ( inNode );
308 for ( int z = 1; z < nbSegments; ++z )
310 double r = myLayerPositions[ z - 1 ];
311 gp_XYZ p = ( 1 - r ) * p1 + r * p2;
312 SMDS_MeshNode* n = meshDS->AddNode( p.X(), p.Y(), p.Z() );
313 meshDS->SetNodeInVolume( n, shapeID );
320 //================================================================================
321 //================================================================================
323 * \brief Class computing layers distribution using data of
324 * StdMeshers_LayerDistribution hypothesis
326 //================================================================================
327 //================================================================================
329 class TNodeDistributor: public StdMeshers_Regular_1D
331 list <const SMESHDS_Hypothesis *> myUsedHyps;
333 // -----------------------------------------------------------------------------
334 static TNodeDistributor* GetDistributor(SMESH_Mesh& aMesh)
336 const int myID = -1000;
337 TNodeDistributor* myHyp = dynamic_cast<TNodeDistributor*>( aMesh.GetHypothesis( myID ));
339 myHyp = new TNodeDistributor( myID, aMesh.GetGen() );
342 // -----------------------------------------------------------------------------
343 bool Compute( vector< double > & positions,
347 const StdMeshers_LayerDistribution* hyp)
349 double len = pIn.Distance( pOut );
350 if ( len <= DBL_MIN ) return error("Too close points of inner and outer shells");
352 if ( !hyp || !hyp->GetLayerDistribution() )
353 return error( "Invalid LayerDistribution hypothesis");
355 myUsedHyps.push_back( hyp->GetLayerDistribution() );
357 TopoDS_Edge edge = BRepBuilderAPI_MakeEdge( pIn, pOut );
358 SMESH_Hypothesis::Hypothesis_Status aStatus;
359 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, edge, aStatus ))
360 return error( "StdMeshers_Regular_1D::CheckHypothesis() failed "
361 "with LayerDistribution hypothesis");
363 BRepAdaptor_Curve C3D(edge);
364 double f = C3D.FirstParameter(), l = C3D.LastParameter();
365 list< double > params;
366 if ( !StdMeshers_Regular_1D::computeInternalParameters( aMesh, C3D, len, f, l, params, false ))
367 return error("StdMeshers_Regular_1D failed to compute layers distribution");
370 positions.reserve( params.size() );
371 for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++)
372 positions.push_back( *itU / len );
376 // -----------------------------------------------------------------------------
377 TNodeDistributor( int hypId, SMESH_Gen* gen)
378 : StdMeshers_Regular_1D( hypId, gen)
381 // -----------------------------------------------------------------------------
382 virtual const list <const SMESHDS_Hypothesis *> &
383 GetUsedHypothesis(SMESH_Mesh &, const TopoDS_Shape &, const bool)
387 // -----------------------------------------------------------------------------
390 //================================================================================
392 * \brief Compute positions of nodes between the internal and the external surfaces
393 * \retval bool - is a success
395 //================================================================================
397 bool StdMeshers_RadialPrism_3D::computeLayerPositions(const gp_Pnt& pIn,
402 int nbSegments = myNbLayerHypo->GetNumberOfLayers();
403 myLayerPositions.resize( nbSegments - 1 );
404 for ( int z = 1; z < nbSegments; ++z )
405 myLayerPositions[ z - 1 ] = double( z )/ double( nbSegments );
408 if ( myDistributionHypo ) {
409 SMESH_Mesh * mesh = myHelper->GetMesh();
410 if ( !TNodeDistributor::GetDistributor(*mesh)->Compute( myLayerPositions, pIn, pOut,
411 *mesh, myDistributionHypo ))
413 error( TNodeDistributor::GetDistributor(*mesh)->GetComputeError() );
417 RETURN_BAD_RESULT("Bad hypothesis");
421 //=======================================================================
422 //function : Evaluate
424 //=======================================================================
426 bool StdMeshers_RadialPrism_3D::Evaluate(SMESH_Mesh& aMesh,
427 const TopoDS_Shape& aShape,
428 MapShapeNbElems& aResMap)
431 TopoDS_Solid solid = TopoDS::Solid( aShape );
432 TopoDS_Shell outerShell = BRepClass3d::OuterShell( solid );
433 TopoDS_Shape innerShell;
435 for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
436 if ( !outerShell.IsSame( It.Value() ))
437 innerShell = It.Value();
438 if ( nbShells != 2 ) {
439 std::vector<int> aResVec(SMDSEntity_Last);
440 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
441 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
442 aResMap.insert(std::make_pair(sm,aResVec));
443 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
444 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
448 // Associate sub-shapes of the shells
449 ProjectionUtils::TShapeShapeMap shape2ShapeMap;
450 if ( !ProjectionUtils::FindSubShapeAssociation( outerShell, &aMesh,
453 std::vector<int> aResVec(SMDSEntity_Last);
454 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
455 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
456 aResMap.insert(std::make_pair(sm,aResVec));
457 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
458 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
462 // get info for outer shell
463 int nb0d_Out=0, nb2d_3_Out=0, nb2d_4_Out=0;
464 //TopTools_SequenceOfShape FacesOut;
465 for (TopExp_Explorer exp(outerShell, TopAbs_FACE); exp.More(); exp.Next()) {
466 //FacesOut.Append(exp.Current());
467 SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
468 MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
469 std::vector<int> aVec = (*anIt).second;
470 nb0d_Out += aVec[SMDSEntity_Node];
471 nb2d_3_Out += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
472 nb2d_4_Out += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
475 TopTools_MapOfShape tmpMap;
476 for (TopExp_Explorer exp(outerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
477 if( tmpMap.Contains( exp.Current() ) )
479 tmpMap.Add( exp.Current() );
480 SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
481 MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
482 std::vector<int> aVec = (*anIt).second;
483 nb0d_Out += aVec[SMDSEntity_Node];
484 nb1d_Out += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
487 for (TopExp_Explorer exp(outerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
488 if( tmpMap.Contains( exp.Current() ) )
490 tmpMap.Add( exp.Current() );
494 // get info for inner shell
495 int nb0d_In=0, nb2d_3_In=0, nb2d_4_In=0;
496 //TopTools_SequenceOfShape FacesIn;
497 for (TopExp_Explorer exp(innerShell, TopAbs_FACE); exp.More(); exp.Next()) {
498 //FacesIn.Append(exp.Current());
499 SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
500 MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
501 std::vector<int> aVec = (*anIt).second;
502 nb0d_In += aVec[SMDSEntity_Node];
503 nb2d_3_In += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
504 nb2d_4_In += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
508 bool IsQuadratic = false;
510 for (TopExp_Explorer exp(innerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
511 if( tmpMap.Contains( exp.Current() ) )
513 tmpMap.Add( exp.Current() );
514 SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
515 MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
516 std::vector<int> aVec = (*anIt).second;
517 nb0d_In += aVec[SMDSEntity_Node];
518 nb1d_In += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
520 IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
525 for (TopExp_Explorer exp(innerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
526 if( tmpMap.Contains( exp.Current() ) )
528 tmpMap.Add( exp.Current() );
532 bool IsOK = (nb0d_Out==nb0d_In) && (nb1d_Out==nb1d_In) &&
533 (nb2d_3_Out==nb2d_3_In) && (nb2d_4_Out==nb2d_4_In);
535 std::vector<int> aResVec(SMDSEntity_Last);
536 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
537 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
538 aResMap.insert(std::make_pair(sm,aResVec));
539 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
540 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
545 if( myNbLayerHypo ) {
546 nbLayers = myNbLayerHypo->GetNumberOfLayers();
548 if ( myDistributionHypo ) {
549 if ( !myDistributionHypo->GetLayerDistribution() ) {
550 std::vector<int> aResVec(SMDSEntity_Last);
551 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
552 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
553 aResMap.insert(std::make_pair(sm,aResVec));
554 SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
555 smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
558 TopExp_Explorer exp(outerShell, TopAbs_VERTEX);
559 TopoDS_Vertex Vout = TopoDS::Vertex(exp.Current());
560 TopoDS_Vertex Vin = TopoDS::Vertex( shape2ShapeMap(Vout) );
561 if ( myLayerPositions.empty() ) {
562 gp_Pnt pIn = BRep_Tool::Pnt(Vin);
563 gp_Pnt pOut = BRep_Tool::Pnt(Vout);
564 computeLayerPositions( pIn, pOut );
566 nbLayers = myLayerPositions.size() + 1;
569 std::vector<int> aResVec(SMDSEntity_Last);
570 for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
572 aResVec[SMDSEntity_Quad_Penta] = nb2d_3_Out * nbLayers;
573 aResVec[SMDSEntity_Quad_Hexa] = nb2d_4_Out * nbLayers;
574 int nb1d = ( nb2d_3_Out*3 + nb2d_4_Out*4 ) / 2;
575 aResVec[SMDSEntity_Node] = nb0d_Out * ( 2*nbLayers - 1 ) - nb1d * nbLayers;
578 aResVec[SMDSEntity_Node] = nb0d_Out * ( nbLayers - 1 );
579 aResVec[SMDSEntity_Penta] = nb2d_3_Out * nbLayers;
580 aResVec[SMDSEntity_Hexa] = nb2d_4_Out * nbLayers;
582 SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
583 aResMap.insert(std::make_pair(sm,aResVec));
588 //================================================================================
590 * \brief Return true if the algorithm can mesh this shape
591 * \param [in] aShape - shape to check
592 * \param [in] toCheckAll - if true, this check returns OK if all shapes are OK,
593 * else, returns OK if at least one shape is OK
595 //================================================================================
597 bool StdMeshers_RadialPrism_3D::IsApplicable( const TopoDS_Shape & aShape, bool toCheckAll )
599 int nbFoundSolids = 0;
600 for (TopExp_Explorer exp( aShape, TopAbs_SOLID ); exp.More(); exp.Next(), ++nbFoundSolids )
602 TopoDS_Shape shell[2];
604 for ( TopoDS_Iterator It (exp.Current()); It.More(); It.Next() )
607 if ( nbShells > 2 ) {
608 if ( toCheckAll ) return false;
611 shell[ nbShells-1 ] = It.Value();
613 if ( nbShells != 2 ) {
614 if ( toCheckAll ) return false;
618 int nbFaces1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_FACE, 0 );
619 int nbFaces2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_FACE, 0 );
620 if ( nbFaces1 != nbFaces2 ){
621 if( toCheckAll ) return false;
624 int nbEdges1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_EDGE, 0 );
625 int nbEdges2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_EDGE, 0 );
626 if ( nbEdges1 != nbEdges2 ){
627 if( toCheckAll ) return false;
630 int nbVertices1 = SMESH_MesherHelper:: Count( shell[0], TopAbs_VERTEX, 0 );
631 int nbVertices2 = SMESH_MesherHelper:: Count( shell[1], TopAbs_VERTEX, 0 );
632 if ( nbVertices1 != nbVertices2 ){
633 if( toCheckAll ) return false;
636 if ( !toCheckAll ) return true;
638 return ( toCheckAll && nbFoundSolids != 0);