src/MEDMEMCppTest/MEDMEMTest_Grid_fault.cxx

   1 //  Copyright (C) 2007-2008  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 #include "MEDMEMTest.hxx"
  23 #include <cppunit/TestAssert.h>
  24
  25 #include "MEDMEM_define.hxx"
  26 #include "MEDMEM_Grid.hxx"
  27 #include "MEDMEM_Mesh.hxx"
  28 #include "MEDMEM_Med.hxx"
  29 #include "MEDMEM_MedMedDriver.hxx"
  30 #include "MEDMEM_MedMeshDriver.hxx"
  31
  32 #include <sstream>
  33 #include <cmath>
  34
  35 // use this define to enable lines, execution of which leads to Segmentation Fault
  36 //#define ENABLE_FAULTS
  37
  38 // use this define to enable CPPUNIT asserts and fails, showing bugs
  39 //#define ENABLE_FORCED_FAILURES
  40
  41 using namespace std;
  42 using namespace MEDMEM;
  43 using namespace MED_EN;
  44
  45
  46 /*!
  47  *  Check methods (44), defined in MEDMEM_Grid.hxx:
  48  *  class GRID: public MESH {
  49  *   (+) GRID();
  50  *   (+) GRID(const MED_EN::med_grid_type type);
  51  *   (BUG:operator=() not implemented but init() not called) GRID(const GRID &m);
  52  *   (+) GRID(driverTypes driverType, const string & fileName="",const string & meshName="");
  53  *   (+) GRID(const std::vector<std::vector<double> >& xyz_array,
  54  *                const std::vector<std::string>& coord_name,
  55  *                const std::vector<std::string>& coord_unit,
  56  *                const MED_EN::med_grid_type type=MED_EN::MED_CARTESIAN);
  57  *   (NOT IMPLEMENTED) GRID & operator=(const GRID &m);
  58  *   (+) virtual ~GRID();
  59  *   (+) virtual void init();
  60  *   (tested together with getCoordinateptr() as it is called
  61  *    internally from there first of all.
  62  *    Moreover, fillCoordinates should be made private to avoid
  63  *    ambiguity druing in GRID class usage.) void fillCoordinates() const;
  64  *   (tested together with getConnectivityptr()) void fillConnectivity() const;
  65  *   (+) inline void makeUnstructured();//fill coordinates and connectivity of MESH
  66  *   (+) void fillMeshAfterRead();
  67  *   (+) void writeUnstructured(int index=0, const string & driverName = "");
  68  *   (+) void read(int index=0);
  69  *   (+) inline int getNodeNumber(const int i, const int j=0, const int k=0) const;
  70  *   (+) inline int getCellNumber(const int i, const int j=0, const int k=0) const;
  71  *   (+) int getEdgeNumber
  72  *           (const int Axis, const int i, const int j=0, const int k=0) const throw (MEDEXCEPTION);
  73  *   (+) int getFaceNumber
  74  *           (const int Axis, const int i, const int j=0, const int k=0) const throw (MEDEXCEPTION);
  75  *   (+) void getNodePosition(const int Number, int& i, int& j, int& k) const throw (MEDEXCEPTION);
  76  *   (+) void getCellPosition(const int Number, int& i, int& j, int& k) const throw (MEDEXCEPTION);
  77  *   (+) void getEdgePosition
  78  *           (const int Number, int& Axis, int& i, int& j, int& k) const throw (MEDEXCEPTION);
  79  *   (+) void getFacePosition
  80  *           (const int Number, int& Axis, int& i, int& j, int& k) const throw (MEDEXCEPTION);
  81  *   (+) inline MED_EN::med_grid_type getGridType() const;
  82  *   (+) int getArrayLength(const int Axis) const throw (MEDEXCEPTION);
  83  *   (+) const double getArrayValue (const int Axis, const int i) const throw (MEDEXCEPTION);
  84  *   (+) inline const COORDINATE * getCoordinateptr() const;
  85  *   (+) inline const double * getCoordinates(MED_EN::medModeSwitch Mode) const;
  86  *   (+) inline const double getCoordinate(int Number,int Axis) const;
  87  *   (+) inline int getNumberOfTypes(MED_EN::medEntityMesh Entity) const;
  88  *   (+) inline int getNumberOfTypesWithPoly(MED_EN::medEntityMesh Entity) const;
  89  *   (+) inline const MED_EN::medGeometryElement * getTypes(MED_EN::medEntityMesh Entity) const;
  90  *   (+) MED_EN::medGeometryElement * getTypesWithPoly(MED_EN::medEntityMesh Entity) const;
  91  *   (+) inline const CELLMODEL * getCellsTypes(MED_EN::medEntityMesh Entity) const;
  92  *   (+) const int * getGlobalNumberingIndex(MED_EN::medEntityMesh Entity) const;
  93  *   (+) inline int getNumberOfElements
  94  *                (MED_EN::medEntityMesh Entity, MED_EN::medGeometryElement Type) const;
  95  *   (+) inline int getNumberOfElementsWithPoly
  96  *                (MED_EN::medEntityMesh Entity, MED_EN::medGeometryElement Type) const;
  97  *   (+) inline bool existConnectivity
  98  *                (MED_EN::medConnectivity ConnectivityType, MED_EN::medEntityMesh Entity) const;
  99  *   (+) inline MED_EN::medGeometryElement getElementType
 100  *                (MED_EN::medEntityMesh Entity, int Number) const;
 101  *   (+) inline MED_EN::medGeometryElement getElementTypeWithPoly
 102  *                (MED_EN::medEntityMesh Entity, int Number) const;
 103  *   (+) inline void calculateConnectivity(MED_EN::medModeSwitch Mode,
 104  *                                             MED_EN::medConnectivity ConnectivityType,
 105  *                                               MED_EN::medEntityMesh Entity) const;
 106  *   (+) inline const CONNECTIVITY* getConnectivityptr() const;
 107  *   (+) inline const int * getConnectivity
 108  *             (MED_EN::medModeSwitch Mode, MED_EN::medConnectivity ConnectivityType,
 109  *                MED_EN::medEntityMesh Entity, MED_EN::medGeometryElement Type) const;
 110  *   (+) inline const int * getConnectivityIndex(MED_EN::medConnectivity ConnectivityType,
 111  *                                          MED_EN::medEntityMesh Entity) const;
 112  *   (+) inline const int * getReverseConnectivity(MED_EN::medConnectivity ConnectivityType,
 113  *                                            MED_EN::medEntityMesh Entity=MED_EN::MED_CELL) const;
 114  *   (+) inline const int * getReverseConnectivityIndex(MED_EN::medConnectivity ConnectivityType,
 115  *                                                 MED_EN::medEntityMesh Entity=MED_EN::MED_CELL) const;
 116  *   (+) inline void setGridType(MED_EN::med_grid_type gridType);
 117  *  }
 118  */
 119 void MEDMEMTest_testGrid()
 120 {
 121   string datadir  = getenv("MED_ROOT_DIR");
 122   string filename = datadir + "/share/salome/resources/med/test19.med" ;
 123   string tmp_dir  = getenv("TMP") ? getenv("TMP") : "/tmp";
 124   if (tmp_dir == "")
 125     tmp_dir = "/tmp";
 126   string filenameout21 = tmp_dir + "/myGridWrite_grid21.med";
 127
 128   // To remove tmp files from disk
 129   MEDMEMTest_TmpFilesRemover aRemover;
 130   aRemover.Register(filenameout21);
 131
 132   // Remove file in advance to ensure it does not exist at the moment of writing,
 133   // because files are not removed by the MEDMEMTest_TmpFilesRemover in case of
 134   // Segmentation Fault (in previous tests execution).
 135   {
 136     MEDMEMTest_TmpFilesRemover aRemover1;
 137     aRemover1.Register(filenameout21);
 138   }
 139
 140   MED * myMed = new MED() ;
 141   MED_MED_RDONLY_DRIVER myMeshDriver (filename, myMed);
 142   myMeshDriver.open();
 143   myMeshDriver.readFileStruct();
 144   myMeshDriver.close();
 145
 146   int nbMeshes = myMed->getNumberOfMeshes();
 147   CPPUNIT_ASSERT(nbMeshes);
 148
 149   deque<string> mesh_names = myMed->getMeshNames();
 150   CPPUNIT_ASSERT(mesh_names.size() != 0);
 151
 152   //////////////////////////////
 153   // test1 "CARTESIAN GRID"   //
 154   //////////////////////////////
 155   {
 156     MESH* myMesh = myMed->getMesh(mesh_names[0]);
 157     myMesh->read();
 158
 159     CPPUNIT_ASSERT(myMesh != NULL);
 160     CPPUNIT_ASSERT(myMesh->getIsAGrid());
 161
 162     GRID* myGrid = dynamic_cast<GRID*>(myMesh);
 163     CPPUNIT_ASSERT(myGrid);
 164
 165     CPPUNIT_ASSERT_THROW(myGrid->getArrayLength(0), MEDEXCEPTION);
 166
 167     int I, J, K;
 168     CPPUNIT_ASSERT_NO_THROW(I = myGrid->getArrayLength(1));
 169     CPPUNIT_ASSERT_NO_THROW(J = myGrid->getArrayLength(2));
 170     CPPUNIT_ASSERT_NO_THROW(K = myGrid->getArrayLength(3));
 171
 172     CPPUNIT_ASSERT(I);
 173     CPPUNIT_ASSERT(J);
 174
 175     med_grid_type grid_type = myGrid->getGridType();
 176     CPPUNIT_ASSERT_MESSAGE("Wrong grid type", grid_type == MED_CARTESIAN);
 177
 178     const double * coordinates = myGrid->getCoordinates(MED_FULL_INTERLACE);
 179     int SpaceDimension = myGrid->getSpaceDimension();
 180     for (int axe = 0; axe < SpaceDimension; axe++) {
 181       for (int num = 0; num < myGrid->getNumberOfNodes(); num++) {
 182         double coordinate;
 183         CPPUNIT_ASSERT_NO_THROW(coordinate = myGrid->getCoordinate(num + 1, axe + 1));
 184         //cout << "coordinate = " << coordinate << endl;
 185         CPPUNIT_ASSERT(fabs(coordinate - coordinates[(num * SpaceDimension)+axe]) < 0.001);
 186       }
 187     }
 188
 189     int nbTypesCell = myGrid->getNumberOfTypes(MED_CELL);
 190     CPPUNIT_ASSERT(nbTypesCell == 1);
 191
 192     const medGeometryElement* types;
 193     CPPUNIT_ASSERT_NO_THROW(types = myGrid->getTypes(MED_CELL));
 194 //#ifdef ENABLE_FORCED_FAILURES
 195     // Compilation warning about GRID::getTypes():
 196     // "inline function
 197     // `virtual const MED_EN::medGeometryElement* MEDMEM::GRID::getTypes(MED_EN::medEntityMesh) const'
 198     // used but never defined".
 199     // In MEDMEM_Grid.hxx:
 200     // inline const MED_EN::medGeometryElement * getTypes(MED_EN::medEntityMesh Entity) const;
 201     // But implemented in MEDMEM_Grid.cxx:
 202     //        const MED_EN::medGeometryElement * GRID::getTypes(MED_EN::medEntityMesh entity) const
 203 //    CPPUNIT_FAIL("Problem with GRID::getTypes() method implementation.");
 204 //#endif
 205     CPPUNIT_ASSERT(types[0] == MED_QUAD4);
 206
 207     int nbElem = 0;
 208     CPPUNIT_ASSERT_NO_THROW(nbElem = myGrid->getNumberOfElements(MED_CELL,types[0]));
 209     CPPUNIT_ASSERT(nbElem);
 210
 211     int nbNodes = myGrid->getNumberOfNodes();
 212     CPPUNIT_ASSERT(nbNodes);
 213
 214     int ijkNode[3];
 215     int NodeNumber;
 216     for (int nbNode = 1; nbNode <= nbNodes; nbNode++) {
 217       CPPUNIT_ASSERT_NO_THROW(myGrid->getNodePosition(nbNode, ijkNode[0], ijkNode[1], ijkNode[2]));
 218       CPPUNIT_ASSERT_NO_THROW(NodeNumber = myGrid->getNodeNumber(ijkNode[0], ijkNode[1], ijkNode[2]));
 219       CPPUNIT_ASSERT(NodeNumber == nbNode);
 220     }
 221
 222     int nbCells = myGrid->getNumberOfElements(MED_CELL, MED_ALL_ELEMENTS);
 223     CPPUNIT_ASSERT(nbCells);
 224
 225     int ijkCell[3];
 226     int CellNumber;
 227     for (int nbCell = 1; nbCell <= nbCells; nbCell++) {
 228       CPPUNIT_ASSERT_NO_THROW(myGrid->getCellPosition(nbCell, ijkCell[0], ijkCell[1], ijkCell[2]));
 229       CPPUNIT_ASSERT_NO_THROW(CellNumber = myGrid->getCellNumber(ijkCell[0], ijkCell[1], ijkCell[2]));
 230       CPPUNIT_ASSERT(CellNumber == nbCell);
 231     }
 232
 233     int nbEdges = myGrid->getNumberOfElements(MED_EDGE, MED_ALL_ELEMENTS);
 234     CPPUNIT_ASSERT(nbEdges);
 235
 236     int ijkAxisEdge[4];
 237     int EdgeNumber;
 238     for (int nbEdge = 1; nbEdge <= nbEdges; nbEdge++) {
 239       CPPUNIT_ASSERT_NO_THROW(myGrid->getEdgePosition(nbEdge, ijkAxisEdge[0], ijkAxisEdge[1],
 240                                                       ijkAxisEdge[2], ijkAxisEdge[3]));
 241       CPPUNIT_ASSERT_NO_THROW(EdgeNumber = myGrid->getEdgeNumber(ijkAxisEdge[0], ijkAxisEdge[1],
 242                                                                  ijkAxisEdge[2], ijkAxisEdge[3]));
 243       CPPUNIT_ASSERT(EdgeNumber == nbEdge);
 244     }
 245
 246     int nbFaces = myGrid->getNumberOfElements(MED_FACE, MED_ALL_ELEMENTS);
 247     CPPUNIT_ASSERT(nbFaces == 0);
 248
 249 //#ifdef ENABLE_FORCED_FAILURES
 250     CPPUNIT_FAIL("ERROR: nbFaces == 0, but getFacePosition(AnyNumber, ...) - return position(i,j,k)");
 251 //#endif
 252     int ijkAxisFace[4];
 253     CPPUNIT_ASSERT_NO_THROW(myGrid->getFacePosition(6, ijkAxisFace[0], ijkAxisFace[1],
 254                                                     ijkAxisFace[2], ijkAxisFace[3]));
 255     CPPUNIT_ASSERT(ijkAxisFace[0]);
 256     CPPUNIT_ASSERT(ijkAxisFace[1]);
 257     CPPUNIT_ASSERT(ijkAxisFace[2]);
 258
 259     /*int FaceNumber;
 260     for(int nbFace = 1; nbFace <= nbFaces; nbFace++)
 261     {
 262       CPPUNIT_ASSERT_NO_THROW( myGrid->getFacePosition(nbFace, ijkAxisFace[0], ijkAxisFace[1],
 263                                                        ijkAxisFace[2], ijkAxisFace[3]));
 264       CPPUNIT_ASSERT_NO_THROW( EdgeNumber = myGrid->getEdgeNumber(ijkAxisFace[0], ijkAxisFace[1],
 265                                                                   ijkAxisFace[2], ijkAxisFace[3]));
 266       CPPUNIT_ASSERT(FaceNumber == nbFace);
 267     }*/
 268
 269     bool existConnect = false;
 270     CPPUNIT_ASSERT_NO_THROW(existConnect = myGrid->existConnectivity(MED_NODAL, MED_CELL));
 271     if (!existConnect) {
 272       CPPUNIT_ASSERT_NO_THROW(myGrid->calculateConnectivity(MED_FULL_INTERLACE, MED_NODAL, MED_CELL));
 273       CPPUNIT_ASSERT(myGrid->existConnectivity(MED_NODAL, MED_CELL));
 274     }
 275
 276     const int* Connectivity;
 277     const int* connectivity_index;
 278     CPPUNIT_ASSERT_NO_THROW(Connectivity = myGrid->getConnectivity(MED_FULL_INTERLACE, MED_NODAL,
 279                                                                    MED_CELL, types[0]));
 280     CPPUNIT_ASSERT_NO_THROW(connectivity_index = myGrid->getConnectivityIndex(MED_NODAL, MED_CELL));
 281     cout << "Nodal connectivity" << endl;
 282     for (int j = 0; j < nbElem; j++) {
 283       cout << "Element "<< j+1 << " : ";
 284       for (int k = connectivity_index[j]; k < connectivity_index[j+1]; k++)
 285         cout << Connectivity[k-1] << " ";
 286       cout << endl;
 287     }
 288
 289     const int * ReverseNodalConnectivity;
 290     const int * ReverseConnectivityIndex;
 291     CPPUNIT_ASSERT_NO_THROW(ReverseNodalConnectivity = myGrid->getReverseConnectivity(MED_NODAL));
 292     CPPUNIT_ASSERT_NO_THROW(ReverseConnectivityIndex = myGrid->getReverseConnectivityIndex(MED_NODAL));
 293     for (int i = 0; i < nbNodes; i++) {
 294       cout << "Node "<< i+1 << " : ";
 295       for (int j = ReverseConnectivityIndex[i]; j < ReverseConnectivityIndex[i+1]; j++)
 296         cout << ReverseNodalConnectivity[j-1] << " ";
 297       cout << endl;
 298     }
 299
 300     const int* myGlobalNbIdx;
 301     CPPUNIT_ASSERT_NO_THROW(myGlobalNbIdx = myGrid->getGlobalNumberingIndex(MED_CELL));
 302     for (int i = 0; i <= nbTypesCell; i++) {
 303       if (i == nbTypesCell) {
 304         CPPUNIT_ASSERT_THROW(myGrid->getElementType(MED_CELL, myGlobalNbIdx[i]), MEDEXCEPTION);
 305         break;
 306       }
 307       medGeometryElement aElem, geomPolyElem;
 308       CPPUNIT_ASSERT_NO_THROW(aElem = myGrid->getElementType(MED_CELL, myGlobalNbIdx[i]));
 309       CPPUNIT_ASSERT_NO_THROW(geomPolyElem = myGrid->getElementTypeWithPoly(MED_CELL, myGlobalNbIdx[i]));
 310       CPPUNIT_ASSERT(types[0] == aElem);
 311       CPPUNIT_ASSERT(geomPolyElem == aElem);
 312     }
 313
 314     CPPUNIT_ASSERT_NO_THROW(existConnect = myGrid->existConnectivity(MED_DESCENDING, MED_CELL));
 315     if (!existConnect) {
 316       CPPUNIT_ASSERT_NO_THROW(myGrid->calculateConnectivity(MED_FULL_INTERLACE, MED_DESCENDING, MED_CELL));
 317       CPPUNIT_ASSERT(myGrid->existConnectivity(MED_DESCENDING, MED_CELL));
 318     }
 319
 320     const int* ConnectivityDes;
 321     const int* connectivity_index_des;
 322     CPPUNIT_ASSERT_NO_THROW(ConnectivityDes = myGrid->getConnectivity(MED_FULL_INTERLACE, MED_DESCENDING,
 323                                                                       MED_CELL, MED_ALL_ELEMENTS));
 324     CPPUNIT_ASSERT_NO_THROW(connectivity_index_des =
 325                             myGrid->getConnectivityIndex(MED_DESCENDING, MED_CELL));
 326     cout<<"Descending connectivity"<<endl;
 327     for (int j = 0; j < nbElem; j++) {
 328       cout << "Element "<< j+1 << " : ";
 329       for (int k = connectivity_index_des[j]; k < connectivity_index_des[j+1]; k++)
 330         cout << ConnectivityDes[k-1] << " ";
 331       cout << endl;
 332     }
 333
 334     const int * ReverseDesConnectivity;
 335     const int * ReverseConnectivityIndexDes;
 336     CPPUNIT_ASSERT_NO_THROW(ReverseDesConnectivity = myGrid->getReverseConnectivity(MED_DESCENDING));
 337     CPPUNIT_ASSERT_NO_THROW(ReverseConnectivityIndexDes =
 338                             myGrid->getReverseConnectivityIndex(MED_DESCENDING));
 339     for (int i = 0; i < nbNodes; i++) {
 340       cout << "Node "<< i+1 << " : ";
 341       for (int j = ReverseConnectivityIndexDes[i]; j < ReverseConnectivityIndexDes[i+1]; j++)
 342         cout << ReverseDesConnectivity[j-1] << " ";
 343       cout << endl;
 344     }
 345
 346     //TEST POLY
 347     {
 348       int nbPolytypes;
 349       //test getNumberOfTypesWithPoly() - a grid has one type
 350       CPPUNIT_ASSERT_NO_THROW(nbPolytypes = myGrid->getNumberOfTypesWithPoly(MED_CELL));
 351       CPPUNIT_ASSERT(nbPolytypes == 1 );
 352
 353       const MED_EN::medGeometryElement * PolyTypes, *Types;
 354       CPPUNIT_ASSERT_NO_THROW(PolyTypes = myGrid->getTypesWithPoly(MED_CELL));
 355       CPPUNIT_ASSERT_NO_THROW(Types = myGrid->getTypes(MED_CELL));
 356       CPPUNIT_ASSERT_EQUAL(PolyTypes[nbPolytypes-1],Types[nbPolytypes-1]);
 357
 358       for (int t = 0; t < nbPolytypes; t++) {
 359         int nbElPoly, nbEl;
 360         CPPUNIT_ASSERT_NO_THROW(nbElPoly = myGrid->getNumberOfElementsWithPoly(MED_CELL, PolyTypes[t]));
 361         CPPUNIT_ASSERT_NO_THROW(nbEl = myGrid->getNumberOfElements(MED_CELL, PolyTypes[t]));
 362         CPPUNIT_ASSERT(nbElPoly == nbEl);
 363       }
 364     }
 365   }
 366
 367   //////////////////////////////
 368   // test2 "MED_BODY_FITTED"  //
 369   //////////////////////////////
 370   {
 371     MESH* myMesh1 = myMed->getMesh(mesh_names[1]);
 372     myMesh1->read();
 373
 374     CPPUNIT_ASSERT(myMesh1 != NULL);
 375     CPPUNIT_ASSERT(myMesh1->getIsAGrid());
 376
 377     GRID* myGrid1 = dynamic_cast<GRID*>(myMesh1);
 378     CPPUNIT_ASSERT(myGrid1);
 379
 380     int I, J, K;
 381     CPPUNIT_ASSERT_NO_THROW(I = myGrid1->getArrayLength(1));
 382     CPPUNIT_ASSERT_NO_THROW(J = myGrid1->getArrayLength(2));
 383     CPPUNIT_ASSERT_NO_THROW(K = myGrid1->getArrayLength(3));
 384
 385     CPPUNIT_ASSERT(I == 2);
 386     CPPUNIT_ASSERT(J == 2);
 387
 388     med_grid_type grid_type = myGrid1->getGridType();
 389     CPPUNIT_ASSERT_MESSAGE("Wrong grid type", grid_type == MED_BODY_FITTED);
 390
 391     int nbTypesCell = myGrid1->getNumberOfTypes(MED_CELL);
 392     CPPUNIT_ASSERT(nbTypesCell == 1);
 393
 394     const medGeometryElement* types1;
 395     CPPUNIT_ASSERT_NO_THROW( types1 = myGrid1->getTypes(MED_CELL) );
 396     CPPUNIT_ASSERT( types1[0] == MED_QUAD4);
 397
 398     int nbElem;
 399     CPPUNIT_ASSERT_NO_THROW(nbElem = myGrid1->getNumberOfElements(MED_CELL, types1[0]));
 400     CPPUNIT_ASSERT(nbElem);
 401
 402     const int* BodyConnectivity;
 403     const int* body_connectivity_index;
 404     int ijkNodeBody[3];
 405     CPPUNIT_ASSERT_NO_THROW(BodyConnectivity = myGrid1->getConnectivity(MED_FULL_INTERLACE, MED_NODAL,
 406                                                                         MED_CELL, types1[0]));
 407     CPPUNIT_ASSERT_NO_THROW(body_connectivity_index = myGrid1->getConnectivityIndex(MED_NODAL, MED_CELL));
 408     cout<<"Nodal connectivity"<<endl;
 409     for (int j = 0; j < nbElem; j++) {
 410       cout << "Element "<< j+1 << " : ";
 411       for (int k = body_connectivity_index[j]; k < body_connectivity_index[j+1]; k++){
 412         CPPUNIT_ASSERT_NO_THROW(myGrid1->getNodePosition(BodyConnectivity[k-1], ijkNodeBody[0],
 413                                                          ijkNodeBody[1], ijkNodeBody[2]));
 414         cout << BodyConnectivity[k-1] << " ";
 415       }
 416       cout << endl;
 417     }
 418   }
 419
 420   ///////////////////////////////////////////////////
 421   // test3 "maa1" which in fact is not a pure GRID //
 422   ///////////////////////////////////////////////////
 423   {
 424     MESH* myMesh2 = NULL;
 425
 426     myMesh2 = myMed->getMesh(mesh_names[2]);
 427     myMesh2->read();
 428
 429     CPPUNIT_ASSERT(myMesh2 != NULL);
 430     CPPUNIT_ASSERT(!(myMesh2->getIsAGrid()));
 431   }
 432
 433   delete myMed;
 434
 435   ////////////////////////////
 436   // test4 create new GRID  //
 437   ////////////////////////////
 438
 439   // Default constructor and destructor
 440   {
 441     GRID* myGrid2 = new GRID();
 442     CPPUNIT_ASSERT(myGrid2->getIsAGrid());
 443     CPPUNIT_ASSERT(myGrid2->getGridType() == MED_CARTESIAN);
 444     CPPUNIT_ASSERT(!myGrid2->getArrayLength(1));
 445     CPPUNIT_ASSERT(!myGrid2->getArrayLength(2));
 446     CPPUNIT_ASSERT(!myGrid2->getArrayLength(3));
 447     delete myGrid2;
 448   }
 449
 450   // Constructor with grid type, setGridType()
 451   {
 452     GRID* myGrid2 = new GRID(MED_POLAR);
 453     CPPUNIT_ASSERT(myGrid2->getGridType() == MED_POLAR);
 454     myGrid2->setGridType(MED_CARTESIAN);
 455     CPPUNIT_ASSERT(myGrid2->getGridType() == MED_CARTESIAN);
 456     delete myGrid2;
 457   }
 458
 459   // Constructor with coordinate values, getArrayValue(), init()
 460   {
 461     vector<vector<double> > xyz;
 462     const int nbCoords = 3;
 463     xyz.resize(nbCoords);
 464     for ( int i = 0; i < nbCoords; i++ )
 465     {
 466       xyz[i].resize(i + 2);
 467       for ( int j = 0; j < i + 2; j++ )
 468         xyz[i][j] = j;
 469     }
 470     vector<string> Coord_Names;
 471     Coord_Names.resize(nbCoords);
 472     Coord_Names[0] = "X";
 473     Coord_Names[1] = "Y";
 474     Coord_Names[2] = "Z";
 475
 476     vector<string> Coord_Units;
 477     Coord_Units.resize(nbCoords);
 478     for(int i = 0; i < 3; i++)
 479       Coord_Units[i] = "cm";
 480
 481     GRID* myGrid2;
 482
 483     try{
 484       myGrid2 = new GRID(xyz, Coord_Names, Coord_Units, MED_CARTESIAN);
 485     }
 486     catch (const std::exception &e)
 487     {
 488       CPPUNIT_FAIL(e.what());
 489     }
 490     catch (...)
 491     {
 492       CPPUNIT_FAIL("Unknown exception");
 493     }
 494
 495     // testing getCoordinateptr() and fillCoordinates()
 496     // We fill a map of all possible coordinate triples.
 497     // After iteration through all coordinates, this map should contain only "true" as data.
 498     // "true" in some map element during iteration means duplicated node position.
 499     // "false" in some map element after iteration means empty node position.
 500     map<int, bool> found;
 501     for ( int i1 = 0; i1 < xyz[0].size(); i1++ )
 502       for ( int i2 = 0; i2 < xyz[1].size(); i2++ )
 503         for ( int i3 = 0; i3 < xyz[2].size(); i3++ )
 504           found[int(xyz[0][i1] * 100 + xyz[1][i2] * 10 + xyz[2][i3])] = false;
 505
 506     COORDINATE* coords = (COORDINATE*)myGrid2->getCoordinateptr();
 507     CPPUNIT_ASSERT(coords);
 508     for (int num = 0; num < myGrid2->getNumberOfNodes(); num++) {
 509       int x = int(coords->getCoordinate(num + 1, 1));
 510       int y = int(coords->getCoordinate(num + 1, 2));
 511       int z = int(coords->getCoordinate(num + 1, 3));
 512       CPPUNIT_ASSERT(!found[x * 100 + y * 10 + z]);
 513       found[x * 100 + y * 10 + z] = true;
 514     }
 515
 516     for ( map<int, bool>::iterator it = found.begin(); it != found.end(); it++ )
 517       CPPUNIT_ASSERT((*it).second);
 518
 519     // Testing fillConnectivity() and getConnectivityptr()
 520     // Basic testing: presence of connectivity arrays, element types and number of elements
 521     CONNECTIVITY* conn = (CONNECTIVITY*)myGrid2->getConnectivityptr();
 522     CPPUNIT_ASSERT(conn);
 523     bool hasFaces = myGrid2->getArrayLength(3), hasEdges = myGrid2->getArrayLength(2);
 524     medGeometryElement aCellGeometry;
 525     if (hasFaces)      aCellGeometry = MED_HEXA8;
 526     else if (hasEdges) aCellGeometry = MED_QUAD4;
 527     else               aCellGeometry = MED_SEG2;
 528     CPPUNIT_ASSERT(conn->getElementType(MED_CELL, 1) == aCellGeometry);
 529     CPPUNIT_ASSERT(conn->existConnectivity(MED_NODAL,      MED_CELL));
 530     CPPUNIT_ASSERT(conn->existConnectivity(MED_DESCENDING, MED_CELL));
 531     //test getCellsTypes
 532     CELLMODEL* cellmodel = (CELLMODEL*)myGrid2->getCellsTypes(MED_CELL);
 533     CPPUNIT_ASSERT(cellmodel);
 534
 535     int nbCells, nbFaces, nbEdges;
 536
 537     int iLen     = myGrid2->getArrayLength(1);
 538     int jLen     = myGrid2->getArrayLength(2);
 539     int kLen     = myGrid2->getArrayLength(3);
 540     int iLenMin1 = myGrid2->getArrayLength(1)-1;
 541     int jLenMin1 = myGrid2->getArrayLength(2)-1;
 542     int kLenMin1 = myGrid2->getArrayLength(3)-1;
 543     const int* aCellCount = conn->getGlobalNumberingIndex(MED_CELL);
 544     nbCells = iLenMin1 * jLenMin1 * kLenMin1;
 545     CPPUNIT_ASSERT(aCellCount[1] - 1 == nbCells);
 546
 547     if (hasFaces){
 548       CPPUNIT_ASSERT(conn->getElementType(MED_FACE, 1) == MED_QUAD4);
 549       nbFaces  = iLen * jLenMin1 * kLenMin1;
 550       nbFaces += jLen * kLenMin1 * iLenMin1;
 551       nbFaces += kLen * iLenMin1 * jLenMin1;
 552       const int* aFaceCount = conn->getGlobalNumberingIndex(MED_FACE);
 553       CPPUNIT_ASSERT(aFaceCount[1] - 1 == nbFaces);
 554       CPPUNIT_ASSERT(conn->existConnectivity(MED_NODAL, MED_FACE));
 555       //test getCellsTypes
 556       CELLMODEL* cellmodelF = (CELLMODEL*)myGrid2->getCellsTypes(MED_FACE);
 557       CPPUNIT_ASSERT(cellmodelF);
 558     }
 559     if (hasEdges){
 560       CPPUNIT_ASSERT(conn->getElementType(MED_EDGE, 1) == MED_SEG2);
 561       if (kLen) { // 3d grid
 562         nbEdges  = iLenMin1 * jLen * kLen;
 563         nbEdges += jLenMin1 * kLen * iLen;
 564         nbEdges += kLenMin1 * iLen * jLen;
 565       }
 566       else if (jLen) { // 2d
 567         nbEdges  = iLenMin1 * jLen;
 568         nbEdges += jLenMin1 * iLen;
 569       }
 570       const int* anEdgeCount = conn->getGlobalNumberingIndex(MED_EDGE);
 571       CPPUNIT_ASSERT(anEdgeCount[1] - 1 == nbEdges);
 572       CPPUNIT_ASSERT(conn->existConnectivity(MED_NODAL, MED_EDGE));
 573       //test getCellsTypes
 574       CELLMODEL* cellmodelE = (CELLMODEL*)myGrid2->getCellsTypes(MED_EDGE);
 575       CPPUNIT_ASSERT(cellmodelE);
 576
 577     }
 578
 579     // Testing getArrayValue()
 580     for ( int ii = 1; ii <= nbCoords; ii++ )
 581       for ( int jj = 0; jj < ii + 1; jj++ )
 582         CPPUNIT_ASSERT(myGrid2->getArrayValue(ii, jj) == xyz[ii - 1][jj]);
 583
 584     CPPUNIT_ASSERT_THROW(myGrid2->getArrayValue(nbCoords + 1, 0), MEDEXCEPTION);
 585     CPPUNIT_ASSERT_THROW(myGrid2->getArrayValue(1, myGrid2->getArrayLength(1) + 1), MEDEXCEPTION);
 586     myGrid2->setGridType(MED_POLAR);
 587
 588     //testing read/write functions
 589
 590     // add new driver
 591     int idGridV21;
 592
 593     try
 594     {
 595       idGridV21 = myGrid2->addDriver(MED_DRIVER,filenameout21);
 596     }
 597     catch(MEDEXCEPTION &e)
 598     {
 599       CPPUNIT_FAIL(e.what());
 600     }
 601     catch( ... )
 602     {
 603       CPPUNIT_FAIL("Unknown exception");
 604     }
 605
 606     // write this driver to file as an unstructured mesh
 607     CPPUNIT_ASSERT_NO_THROW(myGrid2->writeUnstructured(idGridV21));
 608
 609     GRID* myGrid3 = new GRID();
 610     // add new driver for myGrid3
 611     int driver;
 612     CPPUNIT_ASSERT_NO_THROW(driver = myGrid3->addDriver(MED_DRIVER, filenameout21));
 613
 614 //#ifdef ENABLE_FORCED_FAILURES
 615     // ? (BUG) ? "The mesh dimension |-1| seems to be incorrect for the mesh : |Default Mesh Name|"
 616     CPPUNIT_ASSERT_NO_THROW(myGrid3->read(driver));
 617
 618     // Testing getArrayValue()
 619     for ( int ii = 1; ii <= nbCoords; ii++ )
 620       for ( int jj = 0; jj < ii + 1; jj++ )
 621         CPPUNIT_ASSERT(myGrid3->getArrayValue(ii, jj) == xyz[ii - 1][jj]);
 622
 623     CPPUNIT_ASSERT(myGrid3->getGridType() == MED_POLAR);
 624 //#endif
 625
 626     delete myGrid3;
 627
 628     //test init()
 629     CPPUNIT_ASSERT_NO_THROW(myGrid2->init());
 630     CPPUNIT_ASSERT(myGrid2->getGridType() == MED_CARTESIAN);
 631     CPPUNIT_ASSERT(myGrid2->getArrayLength(1) == 0);
 632     CPPUNIT_ASSERT(myGrid2->getArrayLength(2) == 0);
 633     CPPUNIT_ASSERT(myGrid2->getArrayLength(3) == 0);
 634 //#ifdef ENABLE_FAULTS
 635     // (BUG) Segmentation Fault
 636     myGrid2->makeUnstructured();
 637 //#endif
 638 //#ifdef ENABLE_FORCED_FAILURES
 639 //    CPPUNIT_FAIL("ERROR:makeUnstructured() - there is no check if grid is empty or not");
 640 //#endif
 641
 642     delete myGrid2;
 643   }
 644
 645 //#ifdef ENABLE_FORCED_FAILURES
 646   {
 647     GRID* myGrid2;
 648     // ? (BUG) ? MED Exception in /dn20/salome/jfa/V3/SRC/MED_SRC/src/MEDMEM/MEDMEM_MedM
 649     //eshDriver21.cxx [430] : MED_MESH_RDONLY_DRIVER21::getCOORDINATE() : The number
 650     //of nodes |0| seems to be incorrect for the mesh : |bodyfitted|
 651     CPPUNIT_ASSERT_NO_THROW(myGrid2 = new GRID(MED_DRIVER, filename, mesh_names[1]));
 652
 653     // Check if something has been read - full mesh data testing is above
 654     CPPUNIT_ASSERT(myGrid2->getSpaceDimension());
 655     CPPUNIT_ASSERT(myGrid2->getNumberOfNodes());
 656     CPPUNIT_ASSERT(myGrid2->getNumberOfTypes(MED_CELL) == 1);
 657     const medGeometryElement* types2;
 658     CPPUNIT_ASSERT_NO_THROW(types2 = myGrid2->getTypes(MED_CELL));
 659     int nbElem;
 660     CPPUNIT_ASSERT_NO_THROW(nbElem = myGrid2->getNumberOfElements(MED_CELL,types2[0]));
 661     CPPUNIT_ASSERT(nbElem);
 662     delete myGrid2;
 663   }
 664 //#endif
 665
 666   {
 667     GRID* myGrid4 = new GRID();
 668     filename = datadir + "/share/salome/resources/med/pointe.med";
 669     myGrid4->setName("maa1");
 670     MED_MESH_RDONLY_DRIVER myMeshDriver(filename, myGrid4);
 671     myMeshDriver.setMeshName("maa1");
 672
 673     // add new driver for myGrid4
 674     int driver;
 675     CPPUNIT_ASSERT_NO_THROW(driver = myGrid4->addDriver(myMeshDriver));
 676
 677     // ??? ERROR:myGrid4->fillMeshAfterRead()- this method is incomplete:
 678     // currently it only resets _is_coordinates_filled and _is_connectivity_filled
 679     // flags that leads to grid reconstruction
 680
 681     // MED Exception : MED_MESH_RDONLY_DRIVER21::getGRID() : The number
 682     // of nodes |-1| seems to be incorrect for the mesh : |maa1|
 683     // But this exception seems to be correct reaction on attempt
 684     // to read a grid from a file, which does not contain it.
 685     CPPUNIT_ASSERT_THROW(myGrid4->read(driver), MEDEXCEPTION);
 686     //CPPUNIT_ASSERT_NO_THROW(myGrid4->read(driver));
 687     /*CPPUNIT_ASSERT(myGrid4->getArrayLength(1) == 0);
 688     CPPUNIT_ASSERT(myGrid4->getArrayLength(2) == 0);
 689     CPPUNIT_ASSERT(myGrid4->getArrayLength(3) == 0);
 690     myGrid4->fillMeshAfterRead();
 691     CPPUNIT_ASSERT(myGrid4->getArrayLength(1) != 0);
 692     CPPUNIT_ASSERT(myGrid4->getArrayLength(2) != 0);
 693     CPPUNIT_ASSERT(myGrid4->getArrayLength(3) != 0);*/
 694
 695     delete myGrid4;
 696   }
 697 }
 698
 699 int main (int argc, char** argv)
 700 {
 701   MEDMEMTest_testGrid();
 702 }