src/DriverMED/DriverMED_W_SMESHDS_Mesh.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 //  SMESH DriverMED : driver to read and write 'med' files
  24 //  File   : DriverMED_W_SMESHDS_Mesh.cxx
  25 //  Module : SMESH
  26 //
  27 #include <sstream>
  28
  29 #include "DriverMED_W_SMESHDS_Mesh.h"
  30 #include "DriverMED_Family.h"
  31
  32 #include "SMESHDS_Mesh.hxx"
  33 #include "SMDS_MeshElement.hxx"
  34 #include "SMDS_MeshNode.hxx"
  35 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
  36
  37 #include "utilities.h"
  38
  39 #include "MED_Utilities.hxx"
  40
  41 #define _EDF_NODE_IDS_
  42 //#define _ELEMENTS_BY_DIM_
  43
  44 using namespace std;
  45 using namespace MED;
  46
  47
  48 DriverMED_W_SMESHDS_Mesh::DriverMED_W_SMESHDS_Mesh():
  49   myAllSubMeshes (false),
  50   myDoGroupOfNodes (false),
  51   myDoGroupOfEdges (false),
  52   myDoGroupOfFaces (false),
  53   myDoGroupOfVolumes (false),
  54   myDoGroupOf0DElems(false),
  55   myDoGroupOfBalls(false)
  56 {}
  57
  58 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName,
  59                                        MED::EVersion theId)
  60 {
  61   myMed = CrWrapper(theFileName,theId);
  62   Driver_SMESHDS_Mesh::SetFile(theFileName);
  63 }
  64
  65 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName)
  66 {
  67   return SetFile(theFileName,MED::eV2_2);
  68 }
  69
  70 string DriverMED_W_SMESHDS_Mesh::GetVersionString(const MED::EVersion theVersion, int theNbDigits)
  71 {
  72   TInt majeur, mineur, release;
  73   majeur =  mineur = release = 0;
  74 //   if ( theVersion == eV2_1 )
  75 //     MED::GetVersionRelease<eV2_1>(majeur, mineur, release);
  76 //   else
  77     MED::GetVersionRelease<eV2_2>(majeur, mineur, release);
  78   ostringstream name;
  79   if ( theNbDigits > 0 )
  80     name << majeur;
  81   if ( theNbDigits > 1 )
  82     name << "." << mineur;
  83   if ( theNbDigits > 2 )
  84     name << "." << release;
  85   return name.str();
  86 }
  87
  88 void DriverMED_W_SMESHDS_Mesh::AddGroup(SMESHDS_GroupBase* theGroup)
  89 {
  90   myGroups.push_back(theGroup);
  91 }
  92
  93 void DriverMED_W_SMESHDS_Mesh::AddAllSubMeshes()
  94 {
  95   myAllSubMeshes = true;
  96 }
  97
  98 void DriverMED_W_SMESHDS_Mesh::AddSubMesh(SMESHDS_SubMesh* theSubMesh, int theID)
  99 {
 100   mySubMeshes[theID] = theSubMesh;
 101 }
 102
 103 void DriverMED_W_SMESHDS_Mesh::AddGroupOfNodes()
 104 {
 105   myDoGroupOfNodes = true;
 106 }
 107
 108 void DriverMED_W_SMESHDS_Mesh::AddGroupOfEdges()
 109 {
 110   myDoGroupOfEdges = true;
 111 }
 112
 113 void DriverMED_W_SMESHDS_Mesh::AddGroupOfFaces()
 114 {
 115   myDoGroupOfFaces = true;
 116 }
 117
 118 void DriverMED_W_SMESHDS_Mesh::AddGroupOfVolumes()
 119 {
 120   myDoGroupOfVolumes = true;
 121 }
 122
 123 namespace{
 124   typedef double (SMDS_MeshNode::* TGetCoord)() const;
 125   typedef const char* TName;
 126   typedef const char* TUnit;
 127
 128   // name length in a mesh must be equal to 16 :
 129   //         1234567890123456
 130   TName M = "m               ";
 131   TName X = "x               ";
 132   TName Y = "y               ";
 133   TName Z = "z               ";
 134
 135   TUnit aUnit[3] = {M,M,M};
 136
 137   // 3 dim
 138   TGetCoord aXYZGetCoord[3] = {
 139     &SMDS_MeshNode::X,
 140     &SMDS_MeshNode::Y,
 141     &SMDS_MeshNode::Z
 142   };
 143   TName aXYZName[3] = {X,Y,Z};
 144
 145   // 2 dim
 146   TGetCoord aXYGetCoord[2] = {
 147     &SMDS_MeshNode::X,
 148     &SMDS_MeshNode::Y
 149   };
 150   TName aXYName[2] = {X,Y};
 151
 152   TGetCoord aYZGetCoord[2] = {
 153     &SMDS_MeshNode::Y,
 154     &SMDS_MeshNode::Z
 155   };
 156   TName aYZName[2] = {Y,Z};
 157
 158   TGetCoord aXZGetCoord[2] = {
 159     &SMDS_MeshNode::X,
 160     &SMDS_MeshNode::Z
 161   };
 162   TName aXZName[2] = {X,Z};
 163
 164   // 1 dim
 165   TGetCoord aXGetCoord[1] = {
 166     &SMDS_MeshNode::X
 167   };
 168   TName aXName[1] = {X};
 169
 170   TGetCoord aYGetCoord[1] = {
 171     &SMDS_MeshNode::Y
 172   };
 173   TName aYName[1] = {Y};
 174
 175   TGetCoord aZGetCoord[1] = {
 176     &SMDS_MeshNode::Z
 177   };
 178   TName aZName[1] = {Z};
 179
 180
 181   class TCoordHelper{
 182     SMDS_NodeIteratorPtr myNodeIter;
 183     const SMDS_MeshNode* myCurrentNode;
 184     TGetCoord* myGetCoord;
 185     TName* myName;
 186     TUnit* myUnit;
 187   public:
 188     TCoordHelper(const SMDS_NodeIteratorPtr& theNodeIter,
 189                  TGetCoord* theGetCoord,
 190                  TName* theName,
 191                  TUnit* theUnit = aUnit):
 192       myNodeIter(theNodeIter),
 193       myGetCoord(theGetCoord),
 194       myName(theName),
 195       myUnit(theUnit)
 196     {}
 197     virtual ~TCoordHelper(){}
 198     bool Next(){
 199       return myNodeIter->more() &&
 200         (myCurrentNode = myNodeIter->next());
 201     }
 202     const SMDS_MeshNode* GetNode(){
 203       return myCurrentNode;
 204     }
 205     MED::TIntVector::value_type GetID(){
 206       return myCurrentNode->GetID();
 207     }
 208     MED::TFloatVector::value_type GetCoord(TInt theCoodId){
 209       return (myCurrentNode->*myGetCoord[theCoodId])();
 210     }
 211     MED::TStringVector::value_type GetName(TInt theDimId){
 212       return myName[theDimId];
 213     }
 214     MED::TStringVector::value_type GetUnit(TInt theDimId){
 215       return myUnit[theDimId];
 216     }
 217   };
 218   typedef boost::shared_ptr<TCoordHelper> TCoordHelperPtr;
 219
 220   //-------------------------------------------------------
 221   /*!
 222    * \brief Structure describing element type
 223    */
 224   //-------------------------------------------------------
 225   struct TElemTypeData
 226   {
 227     EEntiteMaillage     _entity;
 228     EGeometrieElement   _geomType;
 229     TInt                _nbElems;
 230     SMDSAbs_ElementType _smdsType;
 231
 232     TElemTypeData (EEntiteMaillage entity, EGeometrieElement geom, TInt nb, SMDSAbs_ElementType type)
 233       : _entity(entity), _geomType(geom), _nbElems( nb ), _smdsType( type ) {}
 234   };
 235
 236
 237   typedef NCollection_DataMap< Standard_Address, int > TElemFamilyMap;
 238
 239   //================================================================================
 240   /*!
 241    * \brief Fills element to famaly ID map for element type.
 242    * Removes all families of anElemType
 243    */
 244   //================================================================================
 245
 246   void fillElemFamilyMap( TElemFamilyMap &            anElemFamMap,
 247                           list<DriverMED_FamilyPtr> & aFamilies,
 248                           const SMDSAbs_ElementType   anElemType)
 249   {
 250     anElemFamMap.Clear();
 251     //anElemFamMap.clear();
 252     list<DriverMED_FamilyPtr>::iterator aFamsIter = aFamilies.begin();
 253     while ( aFamsIter != aFamilies.end() )
 254     {
 255       if ((*aFamsIter)->GetType() != anElemType) {
 256         aFamsIter++;
 257       }
 258       else {
 259         int aFamId = (*aFamsIter)->GetId();
 260         const set<const SMDS_MeshElement *>& anElems = (*aFamsIter)->GetElements();
 261         set<const SMDS_MeshElement *>::const_iterator anElemsIter = anElems.begin();
 262         for (; anElemsIter != anElems.end(); anElemsIter++)
 263         {
 264           anElemFamMap.Bind( (Standard_Address)*anElemsIter, aFamId );
 265           //anElemFamMap[*anElemsIter] = aFamId;
 266         }
 267         // remove a family from the list
 268         aFamilies.erase( aFamsIter++ );
 269       }
 270     }
 271   }
 272
 273   //================================================================================
 274   /*!
 275    * \brief For an element, return family ID found in the map or a default one
 276    */
 277   //================================================================================
 278
 279   int getFamilyId( const TElemFamilyMap &  anElemFamMap,
 280                    const SMDS_MeshElement* anElement,
 281                    const int               aDefaultFamilyId)
 282   {
 283     if ( anElemFamMap.IsBound( (Standard_Address) anElement ))
 284       return anElemFamMap( (Standard_Address) anElement );
 285 //     TElemFamilyMap::iterator elem_famNum = anElemFamMap.find( anElement );
 286 //     if ( elem_famNum != anElemFamMap.end() )
 287 //       return elem_famNum->second;
 288     return aDefaultFamilyId;
 289   }
 290 }
 291
 292 Driver_Mesh::Status DriverMED_W_SMESHDS_Mesh::Perform()
 293 {
 294   Status aResult = DRS_OK;
 295   if (myMesh->hasConstructionEdges() || myMesh->hasConstructionFaces()) {
 296     INFOS("SMDS_MESH with hasConstructionEdges() or hasConstructionFaces() do not supports!!!");
 297     return DRS_FAIL;
 298   }
 299   try {
 300     MESSAGE("Perform - myFile : "<<myFile);
 301
 302     // Creating the MED mesh for corresponding SMDS structure
 303     //-------------------------------------------------------
 304     string aMeshName;
 305     if (myMeshId != -1) {
 306       ostringstream aMeshNameStr;
 307       aMeshNameStr<<myMeshId;
 308       aMeshName = aMeshNameStr.str();
 309     } else {
 310       aMeshName = myMeshName;
 311     }
 312
 313     // Mesh dimension definition
 314     TInt aSpaceDimension;
 315     TCoordHelperPtr aCoordHelperPtr;
 316     {
 317       bool anIsXDimension = false;
 318       bool anIsYDimension = false;
 319       bool anIsZDimension = false;
 320       {
 321         SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator();
 322         double aBounds[6];
 323         if(aNodesIter->more()){
 324           const SMDS_MeshNode* aNode = aNodesIter->next();
 325           aBounds[0] = aBounds[1] = aNode->X();
 326           aBounds[2] = aBounds[3] = aNode->Y();
 327           aBounds[4] = aBounds[5] = aNode->Z();
 328         }
 329         while(aNodesIter->more()){
 330           const SMDS_MeshNode* aNode = aNodesIter->next();
 331           aBounds[0] = min(aBounds[0],aNode->X());
 332           aBounds[1] = max(aBounds[1],aNode->X());
 333
 334           aBounds[2] = min(aBounds[2],aNode->Y());
 335           aBounds[3] = max(aBounds[3],aNode->Y());
 336
 337           aBounds[4] = min(aBounds[4],aNode->Z());
 338           aBounds[5] = max(aBounds[5],aNode->Z());
 339         }
 340
 341         double EPS = 1.0E-7;
 342         anIsXDimension = (aBounds[1] - aBounds[0]) + abs(aBounds[1]) + abs(aBounds[0]) > EPS;
 343         anIsYDimension = (aBounds[3] - aBounds[2]) + abs(aBounds[3]) + abs(aBounds[2]) > EPS;
 344         anIsZDimension = (aBounds[5] - aBounds[4]) + abs(aBounds[5]) + abs(aBounds[4]) > EPS;
 345         aSpaceDimension = anIsXDimension + anIsYDimension + anIsZDimension;
 346         if(!aSpaceDimension)
 347           aSpaceDimension = 3;
 348         // PAL16857(SMESH not conform to the MED convention):
 349         if ( aSpaceDimension == 2 && anIsZDimension ) // 2D only if mesh is in XOY plane
 350           aSpaceDimension = 3;
 351         // PAL18941(a saved study with a mesh belong Z is opened and the mesh is belong X)
 352         if ( aSpaceDimension == 1 && !anIsXDimension ) {// 1D only if mesh is along OX
 353           if ( anIsYDimension ) {
 354             aSpaceDimension = 2;
 355             anIsXDimension = true;
 356           } else {
 357             aSpaceDimension = 3;
 358           }
 359         }
 360       }
 361
 362       SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator(/*idInceasingOrder=*/true);
 363       switch(aSpaceDimension){
 364       case 3:
 365         aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYZGetCoord,aXYZName));
 366         break;
 367       case 2:
 368         if(anIsXDimension && anIsYDimension)
 369           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYGetCoord,aXYName));
 370         if(anIsYDimension && anIsZDimension)
 371           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYZGetCoord,aYZName));
 372         if(anIsXDimension && anIsZDimension)
 373           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXZGetCoord,aXZName));
 374         break;
 375       case 1:
 376         if(anIsXDimension)
 377           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXGetCoord,aXName));
 378         if(anIsYDimension)
 379           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYGetCoord,aYName));
 380         if(anIsZDimension)
 381           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aZGetCoord,aZName));
 382         break;
 383       }
 384     }
 385     TInt aMeshDimension = 0;
 386     if ( myMesh->NbEdges() > 0 )
 387       aMeshDimension = 1;
 388     if ( myMesh->NbFaces() > 0 )
 389       aMeshDimension = 2;
 390     if ( myMesh->NbVolumes() > 0 )
 391       aMeshDimension = 3;
 392
 393     PMeshInfo aMeshInfo = myMed->CrMeshInfo(aMeshDimension,aSpaceDimension,aMeshName);
 394     MESSAGE("Add - aMeshName : "<<aMeshName<<"; "<<aMeshInfo->GetName());
 395     myMed->SetMeshInfo(aMeshInfo);
 396
 397     // Storing SMDS groups and sub-meshes as med families
 398     //----------------------------------------------------
 399     int myNodesDefaultFamilyId      = 0;
 400     int my0DElementsDefaultFamilyId = 0;
 401     int myBallsDefaultFamilyId      = 0;
 402     int myEdgesDefaultFamilyId      = 0;
 403     int myFacesDefaultFamilyId      = 0;
 404     int myVolumesDefaultFamilyId    = 0;
 405     int nbNodes      = myMesh->NbNodes();
 406     int nb0DElements = myMesh->Nb0DElements();
 407     int nbBalls      = myMesh->NbBalls();
 408     int nbEdges      = myMesh->NbEdges();
 409     int nbFaces      = myMesh->NbFaces();
 410     int nbVolumes    = myMesh->NbVolumes();
 411     if (myDoGroupOfNodes && nbNodes) myNodesDefaultFamilyId = REST_NODES_FAMILY;
 412     if (myDoGroupOfEdges && nbEdges) myEdgesDefaultFamilyId = REST_EDGES_FAMILY;
 413     if (myDoGroupOfFaces && nbFaces) myFacesDefaultFamilyId = REST_FACES_FAMILY;
 414     if (myDoGroupOfVolumes && nbVolumes) myVolumesDefaultFamilyId = REST_VOLUMES_FAMILY;
 415     if (myDoGroupOf0DElems && nb0DElements) my0DElementsDefaultFamilyId = REST_0DELEM_FAMILY;
 416     if (myDoGroupOfBalls && nbBalls) myBallsDefaultFamilyId = REST_BALL_FAMILY;
 417
 418     MESSAGE("Perform - aFamilyInfo");
 419     //cout << " DriverMED_Family::MakeFamilies() " << endl;
 420     list<DriverMED_FamilyPtr> aFamilies;
 421     if (myAllSubMeshes) {
 422       aFamilies = DriverMED_Family::MakeFamilies
 423         (myMesh->SubMeshes(), myGroups,
 424          myDoGroupOfNodes   && nbNodes,
 425          myDoGroupOfEdges   && nbEdges,
 426          myDoGroupOfFaces   && nbFaces,
 427          myDoGroupOfVolumes && nbVolumes,
 428          myDoGroupOf0DElems && nb0DElements,
 429          myDoGroupOfBalls   && nbBalls);
 430     } else {
 431       aFamilies = DriverMED_Family::MakeFamilies
 432         (mySubMeshes, myGroups,
 433          myDoGroupOfNodes   && nbNodes,
 434          myDoGroupOfEdges   && nbEdges,
 435          myDoGroupOfFaces   && nbFaces,
 436          myDoGroupOfVolumes && nbVolumes,
 437          myDoGroupOf0DElems && nb0DElements,
 438          myDoGroupOfBalls   && nbBalls);
 439     }
 440     //cout << " myMed->SetFamilyInfo() " << endl;
 441     list<DriverMED_FamilyPtr>::iterator aFamsIter;
 442     for (aFamsIter = aFamilies.begin(); aFamsIter != aFamilies.end(); aFamsIter++)
 443     {
 444       PFamilyInfo aFamilyInfo = (*aFamsIter)->GetFamilyInfo(myMed,aMeshInfo);
 445       myMed->SetFamilyInfo(aFamilyInfo);
 446     }
 447
 448     // Storing SMDS nodes to the MED file for the MED mesh
 449     //----------------------------------------------------
 450 #ifdef _EDF_NODE_IDS_
 451     typedef map<TInt,TInt> TNodeIdMap;
 452     TNodeIdMap aNodeIdMap;
 453 #endif
 454     const EModeSwitch   theMode        = eFULL_INTERLACE;
 455     const ERepere       theSystem      = eCART;
 456     const EBooleen      theIsElemNum   = eVRAI;
 457     const EBooleen      theIsElemNames = eFAUX;
 458     const EConnectivite theConnMode    = eNOD;
 459
 460     TInt aNbNodes = myMesh->NbNodes();
 461     //cout << " myMed->CrNodeInfo() aNbNodes = " << aNbNodes << endl;
 462     PNodeInfo aNodeInfo = myMed->CrNodeInfo(aMeshInfo, aNbNodes,
 463                                             theMode, theSystem, theIsElemNum, theIsElemNames);
 464
 465     //cout << " fillElemFamilyMap( SMDSAbs_Node )" << endl;
 466     // find family numbers for nodes
 467     TElemFamilyMap anElemFamMap;
 468     fillElemFamilyMap( anElemFamMap, aFamilies, SMDSAbs_Node );
 469
 470     for (TInt iNode = 0; aCoordHelperPtr->Next(); iNode++)
 471     {
 472       // coordinates
 473       TCoordSlice aTCoordSlice = aNodeInfo->GetCoordSlice( iNode );
 474       for(TInt iCoord = 0; iCoord < aSpaceDimension; iCoord++){
 475         aTCoordSlice[iCoord] = aCoordHelperPtr->GetCoord(iCoord);
 476       }
 477       // node number
 478       int aNodeID = aCoordHelperPtr->GetID();
 479       aNodeInfo->SetElemNum( iNode, aNodeID );
 480 #ifdef _EDF_NODE_IDS_
 481       aNodeIdMap[aNodeID] = iNode+1;
 482 #endif
 483       // family number
 484       const SMDS_MeshNode* aNode = aCoordHelperPtr->GetNode();
 485       int famNum = getFamilyId( anElemFamMap, aNode, myNodesDefaultFamilyId );
 486       aNodeInfo->SetFamNum( iNode, famNum );
 487     }
 488     anElemFamMap.Clear();
 489
 490     // coordinate names and units
 491     for (TInt iCoord = 0; iCoord < aSpaceDimension; iCoord++) {
 492       aNodeInfo->SetCoordName( iCoord, aCoordHelperPtr->GetName(iCoord));
 493       aNodeInfo->SetCoordUnit( iCoord, aCoordHelperPtr->GetUnit(iCoord));
 494     }
 495
 496     //cout << " SetNodeInfo(aNodeInfo)" << endl;
 497     MESSAGE("Perform - aNodeInfo->GetNbElem() = "<<aNbNodes);
 498     myMed->SetNodeInfo(aNodeInfo);
 499     aNodeInfo.reset(); // free memory used for arrays
 500
 501
 502     // Storing SMDS elements to the MED file for the MED mesh
 503     //-------------------------------------------------------
 504     // Write one element type at once in order to minimize memory usage (PAL19276)
 505
 506     const SMDS_MeshInfo& nbElemInfo = myMesh->GetMeshInfo();
 507
 508     // poly elements are not supported by med-2.1
 509     bool polyTypesSupported = myMed->CrPolygoneInfo(aMeshInfo,eMAILLE,ePOLYGONE,0,0);
 510     TInt nbPolygonNodes = 0, nbPolyhedronNodes = 0, nbPolyhedronFaces = 0;
 511
 512     // collect info on all geom types
 513
 514     list< TElemTypeData > aTElemTypeDatas;
 515
 516     EEntiteMaillage anEntity = eMAILLE;
 517 #ifdef _ELEMENTS_BY_DIM_
 518     anEntity = eNOEUD_ELEMENT;
 519 #endif
 520     aTElemTypeDatas.push_back(TElemTypeData(anEntity,
 521                                             ePOINT1,
 522                                             nbElemInfo.Nb0DElements(),
 523                                             SMDSAbs_0DElement));
 524 #ifdef _ELEMENTS_BY_DIM_
 525     anEntity = eSTRUCT_ELEMENT;
 526 #endif
 527     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 528                                              eBALL,
 529                                              nbElemInfo.NbBalls(),
 530                                              SMDSAbs_Ball));
 531 #ifdef _ELEMENTS_BY_DIM_
 532     anEntity = eARETE;
 533 #endif
 534     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 535                                             eSEG2,
 536                                             nbElemInfo.NbEdges( ORDER_LINEAR ),
 537                                             SMDSAbs_Edge));
 538     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 539                                             eSEG3,
 540                                             nbElemInfo.NbEdges( ORDER_QUADRATIC ),
 541                                             SMDSAbs_Edge));
 542 #ifdef _ELEMENTS_BY_DIM_
 543     anEntity = eFACE;
 544 #endif
 545     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 546                                             eTRIA3,
 547                                             nbElemInfo.NbTriangles( ORDER_LINEAR ),
 548                                             SMDSAbs_Face));
 549     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 550                                             eTRIA6,
 551                                             nbElemInfo.NbTriangles( ORDER_QUADRATIC ),
 552                                             SMDSAbs_Face));
 553     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 554                                             eQUAD4,
 555                                             nbElemInfo.NbQuadrangles( ORDER_LINEAR ),
 556                                             SMDSAbs_Face));
 557     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 558                                              eQUAD8,
 559                                              nbElemInfo.NbQuadrangles( ORDER_QUADRATIC ) -
 560                                              nbElemInfo.NbBiQuadQuadrangles(),
 561                                              SMDSAbs_Face));
 562     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 563                                              eQUAD9,
 564                                              nbElemInfo.NbBiQuadQuadrangles(),
 565                                              SMDSAbs_Face));
 566     if ( polyTypesSupported ) {
 567       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 568                                                ePOLYGONE,
 569                                                nbElemInfo.NbPolygons(),
 570                                                SMDSAbs_Face));
 571       // we need one more loop on poly elements to count nb of their nodes
 572       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 573                                                ePOLYGONE,
 574                                                nbElemInfo.NbPolygons(),
 575                                                SMDSAbs_Face));
 576     }
 577 #ifdef _ELEMENTS_BY_DIM_
 578     anEntity = eMAILLE;
 579 #endif
 580     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 581                                             eTETRA4,
 582                                             nbElemInfo.NbTetras( ORDER_LINEAR ),
 583                                             SMDSAbs_Volume));
 584     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 585                                             eTETRA10,
 586                                             nbElemInfo.NbTetras( ORDER_QUADRATIC ),
 587                                             SMDSAbs_Volume));
 588     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 589                                             ePYRA5,
 590                                             nbElemInfo.NbPyramids( ORDER_LINEAR ),
 591                                             SMDSAbs_Volume));
 592     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 593                                             ePYRA13,
 594                                             nbElemInfo.NbPyramids( ORDER_QUADRATIC ),
 595                                             SMDSAbs_Volume));
 596     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 597                                             ePENTA6,
 598                                             nbElemInfo.NbPrisms( ORDER_LINEAR ),
 599                                             SMDSAbs_Volume));
 600     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 601                                              ePENTA15,
 602                                              nbElemInfo.NbPrisms( ORDER_QUADRATIC ),
 603                                              SMDSAbs_Volume));
 604     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 605                                              eHEXA8,
 606                                              nbElemInfo.NbHexas( ORDER_LINEAR ),
 607                                              SMDSAbs_Volume));
 608     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 609                                              eHEXA20,
 610                                              nbElemInfo.NbHexas( ORDER_QUADRATIC )-
 611                                              nbElemInfo.NbTriQuadHexas(),
 612                                              SMDSAbs_Volume));
 613     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 614                                              eHEXA27,
 615                                              nbElemInfo.NbTriQuadHexas(),
 616                                              SMDSAbs_Volume));
 617     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 618                                              eOCTA12,
 619                                              nbElemInfo.NbHexPrisms(),
 620                                              SMDSAbs_Volume));
 621     if ( polyTypesSupported ) {
 622       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 623                                                ePOLYEDRE,
 624                                                nbElemInfo.NbPolyhedrons(),
 625                                                SMDSAbs_Volume));
 626       // we need one more loop on poly elements to count nb of their nodes
 627       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 628                                                ePOLYEDRE,
 629                                                nbElemInfo.NbPolyhedrons(),
 630                                                SMDSAbs_Volume));
 631     }
 632
 633     vector< bool > isElemFamMapBuilt( SMDSAbs_NbElementTypes, false );
 634
 635     // loop on all geom types of elements
 636
 637     list< TElemTypeData >::iterator aElemTypeData = aTElemTypeDatas.begin();
 638     for ( ; aElemTypeData != aTElemTypeDatas.end(); ++aElemTypeData )
 639     {
 640       if ( aElemTypeData->_nbElems == 0 )
 641         continue;
 642
 643       int defaultFamilyId = 0;
 644       switch ( aElemTypeData->_smdsType ) {
 645       case SMDSAbs_0DElement:
 646         defaultFamilyId = my0DElementsDefaultFamilyId;
 647         break;
 648       case SMDSAbs_Ball:
 649         defaultFamilyId = myBallsDefaultFamilyId;
 650         break;
 651       case SMDSAbs_Edge:
 652         defaultFamilyId = myEdgesDefaultFamilyId;
 653         break;
 654       case SMDSAbs_Face:
 655         defaultFamilyId = myFacesDefaultFamilyId;
 656         break;
 657       case SMDSAbs_Volume:
 658         defaultFamilyId = myVolumesDefaultFamilyId;
 659         break;
 660       default:
 661         continue;
 662       }
 663
 664       // iterator on elements of a current type
 665       SMDS_ElemIteratorPtr elemIterator;
 666       int iElem = 0;
 667
 668       // Treat POLYGONs
 669       // ---------------
 670       if ( aElemTypeData->_geomType == ePOLYGONE )
 671       {
 672         elemIterator = myMesh->elementGeomIterator( SMDSGeom_POLYGON );
 673         if ( nbPolygonNodes == 0 ) {
 674           // Count nb of nodes
 675           while ( elemIterator->more() ) {
 676             const SMDS_MeshElement* anElem = elemIterator->next();
 677             nbPolygonNodes += anElem->NbNodes();
 678             if ( ++iElem == aElemTypeData->_nbElems )
 679               break;
 680           }
 681         }
 682         else {
 683           // Store in med file
 684           PPolygoneInfo aPolygoneInfo = myMed->CrPolygoneInfo(aMeshInfo,
 685                                                               aElemTypeData->_entity,
 686                                                               aElemTypeData->_geomType,
 687                                                               aElemTypeData->_nbElems,
 688                                                               nbPolygonNodes,
 689                                                               theConnMode, theIsElemNum,
 690                                                               theIsElemNames);
 691           TElemNum & index = *(aPolygoneInfo->myIndex.get());
 692           index[0] = 1;
 693
 694           while ( elemIterator->more() )
 695           {
 696             const SMDS_MeshElement* anElem = elemIterator->next();
 697             // index
 698             TInt aNbNodes = anElem->NbNodes();
 699             index[ iElem+1 ] = index[ iElem ] + aNbNodes;
 700
 701             // connectivity
 702             TConnSlice aTConnSlice = aPolygoneInfo->GetConnSlice( iElem );
 703             for(TInt iNode = 0; iNode < aNbNodes; iNode++) {
 704               const SMDS_MeshElement* aNode = anElem->GetNode( iNode );
 705 #ifdef _EDF_NODE_IDS_
 706               aTConnSlice[ iNode ] = aNodeIdMap[aNode->GetID()];
 707 #else
 708               aTConnSlice[ iNode ] = aNode->GetID();
 709 #endif
 710             }
 711             // element number
 712             aPolygoneInfo->SetElemNum( iElem, anElem->GetID() );
 713
 714             // family number
 715             int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 716             aPolygoneInfo->SetFamNum( iElem, famNum );
 717
 718             if ( ++iElem == aPolygoneInfo->GetNbElem() )
 719               break;
 720           }
 721           myMed->SetPolygoneInfo(aPolygoneInfo);
 722         }
 723
 724       }
 725
 726       // Treat POLYEDREs
 727       // ----------------
 728       else if (aElemTypeData->_geomType == ePOLYEDRE )
 729       {
 730         elemIterator = myMesh->elementGeomIterator( SMDSGeom_POLYHEDRA );
 731
 732         if ( nbPolyhedronNodes == 0 ) {
 733           // Count nb of nodes
 734           while ( elemIterator->more() ) {
 735             const SMDS_MeshElement* anElem = elemIterator->next();
 736             const SMDS_VtkVolume *aPolyedre = dynamic_cast<const SMDS_VtkVolume*>(anElem);
 737             if ( !aPolyedre ) continue;
 738             nbPolyhedronNodes += aPolyedre->NbNodes();
 739             nbPolyhedronFaces += aPolyedre->NbFaces();
 740             if ( ++iElem == aElemTypeData->_nbElems )
 741               break;
 742           }
 743         }
 744         else {
 745           // Store in med file
 746           PPolyedreInfo aPolyhInfo = myMed->CrPolyedreInfo(aMeshInfo,
 747                                                            aElemTypeData->_entity,
 748                                                            aElemTypeData->_geomType,
 749                                                            aElemTypeData->_nbElems,
 750                                                            nbPolyhedronFaces+1,
 751                                                            nbPolyhedronNodes,
 752                                                            theConnMode,
 753                                                            theIsElemNum,
 754                                                            theIsElemNames);
 755           TElemNum & index = *(aPolyhInfo->myIndex.get());
 756           TElemNum & faces = *(aPolyhInfo->myFaces.get());
 757           TElemNum & conn  = *(aPolyhInfo->myConn.get());
 758           index[0] = 1;
 759           faces[0] = 1;
 760
 761           TInt iFace = 0, iNode = 0;
 762           while ( elemIterator->more() )
 763           {
 764             const SMDS_MeshElement* anElem = elemIterator->next();
 765             const SMDS_VtkVolume *aPolyedre = dynamic_cast<const SMDS_VtkVolume*>(anElem);
 766             if ( !aPolyedre ) continue;
 767             // index
 768             TInt aNbFaces = aPolyedre->NbFaces();
 769             index[ iElem+1 ] = index[ iElem ] + aNbFaces;
 770
 771             // face index
 772             for (TInt f = 1; f <= aNbFaces; ++f, ++iFace ) {
 773               int aNbFaceNodes = aPolyedre->NbFaceNodes( f );
 774               faces[ iFace+1 ] = faces[ iFace ] + aNbFaceNodes;
 775             }
 776             // connectivity
 777             SMDS_ElemIteratorPtr nodeIt = anElem->nodesIterator();
 778             while ( nodeIt->more() ) {
 779               const SMDS_MeshElement* aNode = nodeIt->next();
 780 #ifdef _EDF_NODE_IDS_
 781               conn[ iNode ] = aNodeIdMap[aNode->GetID()];
 782 #else
 783               conn[ iNode ] = aNode->GetID();
 784 #endif
 785               ++iNode;
 786             }
 787             // element number
 788             aPolyhInfo->SetElemNum( iElem, anElem->GetID() );
 789
 790             // family number
 791             int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 792             aPolyhInfo->SetFamNum( iElem, famNum );
 793
 794             if ( ++iElem == aPolyhInfo->GetNbElem() )
 795               break;
 796           }
 797           myMed->SetPolyedreInfo(aPolyhInfo);
 798         }
 799       } // if (aElemTypeData->_geomType == ePOLYEDRE )
 800
 801       // Treat BALLs
 802       // ----------------
 803       else if (aElemTypeData->_geomType == eBALL )
 804       {
 805         // allocate data arrays
 806         PBallInfo aBallInfo = myMed->CrBallInfo( aMeshInfo,
 807                                                  aElemTypeData->_nbElems );
 808
 809         // build map of family numbers for this type
 810         if ( !isElemFamMapBuilt[ aElemTypeData->_smdsType ])
 811         {
 812           fillElemFamilyMap( anElemFamMap, aFamilies, aElemTypeData->_smdsType );
 813           isElemFamMapBuilt[ aElemTypeData->_smdsType ] = true;
 814         }
 815
 816         elemIterator = myMesh->elementsIterator( SMDSAbs_Ball );
 817         while ( elemIterator->more() )
 818         {
 819           const SMDS_MeshElement* anElem = elemIterator->next();
 820           // connectivity
 821           const SMDS_MeshElement* aNode = anElem->GetNode( 0 );
 822 #ifdef _EDF_NODE_IDS_
 823           (*aBallInfo->myConn)[ iElem ] = aNodeIdMap[aNode->GetID()];
 824 #else
 825           (*aBallInfo->myConn)[ iElem ] = aNode->GetID();
 826 #endif
 827           // element number
 828           aBallInfo->SetElemNum( iElem, anElem->GetID() );
 829
 830           // diameter
 831           aBallInfo->myDiameters[ iElem ] =
 832             static_cast<const SMDS_BallElement*>( anElem )->GetDiameter();
 833
 834           // family number
 835           int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 836           aBallInfo->SetFamNum( iElem, famNum );
 837           ++iElem;
 838         }
 839         // store data in a file
 840         myMed->SetBallInfo(aBallInfo);
 841       }
 842
 843       else
 844       {
 845         // Treat standard types
 846         // ---------------------
 847
 848         // allocate data arrays
 849         PCellInfo aCellInfo = myMed->CrCellInfo( aMeshInfo,
 850                                                  aElemTypeData->_entity,
 851                                                  aElemTypeData->_geomType,
 852                                                  aElemTypeData->_nbElems,
 853                                                  theConnMode,
 854                                                  theIsElemNum,
 855                                                  theIsElemNames);
 856         // build map of family numbers for this type
 857         if ( !isElemFamMapBuilt[ aElemTypeData->_smdsType ])
 858         {
 859           //cout << " fillElemFamilyMap()" << endl;
 860           fillElemFamilyMap( anElemFamMap, aFamilies, aElemTypeData->_smdsType );
 861           isElemFamMapBuilt[ aElemTypeData->_smdsType ] = true;
 862         }
 863
 864         TInt aNbNodes = MED::GetNbNodes(aElemTypeData->_geomType);
 865         elemIterator = myMesh->elementsIterator( aElemTypeData->_smdsType );
 866         while ( elemIterator->more() )
 867         {
 868           const SMDS_MeshElement* anElem = elemIterator->next();
 869           if ( anElem->NbNodes() != aNbNodes || anElem->IsPoly() )
 870             continue; // other geometry
 871
 872           // connectivity
 873           TConnSlice aTConnSlice = aCellInfo->GetConnSlice( iElem );
 874           for (TInt iNode = 0; iNode < aNbNodes; iNode++) {
 875             const SMDS_MeshElement* aNode = anElem->GetNode( iNode );
 876 #ifdef _EDF_NODE_IDS_
 877             aTConnSlice[ iNode ] = aNodeIdMap[aNode->GetID()];
 878 #else
 879             aTConnSlice[ iNode ] = aNode->GetID();
 880 #endif
 881           }
 882           // element number
 883           aCellInfo->SetElemNum( iElem, anElem->GetID() );
 884
 885           // family number
 886           int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 887           aCellInfo->SetFamNum( iElem, famNum );
 888
 889           if ( ++iElem == aCellInfo->GetNbElem() )
 890             break;
 891         }
 892         // store data in a file
 893         myMed->SetCellInfo(aCellInfo);
 894       }
 895
 896     } // loop on geom types
 897
 898
 899   }
 900   catch(const std::exception& exc) {
 901     INFOS("The following exception was caught:\n\t"<<exc.what());
 902     throw;
 903   }
 904   catch(...) {
 905     INFOS("Unknown exception was caught !!!");
 906     throw;
 907   }
 908
 909   myMeshId = -1;
 910   myGroups.clear();
 911   mySubMeshes.clear();
 912   return aResult;
 913 }