src/GHS3DPlugin_GHS3D.cxx

   1 //=============================================================================
   2 // File      : GHS3DPlugin_GHS3D.cxx
   3 // Created   :
   4 // Author    : Edward AGAPOV
   5 // Project   : SALOME
   6 // Copyright : CEA 2003
   7 // $Header$
   8 //=============================================================================
   9 using namespace std;
  10
  11 #include "GHS3DPlugin_GHS3D.hxx"
  12 #include "SMESH_Gen.hxx"
  13 #include "SMESH_Mesh.hxx"
  14
  15 #include "SMDS_MeshElement.hxx"
  16 #include "SMDS_MeshNode.hxx"
  17
  18 #include <TopExp_Explorer.hxx>
  19
  20 #include "utilities.h"
  21
  22 #include <qfile.h>
  23
  24 #ifdef _DEBUG_
  25 #define DUMP(txt) \
  26 //  cout << txt
  27 #else
  28 #define DUMP(txt)
  29 #endif
  30
  31 //=============================================================================
  32 /*!
  33  *
  34  */
  35 //=============================================================================
  36
  37 GHS3DPlugin_GHS3D::GHS3DPlugin_GHS3D(int hypId, int studyId, SMESH_Gen* gen)
  38   : SMESH_3D_Algo(hypId, studyId, gen)
  39 {
  40   MESSAGE("GHS3DPlugin_GHS3D::GHS3DPlugin_GHS3D");
  41   _name = "GHS3D_3D";
  42   _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);// 1 bit /shape type
  43 }
  44
  45 //=============================================================================
  46 /*!
  47  *
  48  */
  49 //=============================================================================
  50
  51 GHS3DPlugin_GHS3D::~GHS3DPlugin_GHS3D()
  52 {
  53   MESSAGE("GHS3DPlugin_GHS3D::~GHS3DPlugin_GHS3D");
  54 }
  55
  56 //=============================================================================
  57 /*!
  58  *
  59  */
  60 //=============================================================================
  61
  62 bool GHS3DPlugin_GHS3D::CheckHypothesis
  63                          (SMESH_Mesh& aMesh,
  64                           const TopoDS_Shape& aShape,
  65                           SMESH_Hypothesis::Hypothesis_Status& aStatus)
  66 {
  67 //  MESSAGE("GHS3DPlugin_GHS3D::CheckHypothesis");
  68   aStatus = SMESH_Hypothesis::HYP_OK;
  69   return true;
  70 }
  71
  72 //=======================================================================
  73 //function : writeFaces
  74 //purpose  :
  75 //=======================================================================
  76
  77 static bool writeFaces (ofstream &            theFile,
  78                         SMESHDS_Mesh *        theMesh,
  79                         const map <int,int> & theSmdsToGhs3dIdMap)
  80 {
  81   // record structure:
  82   //
  83   // NB_ELEMS DUMMY_INT
  84   // Loop from 1 to NB_ELEMS
  85   //   NB_NODES NODE_NB_1 NODE_NB_2 ... (NB_NODES + 1) times: DUMMY_INT
  86
  87   int nbFaces = theMesh->NbFaces();
  88   if ( nbFaces == 0 )
  89     return false;
  90
  91   const char* space    = "  ";
  92   const int   dummyint = 0;
  93
  94   // NB_ELEMS DUMMY_INT
  95   theFile << space << nbFaces << space << dummyint << endl;
  96
  97   // Loop from 1 to NB_ELEMS
  98   SMDS_FaceIteratorPtr it = theMesh->facesIterator();
  99   while ( it->more() )
 100   {
 101     // NB_NODES
 102     const SMDS_MeshElement* elem = it->next();
 103     const int nbNodes = elem->NbNodes();
 104     theFile << space << nbNodes;
 105
 106     // NODE_NB_1 NODE_NB_2 ...
 107     SMDS_ElemIteratorPtr nodeIt = elem->nodesIterator();
 108     while ( nodeIt->more() )
 109     {
 110       // find GHS3D ID
 111       int aSmdsID = nodeIt->next()->GetID();
 112       map<int,int>::const_iterator it = theSmdsToGhs3dIdMap.find( aSmdsID );
 113       ASSERT( it != theSmdsToGhs3dIdMap.end() );
 114       theFile << space << (*it).second;
 115     }
 116
 117     // (NB_NODES + 1) times: DUMMY_INT
 118     for ( int i=0; i<=nbNodes; i++)
 119       theFile << space << dummyint;
 120
 121     theFile << endl;
 122   }
 123
 124   return true;
 125 }
 126
 127 //=======================================================================
 128 //function : writePoints
 129 //purpose  :
 130 //=======================================================================
 131
 132 static bool writePoints (ofstream &                       theFile,
 133                          SMESHDS_Mesh *                   theMesh,
 134                          map <int,int> &                  theSmdsToGhs3dIdMap,
 135                          map <int,const SMDS_MeshNode*> & theGhs3dIdToNodeMap)
 136 {
 137   // record structure:
 138   //
 139   // NB_NODES
 140   // Loop from 1 to NB_NODES
 141   //   X Y Z DUMMY_INT
 142
 143   int nbNodes = theMesh->NbNodes();
 144   if ( nbNodes == 0 )
 145     return false;
 146
 147   const char* space    = "  ";
 148   const int   dummyint = 0;
 149
 150   // NB_NODES
 151   theFile << space << nbNodes << endl;
 152
 153   // Loop from 1 to NB_NODES
 154   int aGhs3dID = 1;
 155   SMDS_NodeIteratorPtr it = theMesh->nodesIterator();
 156   while ( it->more() )
 157   {
 158     const SMDS_MeshNode* node = it->next();
 159     theSmdsToGhs3dIdMap.insert( map <int,int>::value_type( node->GetID(), aGhs3dID ));
 160     theGhs3dIdToNodeMap.insert (map <int,const SMDS_MeshNode*>::value_type( aGhs3dID, node ));
 161     aGhs3dID++;
 162
 163     // X Y Z DUMMY_INT
 164     theFile
 165       << space << node->X()
 166       << space << node->Y()
 167       << space << node->Z()
 168       << space << dummyint;
 169
 170     theFile << endl;
 171   }
 172
 173   return true;
 174 }
 175
 176 //=======================================================================
 177 //function : getInt
 178 //purpose  :
 179 //=======================================================================
 180
 181 static bool getInt( int & theValue, char * & theLine )
 182 {
 183   char *ptr;
 184   theValue = strtol( theLine, &ptr, 10 );
 185   if ( ptr == theLine ||
 186       // there must not be neither '.' nor ',' nor 'E' ...
 187       (*ptr != ' ' && *ptr != '\n' && *ptr != '\0'))
 188     return false;
 189
 190   DUMP( "  " << theValue );
 191   theLine = ptr;
 192   return true;
 193 }
 194
 195 //=======================================================================
 196 //function : getDouble
 197 //purpose  :
 198 //=======================================================================
 199
 200 static bool getDouble( double & theValue, char * & theLine )
 201 {
 202   char *ptr;
 203   theValue = strtod( theLine, &ptr );
 204   if ( ptr == theLine )
 205     return false;
 206
 207   DUMP( "   " << theValue );
 208   theLine = ptr;
 209   return true;
 210 }
 211
 212 //=======================================================================
 213 //function : readLine
 214 //purpose  :
 215 //=======================================================================
 216
 217 #define GHS3DPlugin_BUFLENGTH 256
 218 #define GHS3DPlugin_ReadLine(aPtr,aBuf,aFile,aLineNb) \
 219 {  aPtr = fgets( aBuf, GHS3DPlugin_BUFLENGTH - 2, aFile ); aLineNb++; DUMP(endl); }
 220
 221 //=======================================================================
 222 //function : readResult
 223 //purpose  :
 224 //=======================================================================
 225
 226 static bool readResult(FILE *                          theFile,
 227                        SMESHDS_Mesh *                  theMesh,
 228                        const TopoDS_Shape&             theShape,
 229                        map <int,const SMDS_MeshNode*>& theGhs3dIdToNodeMap)
 230 {
 231   // structure:
 232
 233   // record 1:
 234   //  NB_ELEMENTS NB_NODES NB_INPUT_NODES (14 DUMMY_INT)
 235   // record 2:
 236   //  (NB_ELEMENTS * 4) node nbs
 237   // record 3:
 238   //  (NB_NODES) node XYZ
 239
 240   char aBuffer[ GHS3DPlugin_BUFLENGTH ];
 241   char * aPtr;
 242   int aLineNb = 0;
 243
 244   // get shell to set nodes in
 245   TopExp_Explorer exp( theShape, TopAbs_SHELL );
 246   TopoDS_Shell aShell = TopoDS::Shell( exp.Current() );
 247   if ( aShell.IsNull() )
 248     return false;
 249
 250   // ----------------------------------------
 251   // record 1:
 252   // read nb of generated elements and nodes
 253   // ----------------------------------------
 254   int nbElems = 0 , nbNodes = 0, nbInputNodes = 0;
 255   GHS3DPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
 256   if (!aPtr ||
 257       !getInt( nbElems, aPtr ) ||
 258       !getInt( nbNodes, aPtr ) ||
 259       !getInt( nbInputNodes, aPtr))
 260     return false;
 261
 262   // -------------------------------------------
 263   // record 2:
 264   // read element nodes and create tetrahedrons
 265   // -------------------------------------------
 266   GHS3DPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
 267   for (int iElem = 0; iElem < nbElems; iElem++)
 268   {
 269     // read 4 nodes
 270     const SMDS_MeshNode * node[4];
 271     for (int iNode = 0; iNode < 4; iNode++)
 272     {
 273       // read Ghs3d node ID
 274       int ID = 0;
 275       if (!aPtr || ! getInt ( ID, aPtr ))
 276       {
 277         GHS3DPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
 278         if (!aPtr || ! getInt ( ID, aPtr ))
 279         {
 280           MESSAGE( "Cant read " << (iNode+1) << "-th node on line " << aLineNb );
 281           return false;
 282         }
 283       }
 284       // find/create a node with ID
 285       map <int,const SMDS_MeshNode*>::iterator IdNode = theGhs3dIdToNodeMap.find( ID );
 286       if ( IdNode == theGhs3dIdToNodeMap.end())
 287       {
 288         // ID is not yet in theGhs3dIdToNodeMap
 289         ASSERT ( ID > nbInputNodes ); // it should be a new one
 290         SMDS_MeshNode * aNewNode = theMesh->AddNode( 0.,0.,0. ); // read XYZ later
 291         theMesh->SetNodeInVolume( aNewNode, aShell );
 292         theGhs3dIdToNodeMap.insert
 293           ( map <int,const SMDS_MeshNode*>::value_type( ID, aNewNode ));
 294         node[ iNode ] = aNewNode;
 295       }
 296       else
 297       {
 298         node[ iNode ] = (*IdNode).second;
 299       }
 300     }
 301     // create a tetrahedron
 302     SMDS_MeshElement* aTet = theMesh->AddVolume( node[0], node[1], node[2], node[3] );
 303     theMesh->SetMeshElementOnShape( aTet, theShape );
 304   }
 305
 306   // ------------------------
 307   // record 3:
 308   // read and set nodes' XYZ
 309   // ------------------------
 310   GHS3DPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
 311   for (int iNode = 0; iNode < nbNodes; iNode++)
 312   {
 313     // read 3 coordinates
 314     double coord [3];
 315     for (int iCoord = 0; iCoord < 3; iCoord++)
 316     {
 317       if (!aPtr || ! getDouble ( coord[ iCoord ], aPtr ))
 318       {
 319         GHS3DPlugin_ReadLine( aPtr, aBuffer, theFile, aLineNb );
 320         if (!aPtr || ! getDouble ( coord[ iCoord ], aPtr ))
 321         {
 322           MESSAGE( "Cant read " << (iCoord+1) << "-th node coord on line " << aLineNb );
 323           return false;
 324         }
 325       }
 326     }
 327     // do not move old nodes
 328     int ID = iNode + 1;
 329     if (ID <= nbInputNodes)
 330       continue;
 331     // find a node
 332     map <int,const SMDS_MeshNode*>::iterator IdNode = theGhs3dIdToNodeMap.find( ID );
 333     ASSERT ( IdNode != theGhs3dIdToNodeMap.end());
 334     SMDS_MeshNode* node = const_cast<SMDS_MeshNode*> ( (*IdNode).second );
 335
 336     // set XYZ
 337     theMesh->MoveNode( node, coord[0], coord[1], coord[2] );
 338   }
 339
 340   return nbElems;
 341 }
 342
 343 //=============================================================================
 344 /*!
 345  *Here we are going to use the GHS3D mesher
 346  */
 347 //=============================================================================
 348
 349 bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh&         theMesh,
 350                                 const TopoDS_Shape& theShape)
 351 {
 352   MESSAGE("GHS3DPlugin_GHS3D::Compute");
 353
 354   // working dir
 355   QString aTmpDir ( getenv("SALOME_TMP_DIR") );
 356   if ( !aTmpDir.isEmpty() ) {
 357 #ifdef WIN32
 358     if(aTmpDir.at(aTmpDir.length()-1) != '\\') aTmpDir+='\\';
 359 #else
 360     if(aTmpDir.at(aTmpDir.length()-1) != '/') aTmpDir+='/';
 361 #endif
 362   }
 363   else {
 364 #ifdef WIN32
 365     aTmpDir = "C:\\";
 366 #else
 367     aTmpDir = "/tmp/";
 368 #endif
 369   }
 370   // a unique name helps to avoid access to the same files by eg different users
 371   int aUniqueNb;
 372 #ifdef WIN32
 373   aUniqueNb = GetCurrentProcessId();
 374 #else
 375   aUniqueNb = getpid();
 376 #endif
 377
 378   const QString aGenericName    = (aTmpDir + ( "GHS3D_%1" )).arg( aUniqueNb );
 379   const QString aFacesFileName  = aGenericName + ".faces";  // in faces
 380   const QString aPointsFileName = aGenericName + ".points"; // in points
 381   const QString aResultFileName = aGenericName + ".noboite";// out points and volumes
 382   const QString aBadResFileName = aGenericName + ".boite";  // out bad result
 383   const QString aBbResFileName  = aGenericName + ".bb";     // out vertex stepsize
 384   const QString aErrorFileName  = aGenericName + ".log";    // log
 385
 386   // remove possible old files
 387   QFile( aFacesFileName ).remove();
 388   QFile( aPointsFileName ).remove();
 389   QFile( aResultFileName ).remove();
 390   QFile( aErrorFileName ).remove();
 391
 392
 393   // -----------------
 394   // make input files
 395   // -----------------
 396
 397   ofstream aFacesFile  ( aFacesFileName.latin1()  , ios::out);
 398   ofstream aPointsFile ( aPointsFileName.latin1() , ios::out);
 399   bool Ok =
 400 #ifdef WIN32
 401     aFacesFile->is_open() && aPointsFile->is_open();
 402 #else
 403     aFacesFile.rdbuf()->is_open() && aPointsFile.rdbuf()->is_open();
 404 #endif
 405   if (!Ok)
 406   {
 407     MESSAGE( "Can't write into " << aTmpDir << " directory");
 408     return false;
 409   }
 410   SMESHDS_Mesh* meshDS = theMesh.GetMeshDS();
 411   map <int,int> aSmdsToGhs3dIdMap;
 412   map <int,const SMDS_MeshNode*> aGhs3dIdToNodeMap;
 413
 414   Ok =
 415     (writePoints( aPointsFile, meshDS, aSmdsToGhs3dIdMap, aGhs3dIdToNodeMap ) &&
 416      writeFaces ( aFacesFile, meshDS, aSmdsToGhs3dIdMap ));
 417
 418   aFacesFile.close();
 419   aPointsFile.close();
 420
 421   if ( ! Ok )
 422     return false;
 423
 424   // -----------------
 425   // run ghs3d mesher              WIN32???
 426   // -----------------
 427
 428   QString cmd = "ghs3d "
 429     "-m 1000 "                 // memory: 1000 M
 430       "-f " + aGenericName +   // file to read
 431         " 1>" + aErrorFileName; // dump into file
 432   if (system(cmd.latin1()))
 433   {
 434     MESSAGE ("command failed: " << cmd.latin1() );
 435     return false;
 436   }
 437
 438   // --------------
 439   // read a result
 440   // --------------
 441
 442   FILE * aResultFile = fopen( aResultFileName.latin1(), "r" );
 443   if (!aResultFile)
 444   {
 445     MESSAGE( "GHS3D ERROR: see " << aErrorFileName.latin1() );
 446     return false;
 447   }
 448
 449   Ok = readResult( aResultFile, meshDS, theShape, aGhs3dIdToNodeMap );
 450   fclose(aResultFile);
 451
 452   if ( Ok ) {
 453     QFile( aFacesFileName ).remove();
 454     QFile( aPointsFileName ).remove();
 455     QFile( aResultFileName ).remove();
 456     QFile( aErrorFileName ).remove();
 457   }
 458   // remove other possible files
 459   QFile( aBadResFileName ).remove();
 460   QFile( aBbResFileName ).remove();
 461
 462   return Ok;
 463 }
 464
 465
 466 //=============================================================================
 467 /*!
 468  *
 469  */
 470 //=============================================================================
 471
 472 ostream & GHS3DPlugin_GHS3D::SaveTo(ostream & save)
 473 {
 474   return save;
 475 }
 476
 477 //=============================================================================
 478 /*!
 479  *
 480  */
 481 //=============================================================================
 482
 483 istream & GHS3DPlugin_GHS3D::LoadFrom(istream & load)
 484 {
 485   return load;
 486 }
 487
 488 //=============================================================================
 489 /*!
 490  *
 491  */
 492 //=============================================================================
 493
 494 ostream & operator << (ostream & save, GHS3DPlugin_GHS3D & hyp)
 495 {
 496   return hyp.SaveTo( save );
 497 }
 498
 499 //=============================================================================
 500 /*!
 501  *
 502  */
 503 //=============================================================================
 504
 505 istream & operator >> (istream & load, GHS3DPlugin_GHS3D & hyp)
 506 {
 507   return hyp.LoadFrom( load );
 508 }