src/GHS3DPlugin_GHS3D.cxx

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