1 // Copyright (C) 2007-2023 CEA, EDF
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
20 #include "ParaMEDMEMTest.hxx"
21 #include "MEDLoader.hxx"
23 #include "ParaMEDFileMesh.hxx"
24 #include "MEDFileMesh.hxx"
25 #include "MEDFileField1TS.hxx"
26 #include "TestInterpKernelUtils.hxx"
27 #include "MEDCouplingFieldDouble.hxx"
29 #include <cppunit/TestAssert.h>
36 using namespace MEDCoupling;
39 * Generate a 2D mesh that is supposed to match the part that will be loaded by each proc in testParallelLoad1
41 MEDCouplingUMesh* genLocMesh2D(int rk)
45 MCAuto<MEDCouplingCMesh> msh = MEDCouplingCMesh::New("mesh");
46 MCAuto<DataArrayDouble> dax = DataArrayDouble::New(); dax->alloc(nx+1,1);
47 MCAuto<DataArrayDouble> day = DataArrayDouble::New(); day->alloc(ny+1,1);
48 dax->iota(); day->iota();
51 std::transform(dax->begin(), dax->end(),
53 [nxTot](const int& c){return c/(float)nxTot;});
54 std::transform(day->begin(), day->end(),
56 [nyTot](const int& c){return c/(float)nyTot;});
60 std::transform(dax->begin(), dax->end(),
62 [nxTot](const int& c){return c/(float)nxTot+0.5; });
63 std::transform(day->begin(), day->end(),
65 [nyTot](const int& c){return c/(float)nyTot;});
67 msh->setCoords(dax, day);
68 MCAuto<MEDCouplingUMesh> ret = msh->buildUnstructured();
73 * Generate a 2D mesh that is supposed to match the part that will be loaded by proc0 in testParallelLoad2
75 MEDCouplingUMesh* genLocMeshMultipleTypes1()
77 MCAuto<MEDCouplingUMesh> ret= MEDCouplingUMesh::New("mesh",2);
78 double coords[10] = {0.,1., 0.,2., 1.,2., 0.,3., 1.,3.};
79 DataArrayDouble *myCoords=DataArrayDouble::New();
81 std::copy(coords,coords+10,myCoords->getPointer());
82 ret->setCoords(myCoords);
84 mcIdType conn[7]={0,2,1, 1,2,4,3};
85 ret->allocateCells(2);
86 ret->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
87 ret->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+3);
88 ret->finishInsertingCells();
93 * Generate a 2D mesh that is supposed to match the part that will be loaded by proc1 in testParallelLoad2
95 MEDCouplingUMesh* genLocMeshMultipleTypes2()
97 MCAuto<MEDCouplingUMesh> ret= MEDCouplingUMesh::New("mesh",2);
98 double coords[10] = {0.,0., 1.,0., 0.,1., 1.,1., 1.,2.};
99 DataArrayDouble *myCoords=DataArrayDouble::New();
100 myCoords->alloc(5,2);
101 std::copy(coords,coords+10,myCoords->getPointer());
102 ret->setCoords(myCoords);
104 mcIdType conn[7]={2,3,4, 0,1,3,2};
105 ret->allocateCells(2);
106 ret->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
107 ret->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+3);
108 ret->finishInsertingCells();
113 * Generate a 2D mesh that is supposed to match the part that will be loaded by proc2 in testParallelLoad2
115 MEDCouplingUMesh* genLocMeshMultipleTypes3()
117 MCAuto<MEDCouplingUMesh> ret= MEDCouplingUMesh::New("mesh",2);
118 double coords[16] = {1.,0., 2.,0., 1.,1., 2.,1., 1.,2., 2.,2., 1.,3., 2.,3.};
119 DataArrayDouble *myCoords=DataArrayDouble::New();
120 myCoords->alloc(8,2);
121 std::copy(coords,coords+16,myCoords->getPointer());
122 ret->setCoords(myCoords);
124 mcIdType conn[14]={0,1,3, 0,3,2, 2,3,5,4, 4,5,7,6};
125 ret->allocateCells(4);
126 ret->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
127 ret->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+3);
128 ret->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+6);
129 ret->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+10);
130 ret->finishInsertingCells();
135 * Generate a 2D field that is supposed to match the local field loaded by each proc in testParallelLoad4
137 MEDCouplingFieldDouble *genLocFieldCells(int rank)
139 MCAuto<MEDCouplingUMesh> mesh = genLocMesh2D(rank);
140 MCAuto<MEDCouplingFieldDouble> f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
141 f1->setName("field");
144 MCAuto<DataArrayDouble> array(DataArrayDouble::New());
146 std::vector<double> values;
148 values = { 0., 10., 20., 30., 80., 90., 100., 110.};
150 values = { 40., 50., 60., 70., 120., 130., 140., 150.};
151 std::copy(values.data(),values.data()+8,array->getPointer());
152 array->setInfoOnComponent(0,"");
158 * Generate a 2D field that is supposed to match the local field loaded by each proc in testParallelLoad5
160 MEDCouplingFieldDouble *genLocFieldNodes(int rank)
162 MCAuto<MEDCouplingUMesh> mesh = genLocMesh2D(rank);
163 MCAuto<MEDCouplingFieldDouble> f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME);
164 f1->setName("field");
167 MCAuto<DataArrayDouble> array(DataArrayDouble::New());
169 std::vector<double> values;
171 values= { 0., 10., 20., 30., 40., 50., 100., 110., 120., 130., 140., 150., 200., 210., 220., 230., 240., 250. };
173 values= { 40., 50., 60., 70., 80., 90., 140., 150., 160., 170., 180., 190., 240., 250., 260., 270., 280., 290. };
174 std::copy(values.data(),values.data()+18,array->getPointer());
175 array->setInfoOnComponent(0,"");
181 * Test case to load a simple 2D cartesian mesh in parallel on 2 procs
183 void ParaMEDMEMTest::testParallelLoad1()
187 MPI_Comm_size(MPI_COMM_WORLD,&size);
188 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
193 std::map<INTERP_KERNEL::NormalizedCellType, std::vector<mcIdType>> distrib;
195 distrib = { {INTERP_KERNEL::NORM_QUAD4,{0,1,4,5}/*c++ type of indexing: index starts from zero!*/} };
197 distrib = { {INTERP_KERNEL::NORM_QUAD4,{2,3,6,7}} };
199 std::string filename=INTERP_TEST::getResourceFile("SimpleTest2D.med");
200 MCAuto<MEDFileUMesh> mu = ParaMEDFileUMesh::ParaNew(distrib, MPI_COMM_WORLD, MPI_INFO_NULL, filename, "mesh");
201 MCAuto<MEDCouplingUMesh> mesh = mu->getMeshAtLevel(0);
202 MCAuto<MEDCouplingUMesh> meshRef = genLocMesh2D(rank);
203 CPPUNIT_ASSERT(mesh->isEqual(meshRef,1e-12));
204 MPI_Barrier(MPI_COMM_WORLD);
208 * Test case to load a 2D mesh made of squares and triangles in parallel on 3 procs.
209 * Each proc is going to load a part of the mesh.
211 void ParaMEDMEMTest::testParallelLoad2()
215 MPI_Comm_size(MPI_COMM_WORLD,&size);
216 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
221 std::map<INTERP_KERNEL::NormalizedCellType, std::vector<mcIdType>> distrib;
222 // independant numerotation for each geometric type!
224 distrib = { {INTERP_KERNEL::NORM_TRI3,{3}} , {INTERP_KERNEL::NORM_QUAD4,{2}} };
226 distrib = { {INTERP_KERNEL::NORM_TRI3,{2}} , {INTERP_KERNEL::NORM_QUAD4,{0}} };
228 distrib= { {INTERP_KERNEL::NORM_TRI3,{0,1}} , {INTERP_KERNEL::NORM_QUAD4,{1,3}} };
230 std::string filename=INTERP_TEST::getResourceFile("Test2DMultiGeoType.med");
231 MCAuto<MEDFileUMesh> mu = ParaMEDFileUMesh::ParaNew(distrib, MPI_COMM_WORLD, MPI_INFO_NULL, filename, "mesh");
232 MCAuto<MEDCouplingUMesh> mesh = mu->getMeshAtLevel(0);
233 MEDCouplingUMesh *meshRef;
235 meshRef=genLocMeshMultipleTypes1();
237 meshRef=genLocMeshMultipleTypes2();
239 meshRef=genLocMeshMultipleTypes3();
240 //checking that all 3 procs have correctly loaded their part
241 int equal = (int)mesh->isEqual(meshRef,1e-12);
243 MPI_Allreduce(&equal, &allEqual, 1, MPI_INT,MPI_SUM,MPI_COMM_WORLD);
244 CPPUNIT_ASSERT(allEqual==3);
247 MPI_Barrier(MPI_COMM_WORLD);
251 * Test case to load a 3D box meshed with tetras in parallel on 2 procs
253 void ParaMEDMEMTest::testParallelLoad3()
257 MPI_Comm_size(MPI_COMM_WORLD,&size);
258 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
263 std::map<INTERP_KERNEL::NormalizedCellType, std::vector<mcIdType>> distrib;
266 std::vector<mcIdType> distribCells = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,48,49,50,51,52,53,54,55,56,57,
267 58,59,60,61,62,63,64,65,66,67,68,69,70,71,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,
268 119,144,145,146,147,148,149,150,151,152,153,154,155,156,157,158,159,160,161,162,163,164,165,166,167};
269 distrib = { {INTERP_KERNEL::NORM_TETRA4,distribCells} };
273 std::vector<mcIdType> distribCells = {24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,72,73,74,75,76,77,78,79,80,
274 81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,120,121,122,123,124,125,126,127,128,129,130,131,132,133,134,135,136,137,138,139,140,141,142,
275 143,168,169,170,171,172,173,174,175,176,177,178,179,180,181,182,183,184,185,186,187,188,189,190,191};
276 distrib = { {INTERP_KERNEL::NORM_TETRA4,distribCells} };
279 std::string filename=INTERP_TEST::getResourceFile("SimpleTest3D.med");
280 MCAuto<MEDFileUMesh> mu = ParaMEDFileUMesh::ParaNew(distrib, MPI_COMM_WORLD, MPI_INFO_NULL, filename, "mesh");
281 MCAuto<MEDCouplingUMesh> mesh = mu->getMeshAtLevel(0);
282 CPPUNIT_ASSERT_EQUAL(96,(int)mesh->getNumberOfCells());
284 // checking nodal connectivity
285 double nodalConnec[480] = {14, 1, 7, 18, 24, 14, 7, 6, 18, 24, 14, 6, 0, 18, 24, 14, 0, 1, 18, 24, 14, 1, 0, 19, 24, 14, 0, 2, 19, 24, 14, 2, 3, 19, 24,
286 14, 3, 1, 19, 24, 14, 1, 3, 20, 24, 14, 3, 9, 20, 24, 14, 9, 7, 20, 24, 14, 7, 1, 20, 24, 14, 0, 6, 21, 24, 14, 6, 8, 21, 24, 14, 8, 2, 21, 24,
287 14, 2, 0, 21, 24, 14, 7, 9, 22, 24, 14, 9, 8, 22, 24, 14, 8, 6, 22, 24, 14, 6, 7, 22, 24, 14, 2, 8, 23, 24, 14, 8, 9, 23, 24, 14, 9, 3, 23, 24,
288 14, 3, 2, 23, 24, 14, 3, 9, 25, 31, 14, 9, 8, 25, 31, 14, 8, 2, 25, 31, 14, 2, 3, 25, 31, 14, 3, 2, 26, 31, 14, 2, 4, 26, 31, 14, 4, 5, 26, 31,
289 14, 5, 3, 26, 31, 14, 3, 5, 27, 31, 14, 5, 11, 27, 31, 14, 11, 9, 27, 31, 14, 9, 3, 27, 31, 14, 2, 8, 28, 31, 14, 8, 10, 28, 31, 14, 10, 4, 28, 31,
290 14, 4, 2, 28, 31, 14, 9, 11, 29, 31, 14, 11, 10, 29, 31, 14, 10, 8, 29, 31, 14, 8, 9, 29, 31, 14, 4, 10, 30, 31, 14, 10, 11, 30, 31, 14, 11, 5, 30, 31,
291 14, 5, 4, 30, 31, 14, 7, 13, 32, 38, 14, 13, 12, 32, 38, 14, 12, 6, 32, 38, 14, 6, 7, 32, 38, 14, 7, 6, 33, 38, 14, 6, 8, 33, 38, 14, 8, 9, 33, 38,
292 14, 9, 7, 33, 38, 14, 7, 9, 34, 38, 14, 9, 15, 34, 38, 14, 15, 13, 34, 38, 14, 13, 7, 34, 38, 14, 6, 12, 35, 38, 14, 12, 14, 35, 38, 14, 14, 8, 35, 38,
293 14, 8, 6, 35, 38, 14, 13, 15, 36, 38, 14, 15, 14, 36, 38, 14, 14, 12, 36, 38, 14, 12, 13, 36, 38, 14, 8, 14, 37, 38, 14, 14, 15, 37, 38, 14, 15, 9, 37, 38,
294 14, 9, 8, 37, 38, 14, 9, 15, 39, 45, 14, 15, 14, 39, 45, 14, 14, 8, 39, 45, 14, 8, 9, 39, 45, 14, 9, 8, 40, 45, 14, 8, 10, 40, 45, 14, 10, 11, 40, 45,
295 14, 11, 9, 40, 45, 14, 9, 11, 41, 45, 14, 11, 17, 41, 45, 14, 17, 15, 41, 45, 14, 15, 9, 41, 45, 14, 8, 14, 42, 45, 14, 14, 16, 42, 45, 14, 16, 10, 42, 45,
296 14, 10, 8, 42, 45, 14, 15, 17, 43, 45, 14, 17, 16, 43, 45, 14, 16, 14, 43, 45, 14, 14, 15, 43, 45, 14, 10, 16, 44, 45, 14, 16, 17, 44, 45, 14, 17, 11, 44, 45,
298 const mcIdType *nc=mesh->getNodalConnectivity()->getConstPointer();
299 CPPUNIT_ASSERT_EQUAL(480,(int)mesh->getNodalConnectivity()->getNumberOfTuples());
300 CPPUNIT_ASSERT(std::equal(nodalConnec,nodalConnec+480,nc));
302 double nodalConnecInd[97] = {0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 85, 90, 95, 100, 105, 110, 115, 120, 125, 130, 135, 140, 145, 150, 155,
303 160, 165, 170, 175, 180, 185, 190, 195, 200, 205, 210, 215, 220, 225, 230, 235, 240, 245, 250, 255, 260, 265 ,270, 275, 280, 285, 290, 295, 300, 305, 310, 315, 320,
304 325, 330, 335, 340, 345, 350, 355, 360, 365, 370, 375, 380, 385, 390, 395, 400, 405, 410, 415, 420, 425, 430, 435, 440, 445, 450, 455, 460, 465, 470, 475, 480};
305 const mcIdType *ncIndx=mesh->getNodalConnectivityIndex()->getConstPointer();
306 CPPUNIT_ASSERT_EQUAL(97,(int)mesh->getNodalConnectivityIndex()->getNumberOfTuples());
307 CPPUNIT_ASSERT(std::equal(nodalConnecInd,nodalConnecInd+97,ncIndx));
310 std::vector<double> coords(138);
312 coords = {0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 2.0, 0.0, 2.0, 2.0, 0.0, 0.0, 4.0, 0.0, 2.0, 4.0, 0.0, 0.0, 0.0, 2.0, 2.0, 0.0, 2.0, 0.0, 2.0,
313 2.0, 2.0, 2.0, 2.0, 0.0, 4.0, 2.0, 2.0, 4.0, 2.0, 0.0, 0.0, 4.0, 2.0, 0.0, 4.0, 0.0, 2.0, 4.0, 2.0, 2.0, 4.0, 0.0, 4.0, 4.0, 2.0, 4.0, 4.0, 1.0, 0.0,
314 1.0, 1.0, 1.0, 0.0, 2.0, 1.0, 1.0, 0.0, 1.0, 1.0, 1.0, 1.0, 2.0, 1.0, 2.0, 1.0, 1.0, 1.0, 1.0, 1.0, 2.0, 1.0, 1.0, 3.0, 0.0, 2.0, 3.0, 1.0, 0.0, 3.0,
315 1.0, 1.0, 3.0, 2.0, 1.0, 4.0, 1.0, 1.0, 3.0, 1.0, 1.0, 0.0, 3.0, 1.0, 1.0, 2.0, 2.0, 1.0, 3.0, 0.0, 1.0, 3.0, 1.0, 1.0, 4.0, 1.0, 2.0, 3.0, 1.0, 1.0,
316 3.0, 1.0, 2.0, 3.0, 1.0, 3.0, 2.0, 2.0, 3.0, 3.0, 0.0, 3.0, 3.0, 1.0, 3.0, 4.0, 1.0, 4.0, 3.0, 1.0, 3.0, 3.0 };
318 coords = {2.0, 0.0, 0.0, 4.0, 0.0, 0.0, 2.0, 2.0, 0.0, 4.0, 2.0, 0.0, 2.0, 4.0, 0.0, 4.0, 4.0, 0.0, 2.0, 0.0, 2.0, 4.0, 0.0, 2.0, 2.0, 2.0, 2.0, 4.0,
319 2.0, 2.0, 2.0, 4.0, 2.0, 4.0, 4.0, 2.0, 2.0, 0.0, 4.0, 4.0, 0.0, 4.0, 2.0, 2.0, 4.0, 4.0, 2.0, 4.0, 2.0, 4.0, 4.0, 4.0, 4.0, 4.0, 3.0, 0.0, 1.0, 3.0,
320 1.0, 0.0, 4.0, 1.0, 1.0, 2.0, 1.0, 1.0, 3.0, 1.0, 2.0, 3.0, 2.0, 1.0, 3.0, 1.0, 1.0, 3.0, 2.0, 1.0, 3.0, 3.0, 0.0, 4.0, 3.0, 1.0, 2.0, 3.0, 1.0, 3.0,
321 3.0, 2.0, 3.0, 4.0, 1.0, 3.0, 3.0, 1.0, 3.0, 0.0, 3.0, 3.0, 1.0, 2.0, 4.0, 1.0, 3.0, 2.0, 1.0, 3.0, 3.0, 1.0, 4.0, 3.0, 2.0, 3.0, 3.0, 1.0, 3.0, 3.0,
322 2.0, 3.0, 3.0, 3.0, 2.0, 4.0, 3.0, 3.0, 2.0, 3.0, 3.0, 3.0, 3.0, 4.0, 3.0, 4.0, 3.0, 3.0, 3.0, 3.0 };
323 const double *coo=mesh->getCoords()->getConstPointer();
324 CPPUNIT_ASSERT_EQUAL(46,(int)mesh->getCoords()->getNumberOfTuples());
325 CPPUNIT_ASSERT(std::equal(coords.data(),coords.data()+138,coo));
327 MPI_Barrier(MPI_COMM_WORLD);
331 * Test case to load a field located on cells in parallel on 2 procs.
333 void ParaMEDMEMTest::testParallelLoad4()
337 MPI_Comm_size(MPI_COMM_WORLD,&size);
338 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
343 std::vector<mcIdType> distrib;
345 distrib = {0,1,4,5}; //c++ type of indexing: index starts from zero!
349 std::string filename=INTERP_TEST::getResourceFile("SimpleTest2D.med");
350 MCAuto<MEDFileField1TS> f1TS = ParaMEDFileField1TS::ParaNew(MPI_COMM_WORLD, MPI_INFO_NULL,filename,"fieldOnCells","mesh",distrib,ON_CELLS);
351 MCAuto<MEDCouplingFieldDouble> fieldRef = genLocFieldCells(rank);
352 CPPUNIT_ASSERT(f1TS->getUndergroundDataArray()->isEqual(*fieldRef->getArray(),1e-12));
353 MPI_Barrier(MPI_COMM_WORLD);
357 * Test case to load a field located on nodes in parallel on 2 procs.
359 void ParaMEDMEMTest::testParallelLoad5()
363 MPI_Comm_size(MPI_COMM_WORLD,&size);
364 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
369 std::vector<mcIdType> distrib;
371 distrib = {0,1,2,5,6,7,10,11,12}; //c++ type of indexing: index starts from zero!
373 distrib = {2,3,4,7,8,9,12,13,14};
375 std::string filename=INTERP_TEST::getResourceFile("SimpleTest2D.med");
376 MCAuto<MEDFileField1TS> f1TS = ParaMEDFileField1TS::ParaNew(MPI_COMM_WORLD, MPI_INFO_NULL,filename,"fieldOnNodes","mesh",distrib,ON_NODES);
377 MCAuto<MEDCouplingFieldDouble> fieldRef = genLocFieldNodes(rank);
378 CPPUNIT_ASSERT(f1TS->getUndergroundDataArray()->isEqual(*fieldRef->getArray(),1e-12));
379 MPI_Barrier(MPI_COMM_WORLD);