src/StdMeshers/StdMeshers_Penta_3D.cxx

   1 //  SMESH StdMeshers_Penta_3D implementaion of SMESH idl descriptions
   2 //
   3 //  Copyright (C) 2003  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 //  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   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.
  10 //
  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.
  15 //
  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
  19 //
  20 //  See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
  21 //
  22 //
  23 //
  24 //  File   : StdMeshers_Penta_3D.cxx
  25 //  Module : SMESH
  26
  27 #include "StdMeshers_Penta_3D.hxx"
  28
  29 #include "utilities.h"
  30 #include "Utils_ExceptHandlers.hxx"
  31
  32 #include "SMDS_MeshElement.hxx"
  33 #include "SMDS_VolumeOfNodes.hxx"
  34 #include "SMDS_VolumeTool.hxx"
  35 #include "SMESHDS_SubMesh.hxx"
  36 #include "SMESH_Mesh.hxx"
  37 #include "SMESH_subMesh.hxx"
  38
  39 #include <BRep_Tool.hxx>
  40 #include <TopAbs_ShapeEnum.hxx>
  41 #include <TopExp.hxx>
  42 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
  43 #include <TopTools_IndexedMapOfShape.hxx>
  44 #include <TopTools_IndexedMapOfShape.hxx>
  45 #include <TopTools_ListIteratorOfListOfShape.hxx>
  46 #include <TopTools_ListOfShape.hxx>
  47 #include <TopoDS.hxx>
  48 #include <TopoDS_Edge.hxx>
  49 #include <TopoDS_Shell.hxx>
  50 #include <TopoDS_Vertex.hxx>
  51 #include <gp_Pnt.hxx>
  52
  53 #include <stdio.h>
  54 #include <algorithm>
  55
  56 using namespace std;
  57
  58 typedef map < int, int, less<int> >::iterator   \
  59   StdMeshers_IteratorOfDataMapOfIntegerInteger;
  60
  61 //=======================================================================
  62 //
  63 //           StdMeshers_Penta_3D
  64 //
  65 //=======================================================================
  66 //function : StdMeshers_Penta_3D
  67 //purpose  :
  68 //=======================================================================
  69 StdMeshers_Penta_3D::StdMeshers_Penta_3D()
  70 : myErrorStatus(1)
  71 {
  72   myTol3D=0.1;
  73 }
  74 //=======================================================================
  75 //function : Compute
  76 //purpose  :
  77 //=======================================================================
  78 bool StdMeshers_Penta_3D::Compute(SMESH_Mesh& aMesh,
  79                                   const TopoDS_Shape& aShape)
  80 {
  81   myErrorStatus=0;
  82   //
  83   bool bOK=false;
  84   //
  85   myShape=aShape;
  86   SetMesh(aMesh);
  87   //
  88   CheckData();
  89   if (myErrorStatus){
  90     return bOK;
  91   }
  92   //
  93   MakeBlock();
  94     if (myErrorStatus){
  95     return bOK;
  96   }
  97   //
  98   MakeNodes();
  99   if (myErrorStatus){
 100     return bOK;
 101   }
 102   //
 103   MakeConnectingMap();
 104   //
 105   ClearMeshOnFxy1();
 106   if (myErrorStatus) {
 107     return bOK;
 108   }
 109   //
 110   MakeMeshOnFxy1();
 111   if (myErrorStatus) {
 112     return bOK;
 113   }
 114   //
 115   MakeVolumeMesh();
 116   //
 117   return !bOK;
 118 }
 119 //=======================================================================
 120 //function : MakeNodes
 121 //purpose  :
 122 //=======================================================================
 123 void StdMeshers_Penta_3D::MakeNodes()
 124 {
 125   myErrorStatus=0;
 126   //
 127   const int aNbSIDs=9;
 128   int i, j, k, ij, iNbN, aNodeID, aSize, iErr;
 129   double aX, aY, aZ;
 130   SMESH_Block::TShapeID aSID, aSIDs[aNbSIDs]={
 131     SMESH_Block::ID_V000, SMESH_Block::ID_V100,
 132     SMESH_Block::ID_V110, SMESH_Block::ID_V010,
 133     SMESH_Block::ID_Ex00, SMESH_Block::ID_E1y0,
 134     SMESH_Block::ID_Ex10, SMESH_Block::ID_E0y0,
 135     SMESH_Block::ID_Fxy0
 136   };
 137   //
 138   SMESH_Mesh* pMesh=GetMesh();
 139   //
 140   // 1. Define the sizes of mesh
 141   //
 142   // 1.1 Horizontal size
 143   myJSize=0;
 144   for (i=0; i<aNbSIDs; ++i) {
 145     const TopoDS_Shape& aS=myBlock.Shape(aSIDs[i]);
 146     SMESH_subMesh *aSubMesh = pMesh->GetSubMeshContaining(aS);
 147     ASSERT(aSubMesh);
 148     SMESHDS_SubMesh *aSM=aSubMesh->GetSubMeshDS();
 149     iNbN=aSM->NbNodes();
 150     myJSize+=iNbN;
 151   }
 152   //printf("***  Horizontal: number of nodes summary=%d\n", myJSize);
 153   //
 154   // 1.2 Vertical size
 155   myISize=2;
 156   {
 157     const TopoDS_Shape& aS=myBlock.Shape(SMESH_Block::ID_E00z);
 158     SMESH_subMesh *aSubMesh = pMesh->GetSubMeshContaining(aS);
 159     ASSERT(aSubMesh);
 160     SMESHDS_SubMesh *aSM=aSubMesh->GetSubMeshDS();
 161     iNbN=aSM->NbNodes();
 162     myISize+=iNbN;
 163   }
 164   //printf("***  Vertical: number of nodes on edges and vertices=%d\n", myISize);
 165   //
 166   aSize=myISize*myJSize;
 167   myTNodes.resize(aSize);
 168   //
 169   StdMeshers_TNode aTNode;
 170   gp_XYZ aCoords;
 171   gp_Pnt aP3D;
 172   //
 173   // 2. Fill the repers on base face (Z=0)
 174   i=0; j=0;
 175   // vertices
 176   for (k=0; k<aNbSIDs; ++k) {
 177     aSID=aSIDs[k];
 178     const TopoDS_Shape& aS=myBlock.Shape(aSID);
 179     SMDS_NodeIteratorPtr ite =pMesh->GetSubMeshContaining(aS)->GetSubMeshDS()->GetNodes();
 180     while(ite->more()) {
 181       const SMDS_MeshNode* aNode = ite->next();
 182       aNodeID=aNode->GetID();
 183       //
 184       aTNode.SetNode(aNode);
 185       aTNode.SetShapeSupportID(aSID);
 186       aTNode.SetBaseNodeID(aNodeID);
 187       //
 188       switch (aSID){
 189         case SMESH_Block::ID_V000:
 190           aCoords.SetCoord(0., 0., 0.);
 191           break;
 192         case SMESH_Block::ID_V100:
 193           aCoords.SetCoord(1., 0., 0.);
 194           break;
 195         case SMESH_Block::ID_V110:
 196           aCoords.SetCoord(1., 1., 0.);
 197           break;
 198         case SMESH_Block::ID_V010:
 199           aCoords.SetCoord(0., 1., 0.);
 200           break;
 201         case SMESH_Block::ID_Ex00:
 202         case SMESH_Block::ID_E1y0:
 203         case SMESH_Block::ID_Ex10:
 204         case SMESH_Block::ID_E0y0:
 205         case SMESH_Block::ID_Fxy0:{
 206           aX=aNode->X();
 207           aY=aNode->Y();
 208           aZ=aNode->Z();
 209           aP3D.SetCoord(aX, aY, aZ);
 210           myBlock.ComputeParameters(aP3D, aS, aCoords);
 211           iErr=myBlock.ErrorStatus();
 212           if (iErr) {
 213             MESSAGE("StdMeshers_Penta_3D::MakeNodes()," <<
 214                     "SMESHBlock: ComputeParameters operation failed");
 215             myErrorStatus=101; // SMESHBlock: ComputeParameters operation failed
 216             return;
 217           }
 218         }
 219           break;
 220         default:
 221           break;
 222       }
 223       aTNode.SetNormCoord(aCoords);
 224       ij=i*myJSize+j;
 225       myTNodes[ij]=aTNode;
 226       ++j;
 227     }
 228   }
 229   /*
 230   //DEB
 231   {
 232     int iShapeSupportID, iBaseNodeID;
 233     //
 234     //printf("\n\n*** Base Face\n");
 235     i=0;
 236     for (j=0; j<myJSize; ++j) {
 237       ij=i*myJSize+j;
 238       const StdMeshers_TNode& aTNode=myTNodes[ij];
 239       iShapeSupportID=aTNode.ShapeSupportID();
 240       iBaseNodeID=aTNode.BaseNodeID();
 241       const gp_XYZ& aXYZ=aTNode.NormCoord();
 242       printf("*** j:%d bID#%d iSS:%d { %lf %lf %lf }\n",
 243              j,  iBaseNodeID, iShapeSupportID, aXYZ.X(),  aXYZ.Y(), aXYZ.Z());
 244     }
 245   }
 246   */
 247   //DEB
 248   //return; //zz
 249   //
 250   // 3. Finding of Z-layers
 251   vector<double> aZL(myISize);
 252   vector<double>::iterator aItZL1, aItZL2 ;
 253   //
 254   const TopoDS_Shape& aE00z=myBlock.Shape(SMESH_Block::ID_E00z);
 255   SMDS_NodeIteratorPtr aItaE00z =
 256     pMesh->GetSubMeshContaining(aE00z)->GetSubMeshDS()->GetNodes();
 257   //
 258   aZL[0]=0.;
 259   i=1;
 260   while (aItaE00z->more()) {
 261     const SMDS_MeshNode* aNode=aItaE00z->next();
 262     aX=aNode->X(); aY=aNode->Y(); aZ=aNode->Z();
 263     aP3D.SetCoord(aX, aY, aZ);
 264     myBlock.ComputeParameters(aP3D, aE00z, aCoords);
 265     iErr=myBlock.ErrorStatus();
 266     if (iErr) {
 267       MESSAGE("StdMeshers_Penta_3D::MakeNodes()," <<
 268               "SMESHBlock: ComputeParameters operation failed");
 269       myErrorStatus=101; // SMESHBlock: ComputeParameters operation failed
 270       return;
 271     }
 272     aZL[i]=aCoords.Z();
 273     ++i;
 274   }
 275   aZL[i]=1.;
 276   //
 277   aItZL1=aZL.begin();
 278   aItZL2=aZL.end();
 279   //
 280   // Sorting the layers
 281   sort(aItZL1, aItZL2);
 282   //DEB
 283   /*
 284   printf("** \n\n Layers begin\n");
 285   for(i=0, aItZL=aItZL1; aItZL!=aItZL2; ++aItZL, ++i) {
 286     printf(" #%d : %lf\n", i, *aItZL);
 287   }
 288   printf("** Layers end\n");
 289   */
 290   //DEB
 291   //
 292   //
 293   // 4. Fill the rest repers
 294   bool bIsUpperLayer;
 295   int iBNID;
 296   SMESH_Block::TShapeID aSSID, aBNSSID;
 297   StdMeshers_TNode aTN;
 298   //
 299   for (j=0; j<myJSize; ++j) {
 300     for (i=1; i<myISize; ++i) {
 301       //
 302       // base node info
 303       const StdMeshers_TNode& aBN=myTNodes[j];
 304       aBNSSID=(SMESH_Block::TShapeID)aBN.ShapeSupportID();
 305       iBNID=aBN.BaseNodeID();
 306       const gp_XYZ& aBNXYZ=aBN.NormCoord();
 307       //
 308       // fill current node info
 309       //   -index in aTNodes
 310       ij=i*myJSize+j;
 311       //   -normalized coordinates
 312       aX=aBNXYZ.X();
 313       aY=aBNXYZ.Y();
 314       aZ=aZL[i];
 315       aCoords.SetCoord(aX, aY, aZ);
 316       //
 317       //   suporting shape ID
 318       bIsUpperLayer=(i==(myISize-1));
 319       ShapeSupportID(bIsUpperLayer, aBNSSID, aSSID);
 320       if (myErrorStatus) {
 321         MESSAGE("StdMeshers_Penta_3D::MakeNodes() ");
 322         return;
 323       }
 324       //
 325       aTN.SetShapeSupportID(aSSID);
 326       aTN.SetNormCoord(aCoords);
 327       aTN.SetBaseNodeID(iBNID);
 328       //
 329       if (aSSID!=SMESH_Block::ID_NONE){
 330         // try to find the node
 331         const TopoDS_Shape& aS=myBlock.Shape((int)aSSID);
 332         FindNodeOnShape(aS, aCoords, aTN);
 333       }
 334       else{
 335         // create node and get it id
 336         CreateNode (bIsUpperLayer, aCoords, aTN);
 337       }
 338       if (myErrorStatus) {
 339         MESSAGE("StdMeshers_Penta_3D::MakeNodes() ");
 340         return;
 341       }
 342       //
 343       myTNodes[ij]=aTN;
 344     }
 345   }
 346   //DEB
 347   /*
 348   {
 349     int iSSID, iBNID, aID;
 350     //
 351     for (i=0; i<myISize; ++i) {
 352       printf(" Layer# %d\n", i);
 353       for (j=0; j<myJSize; ++j) {
 354         ij=i*myJSize+j;
 355         const StdMeshers_TNode& aTN=myTNodes[ij];
 356         //const StdMeshers_TNode& aTN=aTNodes[ij];
 357         const gp_XYZ& aXYZ=aTN.NormCoord();
 358         iSSID=aTN.ShapeSupportID();
 359         iBNID=aTN.BaseNodeID();
 360         //
 361         const SMDS_MeshNode* aNode=aTN.Node();
 362         aID=aNode->GetID();
 363         aX=aNode->X();
 364         aY=aNode->Y();
 365         aZ=aNode->Z();
 366         printf("*** j:%d BNID#%d iSSID:%d ID:%d { %lf %lf %lf },  { %lf %lf %lf }\n",
 367                j,  iBNID, iSSID, aID, aXYZ.X(),  aXYZ.Y(), aXYZ.Z(), aX, aY, aZ);
 368       }
 369     }
 370   }
 371   */
 372   //DEB t
 373 }
 374 //=======================================================================
 375 //function : FindNodeOnShape
 376 //purpose  :
 377 //=======================================================================
 378 void StdMeshers_Penta_3D::FindNodeOnShape(const TopoDS_Shape& aS,
 379                                           const gp_XYZ& aParams,
 380                                           StdMeshers_TNode& aTN)
 381 {
 382   myErrorStatus=0;
 383   //
 384   double aX, aY, aZ, aD, aTol2, minD;
 385   gp_Pnt aP1, aP2;
 386   //
 387   SMESH_Mesh* pMesh=GetMesh();
 388   aTol2=myTol3D*myTol3D;
 389   minD = 1.e100;
 390   SMDS_MeshNode* pNode=NULL;
 391   //
 392   myBlock.Point(aParams, aS, aP1);
 393   //
 394   SMDS_NodeIteratorPtr ite=
 395     pMesh->GetSubMeshContaining(aS)->GetSubMeshDS()->GetNodes();
 396   while(ite->more()) {
 397     const SMDS_MeshNode* aNode = ite->next();
 398     aX=aNode->X();
 399     aY=aNode->Y();
 400     aZ=aNode->Z();
 401     aP2.SetCoord(aX, aY, aZ);
 402     aD=(double)aP1.SquareDistance(aP2);
 403     //printf("** D=%lf ", aD, aTol2);
 404     if (aD < minD) {
 405       pNode=(SMDS_MeshNode*)aNode;
 406       aTN.SetNode(pNode);
 407       minD = aD;
 408       //printf(" Ok\n");
 409       if (aD<aTol2)
 410         return;
 411     }
 412   }
 413   //
 414   //printf(" KO\n");
 415   //aTN.SetNode(pNode);
 416   //MESSAGE("StdMeshers_Penta_3D::FindNodeOnShape(), can not find the node");
 417   //myErrorStatus=11; // can not find the node;
 418 }
 419
 420 //=======================================================================
 421 //function : MakeVolumeMesh
 422 //purpose  :
 423 //=======================================================================
 424 void StdMeshers_Penta_3D::MakeVolumeMesh()
 425 {
 426   myErrorStatus=0;
 427   //
 428   int i, j, ij, ik, i1, i2, aSSID;
 429   //
 430   TopoDS_Shell aShell;
 431   TopExp_Explorer aExp;
 432   //
 433   SMESH_Mesh*   pMesh =GetMesh();
 434   SMESHDS_Mesh* meshDS=pMesh->GetMeshDS();
 435   //
 436   aExp.Init(myShape, TopAbs_SHELL);
 437   for (; aExp.More(); aExp.Next()){
 438     aShell=TopoDS::Shell(aExp.Current());
 439     break;
 440   }
 441   //
 442   // 1. Set Node In Volume
 443   ik=myISize-1;
 444   for (i=1; i<ik; ++i){
 445     for (j=0; j<myJSize; ++j){
 446       ij=i*myJSize+j;
 447       const StdMeshers_TNode& aTN=myTNodes[ij];
 448       aSSID=aTN.ShapeSupportID();
 449       if (aSSID==SMESH_Block::ID_NONE) {
 450         SMDS_MeshNode* aNode=(SMDS_MeshNode*)aTN.Node();
 451         meshDS->SetNodeInVolume(aNode, aShell);
 452       }
 453     }
 454   }
 455   //
 456   // 2. Make pentahedrons
 457   int aID0, k , aJ[3];
 458   vector<const SMDS_MeshNode*> aN;
 459   //
 460   SMDS_ElemIteratorPtr itf, aItNodes;
 461   //
 462   const TopoDS_Face& aFxy0=
 463     TopoDS::Face(myBlock.Shape(SMESH_Block::ID_Fxy0));
 464   SMESH_subMesh *aSubMesh0 = pMesh->GetSubMeshContaining(aFxy0);
 465   SMESHDS_SubMesh *aSM0=aSubMesh0->GetSubMeshDS();
 466   //
 467   itf=aSM0->GetElements();
 468   while(itf->more()) {
 469     const SMDS_MeshElement* pE0=itf->next();
 470     //
 471     int nbFaceNodes = pE0->NbNodes();
 472     if ( aN.size() < nbFaceNodes * 2 )
 473       aN.resize( nbFaceNodes * 2 );
 474     //
 475     k=0;
 476     aItNodes=pE0->nodesIterator();
 477     while (aItNodes->more()) {
 478       const SMDS_MeshElement* pNode=aItNodes->next();
 479       aID0=pNode->GetID();
 480       aJ[k]=GetIndexOnLayer(aID0);
 481       if (myErrorStatus) {
 482         MESSAGE("StdMeshers_Penta_3D::MakeVolumeMesh");
 483         return;
 484       }
 485       //
 486       ++k;
 487     }
 488     //
 489     bool forward = true;
 490     for (i=0; i<ik; ++i){
 491       i1=i;
 492       i2=i+1;
 493       for(j=0; j<nbFaceNodes; ++j) {
 494         ij=i1*myJSize+aJ[j];
 495         const StdMeshers_TNode& aTN1=myTNodes[ij];
 496         const SMDS_MeshNode* aN1=aTN1.Node();
 497         aN[j]=aN1;
 498         //
 499         ij=i2*myJSize+aJ[j];
 500         const StdMeshers_TNode& aTN2=myTNodes[ij];
 501         const SMDS_MeshNode* aN2=aTN2.Node();
 502         aN[j+nbFaceNodes]=aN2;
 503       }
 504       // check if volume orientation will be ok
 505       if ( i == 0 ) {
 506         SMDS_VolumeTool vTool;
 507         switch ( nbFaceNodes ) {
 508         case 3: {
 509           SMDS_VolumeOfNodes tmpVol (aN[0], aN[1], aN[2],
 510                                      aN[3], aN[4], aN[5]);
 511           vTool.Set( &tmpVol );
 512           break;
 513         }
 514         case 4: {
 515           SMDS_VolumeOfNodes tmpVol(aN[0], aN[1], aN[2], aN[3],
 516                                     aN[4], aN[5], aN[6], aN[7]);
 517           vTool.Set( &tmpVol );
 518           break;
 519         }
 520         default:
 521           continue;
 522         }
 523         forward = vTool.IsForward();
 524       }
 525       // add volume
 526       SMDS_MeshVolume* aV = 0;
 527       switch ( nbFaceNodes ) {
 528       case 3:
 529         if ( forward )
 530           aV = meshDS->AddVolume(aN[0], aN[1], aN[2],
 531                                  aN[3], aN[4], aN[5]);
 532         else
 533           aV = meshDS->AddVolume(aN[0], aN[2], aN[1],
 534                                  aN[3], aN[5], aN[4]);
 535         break;
 536       case 4:
 537         if ( forward )
 538           aV = meshDS->AddVolume(aN[0], aN[1], aN[2], aN[3],
 539                                  aN[4], aN[5], aN[6], aN[7]);
 540         else
 541           aV = meshDS->AddVolume(aN[0], aN[3], aN[2], aN[1],
 542                                  aN[4], aN[7], aN[6], aN[5]);
 543         break;
 544       default:
 545         continue;
 546       }
 547       meshDS->SetMeshElementOnShape(aV, aShell);
 548     }
 549   }
 550 }
 551
 552 //=======================================================================
 553 //function : MakeMeshOnFxy1
 554 //purpose  :
 555 //=======================================================================
 556 void StdMeshers_Penta_3D::MakeMeshOnFxy1()
 557 {
 558   myErrorStatus=0;
 559   //
 560   int aID0, aJ, aLevel, ij, aNbNodes, k;
 561   //
 562   SMDS_NodeIteratorPtr itn;
 563   SMDS_ElemIteratorPtr itf, aItNodes;
 564   SMDSAbs_ElementType aElementType;
 565   //
 566   const TopoDS_Face& aFxy0=
 567     TopoDS::Face(myBlock.Shape(SMESH_Block::ID_Fxy0));
 568   const TopoDS_Face& aFxy1=
 569     TopoDS::Face(myBlock.Shape(SMESH_Block::ID_Fxy1));
 570   //
 571   SMESH_Mesh* pMesh=GetMesh();
 572   SMESHDS_Mesh * meshDS = pMesh->GetMeshDS();
 573   //
 574   SMESH_subMesh *aSubMesh0 = pMesh->GetSubMeshContaining(aFxy0);
 575   SMESHDS_SubMesh *aSM0=aSubMesh0->GetSubMeshDS();
 576   //
 577   // set nodes on aFxy1
 578   aLevel=myISize-1;
 579   itn=aSM0->GetNodes();
 580   aNbNodes=aSM0->NbNodes();
 581   //printf("** aNbNodes=%d\n", aNbNodes);
 582   while(itn->more()) {
 583     const SMDS_MeshNode* aN0=itn->next();
 584     aID0=aN0->GetID();
 585     aJ=GetIndexOnLayer(aID0);
 586     if (myErrorStatus) {
 587       MESSAGE("StdMeshers_Penta_3D::MakeMeshOnFxy1() ");
 588       return;
 589     }
 590     //
 591     ij=aLevel*myJSize+aJ;
 592     const StdMeshers_TNode& aTN1=myTNodes[ij];
 593     SMDS_MeshNode* aN1=(SMDS_MeshNode*)aTN1.Node();
 594     //
 595     meshDS->SetNodeOnFace(aN1, aFxy1);
 596   }
 597   //
 598   // set elements on aFxy1
 599   vector<const SMDS_MeshNode*> aNodes1;
 600   //
 601   itf=aSM0->GetElements();
 602   while(itf->more()) {
 603     const SMDS_MeshElement * pE0=itf->next();
 604     aElementType=pE0->GetType();
 605     if (!aElementType==SMDSAbs_Face) {
 606       continue;
 607     }
 608     aNbNodes=pE0->NbNodes();
 609 //     if (aNbNodes!=3) {
 610 //       continue;
 611 //     }
 612     if ( aNodes1.size() < aNbNodes )
 613       aNodes1.resize( aNbNodes );
 614     //
 615     k=aNbNodes-1; // reverse a face
 616     aItNodes=pE0->nodesIterator();
 617     while (aItNodes->more()) {
 618       const SMDS_MeshElement* pNode=aItNodes->next();
 619       aID0=pNode->GetID();
 620       aJ=GetIndexOnLayer(aID0);
 621       if (myErrorStatus) {
 622         MESSAGE("StdMeshers_Penta_3D::MakeMeshOnFxy1() ");
 623         return;
 624       }
 625       //
 626       ij=aLevel*myJSize+aJ;
 627       const StdMeshers_TNode& aTN1=myTNodes[ij];
 628       const SMDS_MeshNode* aN1=aTN1.Node();
 629       aNodes1[k]=aN1;
 630       --k;
 631     }
 632     SMDS_MeshFace * face = 0;
 633     switch ( aNbNodes ) {
 634     case 3:
 635       face = meshDS->AddFace(aNodes1[0], aNodes1[1], aNodes1[2]);
 636       break;
 637     case 4:
 638       face = meshDS->AddFace(aNodes1[0], aNodes1[1], aNodes1[2], aNodes1[3]);
 639       break;
 640     default:
 641       continue;
 642     }
 643     meshDS->SetMeshElementOnShape(face, aFxy1);
 644   }
 645 }
 646 //=======================================================================
 647 //function : ClearMeshOnFxy1
 648 //purpose  :
 649 //=======================================================================
 650 void StdMeshers_Penta_3D::ClearMeshOnFxy1()
 651 {
 652   myErrorStatus=0;
 653   //
 654   SMESH_subMesh* aSubMesh;
 655   SMESH_Mesh* pMesh=GetMesh();
 656   //
 657   const TopoDS_Shape& aFxy1=myBlock.Shape(SMESH_Block::ID_Fxy1);
 658   aSubMesh = pMesh->GetSubMeshContaining(aFxy1);
 659   if (aSubMesh)
 660     aSubMesh->ComputeStateEngine( SMESH_subMesh::CLEAN );
 661 }
 662
 663 //=======================================================================
 664 //function : GetIndexOnLayer
 665 //purpose  :
 666 //=======================================================================
 667 int StdMeshers_Penta_3D::GetIndexOnLayer(const int aID)
 668 {
 669   myErrorStatus=0;
 670   //
 671   int j=-1;
 672   StdMeshers_IteratorOfDataMapOfIntegerInteger aMapIt;
 673   //
 674   aMapIt=myConnectingMap.find(aID);
 675   if (aMapIt==myConnectingMap.end()) {
 676     myErrorStatus=200;
 677     return j;
 678   }
 679   j=(*aMapIt).second;
 680   return j;
 681 }
 682 //=======================================================================
 683 //function : MakeConnectingMap
 684 //purpose  :
 685 //=======================================================================
 686 void StdMeshers_Penta_3D::MakeConnectingMap()
 687 {
 688   int j, aBNID;
 689   //
 690   for (j=0; j<myJSize; ++j) {
 691     const StdMeshers_TNode& aBN=myTNodes[j];
 692     aBNID=aBN.BaseNodeID();
 693     myConnectingMap[aBNID]=j;
 694   }
 695 }
 696 //=======================================================================
 697 //function : CreateNode
 698 //purpose  :
 699 //=======================================================================
 700 void StdMeshers_Penta_3D::CreateNode(const bool bIsUpperLayer,
 701                                      const gp_XYZ& aParams,
 702                                      StdMeshers_TNode& aTN)
 703 {
 704   myErrorStatus=0;
 705   //
 706   int iErr;
 707   double aX, aY, aZ;
 708   //
 709   gp_Pnt aP;
 710   //
 711   SMDS_MeshNode* pNode=NULL;
 712   aTN.SetNode(pNode);
 713   //
 714   if (bIsUpperLayer) {
 715     // point on face Fxy1
 716     const TopoDS_Shape& aS=myBlock.Shape(SMESH_Block::ID_Fxy1);
 717     myBlock.Point(aParams, aS, aP);
 718   }
 719   else {
 720     // point inside solid
 721     myBlock.Point(aParams, aP);
 722   }
 723   //
 724   iErr=myBlock.ErrorStatus();
 725   if (iErr) {
 726     myErrorStatus=12; // can not find the node point;
 727     return;
 728   }
 729   //
 730   aX=aP.X(); aY=aP.Y(); aZ=aP.Z();
 731   //
 732   SMESH_Mesh* pMesh=GetMesh();
 733   SMESHDS_Mesh* pMeshDS=pMesh->GetMeshDS();
 734   //
 735   pNode = pMeshDS->AddNode(aX, aY, aZ);
 736   aTN.SetNode(pNode);
 737 }
 738 //=======================================================================
 739 //function : ShapeSupportID
 740 //purpose  :
 741 //=======================================================================
 742 void StdMeshers_Penta_3D::ShapeSupportID(const bool bIsUpperLayer,
 743                                          const SMESH_Block::TShapeID aBNSSID,
 744                                          SMESH_Block::TShapeID& aSSID)
 745 {
 746   myErrorStatus=0;
 747   //
 748   switch (aBNSSID) {
 749     case SMESH_Block::ID_V000:
 750       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_V001 : SMESH_Block::ID_E00z;
 751       break;
 752     case SMESH_Block::ID_V100:
 753       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_V101 : SMESH_Block::ID_E10z;
 754       break;
 755     case SMESH_Block::ID_V110:
 756       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_V111 : SMESH_Block::ID_E11z;
 757       break;
 758     case SMESH_Block::ID_V010:
 759       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_V011 : SMESH_Block::ID_E01z;
 760       break;
 761     case SMESH_Block::ID_Ex00:
 762       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_Ex01 : SMESH_Block::ID_Fx0z;
 763       break;
 764     case SMESH_Block::ID_Ex10:
 765       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_Ex11 : SMESH_Block::ID_Fx1z;
 766       break;
 767     case SMESH_Block::ID_E0y0:
 768       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_E0y1 : SMESH_Block::ID_F0yz;
 769       break;
 770     case SMESH_Block::ID_E1y0:
 771       aSSID=(bIsUpperLayer) ?  SMESH_Block::ID_E1y1 : SMESH_Block::ID_F1yz;
 772       break;
 773     case SMESH_Block::ID_Fxy0:
 774       aSSID=SMESH_Block::ID_NONE;//(bIsUpperLayer) ?  Shape_ID_Fxy1 : Shape_ID_NONE;
 775       break;
 776     default:
 777       aSSID=SMESH_Block::ID_NONE;
 778       myErrorStatus=10; // Can not find supporting shape ID
 779       break;
 780   }
 781   return;
 782 }
 783 //=======================================================================
 784 //function : MakeBlock
 785 //purpose  :
 786 //=======================================================================
 787 void StdMeshers_Penta_3D::MakeBlock()
 788 {
 789   myErrorStatus=0;
 790   //
 791   bool bFound;
 792   int i, j, iNbEV, iNbE, iErr, iCnt, iNbNodes, iNbF;
 793   //
 794   TopoDS_Vertex aV000, aV001;
 795   TopoDS_Shape aFTr;
 796   TopTools_IndexedDataMapOfShapeListOfShape aMVES;
 797   TopTools_IndexedMapOfShape aME ,aMEV, aM;
 798   TopTools_ListIteratorOfListOfShape aIt;
 799   //
 800   TopExp::MapShapes(myShape, TopAbs_FACE, aM);
 801   //
 802   // 0. Find triangulated face aFTr
 803   SMDSAbs_ElementType aElementType;
 804   SMESH_Mesh* pMesh=GetMesh();
 805   //
 806   iCnt=0;
 807   iNbF=aM.Extent();
 808   for (i=1; i<=iNbF; ++i) {
 809     const TopoDS_Shape& aF=aM(i);
 810     SMESH_subMesh *aSubMesh = pMesh->GetSubMeshContaining(aF);
 811     ASSERT(aSubMesh);
 812     SMESHDS_SubMesh *aSM=aSubMesh->GetSubMeshDS();
 813     SMDS_ElemIteratorPtr itf=aSM->GetElements();
 814     while(itf->more()) {
 815       const SMDS_MeshElement * pElement=itf->next();
 816       aElementType=pElement->GetType();
 817       if (aElementType==SMDSAbs_Face) {
 818         iNbNodes=pElement->NbNodes();
 819         if (iNbNodes==3) {
 820           aFTr=aF;
 821           ++iCnt;
 822           if (iCnt>1) {
 823             MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 824             myErrorStatus=5; // more than one face has triangulation
 825             return;
 826           }
 827           break; // next face
 828         }
 829       }
 830     }
 831   }
 832   //
 833   // 1. Vetrices V00, V001;
 834   //
 835   TopExp::MapShapes(aFTr, TopAbs_EDGE, aME);
 836   TopExp::MapShapesAndAncestors(myShape, TopAbs_VERTEX, TopAbs_EDGE, aMVES);
 837   //
 838   // 1.1 Base vertex V000
 839   iNbE=aME.Extent();
 840   if (iNbE!=4){
 841     MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 842     myErrorStatus=7; // too few edges are in base face aFTr
 843     return;
 844   }
 845   const TopoDS_Edge& aE1=TopoDS::Edge(aME(1));
 846   aV000=TopExp::FirstVertex(aE1);
 847   //
 848   const TopTools_ListOfShape& aLE=aMVES.FindFromKey(aV000);
 849   aIt.Initialize(aLE);
 850   for (; aIt.More(); aIt.Next()) {
 851     const TopoDS_Shape& aEx=aIt.Value();
 852     aMEV.Add(aEx);
 853   }
 854   iNbEV=aMEV.Extent();
 855   if (iNbEV!=3){
 856     MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 857     myErrorStatus=7; // too few edges meet in base vertex
 858     return;
 859   }
 860   //
 861   // 1.2 Vertex V001
 862   bFound=false;
 863   for (j=1; j<=iNbEV; ++j) {
 864     const TopoDS_Shape& aEx=aMEV(j);
 865     if (!aME.Contains(aEx)) {
 866       TopoDS_Vertex aV[2];
 867       //
 868       const TopoDS_Edge& aE=TopoDS::Edge(aEx);
 869       TopExp::Vertices(aE, aV[0], aV[1]);
 870       for (i=0; i<2; ++i) {
 871         if (!aV[i].IsSame(aV000)) {
 872           aV001=aV[i];
 873           bFound=!bFound;
 874           break;
 875         }
 876       }
 877     }
 878   }
 879   //
 880   if (!bFound) {
 881     MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 882     myErrorStatus=8; // can not find reper V001
 883     return;
 884   }
 885   //DEB
 886   //gp_Pnt aP000, aP001;
 887   //
 888   //aP000=BRep_Tool::Pnt(TopoDS::Vertex(aV000));
 889   //printf("*** aP000 { %lf, %lf, %lf }\n", aP000.X(), aP000.Y(), aP000.Z());
 890   //aP001=BRep_Tool::Pnt(TopoDS::Vertex(aV001));
 891   //printf("*** aP001 { %lf, %lf, %lf }\n", aP001.X(), aP001.Y(), aP001.Z());
 892   //DEB
 893   //
 894   aME.Clear();
 895   TopExp::MapShapes(myShape, TopAbs_SHELL, aME);
 896   iNbE=aME.Extent();
 897   if (iNbE!=1) {
 898     MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 899     myErrorStatus=9; // number of shells in source shape !=1
 900     return;
 901   }
 902   //
 903   // 2. Load Block
 904   const TopoDS_Shell& aShell=TopoDS::Shell(aME(1));
 905   myBlock.Load(aShell, aV000, aV001);
 906   iErr=myBlock.ErrorStatus();
 907   if (iErr) {
 908     MESSAGE("StdMeshers_Penta_3D::MakeBlock() ");
 909     myErrorStatus=100; // SMESHBlock: Load operation failed
 910     return;
 911   }
 912 }
 913 //=======================================================================
 914 //function : CheckData
 915 //purpose  :
 916 //=======================================================================
 917 void StdMeshers_Penta_3D::CheckData()
 918 {
 919   myErrorStatus=0;
 920   //
 921   int i, iNb;
 922   int iNbEx[]={8, 12, 6};
 923   //
 924   TopAbs_ShapeEnum aST;
 925   TopAbs_ShapeEnum aSTEx[]={
 926     TopAbs_VERTEX, TopAbs_EDGE, TopAbs_FACE
 927   };
 928   TopTools_IndexedMapOfShape aM;
 929   //
 930   if (myShape.IsNull()){
 931     MESSAGE("StdMeshers_Penta_3D::CheckData() ");
 932     myErrorStatus=2; // null shape
 933     return;
 934   }
 935   //
 936   aST=myShape.ShapeType();
 937   if (!(aST==TopAbs_SOLID || aST==TopAbs_SHELL)) {
 938     MESSAGE("StdMeshers_Penta_3D::CheckData() ");
 939     myErrorStatus=3; // not compatible type of shape
 940     return;
 941   }
 942   //
 943   for (i=0; i<3; ++i) {
 944     aM.Clear();
 945     TopExp::MapShapes(myShape, aSTEx[i], aM);
 946     iNb=aM.Extent();
 947     if (iNb!=iNbEx[i]){
 948       MESSAGE("StdMeshers_Penta_3D::CheckData() ");
 949       myErrorStatus=4; // number of subshape is not compatible
 950       return;
 951     }
 952   }
 953 }
 954 //////////////////////////////////////////////////////////////////////////
 955 //
 956 //   StdMeshers_SMESHBlock
 957 //
 958 //
 959 #include <TopTools_IndexedMapOfOrientedShape.hxx>
 960 #include <TopoDS_Vertex.hxx>
 961
 962 //=======================================================================
 963 //function : StdMeshers_SMESHBlock
 964 //purpose  :
 965 //=======================================================================
 966 StdMeshers_SMESHBlock::StdMeshers_SMESHBlock()
 967 {
 968   myErrorStatus=1;
 969 }
 970 //=======================================================================
 971 //function : ErrorStatus
 972 //purpose  :
 973 //=======================================================================
 974 int StdMeshers_SMESHBlock::ErrorStatus() const
 975 {
 976   return myErrorStatus;
 977 }
 978 //=======================================================================
 979 //function : Load
 980 //purpose  :
 981 //=======================================================================
 982 void StdMeshers_SMESHBlock::Load(const TopoDS_Shell& theShell)
 983 {
 984
 985   TopoDS_Vertex aV000, aV001;
 986   //
 987   Load(theShell, aV000, aV001);
 988 }
 989 //=======================================================================
 990 //function : Load
 991 //purpose  :
 992 //=======================================================================
 993 void StdMeshers_SMESHBlock::Load(const TopoDS_Shell& theShell,
 994                                  const TopoDS_Vertex& theV000,
 995                                  const TopoDS_Vertex& theV001)
 996 {
 997   myErrorStatus=0;
 998   //
 999   myShell=theShell;
1000   //
1001   bool bOk;
1002   //
1003   myShapeIDMap.Clear();
1004   bOk=myTBlock.LoadBlockShapes(myShell, theV000, theV001, myShapeIDMap);
1005   if (!bOk) {
1006     myErrorStatus=2;
1007     return;
1008   }
1009 }
1010 //=======================================================================
1011 //function : ComputeParameters
1012 //purpose  :
1013 //=======================================================================
1014 void StdMeshers_SMESHBlock::ComputeParameters(const gp_Pnt& thePnt,
1015                                               gp_XYZ& theXYZ)
1016 {
1017   ComputeParameters(thePnt, myShell, theXYZ);
1018 }
1019 //=======================================================================
1020 //function : ComputeParameters
1021 //purpose  :
1022 //=======================================================================
1023 void StdMeshers_SMESHBlock::ComputeParameters(const gp_Pnt& thePnt,
1024                                               const TopoDS_Shape& theShape,
1025                                               gp_XYZ& theXYZ)
1026 {
1027   myErrorStatus=0;
1028   //
1029   int aID;
1030   bool bOk;
1031   //
1032   aID=ShapeID(theShape);
1033   if (myErrorStatus) {
1034     return;
1035   }
1036   bOk=myTBlock.ComputeParameters(thePnt, theXYZ, aID);
1037   if (!bOk) {
1038     myErrorStatus=4; // problems with computation Parameters
1039     return;
1040   }
1041 }
1042 //=======================================================================
1043 //function : Point
1044 //purpose  :
1045 //=======================================================================
1046  void StdMeshers_SMESHBlock::Point(const gp_XYZ& theParams,
1047                                    gp_Pnt& aP3D)
1048 {
1049   TopoDS_Shape aS;
1050   //
1051   Point(theParams, aS, aP3D);
1052 }
1053 //=======================================================================
1054 //function : Point
1055 //purpose  :
1056 //=======================================================================
1057  void StdMeshers_SMESHBlock::Point(const gp_XYZ& theParams,
1058                                    const TopoDS_Shape& theShape,
1059                                    gp_Pnt& aP3D)
1060 {
1061   myErrorStatus=0;
1062   //
1063   int aID;
1064   bool bOk=false;
1065   gp_XYZ aXYZ(99.,99.,99.);
1066   aP3D.SetXYZ(aXYZ);
1067   //
1068   if (theShape.IsNull()) {
1069     bOk=myTBlock.ShellPoint(theParams, aXYZ);
1070   }
1071   //
1072   else {
1073     aID=ShapeID(theShape);
1074     if (myErrorStatus) {
1075       return;
1076     }
1077     //
1078     if (SMESH_Block::IsVertexID(aID)) {
1079       bOk=myTBlock.VertexPoint(aID, aXYZ);
1080     }
1081     else if (SMESH_Block::IsEdgeID(aID)) {
1082       bOk=myTBlock.EdgePoint(aID, theParams, aXYZ);
1083     }
1084     //
1085     else if (SMESH_Block::IsFaceID(aID)) {
1086       bOk=myTBlock.FacePoint(aID, theParams, aXYZ);
1087     }
1088   }
1089   if (!bOk) {
1090     myErrorStatus=4; // problems with point computation
1091     return;
1092   }
1093   aP3D.SetXYZ(aXYZ);
1094 }
1095 //=======================================================================
1096 //function : ShapeID
1097 //purpose  :
1098 //=======================================================================
1099 int StdMeshers_SMESHBlock::ShapeID(const TopoDS_Shape& theShape)
1100 {
1101   myErrorStatus=0;
1102   //
1103   int aID=-1;
1104   TopoDS_Shape aSF, aSR;
1105   //
1106   aSF=theShape;
1107   aSF.Orientation(TopAbs_FORWARD);
1108   aSR=theShape;
1109   aSR.Orientation(TopAbs_REVERSED);
1110   //
1111   if (myShapeIDMap.Contains(aSF)) {
1112     aID=myShapeIDMap.FindIndex(aSF);
1113     return aID;
1114   }
1115   if (myShapeIDMap.Contains(aSR)) {
1116     aID=myShapeIDMap.FindIndex(aSR);
1117     return aID;
1118   }
1119   myErrorStatus=2; // unknown shape;
1120   return aID;
1121 }
1122 //=======================================================================
1123 //function : Shape
1124 //purpose  :
1125 //=======================================================================
1126 const TopoDS_Shape& StdMeshers_SMESHBlock::Shape(const int theID)
1127 {
1128   myErrorStatus=0;
1129   //
1130   int aNb;
1131   //
1132   aNb=myShapeIDMap.Extent();
1133   if (theID<1 || theID>aNb) {
1134     myErrorStatus=3; // ID is out of range
1135     return myEmptyShape;
1136   }
1137   //
1138   const TopoDS_Shape& aS=myShapeIDMap.FindKey(theID);
1139   return aS;
1140 }