src/MEDLoader/SauvWriter.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF 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 // File      : SauvWriter.cxx
  20 // Created   : Wed Aug 24 12:55:55 2011
  21 // Author    : Edward AGAPOV (eap)
  22
  23 #include "SauvWriter.hxx"
  24
  25 #include "InterpKernelException.hxx"
  26 #include "MEDFileMesh.hxx"
  27 #include "MEDFileField.hxx"
  28 #include "MEDFileData.hxx"
  29 #include "CellModel.hxx"
  30
  31 #include <fstream>
  32 #include <sstream>
  33 #include <iostream>
  34 #include <cstdlib>
  35 #include <iomanip>
  36
  37 using namespace ParaMEDMEM;
  38 using namespace SauvUtilities;
  39 using namespace std;
  40
  41 #define INFOS_MED(txt) cout << txt << endl;
  42
  43 namespace
  44 {
  45   const char* zeroI8 = "       0"; // FORMAT(I8)
  46
  47   // ============================================================
  48   // the class writes endl to the file as soon as <limit> fields
  49   // have been written after the last endl
  50   // ============================================================
  51
  52   class TFieldCounter
  53   {
  54     fstream& _file;
  55     int _count, _limit;
  56   public:
  57     TFieldCounter(fstream& f, int limit=0): _file(f), _limit(limit) { init(); }
  58     void init(int limit=0) // init, is done by stop() as well
  59     { if (limit) _limit = limit; _count = 0; }
  60     void operator++(int) // next
  61     { if ( ++_count == _limit ) { _file << endl; init(); }}
  62     void stop() // init() and write endl if there was no endl after the last written field
  63     { if ( _count ) _file << endl; init(); }
  64     ~TFieldCounter() { stop(); }
  65   };
  66
  67   //================================================================================
  68   /*!
  69    * \brief Return a name of a field support on all elements
  70    */
  71   //================================================================================
  72
  73   string noProfileName( INTERP_KERNEL::NormalizedCellType type )
  74   {
  75     return "INTERP_KERNEL::NormalizedCellType_" + SauvUtilities::toString( type );
  76   }
  77
  78   //================================================================================
  79   /*!
  80    * \brief Remove white spaces from the head and tail
  81    */
  82   //================================================================================
  83
  84   string cleanName( const string& theName )
  85   {
  86     string name = theName;
  87     if ( !name.empty() )
  88       {
  89         // cut off leading white spaces
  90         string::size_type firstChar = name.find_first_not_of(" \t");
  91         if (firstChar < name.length())
  92           {
  93             name = name.substr(firstChar);
  94           }
  95         else
  96           {
  97             name = ""; // only whitespaces there - remove them
  98           }
  99         // cut off trailing white spaces
 100         string::size_type lastChar = name.find_last_not_of(" \t");
 101         if (lastChar < name.length())
 102           name = name.substr(0, lastChar + 1);
 103       }
 104     return name;
 105   }
 106
 107   //================================================================================
 108   /*!
 109    * \brief Converts MED long names into SAUVE short ones, returnes a healed long name
 110    */
 111   //================================================================================
 112
 113   string addName (map<string,int>& nameMap,
 114                   map<string,int>& namePrefixesMap,
 115                   const string&    theName,
 116                   const int        index)
 117   {
 118     // Converts names like:
 119     // MED:                       GIBI:
 120     //   TEMPERATURE_FLUIDE   ->    TEMPE001
 121     //   TEMPERATURE_SOLIDE   ->    TEMPE002
 122     //   PRESSION             ->    PRESSION
 123     //   NU                   ->    NU
 124     //   VOLUM001             ->    VOLUM001
 125     //   VOLUMOFOBJECT        ->    VOLUM003
 126     //   VOLUM002             ->    VOLUM002
 127     string healedName = cleanName(theName);
 128     int ind = index;
 129
 130     if (!healedName.empty())
 131       {
 132         string name = healedName;
 133         int len = name.length();
 134         for (int i = 0; i < len; ++i)
 135           name[i] = toupper(name[i]);
 136
 137         bool doResave = false; // only for tracing
 138
 139         // I. Save a short name as it is
 140         if (len <= 8)
 141           {
 142             INFOS_MED("Save <" << theName << "> as <" << name << ">");
 143
 144             map<string,int>::iterator it = nameMap.find(name);
 145             if (it != nameMap.end())
 146               {
 147                 // There is already such name in the map.
 148
 149                 // a. Replace in the map the old pair by the current one
 150                 int old_ind = nameMap[name];
 151                 nameMap[name] = ind;
 152                 // b. Rebuild the old pair (which was in the map,
 153                 //    it seems to be built automatically by step II)
 154                 ind = old_ind;
 155                 // continue with step II
 156                 doResave = true; // only for tracing
 157               }
 158             else
 159               {
 160                 // Save in the map
 161                 nameMap.insert(make_pair(name, ind));
 162
 163                 // Update loc_index for this name (if last free characters represents a number)
 164                 // to avoid conflicts with long names, same in first 5 characters
 165                 if (len == 8)
 166                   {
 167                     int new_loc_index = atoi(name.c_str() + 5);
 168                     if (new_loc_index > 0)
 169                       {
 170                         // prefix
 171                         string str = name.substr(0,5);
 172                         if (namePrefixesMap.find(str) != namePrefixesMap.end())
 173                           {
 174                             int old_loc_index = namePrefixesMap[str];
 175                             if (new_loc_index < old_loc_index) new_loc_index = old_loc_index;
 176                           }
 177                         namePrefixesMap[str] = new_loc_index;
 178                       }
 179                   }
 180                 return healedName;
 181               }
 182           } // if (len <= 8)
 183
 184         // II. Cut long name and add a numeric suffix
 185
 186         // first 5 or less characters of the name
 187         if (len > 5) name = name.substr(0,5);
 188
 189         // numeric suffix
 190         map<string,int>::iterator name2ind = namePrefixesMap.insert( make_pair( name, 0 )).first;
 191         string numSuffix = SauvUtilities::toString( ++(name2ind->second) );
 192
 193         if ( numSuffix.size() + name.size() > 8 )
 194           THROW_IK_EXCEPTION("Can't write not unique name: " << healedName);
 195
 196         if ( numSuffix.size() < 3 )
 197           numSuffix.insert( 0, 3 - numSuffix.size(), '0' );
 198
 199         name += numSuffix;
 200         nameMap.insert(make_pair(name, ind));
 201
 202         if (doResave)
 203           {
 204             INFOS_MED("Resave previous <" << healedName << "> as <" << name << ">");
 205           }
 206         else
 207           {
 208             INFOS_MED("Save <" << theName << "> as <" << name << ">");
 209           }
 210       }
 211     return healedName;
 212   }
 213 }
 214
 215 //================================================================================
 216 /*!
 217  * \brief Creates SauvWriter
 218  */
 219 //================================================================================
 220
 221 SauvWriter* SauvWriter::New()
 222 {
 223   return new SauvWriter;
 224 }
 225
 226 std::size_t SauvWriter::getHeapMemorySizeWithoutChildren() const
 227 {
 228   return 0;
 229 }
 230
 231 std::vector<const BigMemoryObject *> SauvWriter::getDirectChildren() const
 232 {
 233   return std::vector<const BigMemoryObject *>();
 234 }
 235
 236 //================================================================================
 237 /*!
 238  * \brief Fills own DS by MEDFileData
 239  */
 240 //================================================================================
 241
 242 void SauvWriter::setMEDFileDS(const MEDFileData* medData,
 243                               unsigned           meshIndex)
 244 {
 245   if ( !medData) THROW_IK_EXCEPTION("NULL MEDFileData");
 246
 247   MEDFileMeshes * meshes = medData->getMeshes();
 248   MEDFileFields * fields = medData->getFields();
 249   if ( !meshes) THROW_IK_EXCEPTION("No meshes in MEDFileData");
 250
 251   _fileMesh = meshes->getMeshAtPos( meshIndex );
 252   _fileMesh->incrRef();
 253
 254   if ( fields )
 255     for ( int i = 0; i < fields->getNumberOfFields(); ++i )
 256       {
 257         MEDFileAnyTypeFieldMultiTS * fB = fields->getFieldAtPos(i);
 258         MEDFileFieldMultiTS * f = dynamic_cast<MEDFileFieldMultiTS *>(fB);
 259         if(!f)
 260           continue;// fields on int32 not managed
 261         if ( f->getMeshName() == _fileMesh->getName() )
 262           {
 263             vector< vector<TypeOfField> > fTypes = f->getTypesOfFieldAvailable();
 264             if ( fTypes[0].size() == 1 && fTypes[0][0] == ON_NODES )
 265               _nodeFields.push_back( f );
 266             else
 267               _cellFields.push_back( f );
 268           }
 269       }
 270 }
 271
 272 //================================================================================
 273 /*!
 274  * \brief Adds a submesh
 275  */
 276 //================================================================================
 277
 278 SauvWriter::SubMesh* SauvWriter::addSubMesh(const std::string& name, int dimRelExt)
 279 {
 280   if ( _subs.capacity() < _subs.size() + 1 )
 281     THROW_IK_EXCEPTION("SauvWriter: INTERNAL error, wrong evaluation of nb of sub-meshes");
 282   _subs.resize( _subs.size() + 1 );
 283   SubMesh& sm = _subs.back();
 284   sm._name = name;
 285   sm._dimRelExt = dimRelExt;
 286   return &sm;
 287 }
 288 //================================================================================
 289 /*!
 290  * \brief Returns nb of cell types
 291  */
 292 //================================================================================
 293
 294 int SauvWriter::SubMesh::nbTypes() const
 295 {
 296   int nb = 0;
 297   for (int i = 0; i < cellIDsByTypeSize(); ++i )
 298     nb += int( !_cellIDsByType[i].empty() );
 299   return nb;
 300 }
 301
 302 //================================================================================
 303 /*!
 304  * \brief Fill _subs
 305  */
 306 //================================================================================
 307
 308 void SauvWriter::fillSubMeshes( int& nbSauvObjects, map<string,int>& nameNbMap )
 309 {
 310   // evaluate nb of _subs in order to avoid re-allocation of _subs
 311   int nbSubs = 1; // for the very mesh
 312   nbSubs += _fileMesh->getFamilyInfo().size() + 4;  // + 4 zero families (for each dimRelExt)
 313   nbSubs += _fileMesh->getGroupInfo().size();
 314   nbSubs += evaluateNbProfileSubMeshes();
 315   _subs.clear();
 316   _subs.reserve( nbSubs );
 317
 318   fillFamilySubMeshes();
 319   fillGroupSubMeshes();
 320   fillProfileSubMeshes();
 321
 322   // fill names of SubMesh'es and count nb of sauv sub-meshes they will be stored into
 323   nbSauvObjects = 0;
 324   map<string,int> namePrefixMap;
 325   for ( size_t i = 0; i < _subs.size(); ++i )
 326     {
 327       SubMesh& sm = _subs[i];
 328
 329       sm._nbSauvObjects = 0;
 330       if ( sm._subs.empty() )
 331         {
 332           sm._nbSauvObjects = sm.nbTypes();
 333         }
 334       else
 335         {
 336           sm._nbSauvObjects = 1;
 337         }
 338
 339       sm._id = nbSauvObjects+1;
 340       nbSauvObjects += sm._nbSauvObjects;
 341
 342       if ( sm._nbSauvObjects )
 343         sm._name = addName( nameNbMap, namePrefixMap, sm._name, sm._id );
 344
 345       if ( sm._nbSauvObjects && !sm._name.empty() )
 346         {
 347           nameGIBItoMED aMEDName;
 348           aMEDName.gibi_pile = PILE_SOUS_MAILLAGE;
 349           aMEDName.gibi_id   = sm._id;
 350           aMEDName.med_name  = sm._name;
 351           _longNames[ LN_MAIL ].push_back(aMEDName);
 352         }
 353     }
 354 }
 355
 356 //================================================================================
 357 /*!
 358  * \brief fill sub-meshes of families
 359  */
 360 //================================================================================
 361
 362 void SauvWriter::fillFamilySubMeshes()
 363 {
 364   SubMesh* nilSm = (SubMesh*) 0;
 365   std::vector<int> dims = _fileMesh->getNonEmptyLevelsExt();
 366   for ( size_t iDim = 0; iDim < dims.size(); ++iDim )
 367     {
 368       int dimRelExt = dims[ iDim ];
 369       MEDCouplingAutoRefCountObjectPtr< MEDCouplingMesh > mesh = _fileMesh->getGenMeshAtLevel(dimRelExt);
 370       const DataArrayInt * famIds = _fileMesh->getFamilyFieldAtLevel(dimRelExt);
 371       if ( !famIds ) continue;
 372
 373       int curFamID = 0;
 374       SubMesh* curSubMesh = addSubMesh( "", dimRelExt ); // submesh of zero family
 375       _famIDs2Sub[0] = curSubMesh;
 376       int sub0Index = _subs.size()-1;
 377
 378       const int * famID = famIds->begin(), * famIDEnd = famIds->end();
 379       for ( int cellID = 0; famID < famIDEnd; ++famID, cellID++ )
 380         {
 381           if ( *famID != curFamID )
 382             {
 383               curFamID = *famID;
 384               map< int, SubMesh* >::iterator f2s = _famIDs2Sub.insert( make_pair( curFamID, nilSm )).first;
 385               if ( !f2s->second )
 386                 f2s->second = addSubMesh( "", dimRelExt ); // no names for families
 387               curSubMesh = f2s->second;
 388             }
 389           INTERP_KERNEL::NormalizedCellType cellType =
 390             dimRelExt == 1 ? INTERP_KERNEL::NORM_POINT1 : mesh->getTypeOfCell( cellID );
 391           curSubMesh->_cellIDsByType[ cellType ].push_back( cellID );
 392         }
 393
 394       if ( dimRelExt == 1 )
 395         {
 396           // clear submesh of nodal zero family
 397           _famIDs2Sub[0]->_cellIDsByType[ INTERP_KERNEL::NORM_POINT1 ].clear();
 398         }
 399       else if ( dimRelExt == 0 )
 400         {
 401           // make a submesh including all cells
 402           if ( sub0Index == (int)(_subs.size()-1) )
 403             {
 404               _famIDs2Sub[0]->_name = _fileMesh->getName(); // there is the zero family only
 405             }
 406           else
 407             {
 408               curSubMesh = addSubMesh( _fileMesh->getName(), dimRelExt );
 409               if ( _famIDs2Sub[0]->nbTypes() == 0 )
 410                 sub0Index++; // skip an empty zero family
 411               for ( size_t i = sub0Index; i < _subs.size()-1; ++i )
 412                 curSubMesh->_subs.push_back( & _subs[i] );
 413             }
 414         }
 415     }
 416 }
 417
 418 //================================================================================
 419 /*!
 420  * \brief fill sub-meshes of groups
 421  */
 422 //================================================================================
 423
 424 void SauvWriter::fillGroupSubMeshes()
 425 {
 426   const map<string, vector<string> >& grpFams = _fileMesh->getGroupInfo();
 427   map<string, vector<string> >::const_iterator g2ff = grpFams.begin();
 428   for ( ; g2ff != grpFams.end(); ++g2ff )
 429     {
 430       const string&        groupName = g2ff->first;
 431       const vector<string>& famNames = g2ff->second;
 432       if ( famNames.empty() ) continue;
 433       std::vector<SubMesh*> famSubMeshes( famNames.size() );
 434       std::size_t k = 0;
 435       for ( size_t i = 0; i < famNames.size(); ++i )
 436         {
 437           int famID = _fileMesh->getFamilyId( famNames[i].c_str() );
 438           map< int, SubMesh* >::iterator i2f = _famIDs2Sub.find( famID );
 439           if ( i2f != _famIDs2Sub.end() )
 440             {
 441               famSubMeshes[ k ] = i2f->second;
 442               ++k;
 443             }
 444         }
 445       // if a family exists but has no element, no submesh has been found for this family
 446       // => we have to resize famSubMeshes with the number of submeshes stored
 447       if (k != famNames.size())
 448           famSubMeshes.resize(k);
 449       SubMesh* grpSubMesh = addSubMesh( groupName, famSubMeshes[0]->_dimRelExt );
 450       grpSubMesh->_subs.swap( famSubMeshes );
 451     }
 452 }
 453
 454
 455 //================================================================================
 456 /*!
 457  * \brief fill sub-meshes of profiles
 458  */
 459 //================================================================================
 460
 461 void SauvWriter::fillProfileSubMeshes()
 462 {
 463   _profile2Sub.clear();
 464   SubMesh* nilSm = (SubMesh*) 0;
 465   for ( int isOnNodes = 0; isOnNodes < 2; ++isOnNodes )
 466     {
 467       vector< MEDCouplingAutoRefCountObjectPtr< MEDFileFieldMultiTS > >
 468         fields = isOnNodes ? _nodeFields : _cellFields;
 469       for ( size_t i = 0; i < fields.size(); ++i )
 470         {
 471           vector< pair<int,int> > iters = fields[i]->getIterations();
 472
 473           vector<INTERP_KERNEL::NormalizedCellType> types;
 474           vector< vector<TypeOfField> > typesF;
 475           vector< vector<string> > pfls, locs;
 476           fields[i]->getFieldSplitedByType( iters[0].first, iters[0].second,
 477                                             _fileMesh->getName().c_str(), types, typesF, pfls, locs);
 478           int dimRelExt;
 479           for ( size_t iType = 0; iType < types.size(); ++iType )
 480             {
 481               if ( types[iType] == INTERP_KERNEL::NORM_ERROR )
 482                 dimRelExt = 1; // on nodes
 483               else
 484                 dimRelExt = getDimension( types[iType] ) - _fileMesh->getMeshDimension();
 485               for ( size_t iPfl = 0; iPfl < pfls[iType].size(); ++iPfl )
 486                 {
 487                   bool isOnAll = pfls[iType][iPfl].empty();
 488                   if ( isOnAll ) pfls[iType][iPfl] = noProfileName( types[iType] );
 489                   map< string, SubMesh* >::iterator pfl2sm =
 490                     _profile2Sub.insert( make_pair( pfls[iType][iPfl], nilSm )).first;
 491                   if ( !pfl2sm->second )
 492                     {
 493                       SubMesh* sm = pfl2sm->second = addSubMesh( "", dimRelExt ); // no names for profiles
 494                       const DataArrayInt * pfl = isOnAll ? 0 : fields[i]->getProfile( pfls[iType][iPfl].c_str() );
 495                       makeProfileIDs( sm, types[iType], pfl );
 496                     }
 497                 }
 498             }
 499         }
 500     }
 501 }
 502
 503 //================================================================================
 504 /*!
 505  * \brief Return max possible nb of sub-meshes to decsribe field supports
 506  */
 507 //================================================================================
 508
 509 int SauvWriter::evaluateNbProfileSubMeshes() const
 510 {
 511   int nb = 0;
 512   for ( size_t i = 0; i < _nodeFields.size(); ++i )
 513     nb += 1 + _nodeFields[i]->getPflsReallyUsed().size();
 514
 515   for ( size_t i = 0; i < _cellFields.size(); ++i )
 516     {
 517       nb += _cellFields[i]->getPflsReallyUsed().size();
 518
 519       vector< pair<int,int> > iters = _cellFields[i]->getIterations();
 520
 521       vector<INTERP_KERNEL::NormalizedCellType> types;
 522       vector< vector<TypeOfField> > typesF;
 523       vector< vector<string> > pfls, locs;
 524       _cellFields[i]->getFieldSplitedByType( iters[0].first, iters[0].second,
 525                                              _fileMesh->getName().c_str(), types, typesF, pfls, locs);
 526       nb += 2 * types.size(); // x 2 - a type can be on nodes and on cells at the same time
 527     }
 528
 529   return nb;
 530 }
 531
 532 //================================================================================
 533 /*!
 534  * \brief Transorm a profile into ids of mesh elements
 535  */
 536 //================================================================================
 537
 538 void SauvWriter::makeProfileIDs( SubMesh*                          sm,
 539                                  INTERP_KERNEL::NormalizedCellType type,
 540                                  const DataArrayInt*               profile )
 541 {
 542   MEDCouplingAutoRefCountObjectPtr< MEDCouplingMesh >
 543     mesh = _fileMesh->getGenMeshAtLevel(sm->_dimRelExt);
 544   const MEDCouplingUMesh* uMesh = dynamic_cast< const MEDCouplingUMesh* > ((const MEDCouplingMesh*) mesh );
 545
 546   if ( sm->_dimRelExt == 1 ) type = INTERP_KERNEL::NORM_POINT1;
 547   vector< int >& ids = sm->_cellIDsByType[ type ];
 548
 549   if ( sm->_dimRelExt == 1 || !uMesh )
 550     {
 551       // profile on nodes or mesh is CARTESIAN
 552       if ( profile )
 553         {
 554           ids.assign( profile->begin(), profile->end() );
 555         }
 556       else // on all
 557         {
 558           ids.resize( sm->_dimRelExt == 1 ? mesh->getNumberOfNodes() : mesh->getNumberOfCells() );
 559           for ( size_t i = 0; i < ids.size(); ++i )
 560             ids[i]=i;
 561         }
 562     }
 563   else
 564     {
 565       // profile on cells
 566       vector<int> code(3);
 567       code[0] = type;
 568       if ( profile ) // on profile
 569         {
 570           code[1] = profile->getNumberOfTuples();
 571           code[2] = 0;
 572         }
 573       else // on all cells
 574         {
 575           code[1] = mesh->getNumberOfCellsWithType( type );
 576           code[2] = -1;
 577         }
 578       vector<const DataArrayInt *> idsPerType( 1, profile );
 579       MEDCouplingAutoRefCountObjectPtr<DataArrayInt>
 580         resIDs = uMesh->checkTypeConsistencyAndContig( code, idsPerType );
 581       if (( const DataArrayInt *) resIDs )
 582       {
 583         ids.assign( resIDs->begin(), resIDs->end() );
 584       }
 585       else // mesh includes only one type
 586       {
 587         int nbE = code[1];
 588         for ( ids.resize( nbE ); nbE; --nbE )
 589           ids[ nbE-1 ] = nbE-1;
 590       }
 591     }
 592 }
 593
 594 //================================================================================
 595 /*!
 596  * \brief Write its data into the SAUVE file
 597  */
 598 //================================================================================
 599
 600 void SauvWriter::write(const char* fileName)
 601 {
 602   std::fstream fileStream;
 603   fileStream.open( fileName, ios::out);
 604   if
 605 #ifdef WIN32
 606     ( !fileStream || !fileStream.is_open() )
 607 #else
 608     ( !fileStream || !fileStream.rdbuf()->is_open() )
 609 #endif
 610       THROW_IK_EXCEPTION("Can't open the file |"<<fileName<<"|");
 611   _sauvFile = &fileStream;
 612
 613   _subs.clear();
 614   _famIDs2Sub.clear();
 615   _profile2Sub.clear();
 616   _longNames[ LN_MAIL ].clear();
 617   _longNames[ LN_CHAM ].clear();
 618   _longNames[ LN_COMP ].clear();
 619
 620   map<string,int> fldNamePrefixMap;
 621
 622   writeFileHead();
 623   writeSubMeshes();
 624   writeNodes();
 625   writeNodalFields(fldNamePrefixMap);
 626   writeElemFields(fldNamePrefixMap);
 627   writeLongNames();
 628   writeLastRecord();
 629
 630   _sauvFile->close();
 631 }
 632 //================================================================================
 633 /*!
 634  * \brief Writes "ENREGISTREMENT DE TYPE" 4 and 7
 635  */
 636 //================================================================================
 637
 638 void SauvWriter::writeFileHead()
 639 {
 640   MEDCouplingAutoRefCountObjectPtr< MEDCouplingMesh > mesh = _fileMesh->getGenMeshAtLevel(0);
 641
 642   *_sauvFile
 643     << " ENREGISTREMENT DE TYPE   4" << endl
 644     << " NIVEAU  16 NIVEAU ERREUR   0 DIMENSION   " << mesh->getSpaceDimension() <<endl
 645     << " DENSITE 0.00000E+00" << endl
 646     << " ENREGISTREMENT DE TYPE   7" << endl
 647     << " NOMBRE INFO CASTEM2000   8" <<endl
 648     << " IFOUR  -1 NIFOUR   0 IFOMOD  -1 IECHO   1 IIMPI   0 IOSPI   0 ISOTYP   1" << endl
 649     << " NSDPGE     0" << endl;
 650 }
 651
 652 //================================================================================
 653 /*!
 654  * \brief Writes names of objects
 655  */
 656 //================================================================================
 657
 658 void SauvWriter::writeNames( const map<string,int>& nameNbMap )
 659 {
 660   if ( !nameNbMap.empty() )
 661   {
 662     // write names of objects
 663     // * 8001       FORMAT(8(1X,A8))
 664     TFieldCounter fcount( *_sauvFile, 8 );
 665     *_sauvFile << left;
 666     map<string,int>::const_iterator nameNbIt = nameNbMap.begin();
 667     for ( ; nameNbIt != nameNbMap.end(); nameNbIt++, fcount++ )
 668       *_sauvFile << " " << setw(8) << nameNbIt->first;
 669     fcount.stop();
 670     *_sauvFile << right;
 671
 672     // write IDs of named objects in the pile
 673     // *  8000 FORMAT(10I8)
 674     nameNbIt = nameNbMap.begin();
 675     for ( fcount.init(10); nameNbIt != nameNbMap.end(); nameNbIt++, fcount++ )
 676       *_sauvFile << setw(8) << nameNbIt->second;
 677   }
 678 }
 679
 680 //================================================================================
 681 /*!
 682  * \brief Writes "PILE NUMERO   1"
 683  */
 684 //================================================================================
 685
 686 void SauvWriter::writeSubMeshes()
 687 {
 688   int nbSauvObjects;
 689   map<string,int> nameNbMap;
 690   fillSubMeshes( nbSauvObjects, nameNbMap );
 691
 692   // * 800   FORMAT (' ENREGISTREMENT DE TYPE', I4)
 693   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl;
 694   // * 801     FORMAT(' PILE NUMERO',I4,'NBRE OBJETS NOMMES',I8,'NBRE OBJETS',I8)
 695   *_sauvFile << " PILE NUMERO   1NBRE OBJETS NOMMES" << setw(8) << nameNbMap.size() <<
 696     "NBRE OBJETS" << setw(8) << nbSauvObjects <<endl;
 697
 698   writeNames( nameNbMap );
 699
 700   TFieldCounter fcount( *_sauvFile, 10 ); // 10 intergers per line
 701
 702   for ( size_t iSub = 0; iSub < _subs.size(); ++iSub )
 703     {
 704       SubMesh& sm = _subs[iSub];
 705       if ( sm._nbSauvObjects < 1 ) continue;
 706
 707       // The first record of each sub-mesh writes
 708       // - type of cells; zero means a compound object whose the 2nd record enumerates its components
 709       // - number of components of a compound object
 710       // - number of references; each reference means a "pointer" to this sub-mesh
 711       // - number of nodes per cell
 712       // - number of cells
 713
 714       if ( !sm._subs.empty() )
 715         {
 716           writeCompoundSubMesh(iSub);
 717         }
 718       else
 719         {
 720           // write each sub-type as a SAUV sub-mesh
 721           MEDCouplingAutoRefCountObjectPtr< MEDCouplingMesh >
 722             mesh = _fileMesh->getGenMeshAtLevel( sm._dimRelExt );
 723           MEDCouplingAutoRefCountObjectPtr< MEDCouplingUMesh>
 724             umesh = mesh->buildUnstructured();
 725
 726           for ( int iType=0; iType < sm.cellIDsByTypeSize(); ++iType )
 727             {
 728               const vector<int>& cellIDs = sm._cellIDsByType[iType];
 729               if ( cellIDs.empty() ) continue;
 730
 731               INTERP_KERNEL::NormalizedCellType
 732                 cellType = INTERP_KERNEL::NormalizedCellType( iType );
 733               const INTERP_KERNEL::CellModel &
 734                 cell = INTERP_KERNEL::CellModel::GetCellModel( cellType );
 735               int castemType       = SauvUtilities::med2gibiGeom( cellType );
 736               unsigned nbElemNodes = cell.getNumberOfNodes();
 737               unsigned nbElems     = cellIDs.size();
 738
 739               *_sauvFile << setw(8) << castemType
 740                         << zeroI8
 741                         << zeroI8
 742                         << setw(8) << nbElemNodes
 743                         << setw(8) << nbElems << endl;
 744
 745               // write color of each element
 746               // * 8000 FORMAT(10I8)
 747               for ( size_t i = 0; i < nbElems; ++i, fcount++ ) *_sauvFile << zeroI8;
 748               fcount.stop();
 749
 750               // write connectivity
 751               // gibi IDs are in FORTRAN mode while MEDCoupling IDs are in C mode
 752               if ( sm._dimRelExt == 1 ) // nodes
 753                 {
 754                   for ( size_t i = 0; i < nbElems; ++i, fcount++ )
 755                     *_sauvFile << setw(8) << ( cellIDs[i] + 1 );
 756                 }
 757               else
 758                 {
 759                   // indices to transform MED connectivity to GIBI one
 760                   const int * toMedConn = getGibi2MedQuadraticInterlace( cellType );
 761
 762                   vector< int > cellConn( nbElemNodes ), transformedConn( nbElemNodes );
 763                   for ( size_t i = 0; i < nbElems; ++i )
 764                     {
 765                       cellConn.clear();
 766                       umesh->getNodeIdsOfCell( cellIDs[i], cellConn );
 767                       if ( toMedConn )
 768                         {
 769                           for ( unsigned j = 0; j < nbElemNodes; ++j )
 770                             transformedConn[ toMedConn[ j ]] = cellConn[ j ];
 771                           cellConn.swap( transformedConn );
 772                         }
 773                       for ( unsigned j = 0; j < nbElemNodes; ++j, fcount++ )
 774                         *_sauvFile << setw(8) << ( cellConn[j] + 1 );
 775                     }
 776                 }
 777               fcount.stop();
 778
 779             } // loop on cell types
 780         } // not a compound object
 781     } // loop on sub-meshes
 782 }
 783
 784 //================================================================================
 785 /*!
 786  * \brief Writes a sum-mesh composed of other sum-meshes
 787  * This submesh corresponds to a med mesh or group composed of families
 788  */
 789 //================================================================================
 790
 791 void SauvWriter::writeCompoundSubMesh(int iSub)
 792 {
 793   SubMesh& sm = _subs[iSub];
 794   if ( sm._nbSauvObjects < 1 || sm._subs.empty()) return;
 795
 796   vector< int > subIDs;
 797   for ( size_t i = 0; i < sm._subs.size(); ++i ) // loop on sub-meshes of families
 798     for ( int j = 0; j < sm._subs[i]->_nbSauvObjects; ++j )
 799       subIDs.push_back( sm._subs[i]->_id + j );
 800
 801   *_sauvFile << zeroI8
 802              << setw(8) << subIDs.size()
 803              << zeroI8
 804              << zeroI8
 805              << zeroI8 << endl;
 806
 807   TFieldCounter fcount( *_sauvFile, 10 ); // 10 intergers per line
 808   for ( size_t i = 0; i < subIDs.size(); ++i, fcount++ )
 809     *_sauvFile << setw(8) << subIDs[i];
 810 }
 811
 812 //================================================================================
 813 /*!
 814  * \brief Write piles relating to nodes
 815  */
 816 //================================================================================
 817
 818 void SauvWriter::writeNodes()
 819 {
 820   MEDCouplingAutoRefCountObjectPtr< MEDCouplingMesh > mesh = _fileMesh->getGenMeshAtLevel( 1 );
 821   MEDCouplingAutoRefCountObjectPtr< MEDCouplingUMesh > umesh = mesh->buildUnstructured();
 822
 823   // write the index connecting nodes with their coodrinates
 824
 825   const int nbNodes = umesh->getNumberOfNodes();
 826   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl
 827              << " PILE NUMERO  32NBRE OBJETS NOMMES       0NBRE OBJETS" << setw(8) << nbNodes << endl;
 828   *_sauvFile << setw(8) << nbNodes << endl;
 829   //
 830   TFieldCounter fcount( *_sauvFile, 10 );// * 8000 FORMAT(10I8)
 831   for ( int i = 0; i < nbNodes; ++i, fcount++ )
 832     *_sauvFile << setw(8) << i + 1;
 833   fcount.stop();
 834
 835   // write coordinates and density of nodes
 836
 837   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl;
 838   *_sauvFile << " PILE NUMERO  33NBRE OBJETS NOMMES       0NBRE OBJETS       1" << endl;
 839   //
 840   const int dim = umesh->getSpaceDimension();
 841   const int nbValues = nbNodes * ( dim + 1 );
 842   *_sauvFile << setw(8) << nbValues << endl;
 843
 844   // * 8003   FORMAT(1P,3E22.14)
 845   const char* density = "  0.00000000000000E+00";
 846   fcount.init(3);
 847   _sauvFile->precision(14);
 848   _sauvFile->setf( ios_base::scientific, ios_base::floatfield );
 849   _sauvFile->setf( ios_base::uppercase );
 850   MEDCouplingAutoRefCountObjectPtr< DataArrayDouble> coordArray = umesh->getCoordinatesAndOwner();
 851   const double precision = 1.e-99; // PAL12077
 852   for ( int i = 0; i < nbNodes; ++i)
 853   {
 854     for ( int j = 0; j < dim; ++j, fcount++ )
 855       {
 856         double coo = coordArray->getIJ( i, j );
 857         bool  zero = ( -precision < coo && coo < precision );
 858         *_sauvFile << setw(22) << ( zero ? 0.0 : coo );
 859       }
 860     *_sauvFile << density;
 861     fcount++;
 862   }
 863 }
 864
 865 //================================================================================
 866 /*!
 867  * \brief Store correspondence between GIBI (short) and MED (long) names
 868  *
 869  * IMP 0020434: mapping GIBI names to MED names
 870  * Store correspondence between GIBI and MED names as one PILE_STRINGS and one
 871  * PILE_TABLES (in three tables: MED_MAIL, MED_CHAM and MED_COMP)
 872  */
 873 //================================================================================
 874
 875 void SauvWriter::writeLongNames()
 876 {
 877   int nbTables =
 878     3 - _longNames[ LN_MAIL ].empty() - _longNames[ LN_CHAM ].empty() - _longNames[ LN_COMP ].empty();
 879   if (nbTables == 0) return;
 880
 881   // ---------------------
 882   // Write the TABLE pile
 883   // ---------------------
 884
 885   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl
 886         << " PILE NUMERO  10NBRE OBJETS NOMMES" << setw(8) << nbTables
 887         << "NBRE OBJETS" << setw(8) << nbTables << endl;
 888   // table names
 889   if (!_longNames[ LN_MAIL ].empty()) *_sauvFile << " MED_MAIL";
 890   if (!_longNames[ LN_CHAM ].empty()) *_sauvFile << " MED_CHAM";
 891   if (!_longNames[ LN_COMP ].empty()) *_sauvFile << " MED_COMP";
 892   *_sauvFile << endl;
 893   // table indices
 894   for ( int i = 0; i < nbTables; ++i ) *_sauvFile << setw(8) << i+1;
 895   *_sauvFile << endl;
 896
 897   string theWholeString; // concatenated long names
 898   vector<int> theOffsets;
 899   int iStr = 1;
 900   TFieldCounter fcount (*_sauvFile, 10);
 901
 902   for ( int iTbl = 0; iTbl < LN_NB; ++iTbl )
 903     {
 904       vector<nameGIBItoMED>& longNames = _longNames[ iTbl ];
 905       if ( longNames.empty() ) continue;
 906       const bool isComp = ( iTbl == LN_COMP);
 907
 908       // to assure unique MED names
 909       set<string> medUniqueNames;
 910
 911       *_sauvFile << setw(8) << longNames.size()*4 << endl; // Nb of table values
 912
 913       vector<nameGIBItoMED>::iterator itGIBItoMED = longNames.begin();
 914       for (; itGIBItoMED != longNames.end(); itGIBItoMED++, iStr++)
 915         {
 916           // PILE of i-th key (med name)
 917           *_sauvFile << setw(8) << PILE_STRINGS;
 918           fcount++;
 919           // ID of i-th key (med name)
 920           *_sauvFile << setw(8) << iStr;
 921           fcount++;
 922           // PILE of i-th value (gibi name)
 923           *_sauvFile << setw(8) << itGIBItoMED->gibi_pile;
 924           fcount++;
 925           // ID of i-th value (gibi name)
 926           *_sauvFile << setw(8) << ( isComp ? ++iStr : itGIBItoMED->gibi_id );
 927           fcount++;
 928
 929           // add a MED name to the string (while making it be unique for sub-meshes and fields)
 930           string aMedName = itGIBItoMED->med_name;
 931           if ( !isComp )
 932             for (int ind = 1; !medUniqueNames.insert(aMedName).second; ++ind )
 933               aMedName = itGIBItoMED->med_name + "_" + SauvUtilities::toString( ind );
 934           theWholeString += aMedName;
 935
 936           // add an offset
 937           theOffsets.push_back( theWholeString.size() );
 938           if ( isComp )
 939             {
 940               theWholeString += itGIBItoMED->gibi_name;
 941               theOffsets.push_back( theWholeString.size() );
 942             }
 943         }
 944       fcount.stop();
 945     }
 946
 947   // ----------------------
 948   // Write the STRING pile
 949   // ----------------------
 950
 951   const int nbNames = theOffsets.size();
 952   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl
 953         << " PILE NUMERO  27NBRE OBJETS NOMMES" << zeroI8 << "NBRE OBJETS" << setw(8) << nbNames << endl
 954         << setw(8) << theWholeString.length() << setw(8) << nbNames << endl;
 955
 956   // write the whole string
 957   const int fixedLength = 71;
 958   for ( string::size_type aPos = 0; aPos < theWholeString.length(); aPos += fixedLength)
 959     *_sauvFile << setw(72) << theWholeString.substr(aPos, fixedLength) << endl;
 960
 961   // write the offsets
 962   for ( size_t i = 0; i < theOffsets.size(); ++i, fcount++ )
 963     *_sauvFile << setw(8) << theOffsets[i];
 964 }
 965
 966 //================================================================================
 967 /*!
 968  * \brief Write beginning of field record
 969  */
 970 //================================================================================
 971
 972 void SauvWriter::writeFieldNames( const bool                 isNodal,
 973                                   std::map<std::string,int>& fldNamePrefixMap)
 974 {
 975   vector< MEDCouplingAutoRefCountObjectPtr< MEDFileFieldMultiTS > >&
 976     flds = isNodal ? _nodeFields : _cellFields;
 977   map<string,int> nameNbMap;
 978
 979   for ( size_t iF = 0; iF < flds.size(); ++iF )
 980     {
 981       string name = addName( nameNbMap, fldNamePrefixMap, flds[iF]->getName(), iF+1 );
 982       nameGIBItoMED aMEDName;
 983       aMEDName.gibi_pile = isNodal ? PILE_NODES_FIELD : PILE_FIELD;
 984       aMEDName.gibi_id   = iF+1;
 985       aMEDName.med_name  = name;
 986       _longNames[ LN_CHAM ].push_back(aMEDName);
 987     }
 988
 989   *_sauvFile << " ENREGISTREMENT DE TYPE   2" << endl
 990              << ( isNodal ? " PILE NUMERO   2" : " PILE NUMERO  39")
 991              << "NBRE OBJETS NOMMES" << setw(8) << nameNbMap.size()
 992              << "NBRE OBJETS"        << setw(8) << flds.size() << endl;
 993   writeNames( nameNbMap );
 994 }
 995
 996 //================================================================================
 997 /*!
 998  * \brief Make short names of field components
 999  *
1000  * IMP 0020434: mapping GIBI names to MED names
1001  */
1002 //================================================================================
1003
1004 void SauvWriter::makeCompNames(const string&         fieldName,
1005                                const vector<string>& compInfo,
1006                                map<string, string>&  mapMedToGibi)
1007 {
1008   for ( size_t i = 0; i < compInfo.size(); ++i )
1009     mapMedToGibi[compInfo[i]] = cleanName( compInfo[i] );
1010
1011   int compIndex = 1;
1012   map<string, string>::iterator namesIt = mapMedToGibi.begin();
1013   for (; namesIt != mapMedToGibi.end(); namesIt++)
1014     {
1015       string & compGibiName = (*namesIt).second;
1016       if (compGibiName.size() > 4) {
1017         // use new name in form "CXXX", where "XXX" is a number
1018         do
1019           {
1020             compGibiName = SauvUtilities::toString( compIndex++ );
1021             if ( compGibiName.size() < 3 )
1022               compGibiName.insert( 0, 3 - compGibiName.size(), '0' );
1023             compGibiName = "C" + compGibiName;
1024           }
1025         while (mapMedToGibi.count(compGibiName) > 0); // real component name could be CXXX
1026       }
1027
1028       string compMedName = fieldName + "." + namesIt->first;
1029       nameGIBItoMED aMEDName;
1030       aMEDName.med_name  = compMedName;
1031       aMEDName.gibi_pile = PILE_STRINGS;
1032       aMEDName.gibi_name = compGibiName;
1033       _longNames[ LN_COMP ].push_back(aMEDName);
1034     }
1035 }
1036
1037 //================================================================================
1038 /*!
1039  * \brief Writes "PILE NUMERO   2": fields on nodes
1040  */
1041 //================================================================================
1042
1043 void SauvWriter::writeNodalFields(map<string,int>& fldNamePrefixMap)
1044 {
1045   writeFieldNames( /*isNodal=*/true, fldNamePrefixMap );
1046
1047   TFieldCounter fcount (*_sauvFile, 10);
1048
1049   // EXAMPLE ( with no values )
1050
1051   // (1)       4       7       2       1
1052   // (2)     -88       0       3     -89       0       1     -90       0       2     -91
1053   // (2)       0       1
1054   // (3) FX   FY   FZ   FZ   FX   FY   FLX
1055   // (4)       0       0       0       0       0       0       0
1056   // (5)           cree  par  muc pri
1057   // (6)
1058   // (7)       2
1059   for ( size_t iF = 0; iF < _nodeFields.size(); ++iF )
1060     {
1061       // (1) write nb subcomponents, nb components(total)
1062       vector< pair<int,int> >  iters = _nodeFields[iF]->getIterations();
1063       const vector<string>& compInfo = _nodeFields[iF]->getInfo();
1064       const int nbSub = iters.size();
1065       const int nbComp = compInfo.size();
1066       const int totalNbComp = nbSub * nbComp;
1067       *_sauvFile << setw(8) << nbSub
1068                  << setw(8) << totalNbComp
1069                  << setw(8) << -1         // IFOUR
1070                  << setw(8) << 0 << endl; // nb attributes
1071
1072       // (2) for each sub-component (iteration)
1073       // write support, number of values and number of components
1074       fcount.init(10);
1075       vector< int > vals(3);
1076       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1077         {
1078           pair<int,int> it = iters[iIt];
1079
1080           vector<INTERP_KERNEL::NormalizedCellType> types;
1081           vector< vector<TypeOfField> > typesF;
1082           vector< vector<string> > pfls, locs;
1083           vector< vector< std::pair<int,int> > > valsVec;
1084           valsVec=_nodeFields[iF]->getFieldSplitedByType( it.first, it.second, _fileMesh->getName().c_str(),
1085                                                           types, typesF, pfls, locs);
1086           // believe that there can be only one type in a nodal field,
1087           // so do not use a loop on types
1088           if ( pfls[0][0].empty() ) pfls[0][0] = noProfileName( types[0] );
1089           map< string, SubMesh* >::iterator pfl2Sub = _profile2Sub.find( pfls[0][0] );
1090           if ( pfl2Sub == _profile2Sub.end() )
1091             THROW_IK_EXCEPTION( "SauvWriter::writeNodalFields(): no sub-mesh for profile |"
1092                                 << pfls[0][0] << "|");
1093           vals[0] = -pfl2Sub->second->_id;
1094           vals[1] = (valsVec[0][0].second-valsVec[0][0].first);
1095           vals[2] = compInfo.size();
1096           for ( size_t i = 0; i < vals.size(); ++i, fcount++ )
1097             *_sauvFile << setw(8) << vals[i];
1098         }
1099       fcount.stop();
1100
1101       // (3) Write names of components
1102       map<string, string> mapMedToGibi;
1103       makeCompNames( _nodeFields[iF]->getName(), compInfo, mapMedToGibi );
1104       fcount.init(8);
1105       *_sauvFile << left;
1106       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1107         for ( size_t i = 0; i < compInfo.size(); ++i, fcount++ )
1108           *_sauvFile << " "  << setw(4) << mapMedToGibi[compInfo[i]];
1109       *_sauvFile << right;
1110       fcount.stop();
1111
1112       // (4) nb harmonics
1113       fcount.init(10);
1114       for ( size_t i = 0; i < (std::size_t)totalNbComp; ++i, fcount++ )
1115         *_sauvFile << " "  << setw(8) << 0;
1116       fcount.stop();
1117
1118       string description = _nodeFields[iF]->getName();
1119       *_sauvFile << endl;                                         // (5) TYPE
1120       *_sauvFile << setw(72) << description.substr(0,71) << endl; // (6) TITRE
1121       //*_sauvFile << endl;                                         // (7) 0 attributes
1122
1123       // write values of each component
1124       fcount.init( 3 ); // 3 values per a line
1125       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1126         {
1127           pair<int,int> it = iters[iIt];
1128
1129           vector<INTERP_KERNEL::NormalizedCellType> types;
1130           vector< vector<TypeOfField> > typesF;
1131           vector< vector<string> > pfls, locs;
1132           vector< vector< std::pair<int,int> > > valsVec;
1133           valsVec = _nodeFields[iF]->getFieldSplitedByType( it.first, it.second, _fileMesh->getName().c_str(),
1134                                                             types, typesF, pfls, locs);
1135           // believe that there can be only one type in a nodal field,
1136           // so do not perform a loop on types
1137           const DataArrayDouble* valsArray = _nodeFields[iF]->getUndergroundDataArray(it.first, it.second);
1138           for ( size_t j = 0; j < compInfo.size(); ++j )
1139             {
1140               for ( size_t i = valsVec[0][0].first; i < (std::size_t)valsVec[0][0].second; ++i, fcount++ )
1141                 *_sauvFile << setw(22) << valsArray->getIJ( i, j );
1142               fcount.stop();
1143             }
1144         }
1145     } // loop on fiels
1146 }
1147
1148 //================================================================================
1149 /*!
1150  * \brief Writes "PILE NUMERO  39": fields on cells
1151  */
1152 //================================================================================
1153
1154 void SauvWriter::writeElemFields(map<string,int>& fldNamePrefixMap)
1155 {
1156   writeFieldNames( /*isNodal=*/false, fldNamePrefixMap );
1157
1158   TFieldCounter fcount (*_sauvFile, 10);
1159
1160   // REAL EXAMPLE
1161
1162   // (1)        1       2       6      16
1163   // (2)                                                         CARACTERISTIQUES
1164   // (3)      -15  317773       4       0       0       0      -2       0       3
1165   // (4)             317581
1166   // (5)  0
1167   // (6)   317767  317761  317755  317815
1168   // (7)  YOUN     NU       H        SIGY
1169   // (8)  REAL*8            REAL*8            REAL*8            REAL*8
1170   // (9)        1       1       0       0
1171   // (10)  2.00000000000000E+05
1172   // (9)       1       1       0       0
1173   // (10)  3.30000000000000E-01
1174   // (9)       1       1       0       0
1175   // (10)  1.00000000000000E+04
1176   // (9)       6     706       0       0
1177   // (10)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
1178   // (10)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
1179   // (10)  ...
1180
1181   for ( size_t iF = 0; iF < _cellFields.size(); ++iF )
1182     {
1183       // count nb of sub-components
1184       int iSub, nbSub = 0;
1185       vector< pair<int,int> >  iters = _cellFields[iF]->getIterations();
1186       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1187         {
1188           pair<int,int> it = iters[iIt];
1189
1190           vector<INTERP_KERNEL::NormalizedCellType> types;
1191           vector< vector<TypeOfField> > typesF;
1192           vector< vector<string> > pfls, locs;
1193           vector< vector< std::pair<int,int> > > valsVec;
1194           valsVec = _cellFields[iF]->getFieldSplitedByType( it.first, it.second, _fileMesh->getName().c_str(),
1195                                                             types, typesF, pfls, locs);
1196           for ( size_t i = 0; i < valsVec.size(); ++i )
1197             nbSub += valsVec[i].size();
1198         }
1199       // (1) write nb sub-components, title length
1200       *_sauvFile << setw(8) << nbSub
1201                  << setw(8) << -1 // whatever
1202                  << setw(8) << 6  // whatever
1203                  << setw(8) << 72 << endl; // title length
1204       // (2) title
1205       string title = _cellFields[iF]->getName();
1206       *_sauvFile << setw(72) << title.substr(0,71) << endl;
1207       *_sauvFile << setw(72) << " " << endl;
1208
1209       // (3) support, nb components
1210       vector<int> vals(9, 0);
1211       const vector<string>& compInfo = _cellFields[iF]->getInfo();
1212       vals[2] = compInfo.size();
1213       fcount.init(10);
1214       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1215         {
1216           pair<int,int> it = iters[iIt];
1217
1218           vector<INTERP_KERNEL::NormalizedCellType> types;
1219           vector< vector<TypeOfField> > typesF;
1220           vector< vector<string> > pfls, locs;
1221           _cellFields[iF]->getFieldSplitedByType( it.first, it.second, _fileMesh->getName().c_str(),
1222                                                   types, typesF, pfls, locs);
1223           for ( size_t iType = 0; iType < pfls.size(); ++iType )
1224             for ( size_t iP = 0; iP < pfls[iType].size(); ++iP )
1225               {
1226                 if ( pfls[iType][iP].empty() ) pfls[iType][iP] = noProfileName( types[iType] );
1227                 map< string, SubMesh* >::iterator pfl2Sub = _profile2Sub.find( pfls[iType][iP] );
1228                 if ( pfl2Sub == _profile2Sub.end() )
1229                   THROW_IK_EXCEPTION( "SauvWriter::writeElemFields(): no sub-mesh for profile |"
1230                                       << pfls[iType][iP] << "|");
1231                 const int supportID = pfl2Sub->second->_id;
1232                 vals[0] = -supportID;
1233
1234                 for ( size_t i = 0; i < vals.size(); ++i, fcount++ )
1235                   *_sauvFile << setw(8) << vals[ i ];
1236               }
1237         }
1238       fcount.stop();
1239
1240       // (4) dummy strings
1241       for ( fcount.init(4), iSub = 0; iSub < nbSub; ++iSub, fcount++ )
1242         *_sauvFile << "                  ";
1243       fcount.stop();
1244
1245       // (5) dummy strings
1246       for ( fcount.init(8), iSub = 0; iSub < nbSub; ++iSub, fcount++ )
1247         *_sauvFile << "         ";
1248       fcount.stop();
1249
1250       // loop on sub-components of a field, each of which refers to
1251       // a certain support and has its own number of components
1252       for ( std::size_t iIt = 0; iIt < iters.size(); ++iIt )
1253         {
1254           pair<int,int> it = iters[iIt];
1255           writeElemTimeStamp( iF, it.first, it.second );
1256         }
1257     } // loop on cell fields
1258 }
1259
1260 //================================================================================
1261 /*!
1262  * \brief Write one elemental time stamp
1263  */
1264 //================================================================================
1265
1266 void SauvWriter::writeElemTimeStamp(int iF, int iter, int order)
1267 {
1268   // (6)   317767  317761  317755  317815
1269   // (7)  YOUN     NU       H        SIGY
1270   // (8)  REAL*8            REAL*8            REAL*8            REAL*8
1271   // (9)        1       1       0       0
1272   // (10)  2.00000000000000E+05
1273   // (9)       1       1       0       0
1274   // (10)  3.30000000000000E-01
1275   // (9)       1       1       0       0
1276   // (10)  1.00000000000000E+04
1277   // (9)       6     706       0       0
1278   // (10)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
1279   // (10)  1.00000000000000E+02  1.00000000000000E+02  1.00000000000000E+02
1280
1281   TFieldCounter fcount (*_sauvFile, 10);
1282
1283   vector<INTERP_KERNEL::NormalizedCellType> types;
1284   vector< vector<TypeOfField> > typesF;
1285   vector< vector<string> > pfls, locs;
1286   vector< vector< std::pair<int,int> > > valsVec;
1287   valsVec = _cellFields[iF]->getFieldSplitedByType( iter, order, _fileMesh->getName().c_str(),
1288                                                     types, typesF, pfls, locs);
1289   for ( size_t iType = 0; iType < pfls.size(); ++iType )
1290     for ( size_t iP = 0; iP < pfls[iType].size(); ++iP )
1291       {
1292         const vector<string>& compInfo = _cellFields[iF]->getInfo();
1293
1294         // (6) component addresses
1295         int iComp = 0, nbComp = compInfo.size();
1296         for ( fcount.init(10); iComp < nbComp; ++iComp, fcount++ )
1297           *_sauvFile << setw(8) << 777; // a good number
1298         fcount.stop();
1299
1300         // (7) component names
1301         map<string, string> mapMedToGibi;
1302         makeCompNames( _cellFields[iF]->getName(), compInfo, mapMedToGibi );
1303         *_sauvFile << left;
1304         for ( fcount.init(8), iComp = 0; iComp < nbComp; ++iComp, fcount++ )
1305           *_sauvFile << " "  << setw(8) << mapMedToGibi[compInfo[iComp]];
1306         fcount.stop();
1307
1308         // (8) component types
1309         for ( fcount.init(4), iComp = 0; iComp < nbComp; ++iComp, fcount++ )
1310           *_sauvFile << " "  << setw(17) << "REAL*8";
1311         fcount.stop();
1312         *_sauvFile << right;
1313
1314         // (9) nb values per element, nb of elements
1315         int nbPntPerCell = 1;
1316         if ( !locs[iType][iP].empty() )
1317           {
1318             int locID = _cellFields[iF]->getLocalizationId( locs[iType][iP].c_str() );
1319             nbPntPerCell = _cellFields[iF]->getNbOfGaussPtPerCell( locID );
1320           }
1321         else if ( typesF[iType][iP] == ON_GAUSS_NE )
1322           {
1323             nbPntPerCell = INTERP_KERNEL::CellModel::GetCellModel(types[iType]).getNumberOfNodes();
1324           }
1325
1326         // (10) values
1327         const std::pair<int,int>& bgEnd = valsVec[iType][iP];
1328         const DataArrayDouble* valArray = _cellFields[iF]->getUndergroundDataArray(iter, order);
1329         for ( iComp = 0; iComp < nbComp; ++iComp )
1330           {
1331             *_sauvFile << setw(8) << nbPntPerCell
1332                        << setw(8) << (bgEnd.second-bgEnd.first) / nbPntPerCell
1333                        << setw(8) << 0
1334                        << setw(8) << 0
1335                        << endl;
1336             fcount.init(3);
1337             for ( size_t i = bgEnd.first; i < (size_t) bgEnd.second; ++i, fcount++ )
1338               *_sauvFile << setw(22) << valArray->getIJ( i, iComp );
1339             fcount.stop();
1340           }
1341       }
1342 }
1343
1344 //================================================================================
1345 /*!
1346  * \brief Write the last record of the SAUV file
1347  */
1348 //================================================================================
1349
1350 void SauvWriter::writeLastRecord()
1351 {
1352   *_sauvFile << " ENREGISTREMENT DE TYPE   5" << endl;
1353   *_sauvFile << "LABEL AUTOMATIQUE :   1" << endl;
1354 }