1 // Copyright (C) 2007-2016 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
20 #include "ParaMEDMEMTest.hxx"
21 #include <cppunit/TestAssert.h>
23 #include "CommInterface.hxx"
24 #include "ProcessorGroup.hxx"
25 #include "MPIProcessorGroup.hxx"
26 #include "Topology.hxx"
28 #include "MxN_Mapping.hxx"
29 #include "InterpKernelDEC.hxx"
30 #include "ParaMESH.hxx"
31 #include "ParaFIELD.hxx"
32 #include "ComponentTopology.hxx"
33 #include "ICoCoMEDField.hxx"
34 #include "ParaMEDLoader.hxx"
35 #include "MEDLoader.hxx"
36 #include "TestInterpKernelUtils.hxx"
42 // use this define to enable lines, execution of which leads to Segmentation Fault
45 // use this define to enable CPPUNIT asserts and fails, showing bugs
46 #define ENABLE_FORCED_FAILURES
50 using namespace MEDCoupling;
52 void ParaMEDMEMTest::testInterpKernelDEC_2D()
54 testInterpKernelDEC_2D_("P0","P0");
57 void ParaMEDMEMTest::testInterpKernelDEC2_2D()
59 testInterpKernelDEC2_2D_("P0","P0");
62 void ParaMEDMEMTest::testInterpKernelDEC_3D()
64 testInterpKernelDEC_3D_("P0","P0");
67 void ParaMEDMEMTest::testInterpKernelDEC_2DP0P1()
69 //testInterpKernelDEC_2D_("P0","P1");
72 void ParaMEDMEMTest::testInterpKernelDEC_1D()
76 MPI_Comm_size(MPI_COMM_WORLD,&size);
77 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
83 set<int> procs_source;
84 set<int> procs_target;
86 for (int i=0; i<nproc_source; i++)
87 procs_source.insert(i);
88 for (int i=nproc_source; i<size; i++)
89 procs_target.insert(i);
90 self_procs.insert(rank);
92 MEDCoupling::MEDCouplingUMesh *mesh=0;
93 MEDCoupling::ParaMESH *paramesh=0;
94 MEDCoupling::ParaFIELD *parafieldP0=0;
96 MEDCoupling::CommInterface interface;
98 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
99 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
100 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
102 MPI_Barrier(MPI_COMM_WORLD);
103 if(source_group->containsMyRank())
107 double coords[4]={0.3,0.7, 0.9,1.0};
108 int conn[4]={0,1,2,3};
109 mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
110 mesh->allocateCells(2);
111 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
112 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
113 mesh->finishInsertingCells();
114 DataArrayDouble *myCoords=DataArrayDouble::New();
115 myCoords->alloc(4,1);
116 std::copy(coords,coords+4,myCoords->getPointer());
117 mesh->setCoords(myCoords);
122 double coords[2]={0.7,0.9};
124 mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
125 mesh->allocateCells(1);
126 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
127 mesh->finishInsertingCells();
128 DataArrayDouble *myCoords=DataArrayDouble::New();
129 myCoords->alloc(2,1);
130 std::copy(coords,coords+2,myCoords->getPointer());
131 mesh->setCoords(myCoords);
136 double coords[2]={1.,1.12};
138 mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
139 mesh->allocateCells(1);
140 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
141 mesh->finishInsertingCells();
142 DataArrayDouble *myCoords=DataArrayDouble::New();
143 myCoords->alloc(2,1);
144 std::copy(coords,coords+2,myCoords->getPointer());
145 mesh->setCoords(myCoords);
148 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
149 MEDCoupling::ComponentTopology comptopo;
150 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
151 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
152 parafieldP0->getField()->setNature(IntensiveMaximum);
155 valueP0[0]=7.; valueP0[1]=8.;
168 const char targetMeshName[]="target mesh";
171 double coords[2]={0.5,0.75};
173 mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
174 mesh->allocateCells(1);
175 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
176 mesh->finishInsertingCells();
177 DataArrayDouble *myCoords=DataArrayDouble::New();
178 myCoords->alloc(2,1);
179 std::copy(coords,coords+2,myCoords->getPointer());
180 mesh->setCoords(myCoords);
182 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
186 double coords[2]={0.75,1.2};
188 mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
189 mesh->allocateCells(1);
190 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
191 mesh->finishInsertingCells();
192 DataArrayDouble *myCoords=DataArrayDouble::New();
193 myCoords->alloc(2,1);
194 std::copy(coords,coords+2,myCoords->getPointer());
195 mesh->setCoords(myCoords);
197 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
199 MEDCoupling::ComponentTopology comptopo;
200 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
201 parafieldP0->getField()->setNature(IntensiveMaximum);
204 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
205 if (source_group->containsMyRank())
208 dec.attachLocalField(parafieldP0);
210 dec.setForcedRenormalization(false);
213 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
216 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
217 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
221 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
225 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
231 dec.attachLocalField(parafieldP0);
233 dec.setForcedRenormalization(false);
235 const double *res=parafieldP0->getField()->getArray()->getConstPointer();
238 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
239 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
240 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
244 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
245 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
246 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
258 MPI_Barrier(MPI_COMM_WORLD);
261 void ParaMEDMEMTest::testInterpKernelDEC_2DCurve()
265 MPI_Comm_size(MPI_COMM_WORLD,&size);
266 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
270 int nproc_source = 3;
272 set<int> procs_source;
273 set<int> procs_target;
275 for (int i=0; i<nproc_source; i++)
276 procs_source.insert(i);
277 for (int i=nproc_source; i<size; i++)
278 procs_target.insert(i);
279 self_procs.insert(rank);
281 MEDCoupling::MEDCouplingUMesh *mesh=0;
282 MEDCoupling::ParaMESH *paramesh=0;
283 MEDCoupling::ParaFIELD *parafieldP0=0;
285 MEDCoupling::CommInterface interface;
287 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
288 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
289 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
291 MPI_Barrier(MPI_COMM_WORLD);
292 if(source_group->containsMyRank())
296 double coords[8]={0.3,0.3,0.7,0.7, 0.9,0.9,1.0,1.0};
297 int conn[4]={0,1,2,3};
298 mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
299 mesh->allocateCells(2);
300 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
301 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
302 mesh->finishInsertingCells();
303 DataArrayDouble *myCoords=DataArrayDouble::New();
304 myCoords->alloc(4,2);
305 std::copy(coords,coords+8,myCoords->getPointer());
306 mesh->setCoords(myCoords);
311 double coords[4]={0.7,0.7,0.9,0.9};
313 mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
314 mesh->allocateCells(1);
315 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
316 mesh->finishInsertingCells();
317 DataArrayDouble *myCoords=DataArrayDouble::New();
318 myCoords->alloc(2,2);
319 std::copy(coords,coords+4,myCoords->getPointer());
320 mesh->setCoords(myCoords);
325 double coords[4]={1.,1.,1.12,1.12};
327 mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
328 mesh->allocateCells(1);
329 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
330 mesh->finishInsertingCells();
331 DataArrayDouble *myCoords=DataArrayDouble::New();
332 myCoords->alloc(2,2);
333 std::copy(coords,coords+4,myCoords->getPointer());
334 mesh->setCoords(myCoords);
337 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
338 MEDCoupling::ComponentTopology comptopo;
339 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
340 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
341 parafieldP0->getField()->setNature(IntensiveMaximum);
344 valueP0[0]=7.; valueP0[1]=8.;
357 const char targetMeshName[]="target mesh";
360 double coords[4]={0.5,0.5,0.75,0.75};
362 mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
363 mesh->allocateCells(1);
364 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
365 mesh->finishInsertingCells();
366 DataArrayDouble *myCoords=DataArrayDouble::New();
367 myCoords->alloc(2,2);
368 std::copy(coords,coords+4,myCoords->getPointer());
369 mesh->setCoords(myCoords);
371 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
375 double coords[4]={0.75,0.75,1.2,1.2};
377 mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
378 mesh->allocateCells(1);
379 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
380 mesh->finishInsertingCells();
381 DataArrayDouble *myCoords=DataArrayDouble::New();
382 myCoords->alloc(2,2);
383 std::copy(coords,coords+4,myCoords->getPointer());
384 mesh->setCoords(myCoords);
386 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
388 MEDCoupling::ComponentTopology comptopo;
389 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
390 parafieldP0->getField()->setNature(IntensiveMaximum);
393 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
394 if (source_group->containsMyRank())
397 dec.attachLocalField(parafieldP0);
399 dec.setForcedRenormalization(false);
402 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
405 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
406 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
410 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
414 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
420 dec.attachLocalField(parafieldP0);
422 dec.setForcedRenormalization(false);
424 const double *res=parafieldP0->getField()->getArray()->getConstPointer();
427 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
428 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
429 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
433 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
434 CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
435 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
447 MPI_Barrier(MPI_COMM_WORLD);
452 * Check methods defined in InterpKernelDEC.hxx
455 InterpKernelDEC(ProcessorGroup& local_group, ProcessorGroup& distant_group);
456 virtual ~InterpKernelDEC();
462 void ParaMEDMEMTest::testInterpKernelDEC_2D_(const char *srcMeth, const char *targetMeth)
464 std::string srcM(srcMeth);
465 std::string targetM(targetMeth);
468 MPI_Comm_size(MPI_COMM_WORLD,&size);
469 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
471 //the test is meant to run on five processors
472 if (size !=5) return ;
474 int nproc_source = 3;
476 set<int> procs_source;
477 set<int> procs_target;
479 for (int i=0; i<nproc_source; i++)
480 procs_source.insert(i);
481 for (int i=nproc_source; i<size; i++)
482 procs_target.insert(i);
483 self_procs.insert(rank);
485 MEDCoupling::CommInterface interface;
487 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
488 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
489 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
491 //loading the geometry for the source group
493 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
495 MEDCoupling::MEDCouplingUMesh* mesh;
496 MEDCoupling::ParaMESH* paramesh;
497 MEDCoupling::ParaFIELD* parafield;
498 ICoCo::MEDField* icocofield ;
500 string filename_xml1 = "square1_split";
501 string filename_xml2 = "square2_split";
502 //string filename_seq_wr = makeTmpFile("");
503 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
505 // To remove tmp files from disk
506 ParaMEDMEMTest_TmpFilesRemover aRemover;
508 MPI_Barrier(MPI_COMM_WORLD);
509 if (source_group->containsMyRank())
511 string master = filename_xml1;
513 ostringstream strstream;
514 strstream <<master<<rank+1<<".med";
515 string fName = INTERP_TEST::getResourceFile(strstream.str());
516 ostringstream meshname ;
517 meshname<< "Mesh_2_"<< rank+1;
519 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
522 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
524 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
525 MEDCoupling::ComponentTopology comptopo;
528 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
529 parafield->getField()->setNature(IntensiveMaximum);
532 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
535 nb_local=mesh->getNumberOfCells();
537 nb_local=mesh->getNumberOfNodes();
538 // double * value= new double[nb_local];
539 double *value=parafield->getField()->getArray()->getPointer();
540 for(int ielem=0; ielem<nb_local;ielem++)
543 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
544 icocofield=new ICoCo::MEDField(parafield->getField());
545 dec.setMethod(srcMeth);
546 dec.attachLocalField(icocofield);
549 //loading the geometry for the target group
550 if (target_group->containsMyRank())
552 string master= filename_xml2;
553 ostringstream strstream;
554 strstream << master<<(rank-nproc_source+1)<<".med";
555 string fName = INTERP_TEST::getResourceFile(strstream.str());
556 ostringstream meshname ;
557 meshname<< "Mesh_3_"<<rank-nproc_source+1;
558 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
560 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
561 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
562 MEDCoupling::ComponentTopology comptopo;
565 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
566 parafield->getField()->setNature(IntensiveMaximum);
569 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
572 nb_local=mesh->getNumberOfCells();
574 nb_local=mesh->getNumberOfNodes();
575 // double * value= new double[nb_local];
576 double *value=parafield->getField()->getArray()->getPointer();
577 for(int ielem=0; ielem<nb_local;ielem++)
579 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
580 icocofield=new ICoCo::MEDField(parafield->getField());
581 dec.setMethod(targetMeth);
582 dec.attachLocalField(icocofield);
586 //attaching a DEC to the source group
587 double field_before_int;
588 double field_after_int;
590 if (source_group->containsMyRank())
592 field_before_int = parafield->getVolumeIntegral(0,true);
594 cout<<"DEC usage"<<endl;
595 dec.setForcedRenormalization(false);
598 ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
599 if (source_group->myRank()==0)
600 aRemover.Register("./sourcesquareb");
601 ostringstream filename;
602 filename<<"./sourcesquareb_"<<source_group->myRank()+1;
603 aRemover.Register(filename.str().c_str());
604 //WriteField("./sourcesquareb",parafield->getField());
607 cout <<"writing"<<endl;
608 ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
609 if (source_group->myRank()==0)
610 aRemover.Register("./sourcesquare");
611 //WriteField("./sourcesquare",parafield->getField());
614 filename<<"./sourcesquare_"<<source_group->myRank()+1;
615 aRemover.Register(filename.str().c_str());
616 field_after_int = parafield->getVolumeIntegral(0,true);
619 // MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
620 // MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
622 CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
626 //attaching a DEC to the target group
627 if (target_group->containsMyRank())
630 dec.setForcedRenormalization(false);
633 ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
634 //WriteField("./targetsquareb",parafield->getField());
635 if (target_group->myRank()==0)
636 aRemover.Register("./targetsquareb");
637 ostringstream filename;
638 filename<<"./targetsquareb_"<<target_group->myRank()+1;
639 aRemover.Register(filename.str().c_str());
641 ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
642 //WriteField("./targetsquare",parafield->getField());
644 if (target_group->myRank()==0)
645 aRemover.Register("./targetsquareb");
647 filename<<"./targetsquareb_"<<target_group->myRank()+1;
648 aRemover.Register(filename.str().c_str());
649 // double field_before_int, field_after_int;
650 // MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
651 // MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
653 // CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
666 MPI_Barrier(MPI_COMM_WORLD);
667 cout << "end of InterpKernelDEC_2D test"<<endl;
670 void ParaMEDMEMTest::testInterpKernelDEC2_2D_(const char *srcMeth, const char *targetMeth)
672 std::string srcM(srcMeth);
673 std::string targetM(targetMeth);
676 MPI_Comm_size(MPI_COMM_WORLD,&size);
677 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
679 //the test is meant to run on five processors
680 if (size !=5) return ;
682 int nproc_source = 3;
684 set<int> procs_source;
685 set<int> procs_target;
687 for (int i=0; i<nproc_source; i++)
688 procs_source.insert(i);
689 for (int i=nproc_source; i<size; i++)
690 procs_target.insert(i);
691 self_procs.insert(rank);
693 MEDCoupling::CommInterface interface;
695 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
696 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
697 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
699 //loading the geometry for the source group
701 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
703 MEDCoupling::MEDCouplingUMesh* mesh;
704 MEDCoupling::MEDCouplingFieldDouble* mcfield;
706 string filename_xml1 = "square1_split";
707 string filename_xml2 = "square2_split";
709 // To remove tmp files from disk
710 ParaMEDMEMTest_TmpFilesRemover aRemover;
712 MPI_Barrier(MPI_COMM_WORLD);
713 if (source_group->containsMyRank())
715 string master = filename_xml1;
717 ostringstream strstream;
718 strstream <<master<<rank+1<<".med";
719 string fName = INTERP_TEST::getResourceFile(strstream.str());
720 ostringstream meshname ;
721 meshname<< "Mesh_2_"<< rank+1;
723 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
724 MEDCoupling::ComponentTopology comptopo;
727 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
728 mcfield->setMesh(mesh);
729 DataArrayDouble *array=DataArrayDouble::New();
730 array->alloc(mcfield->getNumberOfTuples(),1);
731 mcfield->setArray(array);
733 mcfield->setNature(IntensiveMaximum);
737 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
738 mcfield->setMesh(mesh);
739 DataArrayDouble *array=DataArrayDouble::New();
740 array->alloc(mcfield->getNumberOfTuples(),1);
741 mcfield->setArray(array);
746 nb_local=mesh->getNumberOfCells();
748 nb_local=mesh->getNumberOfNodes();
749 double *value=mcfield->getArray()->getPointer();
750 for(int ielem=0; ielem<nb_local;ielem++)
752 dec.setMethod(srcMeth);
753 dec.attachLocalField(mcfield);
754 dec.attachLocalField(mcfield);
757 //loading the geometry for the target group
758 if (target_group->containsMyRank())
760 string master= filename_xml2;
761 ostringstream strstream;
762 strstream << master<<(rank-nproc_source+1)<<".med";
763 string fName = INTERP_TEST::getResourceFile(strstream.str());
764 ostringstream meshname ;
765 meshname<< "Mesh_3_"<<rank-nproc_source+1;
766 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
767 MEDCoupling::ComponentTopology comptopo;
770 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
771 mcfield->setMesh(mesh);
772 DataArrayDouble *array=DataArrayDouble::New();
773 array->alloc(mcfield->getNumberOfTuples(),1);
774 mcfield->setArray(array);
776 mcfield->setNature(IntensiveMaximum);
780 mcfield = MEDCouplingFieldDouble::New(ON_NODES,NO_TIME);
781 mcfield->setMesh(mesh);
782 DataArrayDouble *array=DataArrayDouble::New();
783 array->alloc(mcfield->getNumberOfTuples(),1);
784 mcfield->setArray(array);
789 nb_local=mesh->getNumberOfCells();
791 nb_local=mesh->getNumberOfNodes();
792 double *value=mcfield->getArray()->getPointer();
793 for(int ielem=0; ielem<nb_local;ielem++)
795 dec.setMethod(targetMeth);
796 dec.attachLocalField(mcfield);
797 dec.attachLocalField(mcfield);
801 //attaching a DEC to the source group
803 if (source_group->containsMyRank())
806 dec.setForcedRenormalization(false);
811 //attaching a DEC to the target group
812 if (target_group->containsMyRank())
815 dec.setForcedRenormalization(false);
825 MPI_Barrier(MPI_COMM_WORLD);
826 cout << "end of InterpKernelDEC2_2D test"<<endl;
829 void ParaMEDMEMTest::testInterpKernelDEC_3D_(const char *srcMeth, const char *targetMeth)
831 std::string srcM(srcMeth);
832 std::string targetM(targetMeth);
835 MPI_Comm_size(MPI_COMM_WORLD,&size);
836 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
838 //the test is meant to run on five processors
839 if (size !=3) return ;
841 int nproc_source = 2;
843 set<int> procs_source;
844 set<int> procs_target;
846 for (int i=0; i<nproc_source; i++)
847 procs_source.insert(i);
848 for (int i=nproc_source; i<size; i++)
849 procs_target.insert(i);
850 self_procs.insert(rank);
852 MEDCoupling::CommInterface interface;
854 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
855 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
856 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
858 //loading the geometry for the source group
860 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
862 MEDCoupling::MEDCouplingUMesh* mesh;
863 MEDCoupling::ParaMESH* paramesh;
864 MEDCoupling::ParaFIELD* parafield;
865 ICoCo::MEDField* icocofield ;
867 char * tmp_dir_c = getenv("TMP");
869 if (tmp_dir_c != NULL)
870 tmp_dir = string(tmp_dir_c);
873 string filename_xml1 = "Mesh3D_10_2d";
874 string filename_xml2 = "Mesh3D_11";
875 //string filename_seq_wr = makeTmpFile("");
876 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
878 // To remove tmp files from disk
879 ParaMEDMEMTest_TmpFilesRemover aRemover;
881 MPI_Barrier(MPI_COMM_WORLD);
882 if (source_group->containsMyRank())
884 string master = filename_xml1;
886 ostringstream strstream;
887 strstream <<master<<rank+1<<".med";
888 std::string fName = INTERP_TEST::getResourceFile(strstream.str());
889 ostringstream meshname ;
890 meshname<< "Mesh_3_"<< rank+1;
892 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
895 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
897 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
898 MEDCoupling::ComponentTopology comptopo;
901 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
902 parafield->getField()->setNature(IntensiveMaximum);
905 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
908 nb_local=mesh->getNumberOfCells();
910 nb_local=mesh->getNumberOfNodes();
911 // double * value= new double[nb_local];
912 double *value=parafield->getField()->getArray()->getPointer();
913 for(int ielem=0; ielem<nb_local;ielem++)
916 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
917 icocofield=new ICoCo::MEDField(parafield->getField());
918 dec.setMethod(srcMeth);
919 dec.attachLocalField(icocofield);
922 //loading the geometry for the target group
923 if (target_group->containsMyRank())
925 string master= filename_xml2;
926 ostringstream strstream;
927 strstream << master << ".med";
928 std::string fName = INTERP_TEST::getResourceFile(strstream.str());
929 ostringstream meshname ;
931 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
933 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
934 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
935 MEDCoupling::ComponentTopology comptopo;
938 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
939 parafield->getField()->setNature(IntensiveMaximum);
942 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
945 nb_local=mesh->getNumberOfCells();
947 nb_local=mesh->getNumberOfNodes();
948 // double * value= new double[nb_local];
949 double *value=parafield->getField()->getArray()->getPointer();
950 for(int ielem=0; ielem<nb_local;ielem++)
952 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
953 icocofield=new ICoCo::MEDField(parafield->getField());
954 dec.setMethod(targetMeth);
955 dec.attachLocalField(icocofield);
957 //attaching a DEC to the source group
958 double field_before_int;
959 double field_after_int;
961 if (source_group->containsMyRank())
963 field_before_int = parafield->getVolumeIntegral(0,true);
965 cout<<"DEC usage"<<endl;
966 dec.setForcedRenormalization(false);
969 ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
970 if (source_group->myRank()==0)
971 aRemover.Register("./sourcesquareb");
972 ostringstream filename;
973 filename<<"./sourcesquareb_"<<source_group->myRank()+1;
974 aRemover.Register(filename.str().c_str());
975 //WriteField("./sourcesquareb",parafield->getField());
978 cout <<"writing"<<endl;
979 ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
980 if (source_group->myRank()==0)
981 aRemover.Register("./sourcesquare");
982 //WriteField("./sourcesquare",parafield->getField());
985 filename<<"./sourcesquare_"<<source_group->myRank()+1;
986 aRemover.Register(filename.str().c_str());
987 field_after_int = parafield->getVolumeIntegral(0,true);
989 CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
993 //attaching a DEC to the target group
994 if (target_group->containsMyRank())
997 dec.setForcedRenormalization(false);
1000 ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
1001 //WriteField("./targetsquareb",parafield->getField());
1002 if (target_group->myRank()==0)
1003 aRemover.Register("./targetsquareb");
1004 ostringstream filename;
1005 filename<<"./targetsquareb_"<<target_group->myRank()+1;
1006 aRemover.Register(filename.str().c_str());
1008 ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
1009 //WriteField("./targetsquare",parafield->getField());
1011 if (target_group->myRank()==0)
1012 aRemover.Register("./targetsquareb");
1014 filename<<"./targetsquareb_"<<target_group->myRank()+1;
1015 aRemover.Register(filename.str().c_str());
1017 delete source_group;
1018 delete target_group;
1026 MPI_Barrier(MPI_COMM_WORLD);
1027 cout << "end of InterpKernelDEC_3D test"<<endl;
1030 //Synchronous tests without interpolation with native mode (AllToAll(v) from lam/MPI:
1031 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpNativeDEC_2D()
1033 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,false,false,false,"P0","P0");
1036 //Synchronous tests without interpolation :
1037 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpDEC_2D()
1039 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,false,"P0","P0");
1042 //Synchronous tests with interpolation :
1043 void ParaMEDMEMTest::testSynchronousEqualInterpKernelDEC_2D()
1045 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,true,"P0","P0");
1047 void ParaMEDMEMTest::testSynchronousFasterSourceInterpKernelDEC_2D()
1049 testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,false,true,"P0","P0");
1051 void ParaMEDMEMTest::testSynchronousSlowerSourceInterpKernelDEC_2D()
1053 testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,false,true,"P0","P0");
1055 void ParaMEDMEMTest::testSynchronousSlowSourceInterpKernelDEC_2D()
1057 testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,false,true,"P0","P0");
1059 void ParaMEDMEMTest::testSynchronousFastSourceInterpKernelDEC_2D()
1061 testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,false,true,"P0","P0");
1064 //Asynchronous tests with interpolation :
1065 void ParaMEDMEMTest::testAsynchronousEqualInterpKernelDEC_2D()
1067 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,true,true,"P0","P0");
1069 void ParaMEDMEMTest::testAsynchronousFasterSourceInterpKernelDEC_2D()
1071 testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,true,true,"P0","P0");
1073 void ParaMEDMEMTest::testAsynchronousSlowerSourceInterpKernelDEC_2D()
1075 testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,true,true,"P0","P0");
1077 void ParaMEDMEMTest::testAsynchronousSlowSourceInterpKernelDEC_2D()
1079 testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,true,true,"P0","P0");
1081 void ParaMEDMEMTest::testAsynchronousFastSourceInterpKernelDEC_2D()
1083 testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,true,true,"P0","P0");
1086 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P0()
1089 const double sourceCoordsAll[2][8]={{0.4,0.5,0.4,1.5,1.6,1.5,1.6,0.5},
1090 {0.3,-0.5,1.6,-0.5,1.6,-1.5,0.3,-1.5}};
1091 const double targetCoordsAll[3][16]={{0.7,1.45,0.7,1.65,0.9,1.65,0.9,1.45, 1.1,1.4,1.1,1.6,1.3,1.6,1.3,1.4},
1092 {0.7,-0.6,0.7,0.7,0.9,0.7,0.9,-0.6, 1.1,-0.7,1.1,0.6,1.3,0.6,1.3,-0.7},
1093 {0.7,-1.55,0.7,-1.35,0.9,-1.35,0.9,-1.55, 1.1,-1.65,1.1,-1.45,1.3,-1.45,1.3,-1.65}};
1094 int conn4All[8]={0,1,2,3,4,5,6,7};
1095 double targetResults[3][2]={{34.,34.},{38.333333333333336,42.666666666666664},{47.,47.}};
1096 double targetResults2[3][2]={{0.28333333333333344,0.56666666666666687},{1.8564102564102569,2.0128205128205132},{1.0846153846153845,0.36153846153846159}};
1097 double targetResults3[3][2]={{3.7777777777777781,7.5555555555555562},{24.511111111111113,26.355555555555558},{14.1,4.7}};
1098 double targetResults4[3][2]={{8.5,17},{8.8461538461538431, 9.8461538461538449},{35.25,11.75}};
1102 MPI_Comm_size(MPI_COMM_WORLD,&size);
1103 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1107 int nproc_source = 2;
1108 set<int> self_procs;
1109 set<int> procs_source;
1110 set<int> procs_target;
1112 for (int i=0; i<nproc_source; i++)
1113 procs_source.insert(i);
1114 for (int i=nproc_source; i<size; i++)
1115 procs_target.insert(i);
1116 self_procs.insert(rank);
1118 MEDCoupling::MEDCouplingUMesh *mesh=0;
1119 MEDCoupling::ParaMESH *paramesh=0;
1120 MEDCoupling::ParaFIELD* parafield=0;
1122 MEDCoupling::CommInterface interface;
1124 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1125 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1126 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1128 MPI_Barrier(MPI_COMM_WORLD);
1129 if(source_group->containsMyRank())
1131 std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
1132 mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
1133 mesh->allocateCells(2);
1134 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
1135 mesh->finishInsertingCells();
1136 DataArrayDouble *myCoords=DataArrayDouble::New();
1137 myCoords->alloc(4,2);
1138 const double *sourceCoords=sourceCoordsAll[rank];
1139 std::copy(sourceCoords,sourceCoords+8,myCoords->getPointer());
1140 mesh->setCoords(myCoords);
1141 myCoords->decrRef();
1142 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1143 MEDCoupling::ComponentTopology comptopo;
1144 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1145 double *value=parafield->getField()->getArray()->getPointer();
1146 value[0]=34+13*((double)rank);
1150 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
1151 mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
1152 mesh->allocateCells(2);
1153 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
1154 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All+4);
1155 mesh->finishInsertingCells();
1156 DataArrayDouble *myCoords=DataArrayDouble::New();
1157 myCoords->alloc(8,2);
1158 const double *targetCoords=targetCoordsAll[rank-nproc_source];
1159 std::copy(targetCoords,targetCoords+16,myCoords->getPointer());
1160 mesh->setCoords(myCoords);
1161 myCoords->decrRef();
1162 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
1163 MEDCoupling::ComponentTopology comptopo;
1164 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1166 //test 1 - Conservative volumic
1167 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1168 parafield->getField()->setNature(IntensiveMaximum);
1169 if (source_group->containsMyRank())
1171 dec.setMethod("P0");
1172 dec.attachLocalField(parafield);
1174 dec.setForcedRenormalization(false);
1179 dec.setMethod("P0");
1180 dec.attachLocalField(parafield);
1182 dec.setForcedRenormalization(false);
1184 const double *res=parafield->getField()->getArray()->getConstPointer();
1185 const double *expected=targetResults[rank-nproc_source];
1186 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1187 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1189 //test 2 - ExtensiveMaximum
1190 MEDCoupling::InterpKernelDEC dec2(*source_group,*target_group);
1191 parafield->getField()->setNature(ExtensiveMaximum);
1192 if (source_group->containsMyRank())
1194 dec2.setMethod("P0");
1195 dec2.attachLocalField(parafield);
1197 dec2.setForcedRenormalization(false);
1202 dec2.setMethod("P0");
1203 dec2.attachLocalField(parafield);
1205 dec2.setForcedRenormalization(false);
1207 const double *res=parafield->getField()->getArray()->getConstPointer();
1208 const double *expected=targetResults2[rank-nproc_source];
1209 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1210 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1212 //test 3 - ExtensiveMaximum with global constraint
1213 MEDCoupling::InterpKernelDEC dec3(*source_group,*target_group);
1214 parafield->getField()->setNature(ExtensiveConservation);
1215 if (source_group->containsMyRank())
1217 dec3.setMethod("P0");
1218 dec3.attachLocalField(parafield);
1220 dec3.setForcedRenormalization(false);
1225 dec3.setMethod("P0");
1226 dec3.attachLocalField(parafield);
1228 dec3.setForcedRenormalization(false);
1230 const double *res=parafield->getField()->getArray()->getConstPointer();
1231 const double *expected=targetResults3[rank-nproc_source];
1232 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1233 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1235 //test 4 - IntensiveConservation
1236 MEDCoupling::InterpKernelDEC dec4(*source_group,*target_group);
1237 parafield->getField()->setNature(IntensiveConservation);
1238 if (source_group->containsMyRank())
1240 dec4.setMethod("P0");
1241 dec4.attachLocalField(parafield);
1243 dec4.setForcedRenormalization(false);
1248 dec4.setMethod("P0");
1249 dec4.attachLocalField(parafield);
1251 dec4.setForcedRenormalization(false);
1253 const double *res=parafield->getField()->getArray()->getConstPointer();
1254 const double *expected=targetResults4[rank-nproc_source];
1255 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1256 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1258 //test 5 - Conservative volumic reversed
1259 MEDCoupling::InterpKernelDEC dec5(*source_group,*target_group);
1260 parafield->getField()->setNature(IntensiveMaximum);
1261 if (source_group->containsMyRank())
1263 dec5.setMethod("P0");
1264 dec5.attachLocalField(parafield);
1266 dec5.setForcedRenormalization(false);
1268 const double *res=parafield->getField()->getArray()->getConstPointer();
1269 CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
1270 const double expected[]={37.8518518518519,43.5333333333333};
1271 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1275 dec5.setMethod("P0");
1276 dec5.attachLocalField(parafield);
1278 dec5.setForcedRenormalization(false);
1279 double *res=parafield->getField()->getArray()->getPointer();
1280 const double *toSet=targetResults[rank-nproc_source];
1285 //test 6 - ExtensiveMaximum reversed
1286 MEDCoupling::InterpKernelDEC dec6(*source_group,*target_group);
1287 parafield->getField()->setNature(ExtensiveMaximum);
1288 if (source_group->containsMyRank())
1290 dec6.setMethod("P0");
1291 dec6.attachLocalField(parafield);
1293 dec6.setForcedRenormalization(false);
1295 const double *res=parafield->getField()->getArray()->getConstPointer();
1296 CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
1297 const double expected[]={0.794600591715977,1.35631163708087};
1298 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1302 dec6.setMethod("P0");
1303 dec6.attachLocalField(parafield);
1305 dec6.setForcedRenormalization(false);
1306 double *res=parafield->getField()->getArray()->getPointer();
1307 const double *toSet=targetResults2[rank-nproc_source];
1312 //test 7 - ExtensiveMaximum with global constraint reversed
1313 MEDCoupling::InterpKernelDEC dec7(*source_group,*target_group);
1314 parafield->getField()->setNature(ExtensiveConservation);
1315 if (source_group->containsMyRank())
1317 dec7.setMethod("P0");
1318 dec7.attachLocalField(parafield);
1320 dec7.setForcedRenormalization(false);
1322 const double *res=parafield->getField()->getArray()->getConstPointer();
1323 CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
1324 const double expected[]={36.4592592592593,44.5407407407407};
1325 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1329 dec7.setMethod("P0");
1330 dec7.attachLocalField(parafield);
1332 dec7.setForcedRenormalization(false);
1333 double *res=parafield->getField()->getArray()->getPointer();
1334 const double *toSet=targetResults3[rank-nproc_source];
1339 //test 8 - ExtensiveMaximum with IntensiveConservation reversed
1340 MEDCoupling::InterpKernelDEC dec8(*source_group,*target_group);
1341 parafield->getField()->setNature(IntensiveConservation);
1342 if (source_group->containsMyRank())
1344 dec8.setMethod("P0");
1345 dec8.attachLocalField(parafield);
1347 dec8.setForcedRenormalization(false);
1349 const double *res=parafield->getField()->getArray()->getConstPointer();
1350 CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
1351 const double expected[]={0.81314102564102553,1.3428994082840233};
1352 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1356 dec8.setMethod("P0");
1357 dec8.attachLocalField(parafield);
1359 dec8.setForcedRenormalization(false);
1360 double *res=parafield->getField()->getArray()->getPointer();
1361 const double *toSet=targetResults4[rank-nproc_source];
1371 delete target_group;
1372 delete source_group;
1374 MPI_Barrier(MPI_COMM_WORLD);
1377 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P1P1P0()
1381 MPI_Comm_size(MPI_COMM_WORLD,&size);
1382 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1386 int nproc_source = 2;
1387 set<int> self_procs;
1388 set<int> procs_source;
1389 set<int> procs_target;
1391 for (int i=0; i<nproc_source; i++)
1392 procs_source.insert(i);
1393 for (int i=nproc_source; i<size; i++)
1394 procs_target.insert(i);
1395 self_procs.insert(rank);
1397 MEDCoupling::MEDCouplingUMesh *mesh=0;
1398 MEDCoupling::ParaMESH *paramesh=0;
1399 MEDCoupling::ParaFIELD *parafieldP0=0,*parafieldP1=0;
1401 MEDCoupling::CommInterface interface;
1403 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1404 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1405 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1407 MPI_Barrier(MPI_COMM_WORLD);
1408 if(source_group->containsMyRank())
1412 double coords[6]={-0.3,-0.3, 0.7,0.7, 0.7,-0.3};
1413 int conn[3]={0,1,2};
1414 //int globalNode[3]={1,2,0};
1415 mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
1416 mesh->allocateCells(1);
1417 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
1418 mesh->finishInsertingCells();
1419 DataArrayDouble *myCoords=DataArrayDouble::New();
1420 myCoords->alloc(3,2);
1421 std::copy(coords,coords+6,myCoords->getPointer());
1422 mesh->setCoords(myCoords);
1423 myCoords->decrRef();
1427 double coords[6]={-0.3,-0.3, -0.3,0.7, 0.7,0.7};
1428 int conn[3]={0,1,2};
1429 //int globalNode[3]={1,3,2};
1430 mesh=MEDCouplingUMesh::New("Source mesh Proc1",2);
1431 mesh->allocateCells(1);
1432 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
1433 mesh->finishInsertingCells();
1434 DataArrayDouble *myCoords=DataArrayDouble::New();
1435 myCoords->alloc(3,2);
1436 std::copy(coords,coords+6,myCoords->getPointer());
1437 mesh->setCoords(myCoords);
1438 myCoords->decrRef();
1440 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1441 MEDCoupling::ComponentTopology comptopo;
1442 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1443 parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
1444 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1445 double *valueP1=parafieldP1->getField()->getArray()->getPointer();
1446 parafieldP0->getField()->setNature(IntensiveMaximum);
1447 parafieldP1->getField()->setNature(IntensiveMaximum);
1451 valueP1[0]=34.; valueP1[1]=77.; valueP1[2]=53.;
1456 valueP1[0]=34.; valueP1[1]=57.; valueP1[2]=77.;
1461 const char targetMeshName[]="target mesh";
1464 double coords[10]={-0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2 };
1465 int conn[7]={0,3,4,1, 1,4,2};
1466 //int globalNode[5]={4,3,0,2,1};
1467 mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
1468 mesh->allocateCells(2);
1469 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1470 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
1471 mesh->finishInsertingCells();
1472 DataArrayDouble *myCoords=DataArrayDouble::New();
1473 myCoords->alloc(5,2);
1474 std::copy(coords,coords+10,myCoords->getPointer());
1475 mesh->setCoords(myCoords);
1476 myCoords->decrRef();
1477 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
1478 DataArrayInt *da=DataArrayInt::New();
1479 const int globalNumberingP2[5]={0,1,2,3,4};
1480 da->useArray(globalNumberingP2,false,CPP_DEALLOC,5,1);
1481 paramesh->setNodeGlobal(da);
1486 double coords[6]={0.2,0.2, 0.7,-0.3, 0.7,0.2};
1487 int conn[3]={0,2,1};
1488 //int globalNode[3]={1,0,5};
1489 mesh=MEDCouplingUMesh::New("Target mesh Proc3",2);
1490 mesh->allocateCells(1);
1491 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
1492 mesh->finishInsertingCells();
1493 DataArrayDouble *myCoords=DataArrayDouble::New();
1494 myCoords->alloc(3,2);
1495 std::copy(coords,coords+6,myCoords->getPointer());
1496 mesh->setCoords(myCoords);
1497 myCoords->decrRef();
1498 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
1499 DataArrayInt *da=DataArrayInt::New();
1500 const int globalNumberingP3[3]={4,2,5};
1501 da->useArray(globalNumberingP3,false,CPP_DEALLOC,3,1);
1502 paramesh->setNodeGlobal(da);
1507 double coords[12]={-0.3,0.2, -0.3,0.7, 0.2,0.7, 0.2,0.2, 0.7,0.7, 0.7,0.2};
1508 int conn[8]={0,1,2,3, 3,2,4,5};
1509 //int globalNode[6]={2,6,7,1,8,5};
1510 mesh=MEDCouplingUMesh::New("Target mesh Proc4",2);
1511 mesh->allocateCells(2);
1512 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1513 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
1514 mesh->finishInsertingCells();
1515 DataArrayDouble *myCoords=DataArrayDouble::New();
1516 myCoords->alloc(6,2);
1517 std::copy(coords,coords+12,myCoords->getPointer());
1518 mesh->setCoords(myCoords);
1519 myCoords->decrRef();
1520 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
1521 DataArrayInt *da=DataArrayInt::New();
1522 const int globalNumberingP4[6]={3,6,7,4,8,5};
1523 da->useArray(globalNumberingP4,false,CPP_DEALLOC,6,1);
1524 paramesh->setNodeGlobal(da);
1527 MEDCoupling::ComponentTopology comptopo;
1528 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1529 parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
1530 parafieldP0->getField()->setNature(IntensiveMaximum);
1531 parafieldP1->getField()->setNature(IntensiveMaximum);
1534 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1535 if (source_group->containsMyRank())
1537 dec.setMethod("P0");
1538 dec.attachLocalField(parafieldP0);
1540 dec.setForcedRenormalization(false);
1543 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1546 CPPUNIT_ASSERT_DOUBLES_EQUAL(34.42857143,valueP0[0],1e-7);
1550 CPPUNIT_ASSERT_DOUBLES_EQUAL(44.,valueP0[0],1e-7);
1555 dec.setMethod("P1");
1556 dec.attachLocalField(parafieldP1);
1558 dec.setForcedRenormalization(false);
1560 const double *res=parafieldP1->getField()->getArray()->getConstPointer();
1563 const double expectP2[5]={39.0, 31.0, 31.0, 47.0, 39.0};
1564 CPPUNIT_ASSERT_EQUAL(5,parafieldP1->getField()->getNumberOfTuples());
1565 CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
1566 for(int kk=0;kk<5;kk++)
1567 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP2[kk],res[kk],1e-12);
1571 const double expectP3[3]={39.0, 31.0, 31.0};
1572 CPPUNIT_ASSERT_EQUAL(3,parafieldP1->getField()->getNumberOfTuples());
1573 CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
1574 for(int kk=0;kk<3;kk++)
1575 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP3[kk],res[kk],1e-12);
1579 const double expectP4[6]={47.0, 47.0, 47.0, 39.0, 39.0, 31.0};
1580 CPPUNIT_ASSERT_EQUAL(6,parafieldP1->getField()->getNumberOfTuples());
1581 CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
1582 for(int kk=0;kk<6;kk++)
1583 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP4[kk],res[kk],1e-12);
1593 delete target_group;
1594 delete source_group;
1596 MPI_Barrier(MPI_COMM_WORLD);
1599 void ParaMEDMEMTest::testInterpKernelDEC2DM1D_P0P0()
1603 MPI_Comm_size(MPI_COMM_WORLD,&size);
1604 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1609 set<int> procs_source;
1610 set<int> procs_target;
1612 for (int i=0; i<nproc_source; i++)
1613 procs_source.insert(i);
1614 for (int i=nproc_source;i<size; i++)
1615 procs_target.insert(i);
1617 MEDCoupling::MEDCouplingUMesh *mesh=0;
1618 MEDCoupling::ParaMESH *paramesh=0;
1619 MEDCoupling::ParaFIELD *parafield=0;
1621 MEDCoupling::CommInterface interface;
1623 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1624 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1626 MPI_Barrier(MPI_COMM_WORLD);
1627 if(source_group->containsMyRank())
1629 double targetCoords[18]={-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 };
1630 mesh=MEDCouplingUMesh::New();
1631 mesh->setMeshDimension(2);
1632 DataArrayDouble *myCoords=DataArrayDouble::New();
1633 myCoords->alloc(9,2);
1634 std::copy(targetCoords,targetCoords+18,myCoords->getPointer());
1635 mesh->setCoords(myCoords);
1636 myCoords->decrRef();
1639 int targetConn[7]={0,3,4,1, 1,4,2};
1640 mesh->allocateCells(2);
1641 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
1642 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+4);
1643 mesh->finishInsertingCells();
1647 int targetConn[11]={4,5,2, 6,7,4,3, 7,8,5,4};
1648 mesh->allocateCells(3);
1649 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1650 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+3);
1651 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+7);
1652 mesh->finishInsertingCells();
1654 MEDCoupling::ComponentTopology comptopo;
1655 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1656 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1657 parafield->getField()->setNature(IntensiveMaximum);
1658 double *vals=parafield->getField()->getArray()->getPointer();
1660 { vals[0]=7.; vals[1]=8.; }
1662 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1666 mesh=MEDCouplingUMesh::New("an example of -1 D mesh",-1);
1667 MEDCoupling::ComponentTopology comptopo;
1668 paramesh=new ParaMESH(mesh,*target_group,"target mesh");
1669 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1670 parafield->getField()->setNature(IntensiveMaximum);
1672 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1673 if(source_group->containsMyRank())
1675 dec.setMethod("P0");
1676 dec.attachLocalField(parafield);
1678 dec.setForcedRenormalization(false);
1681 const double *res=parafield->getField()->getArray()->getConstPointer();
1684 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1685 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1689 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1690 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1691 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
1696 dec.setMethod("P0");
1697 dec.attachLocalField(parafield);
1699 dec.setForcedRenormalization(false);
1701 const double *res=parafield->getField()->getArray()->getConstPointer();
1702 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1705 MEDCoupling::InterpKernelDEC dec2(*source_group,*target_group);
1706 dec2.setMethod("P0");
1707 parafield->getField()->setNature(ExtensiveConservation);
1708 if(source_group->containsMyRank())
1710 double *vals=parafield->getField()->getArray()->getPointer();
1712 { vals[0]=7.; vals[1]=8.; }
1714 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1715 dec2.attachLocalField(parafield);
1719 const double *res=parafield->getField()->getArray()->getConstPointer();
1722 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
1723 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
1727 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
1728 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
1729 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
1734 dec2.attachLocalField(parafield);
1737 const double *res=parafield->getField()->getArray()->getConstPointer();
1738 CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
1742 MEDCoupling::InterpKernelDEC dec3(*source_group,*target_group);
1743 dec3.setMethod("P0");
1744 parafield->getField()->setNature(ExtensiveMaximum);
1745 if(source_group->containsMyRank())
1747 double *vals=parafield->getField()->getArray()->getPointer();
1749 { vals[0]=7.; vals[1]=8.; }
1751 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1752 dec3.attachLocalField(parafield);
1756 const double *res=parafield->getField()->getArray()->getConstPointer();
1759 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
1760 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
1764 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
1765 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
1766 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
1771 dec3.attachLocalField(parafield);
1774 const double *res=parafield->getField()->getArray()->getConstPointer();
1775 CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
1779 MEDCoupling::InterpKernelDEC dec4(*source_group,*target_group);
1780 dec4.setMethod("P0");
1781 parafield->getField()->setNature(IntensiveConservation);
1782 if(source_group->containsMyRank())
1784 double *vals=parafield->getField()->getArray()->getPointer();
1786 { vals[0]=7.; vals[1]=8.; }
1788 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1789 dec4.attachLocalField(parafield);
1793 const double *res=parafield->getField()->getArray()->getConstPointer();
1796 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1797 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1801 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1802 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1803 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
1808 dec4.attachLocalField(parafield);
1811 const double *res=parafield->getField()->getArray()->getConstPointer();
1812 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1818 delete target_group;
1819 delete source_group;
1821 MPI_Barrier(MPI_COMM_WORLD);
1824 void ParaMEDMEMTest::testInterpKernelDECPartialProcs()
1828 MPI_Comm_size(MPI_COMM_WORLD,&size);
1829 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1833 set<int> procs_source;
1834 set<int> procs_target;
1836 procs_source.insert(0);
1837 procs_target.insert(1);
1839 MEDCoupling::MEDCouplingUMesh *mesh=0;
1840 MEDCoupling::ParaMESH *paramesh=0;
1841 MEDCoupling::ParaFIELD *parafield=0;
1843 MEDCoupling::CommInterface interface;
1845 MPI_Barrier(MPI_COMM_WORLD);
1846 double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
1848 int grpIds[2]={0,1};
1849 MPI_Group grp,group_world;
1850 comm.commGroup(MPI_COMM_WORLD,&group_world);
1851 comm.groupIncl(group_world,2,grpIds,&grp);
1852 MPI_Comm partialComm;
1853 comm.commCreate(MPI_COMM_WORLD,grp,&partialComm);
1855 ProcessorGroup* target_group=0;
1856 ProcessorGroup* source_group=0;
1858 MEDCoupling::InterpKernelDEC *dec=0;
1859 if(rank==0 || rank==1)
1861 target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target,partialComm);
1862 source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source,partialComm);
1863 if(source_group->containsMyRank())
1865 mesh=MEDCouplingUMesh::New();
1866 mesh->setMeshDimension(2);
1867 DataArrayDouble *myCoords=DataArrayDouble::New();
1868 myCoords->alloc(4,2);
1869 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
1870 mesh->setCoords(myCoords);
1871 myCoords->decrRef();
1872 int targetConn[4]={0,2,3,1};
1873 mesh->allocateCells(1);
1874 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
1875 mesh->finishInsertingCells();
1876 MEDCoupling::ComponentTopology comptopo;
1877 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1878 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1879 parafield->getField()->setNature(IntensiveMaximum);
1880 double *vals=parafield->getField()->getArray()->getPointer();
1882 dec=new MEDCoupling::InterpKernelDEC(*source_group,*target_group);
1883 dec->attachLocalField(parafield);
1890 mesh=MEDCouplingUMesh::New();
1891 mesh->setMeshDimension(2);
1892 DataArrayDouble *myCoords=DataArrayDouble::New();
1893 myCoords->alloc(4,2);
1894 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
1895 mesh->setCoords(myCoords);
1896 myCoords->decrRef();
1897 int targetConn[6]={0,2,1,2,3,1};
1898 mesh->allocateCells(2);
1899 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1900 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
1901 mesh->finishInsertingCells();
1902 MEDCoupling::ComponentTopology comptopo;
1903 paramesh=new ParaMESH(mesh,*target_group,"target mesh");
1904 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1905 parafield->getField()->setNature(IntensiveMaximum);
1906 dec=new MEDCoupling::InterpKernelDEC(*source_group,*target_group);
1907 dec->attachLocalField(parafield);
1917 delete target_group;
1918 delete source_group;
1920 if(partialComm != MPI_COMM_NULL)
1921 comm.commFree(&partialComm);
1922 comm.groupFree(&grp);
1923 comm.groupFree(&group_world);
1924 MPI_Barrier(MPI_COMM_WORLD);
1928 * This test reproduces bug of Gauthier on 13/9/2010 concerning 3DSurf meshes.
1929 * It is possible to lead to dead lock in InterpKernelDEC when 3DSurfMeshes global bounding boxes intersects whereas cell bounding box intersecting only on one side.
1931 void ParaMEDMEMTest::testInterpKernelDEC3DSurfEmptyBBox()
1935 MPI_Comm_size(MPI_COMM_WORLD,&size);
1936 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1940 int nproc_source = 1;
1941 set<int> self_procs;
1942 set<int> procs_source;
1943 set<int> procs_target;
1945 for (int i=0; i<nproc_source; i++)
1946 procs_source.insert(i);
1947 for (int i=nproc_source; i<size; i++)
1948 procs_target.insert(i);
1949 self_procs.insert(rank);
1951 MEDCoupling::MEDCouplingUMesh *mesh=0;
1952 MEDCoupling::ParaMESH *paramesh=0;
1953 MEDCoupling::ParaFIELD *parafieldP0=0;
1955 MEDCoupling::CommInterface interface;
1957 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1958 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1959 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1961 MPI_Barrier(MPI_COMM_WORLD);
1962 if(source_group->containsMyRank())
1964 double coords[15]={1.,0.,0., 2.,0.,0., 2.,2.,0., 0.,2.,0., 0.5,0.5,1.};
1965 int conn[7]={0,1,2,3,0,3,4};
1966 mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
1967 mesh->allocateCells(2);
1968 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1969 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
1970 mesh->finishInsertingCells();
1971 DataArrayDouble *myCoords=DataArrayDouble::New();
1972 myCoords->alloc(5,3);
1973 std::copy(coords,coords+15,myCoords->getPointer());
1974 mesh->setCoords(myCoords);
1975 myCoords->decrRef();
1977 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1978 MEDCoupling::ComponentTopology comptopo;
1979 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1980 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1981 parafieldP0->getField()->setNature(IntensiveMaximum);
1982 valueP0[0]=7.; valueP0[1]=8.;
1986 const char targetMeshName[]="target mesh";
1989 double coords[12]={0.25,0.25,0.5, 0.,0.25,0.5, 0.,0.,0.5, 0.25,0.,0.5};
1990 int conn[4]={0,1,2,3};
1991 mesh=MEDCouplingUMesh::New("Target mesh Proc1",2);
1992 mesh->allocateCells(1);
1993 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1994 mesh->finishInsertingCells();
1995 DataArrayDouble *myCoords=DataArrayDouble::New();
1996 myCoords->alloc(4,3);
1997 std::copy(coords,coords+12,myCoords->getPointer());
1998 mesh->setCoords(myCoords);
1999 myCoords->decrRef();
2000 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
2004 double coords[12]={0.,0.25,0.5, 0.,0.,0.5, -1.,0.,0.5, -1.,0.25,0.5};
2005 int conn[4]={0,1,2,3};
2006 mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
2007 mesh->allocateCells(1);
2008 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
2009 mesh->finishInsertingCells();
2010 DataArrayDouble *myCoords=DataArrayDouble::New();
2011 myCoords->alloc(4,3);
2012 std::copy(coords,coords+12,myCoords->getPointer());
2013 mesh->setCoords(myCoords);
2014 myCoords->decrRef();
2015 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
2017 MEDCoupling::ComponentTopology comptopo;
2018 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2019 parafieldP0->getField()->setNature(IntensiveMaximum);
2022 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
2023 if (source_group->containsMyRank())
2025 dec.setMethod("P0");
2026 dec.attachLocalField(parafieldP0);
2028 // dec.setForcedRenormalization(false);
2031 // const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
2034 // CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
2035 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
2039 // CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
2043 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
2048 dec.setMethod("P0");
2049 dec.attachLocalField(parafieldP0);
2051 // dec.setForcedRenormalization(false);
2053 // const double *res=parafieldP0->getField()->getArray()->getConstPointer();
2056 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
2057 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
2058 // CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
2062 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
2063 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
2064 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
2073 delete target_group;
2074 delete source_group;
2076 MPI_Barrier(MPI_COMM_WORLD);
2080 * Tests an asynchronous exchange between two codes
2081 * one sends data with dtA as an interval, the max time being tmaxA
2082 * the other one receives with dtB as an interval, the max time being tmaxB
2084 void ParaMEDMEMTest::testAsynchronousInterpKernelDEC_2D(double dtA, double tmaxA,
2085 double dtB, double tmaxB, bool WithPointToPoint, bool Asynchronous,
2086 bool WithInterp, const char *srcMeth, const char *targetMeth)
2088 std::string srcM(srcMeth);
2089 std::string targetM(targetMeth);
2092 MPI_Comm_size(MPI_COMM_WORLD,&size);
2093 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
2095 //the test is meant to run on five processors
2096 if (size !=5) return ;
2098 int nproc_source = 3;
2099 set<int> self_procs;
2100 set<int> procs_source;
2101 set<int> procs_target;
2103 for (int i=0; i<nproc_source; i++)
2104 procs_source.insert(i);
2105 for (int i=nproc_source; i<size; i++)
2106 procs_target.insert(i);
2107 self_procs.insert(rank);
2109 MEDCoupling::CommInterface interface;
2111 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
2112 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
2113 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
2115 //loading the geometry for the source group
2117 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
2119 MEDCoupling::MEDCouplingUMesh* mesh;
2120 MEDCoupling::ParaMESH* paramesh;
2121 MEDCoupling::ParaFIELD* parafield;
2123 ICoCo::MEDField* icocofield ;
2125 char * tmp_dir_c = getenv("TMP");
2127 if (tmp_dir_c != NULL)
2128 tmp_dir = string(tmp_dir_c);
2131 string filename_xml1 = "square1_split";
2132 string filename_xml2 = "square2_split";
2133 //string filename_seq_wr = makeTmpFile("");
2134 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
2136 // To remove tmp files from disk
2137 ParaMEDMEMTest_TmpFilesRemover aRemover;
2139 MPI_Barrier(MPI_COMM_WORLD);
2141 if (source_group->containsMyRank())
2143 string master = filename_xml1;
2145 ostringstream strstream;
2146 strstream <<master<<rank+1<<".med";
2147 string fName = INTERP_TEST::getResourceFile(strstream.str());
2148 ostringstream meshname ;
2149 meshname<< "Mesh_2_"<< rank+1;
2151 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
2153 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
2155 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
2156 MEDCoupling::ComponentTopology comptopo;
2159 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2160 parafield->getField()->setNature(IntensiveMaximum);//InvertIntegral);//IntensiveMaximum);
2163 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
2167 nb_local=mesh->getNumberOfCells();
2169 nb_local=mesh->getNumberOfNodes();
2170 // double * value= new double[nb_local];
2171 double *value=parafield->getField()->getArray()->getPointer();
2172 for(int ielem=0; ielem<nb_local;ielem++)
2175 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
2176 icocofield=new ICoCo::MEDField(parafield->getField());
2178 dec.attachLocalField(icocofield);
2183 //loading the geometry for the target group
2184 if (target_group->containsMyRank())
2186 string master= filename_xml2;
2187 ostringstream strstream;
2188 strstream << master<<(rank-nproc_source+1)<<".med";
2189 string fName = INTERP_TEST::getResourceFile(strstream.str());
2190 ostringstream meshname ;
2191 meshname<< "Mesh_3_"<<rank-nproc_source+1;
2193 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
2195 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
2196 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
2197 MEDCoupling::ComponentTopology comptopo;
2200 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2201 parafield->getField()->setNature(IntensiveMaximum);//InvertIntegral);//IntensiveMaximum);
2204 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
2208 nb_local=mesh->getNumberOfCells();
2210 nb_local=mesh->getNumberOfNodes();
2212 double *value=parafield->getField()->getArray()->getPointer();
2213 for(int ielem=0; ielem<nb_local;ielem++)
2215 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
2216 icocofield=new ICoCo::MEDField(parafield->getField());
2218 dec.attachLocalField(icocofield);
2222 //attaching a DEC to the source group
2224 if (source_group->containsMyRank())
2226 cout<<"DEC usage"<<endl;
2227 dec.setAsynchronous(Asynchronous);
2229 dec.setTimeInterpolationMethod(LinearTimeInterp);
2231 if ( WithPointToPoint ) {
2232 dec.setAllToAllMethod(PointToPoint);
2235 dec.setAllToAllMethod(Native);
2238 dec.setForcedRenormalization(false);
2239 for (double time=0; time<tmaxA+1e-10; time+=dtA)
2241 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2242 << " dtA " << dtA << " tmaxA " << tmaxA << endl ;
2243 if ( time+dtA < tmaxA+1e-7 ) {
2244 dec.sendData( time , dtA );
2247 dec.sendData( time , 0 );
2249 double* value = parafield->getField()->getArray()->getPointer();
2250 int nb_local=parafield->getField()->getMesh()->getNumberOfCells();
2251 for (int i=0; i<nb_local;i++)
2258 //attaching a DEC to the target group
2259 if (target_group->containsMyRank())
2261 cout<<"DEC usage"<<endl;
2262 dec.setAsynchronous(Asynchronous);
2264 dec.setTimeInterpolationMethod(LinearTimeInterp);
2266 if ( WithPointToPoint ) {
2267 dec.setAllToAllMethod(PointToPoint);
2270 dec.setAllToAllMethod(Native);
2273 dec.setForcedRenormalization(false);
2274 vector<double> times;
2275 for (double time=0; time<tmaxB+1e-10; time+=dtB)
2277 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2278 << " dtB " << dtB << " tmaxB " << tmaxB << endl ;
2279 dec.recvData( time );
2280 double vi = parafield->getVolumeIntegral(0,true);
2281 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2282 << " VolumeIntegral " << vi
2283 << " time*10000 " << time*10000 << endl ;
2285 CPPUNIT_ASSERT_DOUBLES_EQUAL(vi,time*10000,0.001);
2290 delete source_group;
2291 delete target_group;
2298 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " MPI_Barrier " << endl ;
2300 if (Asynchronous) MPI_Barrier(MPI_COMM_WORLD);
2301 cout << "end of InterpKernelDEC_2D test"<<endl;