Back to index

salome-smesh  6.5.0
StdMeshers_RadialPrism_3D.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 // This library is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 2.1 of the License.
00010 //
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 // Lesser General Public License for more details.
00015 //
00016 // You should have received a copy of the GNU Lesser General Public
00017 // License along with this library; if not, write to the Free Software
00018 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 
00023 //  SMESH SMESH : implementaion of SMESH idl descriptions
00024 // File      : StdMeshers_RadialPrism_3D.cxx
00025 // Module    : SMESH
00026 // Created   : Fri Oct 20 11:37:07 2006
00027 // Author    : Edward AGAPOV (eap)
00028 //
00029 #include "StdMeshers_RadialPrism_3D.hxx"
00030 
00031 #include "StdMeshers_ProjectionUtils.hxx"
00032 #include "StdMeshers_NumberOfLayers.hxx"
00033 #include "StdMeshers_LayerDistribution.hxx"
00034 #include "StdMeshers_Prism_3D.hxx"
00035 #include "StdMeshers_Regular_1D.hxx"
00036 
00037 #include "SMDS_MeshNode.hxx"
00038 #include "SMESHDS_SubMesh.hxx"
00039 #include "SMESH_Gen.hxx"
00040 #include "SMESH_Mesh.hxx"
00041 #include "SMESH_MesherHelper.hxx"
00042 #include "SMESH_subMesh.hxx"
00043 
00044 #include "utilities.h"
00045 
00046 #include <BRepAdaptor_Curve.hxx>
00047 #include <BRepBuilderAPI_MakeEdge.hxx>
00048 #include <BRepTools.hxx>
00049 #include <BRep_Tool.hxx>
00050 #include <TopExp_Explorer.hxx>
00051 #include <TopoDS.hxx>
00052 #include <TopoDS_Shell.hxx>
00053 #include <TopoDS_Solid.hxx>
00054 #include <TopTools_MapOfShape.hxx>
00055 #include <gp.hxx>
00056 #include <gp_Pnt.hxx>
00057 
00058 
00059 using namespace std;
00060 
00061 #define RETURN_BAD_RESULT(msg) { MESSAGE(")-: Error: " << msg); return false; }
00062 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
00063 
00064 typedef StdMeshers_ProjectionUtils TAssocTool;
00065 
00066 //=======================================================================
00067 //function : StdMeshers_RadialPrism_3D
00068 //purpose  : 
00069 //=======================================================================
00070 
00071 StdMeshers_RadialPrism_3D::StdMeshers_RadialPrism_3D(int hypId, int studyId, SMESH_Gen* gen)
00072   :SMESH_3D_Algo(hypId, studyId, gen)
00073 {
00074   _name = "RadialPrism_3D";
00075   _shapeType = (1 << TopAbs_SOLID);     // 1 bit per shape type
00076 
00077   _compatibleHypothesis.push_back("LayerDistribution");
00078   _compatibleHypothesis.push_back("NumberOfLayers");
00079   myNbLayerHypo = 0;
00080   myDistributionHypo = 0;
00081 }
00082 
00083 //================================================================================
00087 //================================================================================
00088 
00089 StdMeshers_RadialPrism_3D::~StdMeshers_RadialPrism_3D()
00090 {}
00091 
00092 //=======================================================================
00093 //function : CheckHypothesis
00094 //purpose  : 
00095 //=======================================================================
00096 
00097 bool StdMeshers_RadialPrism_3D::CheckHypothesis(SMESH_Mesh&                          aMesh,
00098                                                 const TopoDS_Shape&                  aShape,
00099                                                 SMESH_Hypothesis::Hypothesis_Status& aStatus)
00100 {
00101   // check aShape that must have 2 shells
00102 /*  PAL16229
00103   if ( TAssocTool::Count( aShape, TopAbs_SOLID, 0 ) != 1 ||
00104        TAssocTool::Count( aShape, TopAbs_SHELL, 0 ) != 2 )
00105   {
00106     aStatus = HYP_BAD_GEOMETRY;
00107     return false;
00108   }
00109 */
00110   myNbLayerHypo = 0;
00111   myDistributionHypo = 0;
00112 
00113   list <const SMESHDS_Hypothesis * >::const_iterator itl;
00114 
00115   const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
00116   if ( hyps.size() == 0 )
00117   {
00118     aStatus = SMESH_Hypothesis::HYP_MISSING;
00119     return false;  // can't work with no hypothesis
00120   }
00121 
00122   if ( hyps.size() > 1 )
00123   {
00124     aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
00125     return false;
00126   }
00127 
00128   const SMESHDS_Hypothesis *theHyp = hyps.front();
00129 
00130   string hypName = theHyp->GetName();
00131 
00132   if (hypName == "NumberOfLayers")
00133   {
00134     myNbLayerHypo = static_cast<const StdMeshers_NumberOfLayers *>(theHyp);
00135     aStatus = SMESH_Hypothesis::HYP_OK;
00136     return true;
00137   }
00138   if (hypName == "LayerDistribution")
00139   {
00140     myDistributionHypo = static_cast<const StdMeshers_LayerDistribution *>(theHyp);
00141     aStatus = SMESH_Hypothesis::HYP_OK;
00142     return true;
00143   }
00144   aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
00145   return true;
00146 }
00147 
00148 //=======================================================================
00149 //function : Compute
00150 //purpose  : 
00151 //=======================================================================
00152 
00153 bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape)
00154 {
00155   TopExp_Explorer exp;
00156   SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
00157 
00158   myHelper = new SMESH_MesherHelper( aMesh );
00159   myHelper->IsQuadraticSubMesh( aShape );
00160   // to delete helper at exit from Compute()
00161   std::auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
00162 
00163   // get 2 shells
00164   TopoDS_Solid solid = TopoDS::Solid( aShape );
00165   TopoDS_Shell outerShell = BRepTools::OuterShell( solid );
00166   TopoDS_Shape innerShell;
00167   int nbShells = 0;
00168   for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
00169     if ( !outerShell.IsSame( It.Value() ))
00170       innerShell = It.Value();
00171   if ( nbShells != 2 )
00172     return error(COMPERR_BAD_SHAPE, SMESH_Comment("Must be 2 shells but not ")<<nbShells);
00173 
00174   // ----------------------------------
00175   // Associate sub-shapes of the shells
00176   // ----------------------------------
00177 
00178   TAssocTool::TShapeShapeMap shape2ShapeMap;
00179   if ( !TAssocTool::FindSubShapeAssociation( innerShell, &aMesh,
00180                                              outerShell, &aMesh,
00181                                              shape2ShapeMap) )
00182     return error(COMPERR_BAD_SHAPE,"Topology of inner and outer shells seems different" );
00183 
00184   // ------------------
00185   // Make mesh
00186   // ------------------
00187 
00188   TNode2ColumnMap node2columnMap;
00189   myLayerPositions.clear();
00190 
00191   for ( exp.Init( outerShell, TopAbs_FACE ); exp.More(); exp.Next() )
00192   {
00193     // Corresponding sub-shapes
00194     TopoDS_Face outFace = TopoDS::Face( exp.Current() );
00195     TopoDS_Face inFace;
00196     if ( !shape2ShapeMap.IsBound( outFace, /*isOut=*/true )) {
00197       return error(SMESH_Comment("Corresponding inner face not found for face #" )
00198                    << meshDS->ShapeToIndex( outFace ));
00199     } else {
00200       inFace = TopoDS::Face( shape2ShapeMap( outFace, /*isOut=*/true ));
00201     }
00202 
00203     // Find matching nodes of in and out faces
00204     TNodeNodeMap nodeIn2OutMap;
00205     if ( ! TAssocTool::FindMatchingNodesOnFaces( inFace, &aMesh, outFace, &aMesh,
00206                                                  shape2ShapeMap, nodeIn2OutMap ))
00207       return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Mesh on faces #")
00208                    << meshDS->ShapeToIndex( outFace ) << " and "
00209                    << meshDS->ShapeToIndex( inFace ) << " seems different" );
00210 
00211     // Create volumes
00212 
00213     SMDS_ElemIteratorPtr faceIt = meshDS->MeshElements( inFace )->GetElements();
00214     while ( faceIt->more() ) // loop on faces on inFace
00215     {
00216       const SMDS_MeshElement* face = faceIt->next();
00217       if ( !face || face->GetType() != SMDSAbs_Face )
00218         continue;
00219       int nbNodes = face->NbNodes();
00220       if ( face->IsQuadratic() )
00221         nbNodes /= 2;
00222 
00223       // find node columns for each node
00224       vector< const TNodeColumn* > columns( nbNodes );
00225       for ( int i = 0; i < nbNodes; ++i )
00226       {
00227         const SMDS_MeshNode* nIn = face->GetNode( i );
00228         TNode2ColumnMap::iterator n_col = node2columnMap.find( nIn );
00229         if ( n_col != node2columnMap.end() ) {
00230           columns[ i ] = & n_col->second;
00231         }
00232         else {
00233           TNodeNodeMap::iterator nInOut = nodeIn2OutMap.find( nIn );
00234           if ( nInOut == nodeIn2OutMap.end() )
00235             RETURN_BAD_RESULT("No matching node for "<< nIn->GetID() <<
00236                               " in face "<< face->GetID());
00237           columns[ i ] = makeNodeColumn( node2columnMap, nIn, nInOut->second );
00238         }
00239       }
00240 
00241       StdMeshers_Prism_3D::AddPrisms( columns, myHelper );
00242     }
00243   } // loop on faces of out shell
00244 
00245   return true;
00246 }
00247 
00248 //================================================================================
00256 //================================================================================
00257 
00258 TNodeColumn* StdMeshers_RadialPrism_3D::makeNodeColumn( TNode2ColumnMap&     n2ColMap,
00259                                                         const SMDS_MeshNode* outNode,
00260                                                         const SMDS_MeshNode* inNode)
00261 {
00262   SMESHDS_Mesh * meshDS = myHelper->GetMeshDS();
00263   int shapeID = myHelper->GetSubShapeID();
00264 
00265   if ( myLayerPositions.empty() ) {
00266     gp_Pnt pIn = gpXYZ( inNode ), pOut = gpXYZ( outNode );
00267     computeLayerPositions( pIn, pOut );
00268   }
00269   int nbSegments = myLayerPositions.size() + 1;
00270 
00271   TNode2ColumnMap::iterator n_col =
00272     n2ColMap.insert( make_pair( outNode, TNodeColumn() )).first;
00273   TNodeColumn & column = n_col->second;
00274   column.resize( nbSegments + 1 );
00275   column.front() = outNode;
00276   column.back() =  inNode;
00277 
00278   gp_XYZ p1 = gpXYZ( outNode );
00279   gp_XYZ p2 = gpXYZ( inNode );
00280   for ( int z = 1; z < nbSegments; ++z )
00281   {
00282     double r = myLayerPositions[ z - 1 ];
00283     gp_XYZ p = ( 1 - r ) * p1 + r * p2;
00284     SMDS_MeshNode* n = meshDS->AddNode( p.X(), p.Y(), p.Z() );
00285     meshDS->SetNodeInVolume( n, shapeID );
00286     column[ z ] = n;
00287   }
00288 
00289   return & column;
00290 }
00291 
00292 //================================================================================
00293 //================================================================================
00298 //================================================================================
00299 //================================================================================
00300 
00301 class TNodeDistributor: public StdMeshers_Regular_1D
00302 {
00303   list <const SMESHDS_Hypothesis *> myUsedHyps;
00304 public:
00305   // -----------------------------------------------------------------------------
00306   static TNodeDistributor* GetDistributor(SMESH_Mesh& aMesh)
00307   {
00308     const int myID = -1000;
00309     map < int, SMESH_1D_Algo * > & algoMap = aMesh.GetGen()->_map1D_Algo;
00310     map < int, SMESH_1D_Algo * >::iterator id_algo = algoMap.find( myID );
00311     if ( id_algo == algoMap.end() )
00312       return new TNodeDistributor( myID, 0, aMesh.GetGen() );
00313     return static_cast< TNodeDistributor* >( id_algo->second );
00314   }
00315   // -----------------------------------------------------------------------------
00316   bool Compute( vector< double > &                  positions,
00317                 gp_Pnt                              pIn,
00318                 gp_Pnt                              pOut,
00319                 SMESH_Mesh&                         aMesh,
00320                 const StdMeshers_LayerDistribution* hyp)
00321   {
00322     double len = pIn.Distance( pOut );
00323     if ( len <= DBL_MIN ) return error("Too close points of inner and outer shells");
00324 
00325     if ( !hyp || !hyp->GetLayerDistribution() )
00326       return error( "Invalid LayerDistribution hypothesis");
00327     myUsedHyps.clear();
00328     myUsedHyps.push_back( hyp->GetLayerDistribution() );
00329 
00330     TopoDS_Edge edge = BRepBuilderAPI_MakeEdge( pIn, pOut );
00331     SMESH_Hypothesis::Hypothesis_Status aStatus;
00332     if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, edge, aStatus ))
00333       return error( "StdMeshers_Regular_1D::CheckHypothesis() failed "
00334                     "with LayerDistribution hypothesis");
00335 
00336     BRepAdaptor_Curve C3D(edge);
00337     double f = C3D.FirstParameter(), l = C3D.LastParameter();
00338     list< double > params;
00339     if ( !StdMeshers_Regular_1D::computeInternalParameters( aMesh, C3D, len, f, l, params, false ))
00340       return error("StdMeshers_Regular_1D failed to compute layers distribution");
00341 
00342     positions.clear();
00343     positions.reserve( params.size() );
00344     for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++)
00345       positions.push_back( *itU / len );
00346     return true;
00347   }
00348 protected:
00349   // -----------------------------------------------------------------------------
00350   TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen)
00351     : StdMeshers_Regular_1D( hypId, studyId, gen)
00352   {
00353   }
00354   // -----------------------------------------------------------------------------
00355   virtual const list <const SMESHDS_Hypothesis *> &
00356     GetUsedHypothesis(SMESH_Mesh &, const TopoDS_Shape &, const bool)
00357   {
00358     return myUsedHyps;
00359   }
00360   // -----------------------------------------------------------------------------
00361 };
00362 
00363 //================================================================================
00368 //================================================================================
00369 
00370 bool StdMeshers_RadialPrism_3D::computeLayerPositions(const gp_Pnt& pIn,
00371                                                       const gp_Pnt& pOut)
00372 {
00373   if ( myNbLayerHypo )
00374   {
00375     int nbSegments = myNbLayerHypo->GetNumberOfLayers();
00376     myLayerPositions.resize( nbSegments - 1 );
00377     for ( int z = 1; z < nbSegments; ++z )
00378       myLayerPositions[ z - 1 ] = double( z )/ double( nbSegments );
00379     return true;
00380   }
00381   if ( myDistributionHypo ) {
00382     SMESH_Mesh * mesh = myHelper->GetMesh();
00383     if ( !TNodeDistributor::GetDistributor(*mesh)->Compute( myLayerPositions, pIn, pOut,
00384                                                             *mesh, myDistributionHypo ))
00385     {
00386       error( TNodeDistributor::GetDistributor(*mesh)->GetComputeError() );
00387       return false;
00388     }
00389   }
00390   RETURN_BAD_RESULT("Bad hypothesis");
00391 }
00392 
00393 
00394 //=======================================================================
00395 //function : Evaluate
00396 //purpose  : 
00397 //=======================================================================
00398 
00399 bool StdMeshers_RadialPrism_3D::Evaluate(SMESH_Mesh& aMesh,
00400                                          const TopoDS_Shape& aShape,
00401                                          MapShapeNbElems& aResMap)
00402 {
00403   // get 2 shells
00404   TopoDS_Solid solid = TopoDS::Solid( aShape );
00405   TopoDS_Shell outerShell = BRepTools::OuterShell( solid );
00406   TopoDS_Shape innerShell;
00407   int nbShells = 0;
00408   for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
00409     if ( !outerShell.IsSame( It.Value() ))
00410       innerShell = It.Value();
00411   if ( nbShells != 2 ) {
00412     std::vector<int> aResVec(SMDSEntity_Last);
00413     for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
00414     SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00415     aResMap.insert(std::make_pair(sm,aResVec));
00416     SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
00417     smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
00418     return false;
00419   }
00420 
00421   // Associate sub-shapes of the shells
00422   TAssocTool::TShapeShapeMap shape2ShapeMap;
00423   if ( !TAssocTool::FindSubShapeAssociation( outerShell, &aMesh,
00424                                              innerShell, &aMesh,
00425                                              shape2ShapeMap) ) {
00426     std::vector<int> aResVec(SMDSEntity_Last);
00427     for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
00428     SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00429     aResMap.insert(std::make_pair(sm,aResVec));
00430     SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
00431     smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
00432     return false;
00433   }
00434 
00435   // get info for outer shell
00436   int nb0d_Out=0, nb2d_3_Out=0, nb2d_4_Out=0;
00437   //TopTools_SequenceOfShape FacesOut;
00438   for (TopExp_Explorer exp(outerShell, TopAbs_FACE); exp.More(); exp.Next()) {
00439     //FacesOut.Append(exp.Current());
00440     SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
00441     MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
00442     std::vector<int> aVec = (*anIt).second;
00443     nb0d_Out += aVec[SMDSEntity_Node];
00444     nb2d_3_Out += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
00445     nb2d_4_Out += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
00446   }
00447   int nb1d_Out = 0;
00448   TopTools_MapOfShape tmpMap;
00449   for (TopExp_Explorer exp(outerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
00450     if( tmpMap.Contains( exp.Current() ) )
00451       continue;
00452     tmpMap.Add( exp.Current() );
00453     SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
00454     MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
00455     std::vector<int> aVec = (*anIt).second;
00456     nb0d_Out += aVec[SMDSEntity_Node];
00457     nb1d_Out += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
00458   }
00459   tmpMap.Clear();
00460   for (TopExp_Explorer exp(outerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
00461     if( tmpMap.Contains( exp.Current() ) )
00462       continue;
00463     tmpMap.Add( exp.Current() );
00464     nb0d_Out++;
00465   }
00466 
00467   // get info for inner shell
00468   int nb0d_In=0, nb2d_3_In=0, nb2d_4_In=0;
00469   //TopTools_SequenceOfShape FacesIn;
00470   for (TopExp_Explorer exp(innerShell, TopAbs_FACE); exp.More(); exp.Next()) {
00471     //FacesIn.Append(exp.Current());
00472     SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
00473     MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
00474     std::vector<int> aVec = (*anIt).second;
00475     nb0d_In += aVec[SMDSEntity_Node];
00476     nb2d_3_In += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
00477     nb2d_4_In += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
00478   }
00479   int nb1d_In = 0;
00480   tmpMap.Clear();
00481   bool IsQuadratic = false;
00482   bool IsFirst = true;
00483   for (TopExp_Explorer exp(innerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
00484     if( tmpMap.Contains( exp.Current() ) )
00485       continue;
00486     tmpMap.Add( exp.Current() );
00487     SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
00488     MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
00489     std::vector<int> aVec = (*anIt).second;
00490     nb0d_In += aVec[SMDSEntity_Node];
00491     nb1d_In += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
00492     if(IsFirst) {
00493       IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
00494       IsFirst = false;
00495     }
00496   }
00497   tmpMap.Clear();
00498   for (TopExp_Explorer exp(innerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
00499     if( tmpMap.Contains( exp.Current() ) )
00500       continue;
00501     tmpMap.Add( exp.Current() );
00502     nb0d_In++;
00503   }
00504 
00505   bool IsOK = (nb0d_Out==nb0d_In) && (nb1d_Out==nb1d_In) && 
00506               (nb2d_3_Out==nb2d_3_In) && (nb2d_4_Out==nb2d_4_In);
00507   if(!IsOK) {
00508     std::vector<int> aResVec(SMDSEntity_Last);
00509     for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
00510     SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00511     aResMap.insert(std::make_pair(sm,aResVec));
00512     SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
00513     smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
00514     return false;
00515   }
00516 
00517   int nbLayers = 0;
00518   if( myNbLayerHypo ) {
00519     nbLayers = myNbLayerHypo->GetNumberOfLayers();
00520   }
00521   if ( myDistributionHypo ) {
00522     if ( !myDistributionHypo->GetLayerDistribution() ) {
00523       std::vector<int> aResVec(SMDSEntity_Last);
00524       for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
00525       SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00526       aResMap.insert(std::make_pair(sm,aResVec));
00527       SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
00528       smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
00529       return false;
00530     }
00531     TopExp_Explorer exp(outerShell, TopAbs_VERTEX);
00532     TopoDS_Vertex Vout = TopoDS::Vertex(exp.Current());
00533     TopoDS_Vertex Vin = TopoDS::Vertex( shape2ShapeMap(Vout) );
00534     if ( myLayerPositions.empty() ) {
00535       gp_Pnt pIn = BRep_Tool::Pnt(Vin);
00536       gp_Pnt pOut = BRep_Tool::Pnt(Vout);
00537       computeLayerPositions( pIn, pOut );
00538     }
00539     nbLayers = myLayerPositions.size() + 1;
00540   }
00541 
00542   std::vector<int> aResVec(SMDSEntity_Last);
00543   for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
00544   if(IsQuadratic) {
00545     aResVec[SMDSEntity_Quad_Penta] = nb2d_3_Out * nbLayers;
00546     aResVec[SMDSEntity_Quad_Hexa] = nb2d_4_Out * nbLayers;
00547     int nb1d = ( nb2d_3_Out*3 + nb2d_4_Out*4 ) / 2;
00548     aResVec[SMDSEntity_Node] = nb0d_Out * ( 2*nbLayers - 1 ) - nb1d * nbLayers;
00549   }
00550   else {
00551     aResVec[SMDSEntity_Node] = nb0d_Out * ( nbLayers - 1 );
00552     aResVec[SMDSEntity_Penta] = nb2d_3_Out * nbLayers;
00553     aResVec[SMDSEntity_Hexa] = nb2d_4_Out * nbLayers;
00554   }
00555   SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00556   aResMap.insert(std::make_pair(sm,aResVec));
00557 
00558   return true;
00559 }