1 // Copyright (C) 2007-2015 CEA/DEN, EDF R&D
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay (CEA/DEN)
21 #include "MEDLoaderTest.hxx"
22 #include "MEDLoader.hxx"
23 #include "MEDLoaderBase.hxx"
24 #include "MEDCouplingUMesh.hxx"
25 #include "MEDCouplingFieldDouble.hxx"
26 #include "MEDCouplingMemArray.hxx"
30 #include <unistd.h> // get_current_dir_name()
32 using namespace ParaMEDMEM;
34 void MEDLoaderTest::testMesh1DRW()
36 MEDCouplingUMesh *mesh=build1DMesh_1();
37 mesh->checkCoherency();
38 MEDLoader::WriteUMesh("file1.med",mesh,true);
39 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile("file1.med",mesh->getName().c_str(),0);
40 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
45 void MEDLoaderTest::testMesh2DCurveRW()
47 MEDCouplingUMesh *mesh=build2DCurveMesh_1();
48 mesh->checkCoherency();
49 MEDLoader::WriteUMesh("file2.med",mesh,true);
50 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile("file2.med",mesh->getName().c_str(),0);
51 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
56 void MEDLoaderTest::testMesh2DRW()
58 MEDCouplingUMesh *mesh=build2DMesh_1();
59 mesh->checkCoherency();
60 MEDLoader::WriteUMesh("file3.med",mesh,true);
61 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile("file3.med",mesh->getName().c_str(),0);
62 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
67 void MEDLoaderTest::testMesh3DSurfRW()
69 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
70 mesh->checkCoherency();
71 MEDLoader::WriteUMesh("file4.med",mesh,true);
72 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile("file4.med",mesh->getName().c_str(),0);
73 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
78 void MEDLoaderTest::testMesh3DRW()
80 MEDCouplingUMesh *mesh=build3DMesh_1();
81 mesh->checkCoherency();
82 MEDLoader::WriteUMesh("file5.med",mesh,true);
83 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile("file5.med",mesh->getName().c_str(),0);
84 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
90 * Most basic test : one and only one MEDCoupling field in a new file.
92 void MEDLoaderTest::testFieldRW1()
94 MEDCouplingFieldDouble *f1=buildVecFieldOnCells_1();
95 MEDLoader::WriteField("file6.med",f1,true);
96 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldCell("file6.med",f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),0,1);
97 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
101 f1=buildVecFieldOnNodes_1();
102 MEDLoader::WriteField("file7.med",f1,true);
103 f2=MEDLoader::ReadFieldNode("file7.med",f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,3);
104 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
105 // testing kind message on error of field type.
106 CPPUNIT_ASSERT_THROW(MEDLoader::ReadFieldCell("file7.med",f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,3),INTERP_KERNEL::Exception);
113 * Multi field writing in a same file.
115 void MEDLoaderTest::testFieldRW2()
117 const char fileName[]="file8.med";
118 static const double VAL1=12345.67890314;
119 static const double VAL2=-1111111111111.;
120 MEDCouplingFieldDouble *f1=buildVecFieldOnCells_1();
121 MEDLoader::WriteField(fileName,f1,true);
122 f1->setTime(10.,8,9);
123 double *tmp=f1->getArray()->getPointer();
125 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
126 f1->setTime(10.14,18,19);
128 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
129 //retrieving time steps...
130 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),8,9);
131 f1->setTime(10.,8,9);
133 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
135 f2=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),0,1);
136 MEDCouplingFieldDouble *f3=buildVecFieldOnCells_1();
137 CPPUNIT_ASSERT(f3->isEqual(f2,1e-12,1e-12));
140 f2=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),18,19);
141 f1->setTime(10.14,18,19);
143 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
144 //test of throw on invalid (dt,it)
145 CPPUNIT_ASSERT_THROW(MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),28,19),INTERP_KERNEL::Exception);
149 f1=buildVecFieldOnNodes_1();
150 const char fileName2[]="file9.med";
151 MEDLoader::WriteField(fileName2,f1,true);
152 f1->setTime(110.,108,109);
153 tmp=f1->getArray()->getPointer();
155 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName2,f1);
156 f1->setTime(210.,208,209);
158 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName2,f1);
159 f2=MEDLoader::ReadFieldNode(fileName2,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),108,109);
160 f1->setTime(110.,108,109);
162 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
164 f2=MEDLoader::ReadFieldNode(fileName2,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,3);
165 f3=buildVecFieldOnNodes_1();
166 CPPUNIT_ASSERT(f3->isEqual(f2,1e-12,1e-12));
169 f2=MEDLoader::ReadFieldNode(fileName2,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),208,209);
170 f1->setTime(210.,208,209);
172 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
178 * Multi field in a same file, but this field has several
180 void MEDLoaderTest::testFieldRW3()
182 const char fileName[]="file11.med";
183 static const double VAL1=12345.67890314;
184 static const double VAL2=-1111111111111.;
185 const char name1[]="AField";
186 const char name3[]="AMesh1";
187 MEDCouplingFieldDouble *f1=buildVecFieldOnCells_1();
188 (const_cast<MEDCouplingMesh *>(f1->getMesh()))->setName(name3);
190 f1->setTime(10.,8,9);
191 double *tmp=f1->getArray()->getPointer();
193 MEDLoader::WriteField(fileName,f1,true);
194 f1->setTime(10.14,18,19);
196 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
197 f1->setTime(10.55,28,29);
199 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
200 f1->setTime(10.66,38,39);
202 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
203 f1->setTime(10.77,48,49);
205 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
208 f1=buildVecFieldOnNodes_1();
210 (const_cast<MEDCouplingMesh *>(f1->getMesh()))->setName(name3);
211 f1->setTime(110.,8,9);
212 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
213 f1->setTime(110.,108,109);
214 tmp=f1->getArray()->getPointer();
216 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
217 f1->setTime(210.,208,209);
219 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
221 std::vector< std::pair<int,int> > it1=MEDLoader::GetCellFieldIterations(fileName,name3,name1);
222 CPPUNIT_ASSERT_EQUAL(5,(int)it1.size());
223 CPPUNIT_ASSERT_EQUAL(8,it1[0].first); CPPUNIT_ASSERT_EQUAL(9,it1[0].second);
224 CPPUNIT_ASSERT_EQUAL(18,it1[1].first); CPPUNIT_ASSERT_EQUAL(19,it1[1].second);
225 CPPUNIT_ASSERT_EQUAL(28,it1[2].first); CPPUNIT_ASSERT_EQUAL(29,it1[2].second);
226 CPPUNIT_ASSERT_EQUAL(38,it1[3].first); CPPUNIT_ASSERT_EQUAL(39,it1[3].second);
227 CPPUNIT_ASSERT_EQUAL(48,it1[4].first); CPPUNIT_ASSERT_EQUAL(49,it1[4].second);
228 std::vector< std::pair<int,int> > it3=MEDLoader::GetNodeFieldIterations(fileName,name3,name1);
229 CPPUNIT_ASSERT_EQUAL(3,(int)it3.size());
230 CPPUNIT_ASSERT_EQUAL(8,it3[0].first); CPPUNIT_ASSERT_EQUAL(9,it3[0].second);
231 CPPUNIT_ASSERT_EQUAL(108,it3[1].first); CPPUNIT_ASSERT_EQUAL(109,it3[1].second);
232 CPPUNIT_ASSERT_EQUAL(208,it3[2].first); CPPUNIT_ASSERT_EQUAL(209,it3[2].second);
236 f1=MEDLoader::ReadFieldCell(fileName,name3,0,name1,8,9);
237 CPPUNIT_ASSERT_DOUBLES_EQUAL(VAL1,f1->getArray()->getConstPointer()[0],1e-13);
239 f1=MEDLoader::ReadFieldCell(fileName,name3,0,name1,18,19);
240 CPPUNIT_ASSERT_DOUBLES_EQUAL(VAL2,f1->getArray()->getConstPointer()[0],1e-13);
242 f1=MEDLoader::ReadFieldCell(fileName,name3,0,name1,28,29);
243 CPPUNIT_ASSERT_DOUBLES_EQUAL(3*VAL1,f1->getArray()->getConstPointer()[0],1e-13);
245 f1=MEDLoader::ReadFieldCell(fileName,name3,0,name1,38,39);
246 CPPUNIT_ASSERT_DOUBLES_EQUAL(3*VAL2,f1->getArray()->getConstPointer()[0],1e-13);
248 f1=MEDLoader::ReadFieldCell(fileName,name3,0,name1,48,49);
249 CPPUNIT_ASSERT_DOUBLES_EQUAL(4*VAL2,f1->getArray()->getConstPointer()[0],1e-13);
252 f1=MEDLoader::ReadFieldNode(fileName,name3,0,name1,8,9);
253 CPPUNIT_ASSERT_DOUBLES_EQUAL(71.,f1->getArray()->getConstPointer()[3],1e-13);
255 f1=MEDLoader::ReadFieldNode(fileName,name3,0,name1,108,109);
256 CPPUNIT_ASSERT_DOUBLES_EQUAL(VAL1,f1->getArray()->getConstPointer()[3],1e-13);
258 f1=MEDLoader::ReadFieldNode(fileName,name3,0,name1,208,209);
259 CPPUNIT_ASSERT_DOUBLES_EQUAL(VAL2,f1->getArray()->getConstPointer()[3],1e-13);
263 void MEDLoaderTest::testMultiMeshRW1()
265 const char fileName[]="file10.med";
266 MEDCouplingUMesh *mesh1=build3DMesh_1();
267 const int part1[5]={1,2,4,13,15};
268 MEDCouplingUMesh *mesh2=(MEDCouplingUMesh *)mesh1->buildPartOfMySelf(part1,part1+5,true);
269 mesh2->setName("mesh2");
270 const int part2[4]={3,4,13,14};
271 MEDCouplingUMesh *mesh3=(MEDCouplingUMesh *)mesh1->buildPartOfMySelf(part2,part2+4,true);
272 mesh3->setName("mesh3");
273 MEDCouplingUMesh *mesh4=MEDCouplingUMesh::New();
274 mesh4->setName("mesh4");
275 mesh4->setMeshDimension(3);
276 mesh4->allocateCells(1);
277 int conn[4]={0,11,1,3};
278 mesh4->insertNextCell(INTERP_KERNEL::NORM_TETRA4,4,conn);
279 mesh4->finishInsertingCells();
280 mesh4->setCoords(mesh1->getCoords());
281 std::vector<const MEDCouplingUMesh *> meshes;
282 meshes.push_back(mesh1);
283 meshes.push_back(mesh2);
284 meshes.push_back(mesh3);
285 meshes.push_back(mesh4);
286 const char mnane[]="3DToto";
287 MEDLoader::WriteUMeshesPartition(fileName,mnane,meshes,true);
289 MEDCouplingUMesh *mesh5=MEDLoader::ReadUMeshFromFile(fileName,mnane);
290 mesh1->setName(mnane);
291 const int part3[18]={0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17};
292 MEDCouplingUMesh *mesh6=(MEDCouplingUMesh *)mesh5->buildPartOfMySelf(part3,part3+18,true);
293 mesh6->setName(mnane);
295 CPPUNIT_ASSERT(mesh6->isEqual(mesh1,1e-12));
297 std::vector<std::string> grps=MEDLoader::GetMeshGroupsNames(fileName,mnane);
298 CPPUNIT_ASSERT_EQUAL(4,(int)grps.size());
299 CPPUNIT_ASSERT(std::find(grps.begin(),grps.end(),std::string("mesh2"))!=grps.end());
300 CPPUNIT_ASSERT(std::find(grps.begin(),grps.end(),std::string("mesh3"))!=grps.end());
301 CPPUNIT_ASSERT(std::find(grps.begin(),grps.end(),std::string("mesh4"))!=grps.end());
302 CPPUNIT_ASSERT(std::find(grps.begin(),grps.end(),std::string("3DMesh_1"))!=grps.end());
304 std::vector<std::string> vec;
305 vec.push_back(std::string("mesh2"));
306 MEDCouplingUMesh *mesh2_2=MEDLoader::ReadUMeshFromGroups(fileName,mnane,0,vec);
307 CPPUNIT_ASSERT(mesh2_2->isEqual(mesh2,1e-12));
309 vec.clear(); vec.push_back(std::string("mesh3"));
310 MEDCouplingUMesh *mesh3_2=MEDLoader::ReadUMeshFromGroups(fileName,mnane,0,vec);
311 CPPUNIT_ASSERT(mesh3_2->isEqual(mesh3,1e-12));
313 vec.clear(); vec.push_back(std::string("mesh4"));
314 MEDCouplingUMesh *mesh4_2=MEDLoader::ReadUMeshFromGroups(fileName,mnane,0,vec);
315 CPPUNIT_ASSERT(mesh4_2->isEqual(mesh4,1e-12));
317 vec.clear(); vec.push_back(std::string("3DMesh_1"));
318 MEDCouplingUMesh *mesh1_2=MEDLoader::ReadUMeshFromGroups(fileName,mnane,0,vec);
319 mesh1->setName("3DMesh_1");
320 CPPUNIT_ASSERT(mesh1_2->isEqual(mesh1,1e-12));
323 vec.clear(); vec.push_back(std::string("Family_-3")); vec.push_back(std::string("Family_-5"));
324 mesh2_2=MEDLoader::ReadUMeshFromFamilies(fileName,mnane,0,vec);
325 mesh2_2->setName("mesh2");
326 CPPUNIT_ASSERT(mesh2_2->isEqual(mesh2,1e-12));
329 std::vector<std::string> ret=MEDLoader::GetMeshFamiliesNamesOnGroup(fileName,"3DToto","3DMesh_1");
330 CPPUNIT_ASSERT_EQUAL(4,(int)ret.size());
331 CPPUNIT_ASSERT(ret[0]=="Family_-2");
332 CPPUNIT_ASSERT(ret[1]=="Family_-3");
333 CPPUNIT_ASSERT(ret[2]=="Family_-4");
334 CPPUNIT_ASSERT(ret[3]=="Family_-5");
336 std::vector<std::string> ret1=MEDLoader::GetMeshGroupsNamesOnFamily(fileName,"3DToto","Family_-3");
337 CPPUNIT_ASSERT_EQUAL(2,(int)ret1.size());
338 CPPUNIT_ASSERT(ret1[0]=="3DMesh_1");
339 CPPUNIT_ASSERT(ret1[1]=="mesh2");
347 void MEDLoaderTest::testFieldProfilRW1()
349 const char fileName[]="file12.med";
350 MEDCouplingUMesh *mesh1=build3DMesh_1();
353 DataArrayInt *da=mesh1->mergeNodes(1e-12,b,newNbOfNodes);
355 MEDLoader::WriteUMesh(fileName,mesh1,true);
356 const int part1[5]={1,2,4,13,15};
357 MEDCouplingUMesh *mesh2=(MEDCouplingUMesh *)mesh1->buildPartOfMySelf(part1,part1+5,true);
358 mesh2->setName(mesh1->getName().c_str());//<- important for the test
360 int nbOfCells=mesh2->getNumberOfCells();
361 CPPUNIT_ASSERT_EQUAL(5,nbOfCells);
362 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
363 f1->setName("VectorFieldOnCells");
365 DataArrayDouble *array=DataArrayDouble::New();
366 array->alloc(nbOfCells,2);
369 double *tmp=array->getPointer();
370 const double arr1[10]={71.,171.,10.,110.,20.,120.,30.,130.,40.,140.};
371 std::copy(arr1,arr1+10,tmp);
372 f1->setTime(3.14,2,7);
373 f1->checkCoherency();
375 MEDLoader::WriteField(fileName,f1,false);//<- false important for the test
377 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,7);
378 std::vector<ParaMEDMEM::TypeOfField> types=MEDLoader::GetTypesOfField(fileName,f1->getMesh()->getName().c_str(),f1->getName().c_str());
379 CPPUNIT_ASSERT_EQUAL(1,(int)types.size());
380 CPPUNIT_ASSERT(types[0]==ON_CELLS);
381 f2->checkCoherency();
382 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
391 * Test MED file profiles.
393 void MEDLoaderTest::testFieldNodeProfilRW1()
395 const char fileName[]="file19.med";
396 const char fileName2[]="file20.med";
397 MEDCouplingUMesh *m=build2DMesh_1();
398 int nbOfNodes=m->getNumberOfNodes();
399 MEDLoader::WriteUMesh(fileName,m,true);
400 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
401 f1->setName("VFieldOnNodes");
403 DataArrayDouble *array=DataArrayDouble::New();
404 const double arr1[24]={1.,101.,2.,102.,3.,103.,4.,104.,5.,105.,6.,106.,7.,107.,8.,108.,9.,109.,10.,110.,11.,111.,12.,112.};
405 array->alloc(nbOfNodes,2);
406 std::copy(arr1,arr1+24,array->getPointer());
408 array->setInfoOnComponent(0,"tyty [mm]");
409 array->setInfoOnComponent(1,"uiop [MW]");
411 f1->setTime(3.14,2,7);
412 f1->checkCoherency();
413 const int arr2[2]={1,4};//node ids are 2,4,5,3,6,7
414 MEDCouplingFieldDouble *f2=f1->buildSubPart(arr2,arr2+2);
415 (const_cast<MEDCouplingMesh *>(f2->getMesh()))->setName(f1->getMesh()->getName().c_str());
416 MEDLoader::WriteField(fileName,f2,false);//<- false important for the test
418 MEDCouplingFieldDouble *f3=MEDLoader::ReadFieldNode(fileName,f2->getMesh()->getName().c_str(),0,f2->getName().c_str(),2,7);
419 f3->checkCoherency();
420 CPPUNIT_ASSERT(f3->isEqual(f2,1e-12,1e-12));
423 const int arr3[6]={1,3,0,5,2,4};
424 f2->renumberNodes(arr3);
425 MEDLoader::WriteUMesh(fileName2,m,true);
426 MEDLoader::WriteField(fileName2,f2,false);//<- false important for the test
427 f3=MEDLoader::ReadFieldNode(fileName2,f2->getMesh()->getName().c_str(),0,f2->getName().c_str(),2,7);
428 f3->checkCoherency();
429 CPPUNIT_ASSERT(f3->isEqual(f2,1e-12,1e-12));
437 void MEDLoaderTest::testFieldNodeProfilRW2()
439 const char fileName[]="file23.med";
440 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
441 MEDLoader::WriteUMesh(fileName,mesh,true);
443 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
444 f1->setName("FieldMix");
446 const double arr2[24]={
447 1071.,1171.,1010.,1110.,1020.,1120.,1030.,1130.,1040.,1140.,1050.,1150.,
448 1060.,1160.,1070.,1170.,1080.,1180.,1090.,1190.,1091.,1191.,1092.,1192.
450 DataArrayDouble *array=DataArrayDouble::New();
453 array->setInfoOnComponent(0,"plkj [mm]");
454 array->setInfoOnComponent(1,"pqqqss [mm]");
456 double *tmp=array->getPointer();
457 std::copy(arr2,arr2+24,tmp);
458 f1->setTime(3.17,2,7);
460 const int renumArr[12]={3,7,2,1,5,11,10,0,9,6,8,4};
461 f1->renumberNodes(renumArr);
462 f1->checkCoherency();
463 MEDLoader::WriteField(fileName,f1,false);//<- false important for the test
464 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldNode(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,7);
465 CPPUNIT_ASSERT(f2->isEqual(f1,1e-12,1e-12));
472 void MEDLoaderTest::testFieldGaussRW1()
474 const char fileName[]="file13.med";
475 MEDCouplingFieldDouble *f1=buildVecFieldOnGauss_1();
476 MEDLoader::WriteField(fileName,f1,true);
477 MEDCouplingFieldDouble *f2=MEDLoader::ReadField(ON_GAUSS_PT,fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),1,5);
478 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
483 void MEDLoaderTest::testFieldGaussNERW1()
485 const char fileName[]="file14.med";
486 MEDCouplingFieldDouble *f1=buildVecFieldOnGaussNE_1();
487 MEDLoader::WriteField(fileName,f1,true);
488 std::vector<ParaMEDMEM::TypeOfField> tof(MEDLoader::GetTypesOfField(fileName,"2DMesh_2","MyFieldOnGaussNE"));
489 CPPUNIT_ASSERT_EQUAL(1,(int)tof.size());
490 CPPUNIT_ASSERT(ON_GAUSS_NE==tof[0]);
491 MEDCouplingFieldDouble *f2=MEDLoader::ReadField(ON_GAUSS_NE,fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),1,5);
492 CPPUNIT_ASSERT(f1->isEqual(f2,1e-12,1e-12));
497 void MEDLoaderTest::testLittleStrings1()
499 std::string s("azeeeerrrtty");
500 MEDLoaderBase::zipEqualConsChar(s,3);
501 CPPUNIT_ASSERT(s=="azertty");
504 void MEDLoaderTest::testSplitIntoNameAndUnit1()
506 std::string s(" []");
508 MEDLoaderBase::splitIntoNameAndUnit(s,c,u);
509 CPPUNIT_ASSERT(c.empty());
510 CPPUNIT_ASSERT(u.empty());
512 MEDLoaderBase::strip(s);
513 CPPUNIT_ASSERT(s=="lmmm kki jjj");
515 MEDLoaderBase::strip(s);
516 CPPUNIT_ASSERT(s.empty());
518 MEDLoaderBase::strip(s);
519 CPPUNIT_ASSERT(s.empty());
521 MEDLoaderBase::strip(s);
522 CPPUNIT_ASSERT(s.empty());
524 MEDLoaderBase::strip(s);
525 CPPUNIT_ASSERT(s=="pp");
528 void MEDLoaderTest::testMesh3DSurfShuffleRW()
530 const char fileName[]="file15.med";
531 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
532 const int renumber1[6]={2,5,1,0,3,4};
533 mesh->renumberCells(renumber1,false);
534 mesh->checkCoherency();
535 MEDLoader::WriteUMesh(fileName,mesh,true);
536 MEDCouplingUMesh *mesh_rw=MEDLoader::ReadUMeshFromFile(fileName,mesh->getName().c_str(),0);
537 CPPUNIT_ASSERT(mesh->isEqual(mesh_rw,1e-12));
542 void MEDLoaderTest::testFieldShuffleRW1()
544 const char fileName[]="file16.med";
545 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
546 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
547 f1->setName("FieldOnCellsShuffle");
549 DataArrayDouble *array=DataArrayDouble::New();
553 double *tmp=array->getPointer();
554 const double arr1[12]={71.,171.,10.,110.,20.,120.,30.,130.,40.,140.,50.,150.};
555 std::copy(arr1,arr1+12,tmp);
556 f1->setTime(3.14,2,7);
557 f1->checkCoherency();
559 const int renumber1[6]={2,1,5,0,3,4};
560 f1->renumberCells(renumber1,false);
561 MEDLoader::WriteField(fileName,f1,true);
562 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldCell(fileName,mesh->getName().c_str(),0,f1->getName().c_str(),2,7);
563 CPPUNIT_ASSERT(f2->isEqual(f1,1e-12,1e-12));
571 * Shuffle de cells but no profile. Like pointe.med
573 void MEDLoaderTest::testMultiFieldShuffleRW1()
575 const char fileName[]="file17.med";
576 MEDCouplingUMesh *m=build3DMesh_2();
577 CPPUNIT_ASSERT_EQUAL(20,m->getNumberOfCells());
578 CPPUNIT_ASSERT_EQUAL(45,m->getNumberOfNodes());
579 const int polys[3]={1,4,6};
580 std::vector<int> poly2(polys,polys+3);
581 m->convertToPolyTypes(&poly2[0],&poly2[0]+poly2.size());
582 const int renum[20]={1,3,2,8,9,12,13,16,19,0,4,7,5,15,14,17,10,18,6,11};
583 m->renumberCells(renum,false);
584 m->orientCorrectlyPolyhedrons();
586 MEDLoader::WriteUMesh(fileName,m,true);
587 MEDCouplingFieldDouble *f1Tmp=m->getMeasureField(false);
588 MEDCouplingFieldDouble *f1=f1Tmp->buildNewTimeReprFromThis(ONE_TIME,false);
591 MEDCouplingFieldDouble *f_1=f1->cloneWithMesh(true);
592 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
593 f1->applyFunc("2*x");
594 f1->setTime(0.01,3,4);
595 MEDCouplingFieldDouble *f_2=f1->cloneWithMesh(true);
596 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
597 f1->applyFunc("2*x/3");
598 f1->setTime(0.02,5,6);
599 MEDCouplingFieldDouble *f_3=f1->cloneWithMesh(true);
600 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
603 std::vector<std::pair<int,int> > its;
604 its.push_back(std::pair<int,int>(1,2));
605 its.push_back(std::pair<int,int>(3,4));
606 its.push_back(std::pair<int,int>(5,6));
607 std::vector<MEDCouplingFieldDouble *> fs=MEDLoader::ReadFieldsOnSameMesh(ON_CELLS,fileName,f_1->getMesh()->getName().c_str(),0,f_1->getName().c_str(),its);
608 CPPUNIT_ASSERT_EQUAL(3,(int)fs.size());
609 const MEDCouplingMesh *mm=fs[0]->getMesh();
610 CPPUNIT_ASSERT(fs[0]->isEqual(f_1,1e-12,1e-12));
611 CPPUNIT_ASSERT(fs[1]->isEqual(f_2,1e-12,1e-12));
612 CPPUNIT_ASSERT(fs[2]->isEqual(f_3,1e-12,1e-12));
613 CPPUNIT_ASSERT(mm==fs[1]->getMesh());// <- important for the test
614 CPPUNIT_ASSERT(mm==fs[2]->getMesh());// <- important for the test
615 for(std::vector<MEDCouplingFieldDouble *>::iterator iter=fs.begin();iter!=fs.end();iter++)
625 void MEDLoaderTest::testWriteUMeshesRW1()
627 const char fileName[]="file18.med";
628 MEDCouplingUMesh *m3d=build3DMesh_2();
629 const double pt[3]={0.,0.,-0.3};
630 const double vec[3]={0.,0.,1.};
631 std::vector<int> nodes;
632 m3d->findNodesOnPlane(pt,vec,1e-12,nodes);
633 MEDCouplingUMesh *m2d=(MEDCouplingUMesh *)m3d->buildFacePartOfMySelfNode(&nodes[0],&nodes[0]+nodes.size(),true);
634 const int renumber[5]={1,2,0,4,3};
635 m2d->renumberCells(renumber,false);
636 m2d->setName("ExampleOfMultiDimW");
637 std::vector<const MEDCouplingUMesh *> meshes;
638 meshes.push_back(m2d);
639 meshes.push_back(m3d);
640 MEDLoader::WriteUMeshes(fileName,meshes,true);
641 MEDCouplingUMesh *m3d_bis=MEDLoader::ReadUMeshFromFile(fileName,m2d->getName().c_str(),0);
642 CPPUNIT_ASSERT(!m3d_bis->isEqual(m3d,1e-12));
643 m3d_bis->setName(m3d->getName().c_str());
644 CPPUNIT_ASSERT(m3d_bis->isEqual(m3d,1e-12));
645 MEDCouplingUMesh *m2d_bis=MEDLoader::ReadUMeshFromFile(fileName,m2d->getName().c_str(),-1);//-1 for faces
646 CPPUNIT_ASSERT(m2d_bis->isEqual(m2d,1e-12));
647 // Creation of a field on faces.
648 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
649 f1->setName("FieldOnFacesShuffle");
651 DataArrayDouble *array=DataArrayDouble::New();
652 array->alloc(m2d->getNumberOfCells(),2);
653 array->setInfoOnComponent(0,"plkj [mm]");
654 array->setInfoOnComponent(1,"pqqqss [mm]");
657 double *tmp=array->getPointer();
658 const double arr1[10]={71.,171.,10.,110.,20.,120.,30.,130.,40.,140.};
659 std::copy(arr1,arr1+10,tmp);
660 f1->setTime(3.14,2,7);
661 f1->checkCoherency();
662 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
663 MEDCouplingFieldDouble *f2=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),-1,f1->getName().c_str(),2,7);
664 CPPUNIT_ASSERT(f2->isEqual(f1,1e-12,1e-12));
674 void MEDLoaderTest::testMixCellAndNodesFieldRW1()
676 const char fileName[]="file21.med";
677 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
678 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
679 f1->setName("FieldMix");
681 DataArrayDouble *array=DataArrayDouble::New();
684 array->setInfoOnComponent(0,"plkj [mm]");
685 array->setInfoOnComponent(1,"pqqqss [mm]");
687 double *tmp=array->getPointer();
688 const double arr1[12]={71.,171.,10.,110.,20.,120.,30.,130.,40.,140.,50.,150.};
689 std::copy(arr1,arr1+12,tmp);
690 f1->setTime(3.14,2,7);
691 f1->checkCoherency();
693 MEDCouplingFieldDouble *f2=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
694 f2->setName("FieldMix");
696 array=DataArrayDouble::New();
699 array->setInfoOnComponent(0,"plkj [mm]");
700 array->setInfoOnComponent(1,"pqqqss [mm]");
702 tmp=array->getPointer();
703 const double arr2[24]={
704 1071.,1171.,1010.,1110.,1020.,1120.,1030.,1130.,1040.,1140.,1050.,1150.,
705 1060.,1160.,1070.,1170.,1080.,1180.,1090.,1190.,1091.,1191.,1092.,1192.
707 std::copy(arr2,arr2+24,tmp);
708 f2->setTime(3.14,2,7);
709 f2->checkCoherency();
711 MEDLoader::WriteField(fileName,f1,true);
712 std::vector<ParaMEDMEM::TypeOfField> ts=MEDLoader::GetTypesOfField(fileName,f1->getMesh()->getName().c_str(),f1->getName().c_str());
713 CPPUNIT_ASSERT_EQUAL(1,(int)ts.size());
714 CPPUNIT_ASSERT_EQUAL(ON_CELLS,ts[0]);
715 std::vector<std::string> fs=MEDLoader::GetAllFieldNamesOnMesh(fileName,f1->getMesh()->getName().c_str());
716 CPPUNIT_ASSERT_EQUAL(1,(int)fs.size());
717 CPPUNIT_ASSERT(fs[0]=="FieldMix");
718 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f2);
719 fs=MEDLoader::GetAllFieldNamesOnMesh(fileName,f1->getMesh()->getName().c_str());
720 CPPUNIT_ASSERT_EQUAL(1,(int)fs.size());
721 CPPUNIT_ASSERT(fs[0]=="FieldMix");
723 ts=MEDLoader::GetTypesOfField(fileName,f1->getMesh()->getName().c_str(),f1->getName().c_str());
724 CPPUNIT_ASSERT_EQUAL(2,(int)ts.size());
725 CPPUNIT_ASSERT_EQUAL(ON_NODES,ts[0]);
726 CPPUNIT_ASSERT_EQUAL(ON_CELLS,ts[1]);
728 MEDCouplingFieldDouble *f3=MEDLoader::ReadFieldNode(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,7);
729 CPPUNIT_ASSERT(f3->isEqual(f2,1e-12,1e-12));
731 f3=MEDLoader::ReadFieldCell(fileName,f1->getMesh()->getName().c_str(),0,f1->getName().c_str(),2,7);
732 CPPUNIT_ASSERT(f3->isEqual(f1,1e-12,1e-12));
740 void MEDLoaderTest::testGetAllFieldNamesRW1()
742 const char fileName[]="file22.med";
743 MEDCouplingUMesh *mesh=build2DMesh_2();
744 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
745 f1->setName("Field1");
746 f1->setTime(3.44,5,6);
748 f1->fillFromAnalytic(2,"x+y");
749 MEDLoader::WriteField(fileName,f1,true);
750 f1->setTime(1002.3,7,8);
751 f1->fillFromAnalytic(2,"x+77.*y");
752 MEDLoader::WriteFieldUsingAlreadyWrittenMesh(fileName,f1);
753 f1->setName("Field2");
754 MEDLoader::WriteField(fileName,f1,false);
755 f1->setName("Field3");
756 mesh->setName("2DMesh_2Bis");
757 MEDLoader::WriteField(fileName,f1,false);
759 f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
760 f1->setName("Field8");
761 f1->setTime(8.99,7,9);
763 f1->fillFromAnalytic(3,"3*x+y");
764 MEDLoader::WriteField(fileName,f1,false);
766 std::vector<std::string> fs=MEDLoader::GetAllFieldNames(fileName);
767 CPPUNIT_ASSERT_EQUAL(4,(int)fs.size());
768 CPPUNIT_ASSERT(fs[0]=="Field1");
769 CPPUNIT_ASSERT(fs[1]=="Field2");
770 CPPUNIT_ASSERT(fs[2]=="Field3");
771 CPPUNIT_ASSERT(fs[3]=="Field8");
776 void MEDLoaderTest::testMEDLoaderRead1()
779 using namespace INTERP_KERNEL;
781 string fileName=getResourceFile("pointe.med");
782 vector<string> meshNames=MEDLoader::GetMeshNames(fileName.c_str());
783 CPPUNIT_ASSERT_EQUAL(1,(int)meshNames.size());
784 MEDCouplingUMesh *mesh=MEDLoader::ReadUMeshFromFile(fileName.c_str(),meshNames[0].c_str(),0);
785 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
786 CPPUNIT_ASSERT_EQUAL(3,mesh->getMeshDimension());
787 CPPUNIT_ASSERT_EQUAL(16,mesh->getNumberOfCells());
788 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
789 CPPUNIT_ASSERT_EQUAL(3,(int)mesh->getAllGeoTypes().size());
790 for(int i=0;i<12;i++)
791 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,mesh->getTypeOfCell(i));
792 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,mesh->getTypeOfCell(12));
793 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,mesh->getTypeOfCell(13));
794 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,mesh->getTypeOfCell(14));
795 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,mesh->getTypeOfCell(15));
796 CPPUNIT_ASSERT_EQUAL((std::size_t)90,mesh->getNodalConnectivity()->getNbOfElems());
797 CPPUNIT_ASSERT_EQUAL(701,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+90,0));
798 CPPUNIT_ASSERT_EQUAL(705,std::accumulate(mesh->getNodalConnectivityIndex()->getPointer(),mesh->getNodalConnectivityIndex()->getPointer()+17,0));
799 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+57,0),1e-12);
802 vector<string> families=MEDLoader::GetMeshFamiliesNames(fileName.c_str(),meshNames[0].c_str());
803 CPPUNIT_ASSERT_EQUAL(8,(int)families.size());
804 CPPUNIT_ASSERT(families[2]=="FAMILLE_ELEMENT_3");
806 vector<string> families2;
807 families2.push_back(families[2]);
808 mesh=MEDLoader::ReadUMeshFromFamilies(fileName.c_str(),meshNames[0].c_str(),0,families2);
809 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
810 CPPUNIT_ASSERT_EQUAL(3,mesh->getMeshDimension());
811 CPPUNIT_ASSERT_EQUAL(2,mesh->getNumberOfCells());
812 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
813 CPPUNIT_ASSERT_EQUAL(2,(int)mesh->getAllGeoTypes().size());
814 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,mesh->getTypeOfCell(0));
815 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,mesh->getTypeOfCell(1));
816 CPPUNIT_ASSERT_EQUAL((std::size_t)11,mesh->getNodalConnectivity()->getNbOfElems());
817 CPPUNIT_ASSERT_EQUAL(132,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+11,0));
818 CPPUNIT_ASSERT_EQUAL(16,std::accumulate(mesh->getNodalConnectivityIndex()->getPointer(),mesh->getNodalConnectivityIndex()->getPointer()+3,0));
819 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+57,0),1e-12);
822 vector<string> groups=MEDLoader::GetMeshGroupsNames(fileName.c_str(),meshNames[0].c_str());
823 CPPUNIT_ASSERT_EQUAL(5,(int)groups.size());
824 CPPUNIT_ASSERT(groups[0]=="groupe1");
825 CPPUNIT_ASSERT(groups[1]=="groupe2");
826 CPPUNIT_ASSERT(groups[2]=="groupe3");
827 CPPUNIT_ASSERT(groups[3]=="groupe4");
828 CPPUNIT_ASSERT(groups[4]=="groupe5");
829 vector<string> groups2;
830 groups2.push_back(groups[0]);
831 mesh=MEDLoader::ReadUMeshFromGroups(fileName.c_str(),meshNames[0].c_str(),0,groups2);
832 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
833 CPPUNIT_ASSERT_EQUAL(3,mesh->getMeshDimension());
834 CPPUNIT_ASSERT_EQUAL(7,mesh->getNumberOfCells());
835 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
836 CPPUNIT_ASSERT_EQUAL(2,(int)mesh->getAllGeoTypes().size());
838 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,mesh->getTypeOfCell(i));
839 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,mesh->getTypeOfCell(6));
840 CPPUNIT_ASSERT_EQUAL((std::size_t)36,mesh->getNodalConnectivity()->getNbOfElems());
841 CPPUNIT_ASSERT_EQUAL(254,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+36,0));
842 CPPUNIT_ASSERT_EQUAL(141,std::accumulate(mesh->getNodalConnectivityIndex()->getPointer(),mesh->getNodalConnectivityIndex()->getPointer()+8,0));
843 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+57,0),1e-12);
846 std::vector<std::string> fieldsName=MEDLoader::GetCellFieldNamesOnMesh(fileName.c_str(),meshNames[0].c_str());
847 CPPUNIT_ASSERT_EQUAL(2,(int)fieldsName.size());
848 CPPUNIT_ASSERT(fieldsName[0]=="fieldcelldoublescalar");
849 CPPUNIT_ASSERT(fieldsName[1]=="fieldcelldoublevector");
850 std::vector<std::pair<int,int> > its0=MEDLoader::GetCellFieldIterations(fileName.c_str(),meshNames[0].c_str(),fieldsName[0].c_str());
851 CPPUNIT_ASSERT_EQUAL(1,(int)its0.size());
852 CPPUNIT_ASSERT_EQUAL(-1,its0[0].first);
853 CPPUNIT_ASSERT_EQUAL(-1,its0[0].second);
854 std::vector<std::pair<int,int> > its1=MEDLoader::GetCellFieldIterations(fileName.c_str(),meshNames[0].c_str(),fieldsName[1].c_str());
855 CPPUNIT_ASSERT_EQUAL(1,(int)its1.size());
856 CPPUNIT_ASSERT_EQUAL(-1,its1[0].first);
857 CPPUNIT_ASSERT_EQUAL(-1,its1[0].second);
859 MEDCouplingFieldDouble *field0=MEDLoader::ReadFieldCell(fileName.c_str(),meshNames[0].c_str(),0,fieldsName[0].c_str(),its0[0].first,its0[0].second);
860 field0->checkCoherency();
861 CPPUNIT_ASSERT(field0->getName()==fieldsName[0]);
862 CPPUNIT_ASSERT_EQUAL(1,field0->getNumberOfComponents());
863 CPPUNIT_ASSERT_EQUAL(16,field0->getNumberOfTuples());
864 const double expectedValues[16]={1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,2.,3.,3.,2.};
865 double diffValue[16];
866 std::transform(field0->getArray()->getPointer(),field0->getArray()->getPointer()+16,expectedValues,diffValue,std::minus<double>());
867 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue,diffValue+16),1e-12);
868 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue,diffValue+16),1e-12);
869 const MEDCouplingUMesh *constMesh=dynamic_cast<const MEDCouplingUMesh *>(field0->getMesh());
870 CPPUNIT_ASSERT(constMesh);
871 CPPUNIT_ASSERT_EQUAL(3,constMesh->getSpaceDimension());
872 CPPUNIT_ASSERT_EQUAL(3,constMesh->getMeshDimension());
873 CPPUNIT_ASSERT_EQUAL(16,constMesh->getNumberOfCells());
874 CPPUNIT_ASSERT_EQUAL(19,constMesh->getNumberOfNodes());
875 CPPUNIT_ASSERT_EQUAL(3,(int)constMesh->getAllGeoTypes().size());
876 for(int i=0;i<12;i++)
877 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,constMesh->getTypeOfCell(i));
878 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(12));
879 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(13));
880 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(14));
881 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(15));
882 CPPUNIT_ASSERT_EQUAL((std::size_t)90,constMesh->getNodalConnectivity()->getNbOfElems());
883 CPPUNIT_ASSERT_EQUAL(701,std::accumulate(constMesh->getNodalConnectivity()->getConstPointer(),constMesh->getNodalConnectivity()->getConstPointer()+90,0));
884 CPPUNIT_ASSERT_EQUAL(705,std::accumulate(constMesh->getNodalConnectivityIndex()->getConstPointer(),constMesh->getNodalConnectivityIndex()->getConstPointer()+17,0));
885 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+57,0),1e-12);
888 MEDCouplingFieldDouble *field1=MEDLoader::ReadFieldCell(fileName.c_str(),meshNames[0].c_str(),0,fieldsName[1].c_str(),its1[0].first,its1[0].second);
889 field1->checkCoherency();
890 CPPUNIT_ASSERT(field1->getName()==fieldsName[1]);
891 CPPUNIT_ASSERT_EQUAL(3,field1->getNumberOfComponents());
892 CPPUNIT_ASSERT_EQUAL(16,field1->getNumberOfTuples());
893 const double expectedValues2[48]={1.,0.,1.,1.,0.,1.,1.,0.,1.,2.,1.,0.,2.,1.,0.,2.,1.,0.,3.,0.,1.,3.,0.,1.,3.,0.,1.,4.,1.,0.,4.,1.,0.,4.,1.,0.,5.,0.,0.,6.,1.,1.,6.,0.,0.,5.,1.,1.};
894 double diffValue2[48];
895 std::transform(field1->getArray()->getPointer(),field1->getArray()->getPointer()+48,expectedValues2,diffValue2,std::minus<double>());
896 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue2,diffValue2+48),1e-12);
897 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue2,diffValue2+48),1e-12);
898 constMesh=dynamic_cast<const MEDCouplingUMesh *>(field1->getMesh());
899 CPPUNIT_ASSERT(constMesh);
900 CPPUNIT_ASSERT_EQUAL(3,constMesh->getSpaceDimension());
901 CPPUNIT_ASSERT_EQUAL(3,constMesh->getMeshDimension());
902 CPPUNIT_ASSERT_EQUAL(16,constMesh->getNumberOfCells());
903 CPPUNIT_ASSERT_EQUAL(19,constMesh->getNumberOfNodes());
904 CPPUNIT_ASSERT_EQUAL(3,(int)constMesh->getAllGeoTypes().size());
905 for(int i=0;i<12;i++)
906 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,constMesh->getTypeOfCell(i));
907 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(12));
908 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(13));
909 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(14));
910 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(15));
911 CPPUNIT_ASSERT_EQUAL((std::size_t)90,constMesh->getNodalConnectivity()->getNbOfElems());
912 CPPUNIT_ASSERT_EQUAL(701,std::accumulate(constMesh->getNodalConnectivity()->getConstPointer(),constMesh->getNodalConnectivity()->getConstPointer()+90,0));
913 CPPUNIT_ASSERT_EQUAL(705,std::accumulate(constMesh->getNodalConnectivityIndex()->getConstPointer(),constMesh->getNodalConnectivityIndex()->getConstPointer()+17,0));
914 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+57,0),1e-12);
917 std::vector<std::string> fieldsNameNode=MEDLoader::GetNodeFieldNamesOnMesh(fileName.c_str(),meshNames[0].c_str());
918 CPPUNIT_ASSERT_EQUAL(2,(int)fieldsNameNode.size());
919 CPPUNIT_ASSERT(fieldsNameNode[0]=="fieldnodedouble");
920 CPPUNIT_ASSERT(fieldsNameNode[1]=="fieldnodeint");
921 std::vector<std::pair<int,int> > its0Node=MEDLoader::GetNodeFieldIterations(fileName.c_str(),meshNames[0].c_str(),fieldsNameNode[0].c_str());
922 CPPUNIT_ASSERT_EQUAL(3,(int)its0Node.size());
923 CPPUNIT_ASSERT_EQUAL(-1,its0Node[0].first);
924 CPPUNIT_ASSERT_EQUAL(-1,its0Node[0].second);
925 CPPUNIT_ASSERT_EQUAL(1,its0Node[1].first);
926 CPPUNIT_ASSERT_EQUAL(-1,its0Node[1].second);
927 CPPUNIT_ASSERT_EQUAL(2,its0Node[2].first);
928 CPPUNIT_ASSERT_EQUAL(-1,its0Node[2].second);
929 MEDCouplingFieldDouble *field0Nodes=MEDLoader::ReadFieldNode(fileName.c_str(),meshNames[0].c_str(),0,fieldsNameNode[0].c_str(),its0Node[0].first,its0Node[0].second);
930 field0Nodes->checkCoherency();
931 CPPUNIT_ASSERT(field0Nodes->getName()==fieldsNameNode[0]);
932 CPPUNIT_ASSERT_EQUAL(1,field0Nodes->getNumberOfComponents());
933 CPPUNIT_ASSERT_EQUAL(19,field0Nodes->getNumberOfTuples());
934 const double expectedValues3[19]={1.,1.,1.,2.,2.,2.,3.,3.,3.,4.,4.,4.,5.,5.,5.,6.,6.,6.,7.};
935 double diffValue3[19];
936 std::transform(field0Nodes->getArray()->getPointer(),field0Nodes->getArray()->getPointer()+19,expectedValues3,diffValue3,std::minus<double>());
937 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue3,diffValue3+19),1e-12);
938 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue3,diffValue3+19),1e-12);
939 constMesh=dynamic_cast<const MEDCouplingUMesh *>(field0Nodes->getMesh());
940 CPPUNIT_ASSERT(constMesh);
941 field0Nodes->decrRef();
943 field0Nodes=MEDLoader::ReadFieldNode(fileName.c_str(),meshNames[0].c_str(),0,fieldsNameNode[0].c_str(),its0Node[2].first,its0Node[2].second);
944 field0Nodes->checkCoherency();
945 CPPUNIT_ASSERT(field0Nodes->getName()==fieldsNameNode[0]);
946 CPPUNIT_ASSERT_EQUAL(1,field0Nodes->getNumberOfComponents());
947 CPPUNIT_ASSERT_EQUAL(19,field0Nodes->getNumberOfTuples());
948 const double expectedValues4[19]={1.,2.,2.,2.,3.,3.,3.,4.,4.,4.,5.,5.,5.,6.,6.,6.,7.,7.,7.};
949 std::transform(field0Nodes->getArray()->getPointer(),field0Nodes->getArray()->getPointer()+19,expectedValues4,diffValue3,std::minus<double>());
950 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue3,diffValue3+19),1e-12);
951 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue3,diffValue3+19),1e-12);
952 constMesh=dynamic_cast<const MEDCouplingUMesh *>(field0Nodes->getMesh());
953 CPPUNIT_ASSERT(constMesh);
954 CPPUNIT_ASSERT_EQUAL(3,constMesh->getSpaceDimension());
955 CPPUNIT_ASSERT_EQUAL(3,constMesh->getMeshDimension());
956 CPPUNIT_ASSERT_EQUAL(16,constMesh->getNumberOfCells());
957 CPPUNIT_ASSERT_EQUAL(19,constMesh->getNumberOfNodes());
958 CPPUNIT_ASSERT_EQUAL(3,(int)constMesh->getAllGeoTypes().size());
959 for(int i=0;i<12;i++)
960 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,constMesh->getTypeOfCell(i));
961 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(12));
962 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(13));
963 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(14));
964 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(15));
965 CPPUNIT_ASSERT_EQUAL((std::size_t)90,constMesh->getNodalConnectivity()->getNbOfElems());
966 CPPUNIT_ASSERT_EQUAL(701,std::accumulate(constMesh->getNodalConnectivity()->getConstPointer(),constMesh->getNodalConnectivity()->getConstPointer()+90,0));
967 CPPUNIT_ASSERT_EQUAL(705,std::accumulate(constMesh->getNodalConnectivityIndex()->getConstPointer(),constMesh->getNodalConnectivityIndex()->getConstPointer()+17,0));
968 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+57,0),1e-12);
969 field0Nodes->decrRef();
971 field0Nodes=MEDLoader::ReadFieldNode(fileName.c_str(),meshNames[0].c_str(),0,fieldsNameNode[0].c_str(),its0Node[0].first,its0Node[0].second);
972 field0Nodes->checkCoherency();
973 CPPUNIT_ASSERT(field0Nodes->getName()==fieldsNameNode[0]);
974 CPPUNIT_ASSERT_EQUAL(1,field0Nodes->getNumberOfComponents());
975 CPPUNIT_ASSERT_EQUAL(19,field0Nodes->getNumberOfTuples());
976 const double expectedValues5[19]={1.,1.,1.,2.,2.,2.,3.,3.,3.,4.,4.,4.,5.,5.,5.,6.,6.,6.,7.};
977 std::transform(field0Nodes->getArray()->getPointer(),field0Nodes->getArray()->getPointer()+19,expectedValues5,diffValue3,std::minus<double>());
978 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue3,diffValue3+19),1e-12);
979 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue3,diffValue3+19),1e-12);
980 constMesh=dynamic_cast<const MEDCouplingUMesh *>(field0Nodes->getMesh());
981 CPPUNIT_ASSERT(constMesh);
982 CPPUNIT_ASSERT_EQUAL(3,constMesh->getSpaceDimension());
983 CPPUNIT_ASSERT_EQUAL(3,constMesh->getMeshDimension());
984 CPPUNIT_ASSERT_EQUAL(16,constMesh->getNumberOfCells());
985 CPPUNIT_ASSERT_EQUAL(19,constMesh->getNumberOfNodes());
986 CPPUNIT_ASSERT_EQUAL(3,(int)constMesh->getAllGeoTypes().size());
987 for(int i=0;i<12;i++)
988 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,constMesh->getTypeOfCell(i));
989 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(12));
990 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(13));
991 CPPUNIT_ASSERT_EQUAL(NORM_HEXA8,constMesh->getTypeOfCell(14));
992 CPPUNIT_ASSERT_EQUAL(NORM_PYRA5,constMesh->getTypeOfCell(15));
993 CPPUNIT_ASSERT_EQUAL((std::size_t)90,constMesh->getNodalConnectivity()->getNbOfElems());
994 CPPUNIT_ASSERT_EQUAL(701,std::accumulate(constMesh->getNodalConnectivity()->getConstPointer(),constMesh->getNodalConnectivity()->getConstPointer()+90,0));
995 CPPUNIT_ASSERT_EQUAL(705,std::accumulate(constMesh->getNodalConnectivityIndex()->getConstPointer(),constMesh->getNodalConnectivityIndex()->getConstPointer()+17,0));
996 CPPUNIT_ASSERT_DOUBLES_EQUAL(46.,std::accumulate(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+57,0),1e-12);
997 field0Nodes->decrRef();
1000 void MEDLoaderTest::testMEDLoaderPolygonRead()
1002 using namespace std;
1003 using namespace INTERP_KERNEL;
1005 string fileName=getResourceFile("polygones.med");
1006 vector<string> meshNames=MEDLoader::GetMeshNames(fileName.c_str());
1007 CPPUNIT_ASSERT_EQUAL(1,(int)meshNames.size());
1008 CPPUNIT_ASSERT(meshNames[0]=="Bord");
1009 MEDCouplingUMesh *mesh=MEDLoader::ReadUMeshFromFile(fileName.c_str(),meshNames[0].c_str(),0);
1010 mesh->checkCoherency();
1011 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
1012 CPPUNIT_ASSERT_EQUAL(2,mesh->getMeshDimension());
1013 CPPUNIT_ASSERT_EQUAL(538,mesh->getNumberOfCells());
1014 CPPUNIT_ASSERT_EQUAL(579,mesh->getNumberOfNodes());
1015 CPPUNIT_ASSERT_EQUAL(2,(int)mesh->getAllGeoTypes().size());
1016 for(int i=0;i<514;i++)
1017 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(i));
1018 for(int i=514;i<538;i++)
1019 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,mesh->getTypeOfCell(i));
1020 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+1737,0),1e-12);
1021 const double expectedVals1[12]={1.4851585216522212,-0.5,0.,1.4851585216522212,-0.4,0.,1.4851585216522212,-0.3,0., 1.5741585216522211, -0.5, 0. };
1022 double diffValue1[12];
1023 std::transform(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+12,expectedVals1,diffValue1,std::minus<double>());
1024 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue1,diffValue1+12),1e-12);
1025 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue1,diffValue1+12),1e-12);
1026 CPPUNIT_ASSERT_EQUAL((std::size_t)2768,mesh->getNodalConnectivity()->getNbOfElems());
1027 CPPUNIT_ASSERT_EQUAL(651050,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+2768,0));
1028 CPPUNIT_ASSERT_EQUAL(725943,std::accumulate(mesh->getNodalConnectivityIndex()->getPointer(),mesh->getNodalConnectivityIndex()->getPointer()+539,0));
1031 std::vector<std::string> fieldsName=MEDLoader::GetCellFieldNamesOnMesh(fileName.c_str(),meshNames[0].c_str());
1032 CPPUNIT_ASSERT_EQUAL(3,(int)fieldsName.size());
1033 CPPUNIT_ASSERT(fieldsName[0]=="bord_:_distorsion");
1034 CPPUNIT_ASSERT(fieldsName[1]=="bord_:_familles");
1035 CPPUNIT_ASSERT(fieldsName[2]=="bord_:_non-ortho");
1036 std::vector<std::pair<int,int> > its0=MEDLoader::GetCellFieldIterations(fileName.c_str(),meshNames[0].c_str(),fieldsName[0].c_str());
1037 CPPUNIT_ASSERT_EQUAL(1,(int)its0.size());
1038 MEDCouplingFieldDouble *field=MEDLoader::ReadFieldCell(fileName.c_str(),meshNames[0].c_str(),0,fieldsName[0].c_str(),its0[0].first,its0[0].second);
1039 field->checkCoherency();
1040 CPPUNIT_ASSERT(field->getName()==fieldsName[0]);
1041 CPPUNIT_ASSERT_EQUAL(1,field->getNumberOfComponents());
1042 CPPUNIT_ASSERT_EQUAL(538,field->getNumberOfTuples());
1043 const MEDCouplingUMesh *constMesh=dynamic_cast<const MEDCouplingUMesh *>(field->getMesh());
1044 CPPUNIT_ASSERT(constMesh);
1045 CPPUNIT_ASSERT_EQUAL(3,constMesh->getSpaceDimension());
1046 CPPUNIT_ASSERT_EQUAL(2,constMesh->getMeshDimension());
1047 CPPUNIT_ASSERT_EQUAL(538,constMesh->getNumberOfCells());
1048 CPPUNIT_ASSERT_EQUAL(579,constMesh->getNumberOfNodes());
1049 CPPUNIT_ASSERT_EQUAL(2,(int)constMesh->getAllGeoTypes().size());
1050 for(int i=0;i<514;i++)
1051 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,constMesh->getTypeOfCell(i));
1052 for(int i=514;i<538;i++)
1053 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,constMesh->getTypeOfCell(i));
1054 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,std::accumulate(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+1737,0),1e-12);
1055 std::transform(constMesh->getCoords()->getConstPointer(),constMesh->getCoords()->getConstPointer()+12,expectedVals1,diffValue1,std::minus<double>());
1056 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::max_element(diffValue1,diffValue1+12),1e-12);
1057 CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,*std::min_element(diffValue1,diffValue1+12),1e-12);
1058 CPPUNIT_ASSERT_EQUAL((std::size_t)2768,constMesh->getNodalConnectivity()->getNbOfElems());
1059 CPPUNIT_ASSERT_EQUAL(651050,std::accumulate(constMesh->getNodalConnectivity()->getConstPointer(),constMesh->getNodalConnectivity()->getConstPointer()+2768,0));
1060 CPPUNIT_ASSERT_EQUAL(725943,std::accumulate(constMesh->getNodalConnectivityIndex()->getConstPointer(),constMesh->getNodalConnectivityIndex()->getConstPointer()+539,0));
1061 const double *values=field->getArray()->getPointer();
1062 CPPUNIT_ASSERT_DOUBLES_EQUAL(2.87214203182918,std::accumulate(values,values+538,0.),1e-12);
1066 void MEDLoaderTest::testMEDLoaderPolyhedronRead()
1068 using namespace std;
1069 using namespace INTERP_KERNEL;
1071 string fileName=getResourceFile("poly3D.med");
1072 vector<string> meshNames=MEDLoader::GetMeshNames(fileName.c_str());
1073 CPPUNIT_ASSERT_EQUAL(1,(int)meshNames.size());
1074 CPPUNIT_ASSERT(meshNames[0]=="poly3D");
1075 MEDCouplingUMesh *mesh=MEDLoader::ReadUMeshFromFile(fileName.c_str(),meshNames[0].c_str(),0);
1076 mesh->checkCoherency();
1077 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
1078 CPPUNIT_ASSERT_EQUAL(3,mesh->getMeshDimension());
1079 CPPUNIT_ASSERT_EQUAL(3,mesh->getNumberOfCells());
1080 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
1081 CPPUNIT_ASSERT_EQUAL(2,(int)mesh->getAllGeoTypes().size());
1082 CPPUNIT_ASSERT_EQUAL(NORM_TETRA4,mesh->getTypeOfCell(0));
1083 CPPUNIT_ASSERT_EQUAL(NORM_POLYHED,mesh->getTypeOfCell(1));
1084 CPPUNIT_ASSERT_EQUAL(NORM_POLYHED,mesh->getTypeOfCell(2));
1085 CPPUNIT_ASSERT_EQUAL((std::size_t)98,mesh->getNodalConnectivity()->getNbOfElems());
1086 CPPUNIT_ASSERT_EQUAL(725,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+98,0));
1087 CPPUNIT_ASSERT_DOUBLES_EQUAL(110.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+57,0),1e-12);
1088 CPPUNIT_ASSERT_EQUAL(155,std::accumulate(mesh->getNodalConnectivityIndex()->getPointer(),mesh->getNodalConnectivityIndex()->getPointer()+4,0));
1091 mesh=MEDLoader::ReadUMeshFromFile(fileName.c_str(),meshNames[0].c_str(),-1);
1092 mesh->checkCoherency();
1093 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
1094 CPPUNIT_ASSERT_EQUAL(2,mesh->getMeshDimension());
1095 CPPUNIT_ASSERT_EQUAL(17,mesh->getNumberOfCells());
1096 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
1097 CPPUNIT_ASSERT_EQUAL(3,(int)mesh->getAllGeoTypes().size());
1098 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,mesh->getTypeOfCell(0));
1099 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(1));
1100 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(2));
1101 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(3));
1102 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(4));
1103 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(5));
1104 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(6));
1105 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(7));
1106 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,mesh->getTypeOfCell(8));
1107 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,mesh->getTypeOfCell(9));
1108 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(10));
1109 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(11));
1110 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(12));
1111 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(13));
1112 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(14));
1113 CPPUNIT_ASSERT_EQUAL(NORM_QUAD4,mesh->getTypeOfCell(15));
1114 CPPUNIT_ASSERT_EQUAL(NORM_TRI3,mesh->getTypeOfCell(16));
1115 CPPUNIT_ASSERT_DOUBLES_EQUAL(110.,std::accumulate(mesh->getCoords()->getPointer(),mesh->getCoords()->getPointer()+57,0),1e-12);
1116 CPPUNIT_ASSERT_EQUAL((std::size_t)83,mesh->getNodalConnectivity()->getNbOfElems());
1117 CPPUNIT_ASSERT_EQUAL(619,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+83,0));
1120 vector<string> families=MEDLoader::GetMeshFamiliesNames(fileName.c_str(),meshNames[0].c_str());
1121 CPPUNIT_ASSERT_EQUAL(4,(int)families.size());
1122 CPPUNIT_ASSERT(families[0]=="FAMILLE_FACE_POLYGONS3");
1123 CPPUNIT_ASSERT(families[1]=="FAMILLE_FACE_QUAD41");
1124 CPPUNIT_ASSERT(families[2]=="FAMILLE_FACE_TRIA32");
1125 CPPUNIT_ASSERT(families[3]=="FAMILLE_ZERO");
1126 vector<string> families2;
1127 families2.push_back(families[0]);
1128 mesh=MEDLoader::ReadUMeshFromFamilies(fileName.c_str(),meshNames[0].c_str(),-1,families2);
1129 mesh->checkCoherency();
1130 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
1131 CPPUNIT_ASSERT_EQUAL(2,mesh->getMeshDimension());
1132 CPPUNIT_ASSERT_EQUAL(3,mesh->getNumberOfCells());
1133 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
1134 CPPUNIT_ASSERT_EQUAL(1,(int)mesh->getAllGeoTypes().size());
1135 for(int i=0;i<3;i++)
1136 CPPUNIT_ASSERT_EQUAL(NORM_POLYGON,mesh->getTypeOfCell(i));
1137 CPPUNIT_ASSERT_EQUAL((std::size_t)19,mesh->getNodalConnectivity()->getNbOfElems());
1138 CPPUNIT_ASSERT_EQUAL(117,std::accumulate(mesh->getNodalConnectivity()->getPointer(),mesh->getNodalConnectivity()->getPointer()+19,0));
1141 mesh=MEDLoader::ReadUMeshFromFamilies(fileName.c_str(),meshNames[0].c_str(),0,families2);
1142 CPPUNIT_ASSERT_EQUAL(3,mesh->getSpaceDimension());
1143 CPPUNIT_ASSERT_EQUAL(0,mesh->getNumberOfCells());
1144 CPPUNIT_ASSERT_EQUAL(19,mesh->getNumberOfNodes());
1145 CPPUNIT_ASSERT_EQUAL(3,mesh->getMeshDimension());
1146 CPPUNIT_ASSERT_EQUAL(0,(int)mesh->getAllGeoTypes().size());
1150 std::string MEDLoaderTest::getResourceFile( const std::string& filename ) const
1152 std::string resourceFile = "";
1154 if ( getenv("MEDCOUPLING_ROOT_DIR") ) {
1155 // use MEDCOUPLING_ROOT_DIR env.var
1156 resourceFile = getenv("MEDCOUPLING_ROOT_DIR");
1157 resourceFile += "/share/resources/med/";
1160 resourceFile = get_current_dir_name();
1161 resourceFile += "/../../../resources/";
1164 resourceFile += filename;
1165 return resourceFile;
1169 MEDCouplingUMesh *MEDLoaderTest::build1DMesh_1()
1171 double coords[6]={ 0.0, 0.3, 0.75, 1.0, 1.4, 1.3 };
1172 int conn[9]={ 0,1, 1,2, 2,3 , 3,4,5};
1173 MEDCouplingUMesh *mesh=MEDCouplingUMesh::New();
1174 mesh->setName("1DMesh_1");
1175 mesh->setMeshDimension(1);
1176 mesh->allocateCells(4);
1177 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
1178 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
1179 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+4);
1180 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG3,3,conn+6);
1181 mesh->finishInsertingCells();
1182 DataArrayDouble *myCoords=DataArrayDouble::New();
1183 myCoords->alloc(6,1);
1184 myCoords->setInfoOnComponent(0,"tototototototot [m*m*m*m*m*m*m*m]");
1185 std::copy(coords,coords+6,myCoords->getPointer());
1186 mesh->setCoords(myCoords);
1187 myCoords->decrRef();
1191 MEDCouplingUMesh *MEDLoaderTest::build2DCurveMesh_1()
1193 double coords[12]={ 0.0,0.0, 0.3,0.3, 0.75,0.75, 1.0,1.0, 1.4,1.4, 1.3,1.3 };
1194 int conn[9]={ 0,1, 1,2, 2,3 , 3,4,5};
1195 MEDCouplingUMesh *mesh=MEDCouplingUMesh::New();
1196 mesh->setName("2DCurveMesh_1");
1197 mesh->setMeshDimension(1);
1198 mesh->allocateCells(4);
1199 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
1200 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
1201 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+4);
1202 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG3,3,conn+6);
1203 mesh->finishInsertingCells();
1204 DataArrayDouble *myCoords=DataArrayDouble::New();
1205 myCoords->alloc(6,2);
1206 std::copy(coords,coords+12,myCoords->getPointer());
1207 mesh->setCoords(myCoords);
1208 myCoords->decrRef();
1212 MEDCouplingUMesh *MEDLoaderTest::build2DMesh_1()
1214 double targetCoords[24]={-0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2, 0.7,0.2, -0.3,0.7, 0.2,0.7, 0.7,0.7, -0.05,0.95, 0.2,1.2, 0.45,0.95 };
1215 int targetConn[24]={1,4,2, 4,5,2, 6,10,8,9,11,7, 0,3,4,1, 6,7,4,3, 7,8,5,4};
1216 MEDCouplingUMesh *targetMesh=MEDCouplingUMesh::New();
1217 targetMesh->setMeshDimension(2);
1218 targetMesh->allocateCells(6);
1219 targetMesh->setName("2DMesh_1");
1220 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1221 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
1222 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI6,6,targetConn+6);
1223 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+12);
1224 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+16);
1225 targetMesh->insertNextCell(INTERP_KERNEL::NORM_POLYGON,4,targetConn+20);
1226 targetMesh->finishInsertingCells();
1227 DataArrayDouble *myCoords=DataArrayDouble::New();
1228 myCoords->alloc(12,2);
1229 myCoords->setInfoOnComponent(0,"tototototototot [m]");
1230 myCoords->setInfoOnComponent(1,"energie [kW]");
1231 std::copy(targetCoords,targetCoords+24,myCoords->getPointer());
1232 targetMesh->setCoords(myCoords);
1233 myCoords->decrRef();
1237 MEDCouplingUMesh *MEDLoaderTest::build2DMesh_2()
1239 double targetCoords[24]={
1240 -0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2, 0.7,0.2, -0.3,0.7, 0.2,0.7, 0.7,0.7,
1241 -0.05,0.95, 0.2,1.2, 0.45,0.95
1243 int targetConn[24]={1,4,2, 4,5,2, 6,10,8,9,11,7, 0,3,4,1, 6,7,4,3, 7,8,5,4};
1244 MEDCouplingUMesh *targetMesh=MEDCouplingUMesh::New();
1245 targetMesh->setMeshDimension(2);
1246 targetMesh->allocateCells(5);
1247 targetMesh->setName("2DMesh_2");
1248 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1249 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
1250 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+12);
1251 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+16);
1252 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI6,6,targetConn+6);
1253 targetMesh->finishInsertingCells();
1254 DataArrayDouble *myCoords=DataArrayDouble::New();
1255 myCoords->alloc(12,2);
1256 myCoords->setInfoOnComponent(0,"toto [m]");
1257 myCoords->setInfoOnComponent(1,"energie [kW]");
1258 std::copy(targetCoords,targetCoords+24,myCoords->getPointer());
1259 targetMesh->setCoords(myCoords);
1260 myCoords->decrRef();
1264 MEDCouplingUMesh *MEDLoaderTest::build3DSurfMesh_1()
1266 double targetCoords[36]={
1267 -0.3,-0.3,-0.3, 0.2,-0.3,-0.3, 0.7,-0.3,-0.3, -0.3,0.2,-0.3, 0.2,0.2,-0.3, 0.7,0.2,-0.3, -0.3,0.7,-0.3, 0.2,0.7,-0.3, 0.7,0.7,-0.3
1268 ,-0.05,0.95,-0.3, 0.2,1.2,-0.3, 0.45,0.95,-0.3
1270 int targetConn[24]={1,4,2, 4,5,2, 6,10,8,9,11,7, 0,3,4,1, 6,7,4,3, 7,8,5,4};
1271 MEDCouplingUMesh *targetMesh=MEDCouplingUMesh::New();
1272 targetMesh->setMeshDimension(2);
1273 targetMesh->allocateCells(6);
1274 targetMesh->setName("3DSurfMesh_1");
1275 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1276 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
1277 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+12);
1278 targetMesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+16);
1279 targetMesh->insertNextCell(INTERP_KERNEL::NORM_TRI6,6,targetConn+6);
1280 targetMesh->insertNextCell(INTERP_KERNEL::NORM_POLYGON,4,targetConn+20);
1281 targetMesh->finishInsertingCells();
1282 DataArrayDouble *myCoords=DataArrayDouble::New();
1283 myCoords->alloc(12,3);
1284 myCoords->setInfoOnComponent(0,"toto [m]");
1285 myCoords->setInfoOnComponent(2,"ff [km]");//component 1 is not set for test
1286 std::copy(targetCoords,targetCoords+36,myCoords->getPointer());
1287 targetMesh->setCoords(myCoords);
1288 myCoords->decrRef();
1292 MEDCouplingUMesh *MEDLoaderTest::build3DMesh_1()
1294 double coords[180]={
1295 0.,0.,0., 1.,1.,0., 1.,1.25,0., 0.,1.,0., 1.,1.5,0., 2.,0.,0., 2.,1.,0., 1.,2.,0., 0.,2.,0., 3.,1.,0.,
1296 3.,2.,0., 0.,1.,0., 1.,3.,0., 2.,2.,0., 2.,3.,0.,
1297 0.,0.,1., 1.,1.,1., 1.,1.25,1., 0.,1.,1., 1.,1.5,1., 2.,0.,1., 2.,1.,1., 1.,2.,1., 0.,2.,1., 3.,1.,1.,
1298 3.,2.,1., 0.,1.,1., 1.,3.,1., 2.,2.,1., 2.,3.,1.,
1299 0.,0.,2., 1.,1.,2., 1.,1.25,2., 0.,1.,2., 1.,1.5,2., 2.,0.,2., 2.,1.,2., 1.,2.,2., 0.,2.,2., 3.,1.,2.,
1300 3.,2.,2., 0.,1.,2., 1.,3.,2., 2.,2.,2., 2.,3.,2.,
1301 0.,0.,3., 1.,1.,3., 1.,1.25,3., 0.,1.,3., 1.,1.5,3., 2.,0.,3., 2.,1.,3., 1.,2.,3., 0.,2.,3., 3.,1.,3.,
1302 3.,2.,3., 0.,1.,3., 1.,3.,3., 2.,2.,3., 2.,3.,3.};
1306 0,11,1,3,15,26,16,18, 1,2,4,7,13,6,-1,1,16,21,6,-1,6,21,28,13,-1,13,7,22,28,-1,7,4,19,22,-1,4,2,17,19,-1,2,1,16,17,-1,16,21,28,22,19,17,
1307 1,6,5,3,16,21,20,18, 13,10,9,6,28,25,24,21,
1308 11,8,7,4,2,1,-1,11,26,16,1,-1,1,16,17,2,-1,2,17,19,4,-1,4,19,22,7,-1,7,8,23,22,-1,8,11,26,23,-1,26,16,17,19,22,23,
1309 7,12,14,13,22,27,29,28,
1311 15,26,16,18,30,41,31,33, 16,17,19,22,28,21,-1,16,31,36,21,-1,21,36,43,28,-1,28,22,37,43,-1,22,19,34,37,-1,19,17,32,34,-1,17,16,31,32,-1,31,36,43,37,34,32,
1312 16,21,20,18,31,36,35,33, 28,25,24,21,43,40,39,36,
1313 26,23,22,19,17,16,-1,26,41,31,16,-1,16,31,32,17,-1,17,32,34,19,-1,19,34,37,22,-1,22,23,38,37,-1,23,26,41,38,-1,41,31,32,34,37,38,
1314 22,27,29,28,37,42,44,43,
1316 30,41,31,33,45,56,46,48, 31,32,34,37,43,36,-1,31,46,51,36,-1,36,51,58,43,-1,43,37,52,58,-1,37,34,49,52,-1,34,32,47,49,-1,32,31,46,47,-1,46,51,58,52,49,47,
1317 31,36,35,33,46,51,50,48, 43,40,39,36,58,55,54,51,
1318 41,38,37,34,32,31,-1,41,56,46,31,-1,31,46,47,32,-1,32,47,49,34,-1,34,49,52,37,-1,37,38,53,52,-1,38,41,56,53,-1,56,46,47,49,52,53,
1319 37,42,44,43,52,57,59,58
1322 MEDCouplingUMesh *ret=MEDCouplingUMesh::New();
1323 ret->setName("3DMesh_1");
1324 ret->setMeshDimension(3);
1325 ret->allocateCells(18);
1327 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn);
1328 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+51);
1329 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+59);
1330 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+110);
1332 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+118);
1333 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+169);
1334 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+177);
1335 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+228);
1337 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+236);
1338 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+287);
1339 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+295);
1340 ret->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,conn+346);
1342 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+8);
1343 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+67);
1344 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+126);
1345 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+185);
1346 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+244);
1347 ret->insertNextCell(INTERP_KERNEL::NORM_POLYHED,43,conn+303);
1349 ret->finishInsertingCells();
1350 DataArrayDouble *myCoords=DataArrayDouble::New();
1351 myCoords->alloc(60,3);
1352 myCoords->setInfoOnComponent(0,"titi [m]");
1353 myCoords->setInfoOnComponent(1,"density power [MW/m^3]");
1354 myCoords->setInfoOnComponent(2,"t [kW]");
1355 std::copy(coords,coords+180,myCoords->getPointer());
1356 ret->setCoords(myCoords);
1357 myCoords->decrRef();
1361 MEDCouplingUMesh *MEDLoaderTest::build3DMesh_2()
1363 MEDCouplingUMesh *m3dsurfBase=build3DSurfMesh_1();
1364 int numbers[5]={0,1,2,3,5};
1365 MEDCouplingUMesh *m3dsurf=(MEDCouplingUMesh *)m3dsurfBase->buildPartOfMySelf(numbers,numbers+5,false);
1366 m3dsurfBase->decrRef();
1367 MEDCouplingUMesh *m1dBase=build1DMesh_1();
1368 int numbers2[4]={0,1,2,3};
1369 MEDCouplingUMesh *m1d=(MEDCouplingUMesh *)m1dBase->buildPartOfMySelf(numbers2,numbers2+4,false);
1371 m1d->changeSpaceDimension(3);
1372 const double vec[3]={0.,1.,0.};
1373 const double pt[3]={0.,0.,0.};
1374 m1d->rotate(pt,vec,-M_PI/2.);
1375 MEDCouplingUMesh *ret=m3dsurf->buildExtrudedMesh(m1d,0);
1381 MEDCouplingFieldDouble *MEDLoaderTest::buildVecFieldOnCells_1()
1383 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
1384 int nbOfCells=mesh->getNumberOfCells();
1385 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
1386 f1->setName("VectorFieldOnCells");
1388 DataArrayDouble *array=DataArrayDouble::New();
1389 array->alloc(nbOfCells,3);
1390 array->setInfoOnComponent(0,"power [MW/m^3]");
1391 array->setInfoOnComponent(1,"density [g/cm^3]");
1392 array->setInfoOnComponent(2,"temperature [K]");
1393 f1->setArray(array);
1395 double *tmp=array->getPointer();
1396 const double arr1[18]={0.,10.,20.,1.,11.,21.,2.,12.,22.,3.,13.,23.,4.,14.,24.,5.,15.,25.};
1397 std::copy(arr1,arr1+18,tmp);
1398 f1->setTime(2.,0,1);
1399 f1->checkCoherency();
1404 MEDCouplingFieldDouble *MEDLoaderTest::buildVecFieldOnNodes_1()
1406 MEDCouplingUMesh *mesh=build3DSurfMesh_1();
1407 int nbOfNodes=mesh->getNumberOfNodes();
1408 MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
1409 f1->setName("VectorFieldOnNodes");
1411 DataArrayDouble *array=DataArrayDouble::New();
1412 array->alloc(nbOfNodes,3);
1413 f1->setArray(array);
1414 array->setInfoOnComponent(0,"power [MW/m^3]");
1415 array->setInfoOnComponent(1,"density [g/cm^3]");
1416 array->setInfoOnComponent(2,"temperature [K]");
1418 double *tmp=array->getPointer();
1419 const double arr1[36]={
1420 70.,80.,90.,71.,81.,91.,72.,82.,92.,73.,83.,93.,74.,84.,94.,75.,85.,95.,
1421 1000.,10010.,10020.,1001.,10011.,10021.,1002.,10012.,10022.,1003.,10013.,10023.,1004.,10014.,10024.,1005.,10015.,10025.,
1423 std::copy(arr1,arr1+36,tmp);
1424 f1->setTime(2.12,2,3);
1425 f1->checkCoherency();
1430 MEDCouplingFieldDouble *MEDLoaderTest::buildVecFieldOnGauss_1()
1432 const double _a=0.446948490915965;
1433 const double _b=0.091576213509771;
1434 const double _p1=0.11169079483905;
1435 const double _p2=0.0549758718227661;
1436 const double refCoo1[6]={ 0.,0., 1.,0., 0.,1. };
1437 const double gsCoo1[12]={ 2*_b-1, 1-4*_b, 2*_b-1, 2.07*_b-1, 1-4*_b,
1438 2*_b-1, 1-4*_a, 2*_a-1, 2*_a-1, 1-4*_a, 2*_a-1, 2*_a-1 };
1439 const double wg1[6]={ 4*_p2, 4*_p2, 4*_p2, 4*_p1, 4*_p1, 4*_p1 };
1440 std::vector<double> _refCoo1(refCoo1,refCoo1+6);
1441 std::vector<double> _gsCoo1(gsCoo1,gsCoo1+12);
1442 std::vector<double> _wg1(wg1,wg1+6);
1443 MEDCouplingUMesh *m=build2DMesh_2();
1444 MEDCouplingFieldDouble *f=MEDCouplingFieldDouble::New(ON_GAUSS_PT,ONE_TIME);
1445 f->setTime(3.14,1,5);
1447 f->setGaussLocalizationOnType(INTERP_KERNEL::NORM_TRI3,_refCoo1,_gsCoo1,_wg1);
1448 const double refCoo2[12]={-1.0,1.0, -1.0,-1.0, 1.0,-1.0, -1.0,0.0, 0.0,-1.0, 0.0,0.0 };
1449 std::vector<double> _refCoo2(refCoo2,refCoo2+12);
1450 std::vector<double> _gsCoo2(_gsCoo1);
1451 std::vector<double> _wg2(_wg1);
1452 _gsCoo2.resize(6); _wg2.resize(3);
1453 const double refCoo3[8]={ 0.,0., 1.,0., 1.,1., 0.,1. };
1454 std::vector<double> _refCoo3(refCoo3,refCoo3+8);
1455 _gsCoo1.resize(4); _wg1.resize(2);
1456 f->setGaussLocalizationOnType(INTERP_KERNEL::NORM_QUAD4,_refCoo3,_gsCoo1,_wg1);
1457 f->setGaussLocalizationOnType(INTERP_KERNEL::NORM_TRI6,_refCoo2,_gsCoo2,_wg2);
1458 DataArrayDouble *array=DataArrayDouble::New();
1460 double *ptr=array->getPointer();
1461 for(int i=0;i<19*2;i++)
1462 ptr[i]=(double)(i+7);
1464 f->setName("MyFirstFieldOnGaussPoint");
1465 array->setInfoOnComponent(0,"power [MW/m^3]");
1466 array->setInfoOnComponent(1,"density");
1468 f->checkCoherency();
1473 MEDCouplingFieldDouble *MEDLoaderTest::buildVecFieldOnGaussNE_1()
1475 MEDCouplingUMesh *m=build2DMesh_2();
1476 MEDCouplingFieldDouble *f=MEDCouplingFieldDouble::New(ON_GAUSS_NE,ONE_TIME);
1477 f->setTime(3.14,1,5);
1479 DataArrayDouble *array=DataArrayDouble::New();
1481 double *ptr=array->getPointer();
1482 for(int i=0;i<20*2;i++)
1483 ptr[i]=(double)(i+8);
1485 array->setInfoOnComponent(0,"power [W]");
1486 array->setInfoOnComponent(1,"temperature");
1487 f->setName("MyFieldOnGaussNE");
1489 f->checkCoherency();