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