Salome HOME
PAL16834 (smesh Prism don't work with NETGEN_2D algorithm)
[modules/smesh.git] / src / DriverMED / DriverMED_W_SMESHDS_Mesh.cxx
1 //  SMESH DriverMED : driver to read and write 'med' files
2 //
3 //  Copyright (C) 2003  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 //
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_W_SMDS_Mesh.h"
31 #include "DriverMED_Family.h"
32
33 #include "SMESHDS_Mesh.hxx"
34 #include "SMDS_MeshElement.hxx"
35 #include "SMDS_MeshNode.hxx"
36 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
37
38 #include "utilities.h"
39
40 #include "MED_Utilities.hxx"
41
42 #define _EDF_NODE_IDS_
43 //#define _ELEMENTS_BY_DIM_
44
45 using namespace std;
46 using namespace MED;
47
48
49 DriverMED_W_SMESHDS_Mesh::DriverMED_W_SMESHDS_Mesh():
50   myAllSubMeshes (false),
51   myDoGroupOfNodes (false),
52   myDoGroupOfEdges (false),
53   myDoGroupOfFaces (false),
54   myDoGroupOfVolumes (false)
55 {}
56
57 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName, 
58                                        MED::EVersion theId)
59 {
60   myMed = CrWrapper(theFileName,theId);
61   Driver_SMESHDS_Mesh::SetFile(theFileName);
62 }
63
64 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName)
65 {
66   return SetFile(theFileName,MED::eV2_2);
67 }
68
69 void DriverMED_W_SMESHDS_Mesh::SetMeshName(const std::string& theMeshName)
70 {
71   myMeshName = theMeshName;
72 }
73
74 void DriverMED_W_SMESHDS_Mesh::AddGroup(SMESHDS_GroupBase* theGroup)
75 {
76   myGroups.push_back(theGroup);
77 }
78
79 void DriverMED_W_SMESHDS_Mesh::AddAllSubMeshes()
80 {
81   myAllSubMeshes = true;
82 }
83
84 void DriverMED_W_SMESHDS_Mesh::AddSubMesh(SMESHDS_SubMesh* theSubMesh, int theID)
85 {
86   mySubMeshes[theID] = theSubMesh;
87 }
88
89 void DriverMED_W_SMESHDS_Mesh::AddGroupOfNodes()
90 {
91   myDoGroupOfNodes = true;
92 }
93
94 void DriverMED_W_SMESHDS_Mesh::AddGroupOfEdges()
95 {
96   myDoGroupOfEdges = true;
97 }
98
99 void DriverMED_W_SMESHDS_Mesh::AddGroupOfFaces()
100 {
101   myDoGroupOfFaces = true;
102 }
103
104 void DriverMED_W_SMESHDS_Mesh::AddGroupOfVolumes()
105 {
106   myDoGroupOfVolumes = true;
107 }
108
109 namespace{
110   typedef double (SMDS_MeshNode::* TGetCoord)() const;
111   typedef const char* TName;
112   typedef const char* TUnit;
113
114   // name length in a mesh must be equal to 16 :
115   //         1234567890123456
116   TName M = "m               ";
117   TName X = "x               ";
118   TName Y = "y               ";
119   TName Z = "z               ";
120
121   TUnit aUnit[3] = {M,M,M};
122
123   // 3 dim
124   TGetCoord aXYZGetCoord[3] = {
125     &SMDS_MeshNode::X, 
126     &SMDS_MeshNode::Y, 
127     &SMDS_MeshNode::Z
128   };
129   TName aXYZName[3] = {X,Y,Z};
130   
131   // 2 dim
132   TGetCoord aXYGetCoord[2] = {
133     &SMDS_MeshNode::X, 
134     &SMDS_MeshNode::Y
135   };
136   TName aXYName[2] = {X,Y};
137
138   TGetCoord aYZGetCoord[2] = {
139     &SMDS_MeshNode::Y, 
140     &SMDS_MeshNode::Z
141   };
142   TName aYZName[2] = {Y,Z};
143
144   TGetCoord aXZGetCoord[2] = {
145     &SMDS_MeshNode::X, 
146     &SMDS_MeshNode::Z
147   };
148   TName aXZName[2] = {X,Z};
149
150   // 1 dim
151   TGetCoord aXGetCoord[1] = {
152     &SMDS_MeshNode::X
153   };
154   TName aXName[1] = {X};
155
156   TGetCoord aYGetCoord[1] = {
157     &SMDS_MeshNode::Y
158   };
159   TName aYName[1] = {Y};
160
161   TGetCoord aZGetCoord[1] = {
162     &SMDS_MeshNode::Z
163   };
164   TName aZName[1] = {Z};
165
166
167   class TCoordHelper{
168     SMDS_NodeIteratorPtr myNodeIter;
169     const SMDS_MeshNode* myCurrentNode;
170     TGetCoord* myGetCoord;
171     TName* myName;
172     TUnit* myUnit;
173   public:
174     TCoordHelper(const SMDS_NodeIteratorPtr& theNodeIter,
175                  TGetCoord* theGetCoord,
176                  TName* theName,
177                  TUnit* theUnit = aUnit):
178       myNodeIter(theNodeIter),
179       myGetCoord(theGetCoord),
180       myName(theName),
181       myUnit(theUnit)
182     {}
183     virtual ~TCoordHelper(){}
184     bool Next(){ 
185       return myNodeIter->more() && 
186         (myCurrentNode = myNodeIter->next());
187     }
188     const SMDS_MeshNode* GetNode(){
189       return myCurrentNode;
190     }
191     MED::TIntVector::value_type GetID(){
192       return myCurrentNode->GetID();
193     }
194     MED::TFloatVector::value_type GetCoord(TInt theCoodId){
195       return (myCurrentNode->*myGetCoord[theCoodId])();
196     }
197     MED::TStringVector::value_type GetName(TInt theDimId){
198       return myName[theDimId];
199     }
200     MED::TStringVector::value_type GetUnit(TInt theDimId){
201       return myUnit[theDimId];
202     }
203   };
204   typedef boost::shared_ptr<TCoordHelper> TCoordHelperPtr;
205   
206 }
207
208   
209 Driver_Mesh::Status DriverMED_W_SMESHDS_Mesh::Perform()
210 {
211   Status aResult = DRS_OK;
212   if (myMesh->hasConstructionEdges() || myMesh->hasConstructionFaces()) {
213     INFOS("SMDS_MESH with hasConstructionEdges() or hasConstructionFaces() do not supports!!!");
214     return DRS_FAIL;
215   }
216   try{
217     MESSAGE("Perform - myFile : "<<myFile);
218
219     // Creating the MED mesh for corresponding SMDS structure
220     //-------------------------------------------------------
221     string aMeshName;
222     if (myMeshId != -1) {
223       ostringstream aMeshNameStr;
224       aMeshNameStr<<myMeshId;
225       aMeshName = aMeshNameStr.str();
226     } else {
227       aMeshName = myMeshName;
228     }
229
230     // Mesh dimension definition
231     TInt aMeshDimension;
232     TCoordHelperPtr aCoordHelperPtr;
233     {  
234       bool anIsXDimension = false;
235       bool anIsYDimension = false;
236       bool anIsZDimension = false;
237       {
238         SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator();
239         double aBounds[6];
240         if(aNodesIter->more()){
241           const SMDS_MeshNode* aNode = aNodesIter->next();
242           aBounds[0] = aBounds[1] = aNode->X();
243           aBounds[2] = aBounds[3] = aNode->Y();
244           aBounds[4] = aBounds[5] = aNode->Z();
245         }
246         while(aNodesIter->more()){
247           const SMDS_MeshNode* aNode = aNodesIter->next();
248           aBounds[0] = min(aBounds[0],aNode->X());
249           aBounds[1] = max(aBounds[1],aNode->X());
250           
251           aBounds[2] = min(aBounds[2],aNode->Y());
252           aBounds[3] = max(aBounds[3],aNode->Y());
253           
254           aBounds[4] = min(aBounds[4],aNode->Z());
255           aBounds[5] = max(aBounds[5],aNode->Z());
256         }
257
258         double EPS = 1.0E-7;
259         anIsXDimension = (aBounds[1] - aBounds[0]) + abs(aBounds[1]) + abs(aBounds[0]) > EPS;
260         anIsYDimension = (aBounds[3] - aBounds[2]) + abs(aBounds[3]) + abs(aBounds[2]) > EPS;
261         anIsZDimension = (aBounds[5] - aBounds[4]) + abs(aBounds[5]) + abs(aBounds[4]) > EPS;
262         aMeshDimension = anIsXDimension + anIsYDimension + anIsZDimension;
263         if(!aMeshDimension)
264           aMeshDimension = 3;
265       }
266
267       SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator();
268       switch(aMeshDimension){
269       case 3:
270         aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYZGetCoord,aXYZName));
271         break;
272       case 2:
273         if(anIsXDimension && anIsYDimension)
274           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYGetCoord,aXYName));
275         if(anIsYDimension && anIsZDimension)
276           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYZGetCoord,aYZName));
277         if(anIsXDimension && anIsZDimension)
278           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXZGetCoord,aXZName));
279         break;
280       case 1:
281         if(anIsXDimension)
282           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXGetCoord,aXName));
283         if(anIsYDimension)
284           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYGetCoord,aYName));
285         if(anIsZDimension)
286           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aZGetCoord,aZName));
287         break;
288       }
289     }
290
291     
292     PMeshInfo aMeshInfo = myMed->CrMeshInfo(aMeshDimension,aMeshName);
293     MESSAGE("Add - aMeshName : "<<aMeshName<<"; "<<aMeshInfo->GetName());
294     myMed->SetMeshInfo(aMeshInfo);
295
296     // Storing SMDS groups and sub-meshes
297     //-----------------------------------
298     int myNodesDefaultFamilyId   = 0;
299     int myEdgesDefaultFamilyId   = 0;
300     int myFacesDefaultFamilyId   = 0;
301     int myVolumesDefaultFamilyId = 0;
302     int nbNodes   = myMesh->NbNodes();
303     int nbEdges   = myMesh->NbEdges();
304     int nbFaces   = myMesh->NbFaces();
305     int nbVolumes = myMesh->NbVolumes();
306     if (myDoGroupOfNodes && nbNodes)
307       myNodesDefaultFamilyId = REST_NODES_FAMILY;
308     if (myDoGroupOfEdges && nbEdges)
309       myEdgesDefaultFamilyId = REST_EDGES_FAMILY;
310     if (myDoGroupOfFaces && nbFaces)
311       myFacesDefaultFamilyId = REST_FACES_FAMILY;
312     if (myDoGroupOfVolumes && nbVolumes)
313       myVolumesDefaultFamilyId = REST_VOLUMES_FAMILY;
314
315     MESSAGE("Perform - aFamilyInfo");
316     map<const SMDS_MeshElement *, int> anElemFamMap;
317     list<DriverMED_FamilyPtr> aFamilies;
318     if (myAllSubMeshes) {
319       aFamilies = DriverMED_Family::MakeFamilies
320         (myMesh->SubMeshes(), myGroups,
321          myDoGroupOfNodes   && nbNodes,
322          myDoGroupOfEdges   && nbEdges,
323          myDoGroupOfFaces   && nbFaces,
324          myDoGroupOfVolumes && nbVolumes);
325     } else {
326       aFamilies = DriverMED_Family::MakeFamilies
327         (mySubMeshes, myGroups,
328          myDoGroupOfNodes   && nbNodes,
329          myDoGroupOfEdges   && nbEdges,
330          myDoGroupOfFaces   && nbFaces,
331          myDoGroupOfVolumes && nbVolumes);
332     }
333     list<DriverMED_FamilyPtr>::iterator aFamsIter = aFamilies.begin();
334
335     for (; aFamsIter != aFamilies.end(); aFamsIter++)
336     {
337       PFamilyInfo aFamilyInfo = (*aFamsIter)->GetFamilyInfo(myMed,aMeshInfo);
338       myMed->SetFamilyInfo(aFamilyInfo);
339       int aFamId = (*aFamsIter)->GetId();
340
341       const set<const SMDS_MeshElement *>& anElems = (*aFamsIter)->GetElements();
342       set<const SMDS_MeshElement *>::const_iterator anElemsIter = anElems.begin();
343       for (; anElemsIter != anElems.end(); anElemsIter++)
344       {
345         anElemFamMap[*anElemsIter] = aFamId;
346       }
347     }
348
349     // Storing SMDS nodes to the MED file for the MED mesh
350     //----------------------------------------------------
351 #ifdef _EDF_NODE_IDS_
352     typedef map<TInt,TInt> TNodeIdMap;
353     TNodeIdMap aNodeIdMap;
354 #endif
355     TInt aNbElems = myMesh->NbNodes();
356     MED::TIntVector anElemNums(aNbElems);
357     MED::TIntVector aFamilyNums(aNbElems);
358     MED::TFloatVector aCoordinates(aNbElems*aMeshDimension);
359     for(TInt iNode = 0, aStartId = 0; aCoordHelperPtr->Next(); iNode++, aStartId += aMeshDimension){
360       for(TInt iCoord = 0; iCoord < aMeshDimension; iCoord++){
361         aCoordinates[aStartId+iCoord] = aCoordHelperPtr->GetCoord(iCoord);
362       }
363       int aNodeID = aCoordHelperPtr->GetID();
364       anElemNums[iNode] = aNodeID;
365 #ifdef _EDF_NODE_IDS_
366       aNodeIdMap[aNodeID] = iNode+1;
367 #endif
368       const SMDS_MeshNode* aNode = aCoordHelperPtr->GetNode();
369       if (anElemFamMap.find(aNode) != anElemFamMap.end())
370         aFamilyNums[iNode] = anElemFamMap[aNode];
371       else
372         aFamilyNums[iNode] = myNodesDefaultFamilyId;
373     }
374
375     MED::TStringVector aCoordNames(aMeshDimension);
376     MED::TStringVector aCoordUnits(aMeshDimension);
377     for(TInt iCoord = 0; iCoord < aMeshDimension; iCoord++){
378       aCoordNames[iCoord] = aCoordHelperPtr->GetName(iCoord);
379       aCoordUnits[iCoord] = aCoordHelperPtr->GetUnit(iCoord);
380     }
381
382     const ERepere SMDS_COORDINATE_SYSTEM = eCART;
383
384     PNodeInfo aNodeInfo = myMed->CrNodeInfo(aMeshInfo,
385                                             aCoordinates,
386                                             eFULL_INTERLACE,
387                                             SMDS_COORDINATE_SYSTEM,
388                                             aCoordNames,
389                                             aCoordUnits,
390                                             aFamilyNums,
391                                             anElemNums);
392     MESSAGE("Perform - aNodeInfo->GetNbElem() = "<<aNbElems);
393     myMed->SetNodeInfo(aNodeInfo);
394
395
396     // Storing others SMDS elements to the MED file for the MED mesh
397     //--------------------------------------------------------------
398     EEntiteMaillage SMDS_MED_ENTITY = eMAILLE;
399     const EConnectivite SMDS_MED_CONNECTIVITY = eNOD;
400
401     // Storing SMDS Edges
402     if(TInt aNbElems = myMesh->NbEdges()){
403 #ifdef _ELEMENTS_BY_DIM_
404       SMDS_MED_ENTITY = eARETE;
405 #endif
406       // count edges of diff types
407       int aNbSeg3 = 0, aNbSeg2 = 0;
408       SMDS_EdgeIteratorPtr anIter = myMesh->edgesIterator();
409       while ( anIter->more() )
410         if ( anIter->next()->NbNodes() == 3 )
411           ++aNbSeg3;
412       aNbSeg2 = aNbElems - aNbSeg3;
413
414       TInt aNbSeg2Conn = MED::GetNbNodes(eSEG2);
415       MED::TIntVector aSeg2ElemNums, aSeg2FamilyNums, aSeg2Conn;
416       aSeg2ElemNums  .reserve( aNbSeg2 );
417       aSeg2FamilyNums.reserve( aNbSeg2 );
418       aSeg2Conn      .reserve( aNbSeg2*aNbSeg2Conn );
419
420       TInt aNbSeg3Conn = MED::GetNbNodes(eSEG3);
421       MED::TIntVector aSeg3ElemNums, aSeg3FamilyNums, aSeg3Conn;
422       aSeg3ElemNums  .reserve( aNbSeg3 );
423       aSeg3FamilyNums.reserve( aNbSeg3 );
424       aSeg3Conn      .reserve( aNbSeg3*aNbSeg3Conn );
425
426       anIter = myMesh->edgesIterator();
427       while ( anIter->more() ) {
428         const SMDS_MeshEdge* anElem = anIter->next();
429         TInt aNbNodes = anElem->NbNodes();
430
431         TInt aNbConnectivity;
432         MED::TIntVector* anElemNums;
433         MED::TIntVector* aFamilyNums;
434         MED::TIntVector* aConnectivity;
435         switch(aNbNodes){
436         case 2:
437           aNbConnectivity = aNbSeg2Conn;
438           anElemNums      = &aSeg2ElemNums;
439           aFamilyNums     = &aSeg2FamilyNums;
440           aConnectivity   = &aSeg2Conn;
441           break;
442         case 3:
443           aNbConnectivity = aNbSeg3Conn;
444           anElemNums      = &aSeg3ElemNums;
445           aFamilyNums     = &aSeg3FamilyNums;
446           aConnectivity   = &aSeg3Conn;
447           break;
448         default:
449           break;
450         }
451
452         for(TInt iNode = 0; iNode < aNbNodes; iNode++) {
453           const SMDS_MeshElement* aNode = anElem->GetNode( iNode );
454 #ifdef _EDF_NODE_IDS_
455           aConnectivity->push_back( aNodeIdMap[aNode->GetID()] );
456 #else
457           aConnectivity->push_back( aNode->GetID() );
458 #endif
459         }
460
461         anElemNums->push_back(anElem->GetID());
462
463         map<const SMDS_MeshElement*,int>::iterator edge_fam = anElemFamMap.find( anElem );
464         if ( edge_fam != anElemFamMap.end() )
465           aFamilyNums->push_back( edge_fam->second );
466         else
467           aFamilyNums->push_back( myEdgesDefaultFamilyId );
468       }
469       
470       if ( aNbSeg2 ) {
471         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
472                                                 SMDS_MED_ENTITY,
473                                                 eSEG2,
474                                                 aSeg2Conn,
475                                                 SMDS_MED_CONNECTIVITY,
476                                                 aSeg2FamilyNums,
477                                                 aSeg2ElemNums);
478         myMed->SetCellInfo(aCellInfo);
479       }
480       if ( aNbSeg3 ) {
481         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
482                                                 SMDS_MED_ENTITY,
483                                                 eSEG3,
484                                                 aSeg3Conn,
485                                                 SMDS_MED_CONNECTIVITY,
486                                                 aSeg3FamilyNums,
487                                                 aSeg3ElemNums);
488         myMed->SetCellInfo(aCellInfo);
489       }
490     }
491
492     // Storing SMDS Faces
493     if(TInt aNbElems = myMesh->NbFaces()){
494       SMDS_FaceIteratorPtr anIter = myMesh->facesIterator();
495 #ifdef _ELEMENTS_BY_DIM_
496       SMDS_MED_ENTITY = eFACE;
497 #endif
498       TInt aNbTriaConn = MED::GetNbNodes(eTRIA3);
499       MED::TIntVector anTriaElemNums; 
500       anTriaElemNums.reserve(aNbElems);
501       MED::TIntVector aTriaFamilyNums;
502       aTriaFamilyNums.reserve(aNbElems);
503       MED::TIntVector aTriaConn;
504       aTriaConn.reserve(aNbElems*aNbTriaConn);
505
506       TInt aNbTria6Conn = MED::GetNbNodes(eTRIA6);
507       MED::TIntVector anTria6ElemNums; 
508       anTria6ElemNums.reserve(aNbElems);
509       MED::TIntVector aTria6FamilyNums;
510       aTria6FamilyNums.reserve(aNbElems);
511       MED::TIntVector aTria6Conn;
512       aTria6Conn.reserve(aNbElems*aNbTria6Conn);
513
514       TInt aNbQuadConn = MED::GetNbNodes(eQUAD4);
515       MED::TIntVector aQuadElemNums;
516       aQuadElemNums.reserve(aNbElems);
517       MED::TIntVector aQuadFamilyNums;
518       aQuadFamilyNums.reserve(aNbElems);
519       MED::TIntVector aQuadConn;
520       aQuadConn.reserve(aNbElems*aNbQuadConn);
521
522       TInt aNbQuad8Conn = MED::GetNbNodes(eQUAD8);
523       MED::TIntVector aQuad8ElemNums;
524       aQuad8ElemNums.reserve(aNbElems);
525       MED::TIntVector aQuad8FamilyNums;
526       aQuad8FamilyNums.reserve(aNbElems);
527       MED::TIntVector aQuad8Conn;
528       aQuad8Conn.reserve(aNbElems*aNbQuad8Conn);
529
530       MED::TIntVector aPolygoneElemNums;
531       aPolygoneElemNums.reserve(aNbElems);
532       MED::TIntVector aPolygoneInds;
533       aPolygoneInds.reserve(aNbElems + 1);
534       aPolygoneInds.push_back(1); // reference on the first element in the connectivities
535       MED::TIntVector aPolygoneFamilyNums;
536       aPolygoneFamilyNums.reserve(aNbElems);
537       MED::TIntVector aPolygoneConn;
538       aPolygoneConn.reserve(aNbElems*aNbQuadConn);
539
540       for(TInt iElem = 0; iElem < aNbElems && anIter->more(); iElem++){
541         const SMDS_MeshFace* anElem = anIter->next();
542         TInt aNbNodes = anElem->NbNodes();
543         SMDS_ElemIteratorPtr aNodesIter = anElem->nodesIterator();
544         TInt aNbConnectivity;
545         MED::TIntVector* anElemNums;
546         MED::TIntVector* aFamilyNums;
547         MED::TIntVector* aConnectivity;
548         if (anElem->IsPoly()) {
549           aNbConnectivity = aNbNodes;
550           anElemNums = &aPolygoneElemNums;
551           aFamilyNums = &aPolygoneFamilyNums;
552           aConnectivity = &aPolygoneConn;
553         }
554         else {
555           switch(aNbNodes){
556           case 3:
557             aNbConnectivity = aNbTriaConn;
558             anElemNums = &anTriaElemNums;
559             aFamilyNums = &aTriaFamilyNums;
560             aConnectivity = &aTriaConn;
561             break;
562           case 4:
563             aNbConnectivity = aNbQuadConn;
564             anElemNums = &aQuadElemNums;
565             aFamilyNums = &aQuadFamilyNums;
566             aConnectivity = &aQuadConn;
567             break;
568           case 6:
569             aNbConnectivity = aNbTria6Conn;
570             anElemNums = &anTria6ElemNums;
571             aFamilyNums = &aTria6FamilyNums;
572             aConnectivity = &aTria6Conn;
573             break;
574           case 8:
575             aNbConnectivity = aNbQuad8Conn;
576             anElemNums = &aQuad8ElemNums;
577             aFamilyNums = &aQuad8FamilyNums;
578             aConnectivity = &aQuad8Conn;
579             break;
580           default:
581             break;
582           }
583         }
584         MED::TIntVector aVector(aNbNodes);
585         for(TInt iNode = 0; aNodesIter->more(); iNode++){
586           const SMDS_MeshElement* aNode = aNodesIter->next();
587 #ifdef _EDF_NODE_IDS_
588           aVector[iNode] = aNodeIdMap[aNode->GetID()];
589 #else
590           aVector[iNode] = aNode->GetID();
591 #endif
592         }
593
594         TInt aSize = aConnectivity->size();
595         aConnectivity->resize(aSize+aNbConnectivity);
596         // There is some differences between SMDS and MED in cells mapping
597         switch(aNbNodes){
598         case 4:
599           (*aConnectivity)[aSize+0] = aVector[0];
600           (*aConnectivity)[aSize+1] = aVector[1];
601           (*aConnectivity)[aSize+2] = aVector[3];  
602           (*aConnectivity)[aSize+3] = aVector[2];  
603         default:
604           for(TInt iNode = 0; iNode < aNbNodes; iNode++) 
605             (*aConnectivity)[aSize+iNode] = aVector[iNode];
606         }
607
608         if (anElem->IsPoly()) {
609           // fill indices for polygonal element
610           TInt aPrevPos = aPolygoneInds.back();
611           aPolygoneInds.push_back(aPrevPos + aNbNodes);
612         }
613
614         anElemNums->push_back(anElem->GetID());
615
616         if (anElemFamMap.find(anElem) != anElemFamMap.end())
617           aFamilyNums->push_back(anElemFamMap[anElem]);
618         else
619           aFamilyNums->push_back(myFacesDefaultFamilyId);
620       }
621       if(TInt aNbElems = anTriaElemNums.size()){
622         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
623                                                 SMDS_MED_ENTITY,
624                                                 eTRIA3,
625                                                 aTriaConn,
626                                                 SMDS_MED_CONNECTIVITY,
627                                                 aTriaFamilyNums,
628                                                 anTriaElemNums);
629         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eTRIA3<<"; aNbElems = "<<aNbElems);
630         myMed->SetCellInfo(aCellInfo);
631       }
632       if(TInt aNbElems = aQuadElemNums.size()){
633         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
634                                                 SMDS_MED_ENTITY,
635                                                 eQUAD4,
636                                                 aQuadConn,
637                                                 SMDS_MED_CONNECTIVITY,
638                                                 aQuadFamilyNums,
639                                                 aQuadElemNums);
640         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eQUAD4<<"; aNbElems = "<<aNbElems);
641         myMed->SetCellInfo(aCellInfo);
642       }
643       if(TInt aNbElems = anTria6ElemNums.size()){
644         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
645                                                 SMDS_MED_ENTITY,
646                                                 eTRIA6,
647                                                 aTria6Conn,
648                                                 SMDS_MED_CONNECTIVITY,
649                                                 aTria6FamilyNums,
650                                                 anTria6ElemNums);
651         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eTRIA6<<"; aNbElems = "<<aNbElems);
652         myMed->SetCellInfo(aCellInfo);
653       }
654       if(TInt aNbElems = aQuad8ElemNums.size()){
655         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
656                                                 SMDS_MED_ENTITY,
657                                                 eQUAD8,
658                                                 aQuad8Conn,
659                                                 SMDS_MED_CONNECTIVITY,
660                                                 aQuad8FamilyNums,
661                                                 aQuad8ElemNums);
662         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eQUAD8<<"; aNbElems = "<<aNbElems);
663         myMed->SetCellInfo(aCellInfo);
664       }
665       if(TInt aNbElems = aPolygoneElemNums.size()){
666         // add one element in connectivities,
667         // referenced by the last element in indices
668         aPolygoneConn.push_back(0);
669
670         PPolygoneInfo aCellInfo = myMed->CrPolygoneInfo(aMeshInfo,
671                                                         SMDS_MED_ENTITY,
672                                                         ePOLYGONE,
673                                                         aPolygoneInds,
674                                                         aPolygoneConn,
675                                                         SMDS_MED_CONNECTIVITY,
676                                                         aPolygoneFamilyNums,
677                                                         aPolygoneElemNums);
678         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePOLYGONE<<"; aNbElems = "<<aNbElems);
679         myMed->SetPolygoneInfo(aCellInfo);
680       }
681     }
682
683     // Storing SMDS Volumes
684     if(TInt aNbElems = myMesh->NbVolumes()){
685       SMDS_VolumeIteratorPtr anIter = myMesh->volumesIterator();
686 #ifdef _ELEMENTS_BY_DIM_
687       SMDS_MED_ENTITY = eMAILLE;
688 #endif
689       TInt aNbTetraConn = MED::GetNbNodes(eTETRA4);
690       MED::TIntVector anTetraElemNums; 
691       anTetraElemNums.reserve(aNbElems);
692       MED::TIntVector aTetraFamilyNums;
693       aTetraFamilyNums.reserve(aNbElems);
694       MED::TIntVector aTetraConn;
695       aTetraConn.reserve(aNbElems*aNbTetraConn);
696
697       TInt aNbPyraConn = MED::GetNbNodes(ePYRA5);
698       MED::TIntVector anPyraElemNums; 
699       anPyraElemNums.reserve(aNbElems);
700       MED::TIntVector aPyraFamilyNums;
701       aPyraFamilyNums.reserve(aNbElems);
702       MED::TIntVector aPyraConn;
703       aPyraConn.reserve(aNbElems*aNbPyraConn);
704
705       TInt aNbPentaConn = MED::GetNbNodes(ePENTA6);
706       MED::TIntVector anPentaElemNums; 
707       anPentaElemNums.reserve(aNbElems);
708       MED::TIntVector aPentaFamilyNums;
709       aPentaFamilyNums.reserve(aNbElems);
710       MED::TIntVector aPentaConn;
711       aPentaConn.reserve(aNbElems*aNbPentaConn);
712
713       TInt aNbHexaConn = MED::GetNbNodes(eHEXA8);
714       MED::TIntVector aHexaElemNums;
715       aHexaElemNums.reserve(aNbElems);
716       MED::TIntVector aHexaFamilyNums;
717       aHexaFamilyNums.reserve(aNbElems);
718       MED::TIntVector aHexaConn;
719       aHexaConn.reserve(aNbElems*aNbHexaConn);
720
721       TInt aNbTetra10Conn = MED::GetNbNodes(eTETRA10);
722       MED::TIntVector anTetra10ElemNums; 
723       anTetra10ElemNums.reserve(aNbElems);
724       MED::TIntVector aTetra10FamilyNums;
725       aTetra10FamilyNums.reserve(aNbElems);
726       MED::TIntVector aTetra10Conn;
727       aTetra10Conn.reserve(aNbElems*aNbTetra10Conn);
728
729       TInt aNbPyra13Conn = MED::GetNbNodes(ePYRA13);
730       MED::TIntVector anPyra13ElemNums; 
731       anPyra13ElemNums.reserve(aNbElems);
732       MED::TIntVector aPyra13FamilyNums;
733       aPyra13FamilyNums.reserve(aNbElems);
734       MED::TIntVector aPyra13Conn;
735       aPyra13Conn.reserve(aNbElems*aNbPyra13Conn);
736
737       TInt aNbPenta15Conn = MED::GetNbNodes(ePENTA15);
738       MED::TIntVector anPenta15ElemNums; 
739       anPenta15ElemNums.reserve(aNbElems);
740       MED::TIntVector aPenta15FamilyNums;
741       aPenta15FamilyNums.reserve(aNbElems);
742       MED::TIntVector aPenta15Conn;
743       aPenta15Conn.reserve(aNbElems*aNbPenta15Conn);
744
745       TInt aNbHexa20Conn = MED::GetNbNodes(eHEXA20);
746       MED::TIntVector aHexa20ElemNums;
747       aHexa20ElemNums.reserve(aNbElems);
748       MED::TIntVector aHexa20FamilyNums;
749       aHexa20FamilyNums.reserve(aNbElems);
750       MED::TIntVector aHexa20Conn;
751       aHexa20Conn.reserve(aNbElems*aNbHexa20Conn);
752
753       MED::TIntVector aPolyedreElemNums;
754       aPolyedreElemNums.reserve(aNbElems);
755       MED::TIntVector aPolyedreInds;
756       aPolyedreInds.reserve(aNbElems + 1);
757       aPolyedreInds.push_back(1); // reference on the first element in the faces
758       MED::TIntVector aPolyedreFaces;
759       aPolyedreFaces.reserve(aNbElems + 1);
760       aPolyedreFaces.push_back(1); // reference on the first element in the connectivities
761       MED::TIntVector aPolyedreFamilyNums;
762       aPolyedreFamilyNums.reserve(aNbElems);
763       MED::TIntVector aPolyedreConn;
764       aPolyedreConn.reserve(aNbElems*aNbHexaConn);
765
766       for(TInt iElem = 0; iElem < aNbElems && anIter->more(); iElem++){
767         const SMDS_MeshVolume* anElem = anIter->next();
768
769         MED::TIntVector* anElemNums;
770         MED::TIntVector* aFamilyNums;
771
772         if (anElem->IsPoly()) {
773           const SMDS_PolyhedralVolumeOfNodes* aPolyedre =
774             (const SMDS_PolyhedralVolumeOfNodes*) anElem;
775           if (!aPolyedre) {
776             MESSAGE("Warning: bad volumic element");
777             continue;
778           }
779
780           anElemNums = &aPolyedreElemNums;
781           aFamilyNums = &aPolyedreFamilyNums;
782
783           TInt aNodeId, aNbFaces = aPolyedre->NbFaces();
784           for (int iface = 1; iface <= aNbFaces; iface++) {
785             int aNbFaceNodes = aPolyedre->NbFaceNodes(iface);
786             for (int inode = 1; inode <= aNbFaceNodes; inode++) {
787               aNodeId = aPolyedre->GetFaceNode(iface, inode)->GetID();
788 #ifdef _EDF_NODE_IDS_
789               aPolyedreConn.push_back(aNodeIdMap[aNodeId]);
790 #else
791               aPolyedreConn.push_back(aNodeId);
792 #endif
793             }
794             TInt aPrevPos = aPolyedreFaces.back();
795             aPolyedreFaces.push_back(aPrevPos + aNbFaceNodes);
796           }
797           TInt aPrevPos = aPolyedreInds.back();
798           aPolyedreInds.push_back(aPrevPos + aNbFaces);
799
800         }
801         else {
802           TInt aNbNodes = anElem->NbNodes();
803           SMDS_ElemIteratorPtr aNodesIter = anElem->nodesIterator();
804           TInt aNbConnectivity;
805           MED::TIntVector* aConnectivity;
806           switch(aNbNodes){
807           case 4:
808             aNbConnectivity = aNbTetraConn;
809             anElemNums = &anTetraElemNums;
810             aFamilyNums = &aTetraFamilyNums;
811             aConnectivity = &aTetraConn;
812             break;
813           case 5:
814             aNbConnectivity = aNbPyraConn;
815             anElemNums = &anPyraElemNums;
816             aFamilyNums = &aPyraFamilyNums;
817             aConnectivity = &aPyraConn;
818             break;
819           case 6:
820             aNbConnectivity = aNbPentaConn;
821             anElemNums = &anPentaElemNums;
822             aFamilyNums = &aPentaFamilyNums;
823             aConnectivity = &aPentaConn;
824             break;
825           case 8:
826             aNbConnectivity = aNbHexaConn;
827             anElemNums = &aHexaElemNums;
828             aFamilyNums = &aHexaFamilyNums;
829             aConnectivity = &aHexaConn;
830             break;
831           case 10:
832             aNbConnectivity = aNbTetra10Conn;
833             anElemNums = &anTetra10ElemNums;
834             aFamilyNums = &aTetra10FamilyNums;
835             aConnectivity = &aTetra10Conn;
836             break;
837           case 13:
838             aNbConnectivity = aNbPyra13Conn;
839             anElemNums = &anPyra13ElemNums;
840             aFamilyNums = &aPyra13FamilyNums;
841             aConnectivity = &aPyra13Conn;
842             break;
843           case 15:
844             aNbConnectivity = aNbPenta15Conn;
845             anElemNums = &anPenta15ElemNums;
846             aFamilyNums = &aPenta15FamilyNums;
847             aConnectivity = &aPenta15Conn;
848             break;
849           case 20:
850             aNbConnectivity = aNbHexa20Conn;
851             anElemNums = &aHexa20ElemNums;
852             aFamilyNums = &aHexa20FamilyNums;
853             aConnectivity = &aHexa20Conn;
854           }
855
856           TInt aSize = aConnectivity->size();
857           aConnectivity->resize(aSize + aNbConnectivity);
858
859           MED::TIntVector aVector(aNbNodes);
860           for(TInt iNode = 0; aNodesIter->more(); iNode++){
861             const SMDS_MeshElement* aNode = aNodesIter->next();
862 #ifdef _EDF_NODE_IDS_
863             aVector[iNode] = aNodeIdMap[aNode->GetID()];
864 #else
865             aVector[iNode] = aNode->GetID();
866 #endif
867           }
868           // There is some difference between SMDS and MED in cells mapping
869           switch(aNbNodes){
870           case 5:
871             (*aConnectivity)[aSize+0] = aVector[0];
872             (*aConnectivity)[aSize+1] = aVector[3];
873             (*aConnectivity)[aSize+2] = aVector[2];  
874             (*aConnectivity)[aSize+3] = aVector[1];  
875             (*aConnectivity)[aSize+4] = aVector[4];  
876           default:
877             for(TInt iNode = 0; iNode < aNbNodes; iNode++) 
878               (*aConnectivity)[aSize+iNode] = aVector[iNode];
879           }
880         }
881
882         anElemNums->push_back(anElem->GetID());
883
884         if (anElemFamMap.find(anElem) != anElemFamMap.end())
885           aFamilyNums->push_back(anElemFamMap[anElem]);
886         else
887           aFamilyNums->push_back(myVolumesDefaultFamilyId);
888       }
889
890       if(TInt aNbElems = anTetraElemNums.size()){
891         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
892                                                 SMDS_MED_ENTITY,
893                                                 eTETRA4,
894                                                 aTetraConn,
895                                                 SMDS_MED_CONNECTIVITY,
896                                                 aTetraFamilyNums,
897                                                 anTetraElemNums);
898         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eTETRA4<<"; aNbElems = "<<aNbElems);
899         myMed->SetCellInfo(aCellInfo);
900       }
901       if(TInt aNbElems = anPyraElemNums.size()){
902         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
903                                                 SMDS_MED_ENTITY,
904                                                 ePYRA5,
905                                                 aPyraConn,
906                                                 SMDS_MED_CONNECTIVITY,
907                                                 aPyraFamilyNums,
908                                                 anPyraElemNums);
909         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePYRA5<<"; aNbElems = "<<aNbElems);
910         myMed->SetCellInfo(aCellInfo);
911       }
912       if(TInt aNbElems = anPentaElemNums.size()){
913         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
914                                                 SMDS_MED_ENTITY,
915                                                 ePENTA6,
916                                                 aPentaConn,
917                                                 SMDS_MED_CONNECTIVITY,
918                                                 aPentaFamilyNums,
919                                                 anPentaElemNums);
920         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePENTA6<<"; aNbElems = "<<aNbElems);
921         myMed->SetCellInfo(aCellInfo);
922       }
923       if(TInt aNbElems = aHexaElemNums.size()){
924         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
925                                                 SMDS_MED_ENTITY,
926                                                 eHEXA8,
927                                                 aHexaConn,
928                                                 SMDS_MED_CONNECTIVITY,
929                                                 aHexaFamilyNums,
930                                                 aHexaElemNums);
931         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eHEXA8<<"; aNbElems = "<<aNbElems);
932         myMed->SetCellInfo(aCellInfo);
933       }
934       if(TInt aNbElems = anTetra10ElemNums.size()){
935         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
936                                                 SMDS_MED_ENTITY,
937                                                 eTETRA10,
938                                                 aTetra10Conn,
939                                                 SMDS_MED_CONNECTIVITY,
940                                                 aTetra10FamilyNums,
941                                                 anTetra10ElemNums);
942         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eTETRA10<<"; aNbElems = "<<aNbElems);
943         myMed->SetCellInfo(aCellInfo);
944       }
945       if(TInt aNbElems = anPyra13ElemNums.size()){
946         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
947                                                 SMDS_MED_ENTITY,
948                                                 ePYRA13,
949                                                 aPyra13Conn,
950                                                 SMDS_MED_CONNECTIVITY,
951                                                 aPyra13FamilyNums,
952                                                 anPyra13ElemNums);
953         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePYRA13<<"; aNbElems = "<<aNbElems);
954         myMed->SetCellInfo(aCellInfo);
955       }
956       if(TInt aNbElems = anPenta15ElemNums.size()){
957         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
958                                                 SMDS_MED_ENTITY,
959                                                 ePENTA15,
960                                                 aPenta15Conn,
961                                                 SMDS_MED_CONNECTIVITY,
962                                                 aPenta15FamilyNums,
963                                                 anPenta15ElemNums);
964         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePENTA15<<"; aNbElems = "<<aNbElems);
965         myMed->SetCellInfo(aCellInfo);
966       }
967       if(TInt aNbElems = aHexa20ElemNums.size()){
968         PCellInfo aCellInfo = myMed->CrCellInfo(aMeshInfo,
969                                                 SMDS_MED_ENTITY,
970                                                 eHEXA20,
971                                                 aHexa20Conn,
972                                                 SMDS_MED_CONNECTIVITY,
973                                                 aHexa20FamilyNums,
974                                                 aHexa20ElemNums);
975         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<eHEXA20<<"; aNbElems = "<<aNbElems);
976         myMed->SetCellInfo(aCellInfo);
977       }
978
979       if(TInt aNbElems = aPolyedreElemNums.size()){
980         // add one element in connectivities,
981         // referenced by the last element in faces
982         aPolyedreConn.push_back(0);
983
984         PPolyedreInfo aCellInfo = myMed->CrPolyedreInfo(aMeshInfo,
985                                                         SMDS_MED_ENTITY,
986                                                         ePOLYEDRE,
987                                                         aPolyedreInds,
988                                                         aPolyedreFaces,
989                                                         aPolyedreConn,
990                                                         SMDS_MED_CONNECTIVITY,
991                                                         aPolyedreFamilyNums,
992                                                         aPolyedreElemNums);
993         MESSAGE("Perform - anEntity = "<<SMDS_MED_ENTITY<<"; aGeom = "<<ePOLYEDRE<<"; aNbElems = "<<aNbElems);
994         myMed->SetPolyedreInfo(aCellInfo);
995       }
996     }
997   }
998   catch(const std::exception& exc) {
999     INFOS("Follow exception was cought:\n\t"<<exc.what());
1000     throw;
1001   }
1002   catch(...) {
1003     INFOS("Unknown exception was cought !!!");
1004     throw;
1005   }
1006
1007   myMeshId = -1;
1008   myGroups.clear();
1009   mySubMeshes.clear();
1010   return aResult;
1011 }