src/MEDLoader/Test/SauvLoaderTest.cxx

   1 // Copyright (C) 2007-2015  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "SauvLoaderTest.hxx"
  21
  22 #include "SauvReader.hxx"
  23 #include "SauvWriter.hxx"
  24 #include "MEDFileData.hxx"
  25 #include "MEDCouplingFieldDouble.hxx"
  26 #include "MEDCouplingMemArray.hxx"
  27
  28 #ifdef WIN32
  29 #include <windows.h>
  30 #include <direct.h>
  31 #define getcwd _getcwd
  32 #else
  33 # include <unistd.h>
  34 #endif
  35
  36 #include <vector>
  37 #include <string>
  38
  39 using namespace ParaMEDMEM;
  40
  41 void SauvLoaderTest::testSauv2Med()
  42 {
  43   // read a file containing all types of readable piles
  44   std::string file = getResourceFile("allPillesTest.sauv");
  45   MEDCouplingAutoRefCountObjectPtr<SauvReader> sr=SauvReader::New(file.c_str());
  46   MEDCouplingAutoRefCountObjectPtr<MEDFileData> d2=sr->loadInMEDFileDS();
  47   // write MED
  48   d2->write("allPillesTest.med",0);
  49   // check
  50   CPPUNIT_ASSERT_EQUAL(1,d2->getNumberOfMeshes());
  51   CPPUNIT_ASSERT_EQUAL(8+97,d2->getNumberOfFields());
  52   MEDFileMesh * m = d2->getMeshes()->getMeshAtPos(0);
  53   CPPUNIT_ASSERT_EQUAL(17,int(m->getGroupsNames().size()));
  54 }
  55
  56 void SauvLoaderTest::testMed2SauvOnAMeshWithVoidFamily()
  57 {
  58   // Create a mesh with 2 quads.
  59   const int spaceDim = 2;
  60   const int nbOfNodes = 6;
  61   double coords[nbOfNodes*spaceDim] = {0,0, 1,0, 1,1, 0,1, 2,0, 2,1};
  62   int conn[8]={0,1,2,3, 1,4,5,2};
  63   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh2d=MEDCouplingUMesh::New("Mesh",spaceDim);
  64   mesh2d->allocateCells(2);
  65   mesh2d->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
  66   mesh2d->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
  67   mesh2d->finishInsertingCells();
  68   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> myCoords=DataArrayDouble::New();
  69   myCoords->alloc(nbOfNodes,spaceDim);
  70   std::copy(coords,coords+nbOfNodes*spaceDim,myCoords->getPointer());
  71   mesh2d->setCoords(myCoords);
  72
  73   // create a MedFileUMesh
  74   MEDCouplingAutoRefCountObjectPtr<MEDFileUMesh> m= MEDFileUMesh::New();
  75   m->setMeshAtLevel(0,mesh2d);
  76
  77   // Create families and groups
  78
  79   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> fam = DataArrayInt::New();
  80   fam->alloc(2,1);
  81   int elemsFams[2] = {-2,-3};
  82   std::copy(elemsFams,elemsFams+2,fam->getPointer());
  83
  84   m->setFamilyFieldArr(0,fam);
  85
  86   std::map<std::string,int> theFamilies;
  87   theFamilies["FAM_-1"]=-1;
  88   theFamilies["FAM_-2"]=-2;
  89   theFamilies["FAM_-3"]=-3;
  90
  91   std::map<std::string, std::vector<std::string> > theGroups;
  92   theGroups["Group1"].push_back("FAM_-2");
  93   theGroups["Group2"].push_back("FAM_-3");
  94   theGroups["Grouptot"].push_back("FAM_-1");
  95   theGroups["Grouptot"].push_back("FAM_-2");
  96   theGroups["Grouptot"].push_back("FAM_-3");
  97   m->setFamilyInfo(theFamilies);
  98   m->setGroupInfo(theGroups);
  99
 100   // write to MED for visual check
 101   //const char* medFile = "mesh_with_void_family.med";
 102   //m->write(medFile, 2);
 103
 104   // write to SAUV
 105   const char* sauvFile = "mesh_with_void_family.sauv";
 106   MEDCouplingAutoRefCountObjectPtr<MEDFileData> medData = MEDFileData::New();
 107   MEDCouplingAutoRefCountObjectPtr<MEDFileMeshes> medMeshes = MEDFileMeshes::New();
 108   MEDCouplingAutoRefCountObjectPtr<SauvWriter> sw=SauvWriter::New();
 109   medMeshes->setMeshAtPos(0, m);
 110   medData->setMeshes(medMeshes);
 111   sw->setMEDFileDS(medData);
 112   sw->write(sauvFile);
 113
 114   // read SAUV and check groups
 115   MEDCouplingAutoRefCountObjectPtr<SauvReader> sr=SauvReader::New(sauvFile);
 116   MEDCouplingAutoRefCountObjectPtr<MEDFileData> d2=sr->loadInMEDFileDS();
 117   MEDFileUMesh* m2 = static_cast<MEDFileUMesh*>( d2->getMeshes()->getMeshAtPos(0) );
 118   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> group1 = m2->getGroup(0, "Group1");
 119   CPPUNIT_ASSERT_EQUAL(1,(int)group1->getNumberOfCells());
 120   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> group2 = m2->getGroup(0, "Group2");
 121   CPPUNIT_ASSERT_EQUAL(1,(int)group2->getNumberOfCells());
 122   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> grptot = m2->getGroup(0, "Grouptot");
 123   CPPUNIT_ASSERT_EQUAL(2,(int)grptot->getNumberOfCells());
 124 }
 125
 126 void SauvLoaderTest::testSauv2MedOnA3SubsField()
 127 {
 128   // read SAUV
 129   std::string sauvFile = getResourceFile("portico_3subs.sauv");
 130   MEDCouplingAutoRefCountObjectPtr<SauvReader> sr=SauvReader::New(sauvFile.c_str());
 131   MEDCouplingAutoRefCountObjectPtr<MEDFileData> d2=sr->loadInMEDFileDS();
 132   // check mesh
 133   MEDFileUMesh* m2 = static_cast<MEDFileUMesh*>(d2->getMeshes()->getMeshAtPos(0));
 134   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh1d = m2->getMeshAtLevel(0);
 135   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> length1dField = mesh1d->getMeasureField(0);
 136   std::cout << "Length of 1d elements: " << length1dField->accumulate(0) << std::endl;
 137   CPPUNIT_ASSERT_DOUBLES_EQUAL(3, length1dField->accumulate(0), 1e-12);
 138   // check field
 139   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> field =
 140     dynamic_cast<MEDFileFieldMultiTS *>(d2->getFields()->getFieldWithName("CHAM1D"));
 141   std::cout << "Number of components in field: " << field->getInfo().size() << std::endl;
 142   CPPUNIT_ASSERT_EQUAL(6,(int)field->getInfo().size());
 143   std::vector< std::pair<int,int> > timesteps = field->getIterations();
 144
 145   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> field1d =
 146     field->getFieldOnMeshAtLevel(ON_GAUSS_NE, timesteps[0].first, timesteps[0].second, 0, m2);
 147
 148   // Check first component of the field
 149   // 2 gauss points per element => 12 values
 150   double values[12] = {
 151       -7.687500000000e-03,
 152       -7.687500000000e-03,
 153       -4.562500000000e-03,
 154       -4.562500000000e-03,
 155       -8.208333333333e-03,
 156       -8.208333333333e-03,
 157       -6.125000000000e-03,
 158       -6.125000000000e-03,
 159       -4.041666666666e-03,
 160       -4.041666666666e-03,
 161       -6.111413346910e-07,
 162       -6.111413346910e-07};
 163
 164   for (int i=0; i < field1d->getNumberOfTuples(); i++)
 165   {
 166     CPPUNIT_ASSERT_DOUBLES_EQUAL( values[i], field1d->getIJ(i, 0), 1e-12 );
 167   }
 168 }
 169
 170 void SauvLoaderTest::testMed2Sauv()
 171 {
 172   // read pointe.med
 173   std::string file = getResourceFile("pointe.med");
 174   MEDCouplingAutoRefCountObjectPtr<MEDFileData> pointeMed=MEDFileData::New(file.c_str());
 175
 176   // add 3 faces to pointeMed
 177   MEDFileUMesh* pointeMedMesh = static_cast<MEDFileUMesh*>(pointeMed->getMeshes()->getMeshAtPos(0));
 178   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> pointeM1D = MEDCouplingUMesh::New();
 179   DataArrayDouble     *coords = pointeMedMesh->getCoords();
 180   pointeM1D->setCoords( coords );
 181   pointeM1D->setMeshDimension( 2 );
 182   pointeM1D->allocateCells( 3 );
 183   int conn[]=
 184     {
 185       0,1,2, 0,1,3, 10,11,12,13
 186     };
 187   pointeM1D->insertNextCell( INTERP_KERNEL::NORM_TRI3, 3, conn);
 188   pointeM1D->insertNextCell( INTERP_KERNEL::NORM_TRI3, 3, conn+3);
 189   pointeM1D->insertNextCell( INTERP_KERNEL::NORM_QUAD4, 4, conn+6);
 190   pointeM1D->finishInsertingCells();
 191   pointeMedMesh->setMeshAtLevel( -1, pointeM1D );
 192   pointeMed->getMeshes()->setMeshAtPos( 0, pointeMedMesh );
 193
 194   // add a field on 2 faces to pointeMed
 195   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> ff1=MEDFileFieldMultiTS::New();
 196   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> f1=MEDCouplingFieldDouble::New(ON_GAUSS_NE,ONE_TIME);
 197   f1->setMesh( pointeM1D );
 198   f1->setName("Field on 2 faces");
 199   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> d=DataArrayDouble::New();
 200   d->alloc(3+4,2);
 201   d->setInfoOnComponent(0,"sigX [MPa]");
 202   d->setInfoOnComponent(1,"sigY [GPa]");
 203   double vals[2*(3+4)] =
 204     {
 205       311,312,321,322,331,332,411,412,421,422,431,432,441,442
 206     };
 207   std::copy(vals,vals+d->getNbOfElems(),d->getPointer());
 208   f1->setArray(d);
 209   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> da=DataArrayInt::New();
 210   int ids[] =
 211     {
 212       0,2
 213     };
 214   da->alloc(2,1);
 215   std::copy(ids,ids+da->getNbOfElems(),da->getPointer());
 216   da->setName("sup2");
 217   ff1->appendFieldProfile(f1,pointeMedMesh,-1,da);
 218   pointeMed->getFields()->pushField( ff1 );
 219
 220   // remove "fieldnodeint"
 221   MEDFileFields* pointeFields = pointeMed->getFields();
 222   for ( int i = 0; i < pointeFields->getNumberOfFields(); ++i )
 223     {
 224       MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> ts = pointeFields->getFieldAtPos(i);
 225       if ( std::string("fieldnodeint") == ts->getName())
 226         {
 227           pointeFields->destroyFieldAtPos( i );
 228           break;
 229         }
 230     }
 231   // write pointeMed to SAUV
 232   const char* sauvFile = "pointe.sauv";
 233   MEDCouplingAutoRefCountObjectPtr<SauvWriter> sw=SauvWriter::New();
 234   sw->setMEDFileDS(pointeMed);
 235   sw->write(sauvFile);
 236
 237   // read SAUV and check
 238   MEDCouplingAutoRefCountObjectPtr<SauvReader> sr=SauvReader::New(sauvFile);
 239   MEDCouplingAutoRefCountObjectPtr<MEDFileData> d2=sr->loadInMEDFileDS();
 240   CPPUNIT_ASSERT_EQUAL(1,d2->getNumberOfMeshes());
 241   CPPUNIT_ASSERT_EQUAL(4,d2->getNumberOfFields());
 242   MEDFileUMesh * m = static_cast<MEDFileUMesh*>( d2->getMeshes()->getMeshAtPos(0) );
 243   CPPUNIT_ASSERT_EQUAL(std::string("maa1"),std::string(m->getName() ));
 244   CPPUNIT_ASSERT_EQUAL(3,m->getMeshDimension());
 245   std::vector<std::string > groups = m->getGroupsNames();
 246   CPPUNIT_ASSERT_EQUAL(6,(int)groups.size());
 247   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"groupe1") != groups.end() );
 248   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"groupe2") != groups.end() );
 249   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"groupe3") != groups.end() );
 250   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"groupe4") != groups.end() );
 251   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"groupe5") != groups.end() );
 252   CPPUNIT_ASSERT( std::find(groups.begin(),groups.end(),"maa1") != groups.end() );
 253   CPPUNIT_ASSERT_EQUAL(16,m->getSizeAtLevel(0));
 254   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> um0 = m->getGenMeshAtLevel(0);
 255   CPPUNIT_ASSERT_EQUAL(12, um0->getNumberOfCellsWithType( INTERP_KERNEL::NORM_TETRA4 ));
 256   CPPUNIT_ASSERT_EQUAL(2,  um0->getNumberOfCellsWithType( INTERP_KERNEL::NORM_PYRA5 ));
 257   CPPUNIT_ASSERT_EQUAL(2,  um0->getNumberOfCellsWithType( INTERP_KERNEL::NORM_HEXA8 ));
 258   MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> um1 = m->getGenMeshAtLevel(-1);
 259   CPPUNIT_ASSERT_EQUAL(2, um1->getNumberOfCellsWithType( INTERP_KERNEL::NORM_TRI3 ));
 260   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> pointeUM0 =
 261     static_cast<MEDCouplingUMesh*>( pointeMedMesh->getGenMeshAtLevel(0));
 262   DataArrayDouble *coo = m->getCoords();
 263   DataArrayDouble *pointeCoo = pointeMedMesh->getCoords();
 264   CPPUNIT_ASSERT(coo->isEqualWithoutConsideringStr(*pointeCoo,1e-12));
 265   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol = um0->getMeasureField(0);
 266   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> pointeVol = pointeUM0->getMeasureField(0);
 267   CPPUNIT_ASSERT_DOUBLES_EQUAL( vol->accumulate(0), pointeVol->accumulate(0),1e-12);
 268   // check fields
 269   // fieldnodedouble
 270   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> fieldnodedoubleTS1 =
 271     dynamic_cast<MEDFileFieldMultiTS *>(pointeMed->getFields()->getFieldWithName("fieldnodedouble"));
 272   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> fieldnodedoubleTS2 =
 273     dynamic_cast<MEDFileFieldMultiTS *>(d2->getFields()->getFieldWithName("fieldnodedouble"));
 274   CPPUNIT_ASSERT_EQUAL( fieldnodedoubleTS1->getInfo().size(), fieldnodedoubleTS2->getInfo().size());
 275   for ( size_t i = 0; i < fieldnodedoubleTS1->getInfo().size(); ++i )
 276     CPPUNIT_ASSERT_EQUAL( fieldnodedoubleTS1->getInfo()[i], fieldnodedoubleTS2->getInfo()[i]);
 277   CPPUNIT_ASSERT_EQUAL( fieldnodedoubleTS1->getNumberOfTS(), fieldnodedoubleTS2->getNumberOfTS());
 278   std::vector< std::pair<int,int> > io1 = fieldnodedoubleTS1->getIterations();
 279   std::vector< std::pair<int,int> > io2 = fieldnodedoubleTS2->getIterations();
 280   for ( int i =0; i < fieldnodedoubleTS1->getNumberOfTS(); ++i )
 281     {
 282       MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> fnd1 =
 283         fieldnodedoubleTS1->getFieldOnMeshAtLevel(ON_NODES, io1[i].first,io1[i].second,pointeUM0);
 284       MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> fnd2 =
 285         fieldnodedoubleTS2->getFieldOnMeshAtLevel(ON_NODES, io2[i].first,io2[i].second,um0);
 286       CPPUNIT_ASSERT( fnd1->getArray()->isEqual( *fnd2->getArray(), 1e-12 ));
 287     }
 288   // fieldcelldoublevector
 289   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> fieldcelldoublevectorTS1 =
 290     dynamic_cast<MEDFileFieldMultiTS *>(pointeMed->getFields()->getFieldWithName("fieldcelldoublevector"));
 291   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> fieldcelldoublevectorTS2 =
 292     dynamic_cast<MEDFileFieldMultiTS *>(d2->getFields()->getFieldWithName("fieldcelldoublevector"));
 293   CPPUNIT_ASSERT_EQUAL( fieldcelldoublevectorTS1->getInfo().size(), fieldcelldoublevectorTS2->getInfo().size());
 294   for ( size_t i = 0; i < fieldcelldoublevectorTS1->getInfo().size(); ++i )
 295     CPPUNIT_ASSERT_EQUAL( fieldcelldoublevectorTS1->getInfo()[i], fieldcelldoublevectorTS2->getInfo()[i]);
 296   CPPUNIT_ASSERT_EQUAL( fieldcelldoublevectorTS1->getNumberOfTS(), fieldcelldoublevectorTS2->getNumberOfTS());
 297   io1 = fieldcelldoublevectorTS1->getIterations();
 298   io2 = fieldcelldoublevectorTS2->getIterations();
 299   for ( int i =0; i < fieldcelldoublevectorTS1->getNumberOfTS(); ++i )
 300     {
 301       MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> fnd1 =
 302         fieldcelldoublevectorTS1->getFieldOnMeshAtLevel(ON_CELLS, io1[i].first,io1[i].second,pointeUM0);
 303       MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> fnd2 =
 304         fieldcelldoublevectorTS2->getFieldOnMeshAtLevel(ON_CELLS, io2[i].first,io2[i].second,um0);
 305       CPPUNIT_ASSERT_DOUBLES_EQUAL( fnd1->accumulate(0), fnd2->accumulate(0), 1e-12 );
 306       CPPUNIT_ASSERT_DOUBLES_EQUAL( fnd1->accumulate(1), fnd2->accumulate(1), 1e-12 );
 307       CPPUNIT_ASSERT_DOUBLES_EQUAL( fnd1->accumulate(2), fnd2->accumulate(2), 1e-12 );
 308     }
 309   // "Field on 2 faces"
 310   MEDCouplingAutoRefCountObjectPtr<MEDFileFieldMultiTS> fieldOnFaces =
 311     dynamic_cast<MEDFileFieldMultiTS *>(d2->getFields()->getFieldWithName(f1->getName().c_str()));
 312   io1 = fieldOnFaces->getIterations();
 313   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> fof =
 314     fieldOnFaces->getFieldOnMeshAtLevel(f1->getTypeOfField(),io1[0].first,io1[0].second,um1);
 315   CPPUNIT_ASSERT( d->isEqual( *fof->getArray(), 1e-12 ));
 316 }
 317
 318 void SauvLoaderTest::tearDown()
 319 {
 320   const int nbFilesToRemove = 3;
 321   const char* fileToRemove[nbFilesToRemove] = { "allPillesTest.med", "pointe.sauv", "mesh_with_void_family.sauv" };
 322   for ( int i = 0; i < nbFilesToRemove; ++i )
 323     {
 324 #ifdef WIN32
 325       if (GetFileAttributes(fileToRemove[i]) != INVALID_FILE_ATTRIBUTES)
 326 #else
 327         if (access(fileToRemove[i], F_OK) == 0)
 328 #endif
 329       remove(fileToRemove[i]);
 330   }
 331 }
 332
 333 std::string SauvLoaderTest::getResourceFile( const std::string& filename )
 334 {
 335   std::string resourceFile = "";
 336
 337   if ( getenv("MEDCOUPLING_ROOT_DIR") ) {
 338     // use MEDCOUPLING_ROOT_DIR env.var
 339     resourceFile = getenv("MEDCOUPLING_ROOT_DIR");
 340     resourceFile += "/share/resources/med/";
 341   }
 342   else {
 343     resourceFile = getcwd(NULL, 0);
 344     resourceFile += "/../../../resources/";
 345   }
 346
 347   resourceFile += filename;
 348 #ifdef WIN32
 349   std::string fixedpath = resourceFile;
 350   for ( int i=0; i < fixedpath.length(); ++i )
 351     if (fixedpath[i] == '/')
 352       fixedpath[i] = '\\';
 353   return fixedpath;
 354 #endif
 355   return resourceFile;
 356 }