1 // Copyright (C) 2007-2024 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 <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 "ParaMEDLoader.hxx"
34 #include "ICoCoMEDDoubleField.hxx"
35 #include "MEDLoader.hxx"
36 #include "MEDCouplingUMesh.hxx"
37 #include "TestInterpKernelUtils.hxx"
43 // use this define to enable lines, execution of which leads to Segmentation Fault
46 // use this define to enable CPPUNIT asserts and fails, showing bugs
47 #define ENABLE_FORCED_FAILURES
51 using namespace MEDCoupling;
53 void ParaMEDMEMTest::testInterpKernelDEC_2D()
55 testInterpKernelDEC_2D_("P0","P0");
58 void ParaMEDMEMTest::testInterpKernelDEC2_2D()
60 testInterpKernelDEC2_2D_("P0","P0");
63 void ParaMEDMEMTest::testInterpKernelDEC_3D()
65 testInterpKernelDEC_3D_("P0","P0");
68 void ParaMEDMEMTest::testInterpKernelDEC_2DP0P1()
70 //testInterpKernelDEC_2D_("P0","P1");
73 void ParaMEDMEMTest::testInterpKernelDEC_1D()
77 MPI_Comm_size(MPI_COMM_WORLD,&size);
78 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
84 set<int> procs_source;
85 set<int> procs_target;
87 for (int i=0; i<nproc_source; i++)
88 procs_source.insert(i);
89 for (int i=nproc_source; i<size; i++)
90 procs_target.insert(i);
91 self_procs.insert(rank);
93 MEDCoupling::MEDCouplingUMesh *mesh=0;
94 MEDCoupling::ParaMESH *paramesh=0;
95 MEDCoupling::ParaFIELD *parafieldP0=0;
97 MEDCoupling::CommInterface interface;
99 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
100 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
101 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
103 MPI_Barrier(MPI_COMM_WORLD);
104 if(source_group->containsMyRank())
108 double coords[4]={0.3,0.7, 0.9,1.0};
109 mcIdType conn[4]={0,1,2,3};
110 mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
111 mesh->allocateCells(2);
112 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
113 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
114 mesh->finishInsertingCells();
115 DataArrayDouble *myCoords=DataArrayDouble::New();
116 myCoords->alloc(4,1);
117 std::copy(coords,coords+4,myCoords->getPointer());
118 mesh->setCoords(myCoords);
123 double coords[2]={0.7,0.9};
124 mcIdType conn[2]={0,1};
125 mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
126 mesh->allocateCells(1);
127 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
128 mesh->finishInsertingCells();
129 DataArrayDouble *myCoords=DataArrayDouble::New();
130 myCoords->alloc(2,1);
131 std::copy(coords,coords+2,myCoords->getPointer());
132 mesh->setCoords(myCoords);
137 double coords[2]={1.,1.12};
138 mcIdType conn[2]={0,1};
139 mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
140 mesh->allocateCells(1);
141 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
142 mesh->finishInsertingCells();
143 DataArrayDouble *myCoords=DataArrayDouble::New();
144 myCoords->alloc(2,1);
145 std::copy(coords,coords+2,myCoords->getPointer());
146 mesh->setCoords(myCoords);
149 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
150 MEDCoupling::ComponentTopology comptopo;
151 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
152 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
153 parafieldP0->getField()->setNature(IntensiveMaximum);
156 valueP0[0]=7.; valueP0[1]=8.;
169 const char targetMeshName[]="target mesh";
172 double coords[2]={0.5,0.75};
173 mcIdType conn[2]={0,1};
174 mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
175 mesh->allocateCells(1);
176 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
177 mesh->finishInsertingCells();
178 DataArrayDouble *myCoords=DataArrayDouble::New();
179 myCoords->alloc(2,1);
180 std::copy(coords,coords+2,myCoords->getPointer());
181 mesh->setCoords(myCoords);
183 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
187 double coords[2]={0.75,1.2};
188 mcIdType conn[2]={0,1};
189 mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
190 mesh->allocateCells(1);
191 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
192 mesh->finishInsertingCells();
193 DataArrayDouble *myCoords=DataArrayDouble::New();
194 myCoords->alloc(2,1);
195 std::copy(coords,coords+2,myCoords->getPointer());
196 mesh->setCoords(myCoords);
198 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
200 MEDCoupling::ComponentTopology comptopo;
201 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
202 parafieldP0->getField()->setNature(IntensiveMaximum);
205 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
206 if (source_group->containsMyRank())
209 dec.attachLocalField(parafieldP0);
211 dec.setForcedRenormalization(false);
214 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
217 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
218 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
222 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
226 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
232 dec.attachLocalField(parafieldP0);
234 dec.setForcedRenormalization(false);
236 const double *res=parafieldP0->getField()->getArray()->getConstPointer();
239 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfTuples());
240 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfComponents());
241 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
245 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfTuples());
246 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfComponents());
247 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
259 MPI_Barrier(MPI_COMM_WORLD);
262 void ParaMEDMEMTest::testInterpKernelDEC_2DCurve()
266 MPI_Comm_size(MPI_COMM_WORLD,&size);
267 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
271 int nproc_source = 3;
273 set<int> procs_source;
274 set<int> procs_target;
276 for (int i=0; i<nproc_source; i++)
277 procs_source.insert(i);
278 for (int i=nproc_source; i<size; i++)
279 procs_target.insert(i);
280 self_procs.insert(rank);
282 MEDCoupling::MEDCouplingUMesh *mesh=0;
283 MEDCoupling::ParaMESH *paramesh=0;
284 MEDCoupling::ParaFIELD *parafieldP0=0;
286 MEDCoupling::CommInterface interface;
288 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
289 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
290 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
292 MPI_Barrier(MPI_COMM_WORLD);
293 if(source_group->containsMyRank())
297 double coords[8]={0.3,0.3,0.7,0.7, 0.9,0.9,1.0,1.0};
298 mcIdType conn[4]={0,1,2,3};
299 mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
300 mesh->allocateCells(2);
301 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
302 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
303 mesh->finishInsertingCells();
304 DataArrayDouble *myCoords=DataArrayDouble::New();
305 myCoords->alloc(4,2);
306 std::copy(coords,coords+8,myCoords->getPointer());
307 mesh->setCoords(myCoords);
312 double coords[4]={0.7,0.7,0.9,0.9};
313 mcIdType conn[2]={0,1};
314 mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
315 mesh->allocateCells(1);
316 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
317 mesh->finishInsertingCells();
318 DataArrayDouble *myCoords=DataArrayDouble::New();
319 myCoords->alloc(2,2);
320 std::copy(coords,coords+4,myCoords->getPointer());
321 mesh->setCoords(myCoords);
326 double coords[4]={1.,1.,1.12,1.12};
327 mcIdType conn[2]={0,1};
328 mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
329 mesh->allocateCells(1);
330 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
331 mesh->finishInsertingCells();
332 DataArrayDouble *myCoords=DataArrayDouble::New();
333 myCoords->alloc(2,2);
334 std::copy(coords,coords+4,myCoords->getPointer());
335 mesh->setCoords(myCoords);
338 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
339 MEDCoupling::ComponentTopology comptopo;
340 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
341 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
342 parafieldP0->getField()->setNature(IntensiveMaximum);
345 valueP0[0]=7.; valueP0[1]=8.;
358 const char targetMeshName[]="target mesh";
361 double coords[4]={0.5,0.5,0.75,0.75};
362 mcIdType conn[2]={0,1};
363 mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
364 mesh->allocateCells(1);
365 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
366 mesh->finishInsertingCells();
367 DataArrayDouble *myCoords=DataArrayDouble::New();
368 myCoords->alloc(2,2);
369 std::copy(coords,coords+4,myCoords->getPointer());
370 mesh->setCoords(myCoords);
372 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
376 double coords[4]={0.75,0.75,1.2,1.2};
377 mcIdType conn[2]={0,1};
378 mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
379 mesh->allocateCells(1);
380 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
381 mesh->finishInsertingCells();
382 DataArrayDouble *myCoords=DataArrayDouble::New();
383 myCoords->alloc(2,2);
384 std::copy(coords,coords+4,myCoords->getPointer());
385 mesh->setCoords(myCoords);
387 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
389 MEDCoupling::ComponentTopology comptopo;
390 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
391 parafieldP0->getField()->setNature(IntensiveMaximum);
394 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
395 if (source_group->containsMyRank())
398 dec.attachLocalField(parafieldP0);
400 dec.setForcedRenormalization(false);
403 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
406 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
407 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
411 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
415 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
421 dec.attachLocalField(parafieldP0);
423 dec.setForcedRenormalization(false);
425 const double *res=parafieldP0->getField()->getArray()->getConstPointer();
428 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfTuples());
429 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfComponents());
430 CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
434 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfTuples());
435 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP0->getField()->getNumberOfComponents());
436 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
448 MPI_Barrier(MPI_COMM_WORLD);
453 * Check methods defined in InterpKernelDEC.hxx
456 InterpKernelDEC(ProcessorGroup& local_group, ProcessorGroup& distant_group);
457 virtual ~InterpKernelDEC();
463 void ParaMEDMEMTest::testInterpKernelDEC_2D_(const char *srcMeth, const char *targetMeth)
465 std::string srcM(srcMeth);
466 std::string targetM(targetMeth);
469 MPI_Comm_size(MPI_COMM_WORLD,&size);
470 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
472 //the test is meant to run on five processors
473 if (size !=5) return ;
475 int nproc_source = 3;
477 set<int> procs_source;
478 set<int> procs_target;
480 for (int i=0; i<nproc_source; i++)
481 procs_source.insert(i);
482 for (int i=nproc_source; i<size; i++)
483 procs_target.insert(i);
484 self_procs.insert(rank);
486 MEDCoupling::CommInterface interface;
488 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
489 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
490 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
492 //loading the geometry for the source group
494 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
496 MEDCoupling::MEDCouplingUMesh* mesh = nullptr;
497 MEDCoupling::ParaMESH* paramesh = nullptr;
498 MEDCoupling::ParaFIELD* parafield = nullptr;
499 ICoCo::MEDDoubleField* icocofield = nullptr;
501 string filename_xml1 = "square1_split";
502 string filename_xml2 = "square2_split";
503 //string filename_seq_wr = makeTmpFile("");
504 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
506 // To remove tmp files from disk
507 ParaMEDMEMTest_TmpFilesRemover aRemover;
509 MPI_Barrier(MPI_COMM_WORLD);
510 if (source_group->containsMyRank())
512 string master = filename_xml1;
514 ostringstream strstream;
515 strstream <<master<<rank+1<<".med";
516 string fName = INTERP_TEST::getResourceFile(strstream.str());
517 ostringstream meshname ;
518 meshname<< "Mesh_2_"<< rank+1;
520 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
523 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
525 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
526 MEDCoupling::ComponentTopology comptopo;
529 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
530 parafield->getField()->setNature(IntensiveMaximum);
533 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
536 nb_local=mesh->getNumberOfCells();
538 nb_local=mesh->getNumberOfNodes();
539 // double * value= new double[nb_local];
540 double *value=parafield->getField()->getArray()->getPointer();
541 for(int ielem=0; ielem<nb_local;ielem++)
544 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
545 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
546 dec.setMethod(srcMeth);
547 dec.attachLocalField(icocofield);
550 //loading the geometry for the target group
551 if (target_group->containsMyRank())
553 string master= filename_xml2;
554 ostringstream strstream;
555 strstream << master<<(rank-nproc_source+1)<<".med";
556 string fName = INTERP_TEST::getResourceFile(strstream.str());
557 ostringstream meshname ;
558 meshname<< "Mesh_3_"<<rank-nproc_source+1;
559 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
561 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
562 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
563 MEDCoupling::ComponentTopology comptopo;
566 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
567 parafield->getField()->setNature(IntensiveMaximum);
570 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
573 nb_local=mesh->getNumberOfCells();
575 nb_local=mesh->getNumberOfNodes();
576 // double * value= new double[nb_local];
577 double *value=parafield->getField()->getArray()->getPointer();
578 for(int ielem=0; ielem<nb_local;ielem++)
580 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
581 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
582 dec.setMethod(targetMeth);
583 dec.attachLocalField(icocofield);
587 //attaching a DEC to the source group
588 double field_before_int;
589 double field_after_int;
591 if (source_group->containsMyRank())
593 field_before_int = parafield->getVolumeIntegral(0,true);
595 cout<<"DEC usage"<<endl;
596 dec.setForcedRenormalization(false);
599 ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
600 if (source_group->myRank()==0)
601 aRemover.Register("./sourcesquareb");
602 ostringstream filename;
603 filename<<"./sourcesquareb_"<<source_group->myRank()+1;
604 aRemover.Register(filename.str().c_str());
605 //WriteField("./sourcesquareb",parafield->getField());
608 cout <<"writing"<<endl;
609 ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
610 if (source_group->myRank()==0)
611 aRemover.Register("./sourcesquare");
612 //WriteField("./sourcesquare",parafield->getField());
615 filename<<"./sourcesquare_"<<source_group->myRank()+1;
616 aRemover.Register(filename.str().c_str());
617 field_after_int = parafield->getVolumeIntegral(0,true);
620 // MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
621 // MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
623 CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
627 //attaching a DEC to the target group
628 if (target_group->containsMyRank())
631 dec.setForcedRenormalization(false);
634 ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
635 //WriteField("./targetsquareb",parafield->getField());
636 if (target_group->myRank()==0)
637 aRemover.Register("./targetsquareb");
638 ostringstream filename;
639 filename<<"./targetsquareb_"<<target_group->myRank()+1;
640 aRemover.Register(filename.str().c_str());
642 ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
643 //WriteField("./targetsquare",parafield->getField());
645 if (target_group->myRank()==0)
646 aRemover.Register("./targetsquareb");
648 filename<<"./targetsquareb_"<<target_group->myRank()+1;
649 aRemover.Register(filename.str().c_str());
650 // double field_before_int, field_after_int;
651 // MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
652 // MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
654 // CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
667 MPI_Barrier(MPI_COMM_WORLD);
668 cout << "end of InterpKernelDEC_2D test"<<endl;
671 void ParaMEDMEMTest::testInterpKernelDEC2_2D_(const char *srcMeth, const char *targetMeth)
673 std::string srcM(srcMeth);
674 std::string targetM(targetMeth);
677 MPI_Comm_size(MPI_COMM_WORLD,&size);
678 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
680 //the test is meant to run on five processors
681 if (size !=5) return ;
683 int nproc_source = 3;
685 set<int> procs_source;
686 set<int> procs_target;
688 for (int i=0; i<nproc_source; i++)
689 procs_source.insert(i);
690 for (int i=nproc_source; i<size; i++)
691 procs_target.insert(i);
692 self_procs.insert(rank);
694 MEDCoupling::CommInterface interface;
696 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
697 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
698 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
700 //loading the geometry for the source group
702 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
704 MEDCoupling::MEDCouplingUMesh* mesh = nullptr;
705 MEDCoupling::MEDCouplingFieldDouble* mcfield = nullptr;
707 string filename_xml1 = "square1_split";
708 string filename_xml2 = "square2_split";
710 // To remove tmp files from disk
711 ParaMEDMEMTest_TmpFilesRemover aRemover;
713 MPI_Barrier(MPI_COMM_WORLD);
714 if (source_group->containsMyRank())
716 string master = filename_xml1;
718 ostringstream strstream;
719 strstream <<master<<rank+1<<".med";
720 string fName = INTERP_TEST::getResourceFile(strstream.str());
721 ostringstream meshname ;
722 meshname<< "Mesh_2_"<< rank+1;
724 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
725 MEDCoupling::ComponentTopology comptopo;
728 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
729 mcfield->setMesh(mesh);
730 DataArrayDouble *array=DataArrayDouble::New();
731 array->alloc(mcfield->getNumberOfTuples(),1);
732 mcfield->setArray(array);
734 mcfield->setNature(IntensiveMaximum);
738 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
739 mcfield->setMesh(mesh);
740 DataArrayDouble *array=DataArrayDouble::New();
741 array->alloc(mcfield->getNumberOfTuples(),1);
742 mcfield->setArray(array);
747 nb_local=mesh->getNumberOfCells();
749 nb_local=mesh->getNumberOfNodes();
750 double *value=mcfield->getArray()->getPointer();
751 for(int ielem=0; ielem<nb_local;ielem++)
753 dec.setMethod(srcMeth);
754 dec.attachLocalField(mcfield);
755 dec.attachLocalField(mcfield);
758 //loading the geometry for the target group
759 if (target_group->containsMyRank())
761 string master= filename_xml2;
762 ostringstream strstream;
763 strstream << master<<(rank-nproc_source+1)<<".med";
764 string fName = INTERP_TEST::getResourceFile(strstream.str());
765 ostringstream meshname ;
766 meshname<< "Mesh_3_"<<rank-nproc_source+1;
767 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
768 MEDCoupling::ComponentTopology comptopo;
771 mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
772 mcfield->setMesh(mesh);
773 DataArrayDouble *array=DataArrayDouble::New();
774 array->alloc(mcfield->getNumberOfTuples(),1);
775 mcfield->setArray(array);
777 mcfield->setNature(IntensiveMaximum);
781 mcfield = MEDCouplingFieldDouble::New(ON_NODES,NO_TIME);
782 mcfield->setMesh(mesh);
783 DataArrayDouble *array=DataArrayDouble::New();
784 array->alloc(mcfield->getNumberOfTuples(),1);
785 mcfield->setArray(array);
790 nb_local=mesh->getNumberOfCells();
792 nb_local=mesh->getNumberOfNodes();
793 double *value=mcfield->getArray()->getPointer();
794 for(int ielem=0; ielem<nb_local;ielem++)
796 dec.setMethod(targetMeth);
797 dec.attachLocalField(mcfield);
798 dec.attachLocalField(mcfield);
802 //attaching a DEC to the source group
804 if (source_group->containsMyRank())
807 dec.setForcedRenormalization(false);
812 //attaching a DEC to the target group
813 if (target_group->containsMyRank())
816 dec.setForcedRenormalization(false);
826 MPI_Barrier(MPI_COMM_WORLD);
827 cout << "end of InterpKernelDEC2_2D test"<<endl;
830 void ParaMEDMEMTest::testInterpKernelDEC_3D_(const char *srcMeth, const char *targetMeth)
832 std::string srcM(srcMeth);
833 std::string targetM(targetMeth);
836 MPI_Comm_size(MPI_COMM_WORLD,&size);
837 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
839 //the test is meant to run on five processors
840 if (size !=3) return ;
842 int nproc_source = 2;
844 set<int> procs_source;
845 set<int> procs_target;
847 for (int i=0; i<nproc_source; i++)
848 procs_source.insert(i);
849 for (int i=nproc_source; i<size; i++)
850 procs_target.insert(i);
851 self_procs.insert(rank);
853 MEDCoupling::CommInterface interface;
855 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
856 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
857 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
859 //loading the geometry for the source group
861 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
863 MEDCoupling::MEDCouplingUMesh* mesh = nullptr;
864 MEDCoupling::ParaMESH* paramesh = nullptr;
865 MEDCoupling::ParaFIELD* parafield = nullptr;
866 ICoCo::MEDDoubleField* icocofield = nullptr;
868 char * tmp_dir_c = getenv("TMP");
870 if (tmp_dir_c != NULL)
871 tmp_dir = string(tmp_dir_c);
874 string filename_xml1 = "Mesh3D_10_2d";
875 string filename_xml2 = "Mesh3D_11";
876 //string filename_seq_wr = makeTmpFile("");
877 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
879 // To remove tmp files from disk
880 ParaMEDMEMTest_TmpFilesRemover aRemover;
882 MPI_Barrier(MPI_COMM_WORLD);
883 if (source_group->containsMyRank())
885 string master = filename_xml1;
887 ostringstream strstream;
888 strstream <<master<<rank+1<<".med";
889 std::string fName = INTERP_TEST::getResourceFile(strstream.str());
890 ostringstream meshname ;
891 meshname<< "Mesh_3_"<< rank+1;
893 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
896 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
898 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
899 MEDCoupling::ComponentTopology comptopo;
902 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
903 parafield->getField()->setNature(IntensiveMaximum);
906 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
909 nb_local=mesh->getNumberOfCells();
911 nb_local=mesh->getNumberOfNodes();
912 // double * value= new double[nb_local];
913 double *value=parafield->getField()->getArray()->getPointer();
914 for(int ielem=0; ielem<nb_local;ielem++)
917 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
918 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
919 dec.setMethod(srcMeth);
920 dec.attachLocalField(icocofield);
923 //loading the geometry for the target group
924 if (target_group->containsMyRank())
926 string master= filename_xml2;
927 ostringstream strstream;
928 strstream << master << ".med";
929 std::string fName = INTERP_TEST::getResourceFile(strstream.str());
930 ostringstream meshname ;
932 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
934 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
935 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
936 MEDCoupling::ComponentTopology comptopo;
939 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
940 parafield->getField()->setNature(IntensiveMaximum);
943 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
946 nb_local=mesh->getNumberOfCells();
948 nb_local=mesh->getNumberOfNodes();
949 // double * value= new double[nb_local];
950 double *value=parafield->getField()->getArray()->getPointer();
951 for(int ielem=0; ielem<nb_local;ielem++)
953 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
954 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
955 dec.setMethod(targetMeth);
956 dec.attachLocalField(icocofield);
958 //attaching a DEC to the source group
959 double field_before_int;
960 double field_after_int;
962 if (source_group->containsMyRank())
964 field_before_int = parafield->getVolumeIntegral(0,true);
966 cout<<"DEC usage"<<endl;
967 dec.setForcedRenormalization(false);
970 ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
971 if (source_group->myRank()==0)
972 aRemover.Register("./sourcesquareb");
973 ostringstream filename;
974 filename<<"./sourcesquareb_"<<source_group->myRank()+1;
975 aRemover.Register(filename.str().c_str());
976 //WriteField("./sourcesquareb",parafield->getField());
979 cout <<"writing"<<endl;
980 ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
981 if (source_group->myRank()==0)
982 aRemover.Register("./sourcesquare");
983 //WriteField("./sourcesquare",parafield->getField());
986 filename<<"./sourcesquare_"<<source_group->myRank()+1;
987 aRemover.Register(filename.str().c_str());
988 field_after_int = parafield->getVolumeIntegral(0,true);
990 CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
994 //attaching a DEC to the target group
995 if (target_group->containsMyRank())
998 dec.setForcedRenormalization(false);
1001 ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
1002 //WriteField("./targetsquareb",parafield->getField());
1003 if (target_group->myRank()==0)
1004 aRemover.Register("./targetsquareb");
1005 ostringstream filename;
1006 filename<<"./targetsquareb_"<<target_group->myRank()+1;
1007 aRemover.Register(filename.str().c_str());
1009 ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
1010 //WriteField("./targetsquare",parafield->getField());
1012 if (target_group->myRank()==0)
1013 aRemover.Register("./targetsquareb");
1015 filename<<"./targetsquareb_"<<target_group->myRank()+1;
1016 aRemover.Register(filename.str().c_str());
1018 delete source_group;
1019 delete target_group;
1027 MPI_Barrier(MPI_COMM_WORLD);
1028 cout << "end of InterpKernelDEC_3D test"<<endl;
1031 //Synchronous tests without interpolation with native mode (AllToAll(v) from lam/MPI:
1032 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpNativeDEC_2D()
1034 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,false,false,false,"P0","P0");
1037 //Synchronous tests without interpolation :
1038 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpDEC_2D()
1040 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,false,"P0","P0");
1043 //Synchronous tests with interpolation :
1044 void ParaMEDMEMTest::testSynchronousEqualInterpKernelDEC_2D()
1046 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,true,"P0","P0");
1048 void ParaMEDMEMTest::testSynchronousFasterSourceInterpKernelDEC_2D()
1050 testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,false,true,"P0","P0");
1052 void ParaMEDMEMTest::testSynchronousSlowerSourceInterpKernelDEC_2D()
1054 testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,false,true,"P0","P0");
1056 void ParaMEDMEMTest::testSynchronousSlowSourceInterpKernelDEC_2D()
1058 testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,false,true,"P0","P0");
1060 void ParaMEDMEMTest::testSynchronousFastSourceInterpKernelDEC_2D()
1062 testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,false,true,"P0","P0");
1065 //Asynchronous tests with interpolation :
1066 void ParaMEDMEMTest::testAsynchronousEqualInterpKernelDEC_2D()
1068 testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,true,true,"P0","P0");
1070 void ParaMEDMEMTest::testAsynchronousFasterSourceInterpKernelDEC_2D()
1072 testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,true,true,"P0","P0");
1074 void ParaMEDMEMTest::testAsynchronousSlowerSourceInterpKernelDEC_2D()
1076 testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,true,true,"P0","P0");
1078 void ParaMEDMEMTest::testAsynchronousSlowSourceInterpKernelDEC_2D()
1080 testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,true,true,"P0","P0");
1082 void ParaMEDMEMTest::testAsynchronousFastSourceInterpKernelDEC_2D()
1084 testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,true,true,"P0","P0");
1087 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P0()
1090 const double sourceCoordsAll[2][8]={{0.4,0.5,0.4,1.5,1.6,1.5,1.6,0.5},
1091 {0.3,-0.5,1.6,-0.5,1.6,-1.5,0.3,-1.5}};
1092 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},
1093 {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},
1094 {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}};
1095 mcIdType conn4All[8]={0,1,2,3,4,5,6,7};
1096 double targetResults[3][2]={{34.,34.},{38.333333333333336,42.666666666666664},{47.,47.}};
1097 double targetResults2[3][2]={{0.28333333333333344,0.56666666666666687},{1.8564102564102569,2.0128205128205132},{1.0846153846153845,0.36153846153846159}};
1098 double targetResults3[3][2]={{3.7777777777777781,7.5555555555555562},{24.511111111111113,26.355555555555558},{14.1,4.7}};
1099 double targetResults4[3][2]={{8.5,17},{8.8461538461538431, 9.8461538461538449},{35.25,11.75}};
1103 MPI_Comm_size(MPI_COMM_WORLD,&size);
1104 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1108 int nproc_source = 2;
1109 set<int> self_procs;
1110 set<int> procs_source;
1111 set<int> procs_target;
1113 for (int i=0; i<nproc_source; i++)
1114 procs_source.insert(i);
1115 for (int i=nproc_source; i<size; i++)
1116 procs_target.insert(i);
1117 self_procs.insert(rank);
1119 MEDCoupling::MEDCouplingUMesh *mesh=0;
1120 MEDCoupling::ParaMESH *paramesh=0;
1121 MEDCoupling::ParaFIELD* parafield=0;
1123 MEDCoupling::CommInterface interface;
1125 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1126 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1127 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1129 MPI_Barrier(MPI_COMM_WORLD);
1130 if(source_group->containsMyRank())
1132 std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
1133 mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
1134 mesh->allocateCells(2);
1135 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
1136 mesh->finishInsertingCells();
1137 DataArrayDouble *myCoords=DataArrayDouble::New();
1138 myCoords->alloc(4,2);
1139 const double *sourceCoords=sourceCoordsAll[rank];
1140 std::copy(sourceCoords,sourceCoords+8,myCoords->getPointer());
1141 mesh->setCoords(myCoords);
1142 myCoords->decrRef();
1143 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1144 MEDCoupling::ComponentTopology comptopo;
1145 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1146 double *value=parafield->getField()->getArray()->getPointer();
1147 value[0]=34+13*((double)rank);
1151 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
1152 mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
1153 mesh->allocateCells(2);
1154 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
1155 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All+4);
1156 mesh->finishInsertingCells();
1157 DataArrayDouble *myCoords=DataArrayDouble::New();
1158 myCoords->alloc(8,2);
1159 const double *targetCoords=targetCoordsAll[rank-nproc_source];
1160 std::copy(targetCoords,targetCoords+16,myCoords->getPointer());
1161 mesh->setCoords(myCoords);
1162 myCoords->decrRef();
1163 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
1164 MEDCoupling::ComponentTopology comptopo;
1165 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1167 //test 1 - Conservative volumic
1168 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1169 parafield->getField()->setNature(IntensiveMaximum);
1170 if (source_group->containsMyRank())
1172 dec.setMethod("P0");
1173 dec.attachLocalField(parafield);
1175 dec.setForcedRenormalization(false);
1180 dec.setMethod("P0");
1181 dec.attachLocalField(parafield);
1183 dec.setForcedRenormalization(false);
1185 const double *res=parafield->getField()->getArray()->getConstPointer();
1186 const double *expected=targetResults[rank-nproc_source];
1187 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1188 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1190 //test 2 - ExtensiveMaximum
1191 MEDCoupling::InterpKernelDEC dec2(*source_group,*target_group);
1192 parafield->getField()->setNature(ExtensiveMaximum);
1193 if (source_group->containsMyRank())
1195 dec2.setMethod("P0");
1196 dec2.attachLocalField(parafield);
1198 dec2.setForcedRenormalization(false);
1203 dec2.setMethod("P0");
1204 dec2.attachLocalField(parafield);
1206 dec2.setForcedRenormalization(false);
1208 const double *res=parafield->getField()->getArray()->getConstPointer();
1209 const double *expected=targetResults2[rank-nproc_source];
1210 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1211 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1213 //test 3 - ExtensiveMaximum with global constraint
1214 MEDCoupling::InterpKernelDEC dec3(*source_group,*target_group);
1215 parafield->getField()->setNature(ExtensiveConservation);
1216 if (source_group->containsMyRank())
1218 dec3.setMethod("P0");
1219 dec3.attachLocalField(parafield);
1221 dec3.setForcedRenormalization(false);
1226 dec3.setMethod("P0");
1227 dec3.attachLocalField(parafield);
1229 dec3.setForcedRenormalization(false);
1231 const double *res=parafield->getField()->getArray()->getConstPointer();
1232 const double *expected=targetResults3[rank-nproc_source];
1233 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1234 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1236 //test 4 - IntensiveConservation
1237 MEDCoupling::InterpKernelDEC dec4(*source_group,*target_group);
1238 parafield->getField()->setNature(IntensiveConservation);
1239 if (source_group->containsMyRank())
1241 dec4.setMethod("P0");
1242 dec4.attachLocalField(parafield);
1244 dec4.setForcedRenormalization(false);
1249 dec4.setMethod("P0");
1250 dec4.attachLocalField(parafield);
1252 dec4.setForcedRenormalization(false);
1254 const double *res=parafield->getField()->getArray()->getConstPointer();
1255 const double *expected=targetResults4[rank-nproc_source];
1256 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
1257 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
1259 //test 5 - Conservative volumic reversed
1260 MEDCoupling::InterpKernelDEC dec5(*source_group,*target_group);
1261 parafield->getField()->setNature(IntensiveMaximum);
1262 if (source_group->containsMyRank())
1264 dec5.setMethod("P0");
1265 dec5.attachLocalField(parafield);
1267 dec5.setForcedRenormalization(false);
1269 const double *res=parafield->getField()->getArray()->getConstPointer();
1270 CPPUNIT_ASSERT_EQUAL(1,(int)parafield->getField()->getNumberOfTuples());
1271 const double expected[]={37.8518518518519,43.5333333333333};
1272 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1276 dec5.setMethod("P0");
1277 dec5.attachLocalField(parafield);
1279 dec5.setForcedRenormalization(false);
1280 double *res=parafield->getField()->getArray()->getPointer();
1281 const double *toSet=targetResults[rank-nproc_source];
1286 //test 6 - ExtensiveMaximum reversed
1287 MEDCoupling::InterpKernelDEC dec6(*source_group,*target_group);
1288 parafield->getField()->setNature(ExtensiveMaximum);
1289 if (source_group->containsMyRank())
1291 dec6.setMethod("P0");
1292 dec6.attachLocalField(parafield);
1294 dec6.setForcedRenormalization(false);
1296 const double *res=parafield->getField()->getArray()->getConstPointer();
1297 CPPUNIT_ASSERT_EQUAL(1,(int)parafield->getField()->getNumberOfTuples());
1298 const double expected[]={0.794600591715977,1.35631163708087};
1299 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1303 dec6.setMethod("P0");
1304 dec6.attachLocalField(parafield);
1306 dec6.setForcedRenormalization(false);
1307 double *res=parafield->getField()->getArray()->getPointer();
1308 const double *toSet=targetResults2[rank-nproc_source];
1313 //test 7 - ExtensiveMaximum with global constraint reversed
1314 MEDCoupling::InterpKernelDEC dec7(*source_group,*target_group);
1315 parafield->getField()->setNature(ExtensiveConservation);
1316 if (source_group->containsMyRank())
1318 dec7.setMethod("P0");
1319 dec7.attachLocalField(parafield);
1321 dec7.setForcedRenormalization(false);
1323 const double *res=parafield->getField()->getArray()->getConstPointer();
1324 CPPUNIT_ASSERT_EQUAL(1,(int)parafield->getField()->getNumberOfTuples());
1325 const double expected[]={36.4592592592593,44.5407407407407};
1326 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1330 dec7.setMethod("P0");
1331 dec7.attachLocalField(parafield);
1333 dec7.setForcedRenormalization(false);
1334 double *res=parafield->getField()->getArray()->getPointer();
1335 const double *toSet=targetResults3[rank-nproc_source];
1340 //test 8 - ExtensiveMaximum with IntensiveConservation reversed
1341 MEDCoupling::InterpKernelDEC dec8(*source_group,*target_group);
1342 parafield->getField()->setNature(IntensiveConservation);
1343 if (source_group->containsMyRank())
1345 dec8.setMethod("P0");
1346 dec8.attachLocalField(parafield);
1348 dec8.setForcedRenormalization(false);
1350 const double *res=parafield->getField()->getArray()->getConstPointer();
1351 CPPUNIT_ASSERT_EQUAL(1,(int)parafield->getField()->getNumberOfTuples());
1352 const double expected[]={0.81314102564102553,1.3428994082840233};
1353 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
1357 dec8.setMethod("P0");
1358 dec8.attachLocalField(parafield);
1360 dec8.setForcedRenormalization(false);
1361 double *res=parafield->getField()->getArray()->getPointer();
1362 const double *toSet=targetResults4[rank-nproc_source];
1372 delete target_group;
1373 delete source_group;
1375 MPI_Barrier(MPI_COMM_WORLD);
1378 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P1P1P0()
1382 MPI_Comm_size(MPI_COMM_WORLD,&size);
1383 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1387 int nproc_source = 2;
1388 set<int> self_procs;
1389 set<int> procs_source;
1390 set<int> procs_target;
1392 for (int i=0; i<nproc_source; i++)
1393 procs_source.insert(i);
1394 for (int i=nproc_source; i<size; i++)
1395 procs_target.insert(i);
1396 self_procs.insert(rank);
1398 MEDCoupling::MEDCouplingUMesh *mesh=0;
1399 MEDCoupling::ParaMESH *paramesh=0;
1400 MEDCoupling::ParaFIELD *parafieldP0=0,*parafieldP1=0;
1402 MEDCoupling::CommInterface interface;
1404 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1405 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1406 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1408 MPI_Barrier(MPI_COMM_WORLD);
1409 if(source_group->containsMyRank())
1413 double coords[6]={-0.3,-0.3, 0.7,0.7, 0.7,-0.3};
1414 mcIdType conn[3]={0,1,2};
1415 //int globalNode[3]={1,2,0};
1416 mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
1417 mesh->allocateCells(1);
1418 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
1419 mesh->finishInsertingCells();
1420 DataArrayDouble *myCoords=DataArrayDouble::New();
1421 myCoords->alloc(3,2);
1422 std::copy(coords,coords+6,myCoords->getPointer());
1423 mesh->setCoords(myCoords);
1424 myCoords->decrRef();
1428 double coords[6]={-0.3,-0.3, -0.3,0.7, 0.7,0.7};
1429 mcIdType conn[3]={0,1,2};
1430 //int globalNode[3]={1,3,2};
1431 mesh=MEDCouplingUMesh::New("Source mesh Proc1",2);
1432 mesh->allocateCells(1);
1433 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
1434 mesh->finishInsertingCells();
1435 DataArrayDouble *myCoords=DataArrayDouble::New();
1436 myCoords->alloc(3,2);
1437 std::copy(coords,coords+6,myCoords->getPointer());
1438 mesh->setCoords(myCoords);
1439 myCoords->decrRef();
1441 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1442 MEDCoupling::ComponentTopology comptopo;
1443 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1444 parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
1445 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1446 double *valueP1=parafieldP1->getField()->getArray()->getPointer();
1447 parafieldP0->getField()->setNature(IntensiveMaximum);
1448 parafieldP1->getField()->setNature(IntensiveMaximum);
1452 valueP1[0]=34.; valueP1[1]=77.; valueP1[2]=53.;
1457 valueP1[0]=34.; valueP1[1]=57.; valueP1[2]=77.;
1462 const char targetMeshName[]="target mesh";
1465 double coords[10]={-0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2 };
1466 mcIdType conn[7]={0,3,4,1, 1,4,2};
1467 //int globalNode[5]={4,3,0,2,1};
1468 mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
1469 mesh->allocateCells(2);
1470 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1471 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
1472 mesh->finishInsertingCells();
1473 DataArrayDouble *myCoords=DataArrayDouble::New();
1474 myCoords->alloc(5,2);
1475 std::copy(coords,coords+10,myCoords->getPointer());
1476 mesh->setCoords(myCoords);
1477 myCoords->decrRef();
1478 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
1479 std::vector<mcIdType> globalNumberingP2 = {0,1,2,3,4};
1480 MCAuto<DataArrayIdType> da=DataArrayIdType::New(); da->alloc(5,1);
1481 std::copy(globalNumberingP2.begin(), globalNumberingP2.end(), da->rwBegin());
1482 paramesh->setNodeGlobal(da);
1486 double coords[6]={0.2,0.2, 0.7,-0.3, 0.7,0.2};
1487 mcIdType 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 std::vector<mcIdType> globalNumberingP3 = {4,2,5};
1500 MCAuto<DataArrayIdType> da=DataArrayIdType::New(); da->alloc(3,1);
1501 std::copy(globalNumberingP3.begin(), globalNumberingP3.end(), da->rwBegin());
1502 paramesh->setNodeGlobal(da);
1506 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};
1507 mcIdType conn[8]={0,1,2,3, 3,2,4,5};
1508 //int globalNode[6]={2,6,7,1,8,5};
1509 mesh=MEDCouplingUMesh::New("Target mesh Proc4",2);
1510 mesh->allocateCells(2);
1511 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1512 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
1513 mesh->finishInsertingCells();
1514 DataArrayDouble *myCoords=DataArrayDouble::New();
1515 myCoords->alloc(6,2);
1516 std::copy(coords,coords+12,myCoords->getPointer());
1517 mesh->setCoords(myCoords);
1518 myCoords->decrRef();
1519 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
1520 std::vector<mcIdType> globalNumberingP4 = {3,6,7,4,8,5};
1521 MCAuto<DataArrayIdType> da=DataArrayIdType::New(); da->alloc(6,1);
1522 std::copy(globalNumberingP4.begin(), globalNumberingP4.end(), da->rwBegin());
1523 paramesh->setNodeGlobal(da);
1525 MEDCoupling::ComponentTopology comptopo;
1526 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1527 parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
1528 parafieldP0->getField()->setNature(IntensiveMaximum);
1529 parafieldP1->getField()->setNature(IntensiveMaximum);
1532 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1533 if (source_group->containsMyRank())
1535 dec.setMethod("P0");
1536 dec.attachLocalField(parafieldP0);
1538 dec.setForcedRenormalization(false);
1541 const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1544 CPPUNIT_ASSERT_DOUBLES_EQUAL(34.42857143,valueP0[0],1e-7);
1548 CPPUNIT_ASSERT_DOUBLES_EQUAL(44.,valueP0[0],1e-7);
1553 dec.setMethod("P1");
1554 dec.attachLocalField(parafieldP1);
1556 dec.setForcedRenormalization(false);
1558 const double *res=parafieldP1->getField()->getArray()->getConstPointer();
1561 const double expectP2[5]={39.0, 31.0, 31.0, 47.0, 39.0};
1562 CPPUNIT_ASSERT_EQUAL(5,(int)parafieldP1->getField()->getNumberOfTuples());
1563 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP1->getField()->getNumberOfComponents());
1564 for(int kk=0;kk<5;kk++)
1565 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP2[kk],res[kk],1e-12);
1569 const double expectP3[3]={39.0, 31.0, 31.0};
1570 CPPUNIT_ASSERT_EQUAL(3,(int)parafieldP1->getField()->getNumberOfTuples());
1571 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP1->getField()->getNumberOfComponents());
1572 for(int kk=0;kk<3;kk++)
1573 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP3[kk],res[kk],1e-12);
1577 const double expectP4[6]={47.0, 47.0, 47.0, 39.0, 39.0, 31.0};
1578 CPPUNIT_ASSERT_EQUAL(6,(int)parafieldP1->getField()->getNumberOfTuples());
1579 CPPUNIT_ASSERT_EQUAL(1,(int)parafieldP1->getField()->getNumberOfComponents());
1580 for(int kk=0;kk<6;kk++)
1581 CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP4[kk],res[kk],1e-12);
1591 delete target_group;
1592 delete source_group;
1594 MPI_Barrier(MPI_COMM_WORLD);
1597 void ParaMEDMEMTest::testInterpKernelDEC2DM1D_P0P0()
1601 MPI_Comm_size(MPI_COMM_WORLD,&size);
1602 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1607 set<int> procs_source;
1608 set<int> procs_target;
1610 for (int i=0; i<nproc_source; i++)
1611 procs_source.insert(i);
1612 for (int i=nproc_source;i<size; i++)
1613 procs_target.insert(i);
1615 MEDCoupling::MEDCouplingUMesh *mesh=0;
1616 MEDCoupling::ParaMESH *paramesh=0;
1617 MEDCoupling::ParaFIELD *parafield=0;
1619 MEDCoupling::CommInterface interface;
1621 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1622 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1624 MPI_Barrier(MPI_COMM_WORLD);
1625 if(source_group->containsMyRank())
1627 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 };
1628 mesh=MEDCouplingUMesh::New();
1629 mesh->setMeshDimension(2);
1630 DataArrayDouble *myCoords=DataArrayDouble::New();
1631 myCoords->alloc(9,2);
1632 std::copy(targetCoords,targetCoords+18,myCoords->getPointer());
1633 mesh->setCoords(myCoords);
1634 myCoords->decrRef();
1637 mcIdType targetConn[7]={0,3,4,1, 1,4,2};
1638 mesh->allocateCells(2);
1639 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
1640 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+4);
1641 mesh->finishInsertingCells();
1645 mcIdType targetConn[11]={4,5,2, 6,7,4,3, 7,8,5,4};
1646 mesh->allocateCells(3);
1647 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1648 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+3);
1649 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+7);
1650 mesh->finishInsertingCells();
1652 MEDCoupling::ComponentTopology comptopo;
1653 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1654 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1655 parafield->getField()->setNature(IntensiveMaximum);
1656 double *vals=parafield->getField()->getArray()->getPointer();
1658 { vals[0]=7.; vals[1]=8.; }
1660 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1664 mesh=MEDCouplingUMesh::New("an example of -1 D mesh",-1);
1665 MEDCoupling::ComponentTopology comptopo;
1666 paramesh=new ParaMESH(mesh,*target_group,"target mesh");
1667 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1668 parafield->getField()->setNature(IntensiveMaximum);
1670 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
1671 if(source_group->containsMyRank())
1673 dec.setMethod("P0");
1674 dec.attachLocalField(parafield);
1676 dec.setForcedRenormalization(false);
1679 const double *res=parafield->getField()->getArray()->getConstPointer();
1682 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1683 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1687 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1688 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1689 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
1694 dec.setMethod("P0");
1695 dec.attachLocalField(parafield);
1697 dec.setForcedRenormalization(false);
1699 const double *res=parafield->getField()->getArray()->getConstPointer();
1700 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1703 MEDCoupling::InterpKernelDEC dec2(*source_group,*target_group);
1704 dec2.setMethod("P0");
1705 parafield->getField()->setNature(ExtensiveConservation);
1706 if(source_group->containsMyRank())
1708 double *vals=parafield->getField()->getArray()->getPointer();
1710 { vals[0]=7.; vals[1]=8.; }
1712 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1713 dec2.attachLocalField(parafield);
1717 const double *res=parafield->getField()->getArray()->getConstPointer();
1720 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
1721 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
1725 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
1726 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
1727 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
1732 dec2.attachLocalField(parafield);
1735 const double *res=parafield->getField()->getArray()->getConstPointer();
1736 CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
1740 MEDCoupling::InterpKernelDEC dec3(*source_group,*target_group);
1741 dec3.setMethod("P0");
1742 parafield->getField()->setNature(ExtensiveMaximum);
1743 if(source_group->containsMyRank())
1745 double *vals=parafield->getField()->getArray()->getPointer();
1747 { vals[0]=7.; vals[1]=8.; }
1749 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1750 dec3.attachLocalField(parafield);
1754 const double *res=parafield->getField()->getArray()->getConstPointer();
1757 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
1758 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
1762 CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
1763 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
1764 CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
1769 dec3.attachLocalField(parafield);
1772 const double *res=parafield->getField()->getArray()->getConstPointer();
1773 CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
1777 MEDCoupling::InterpKernelDEC dec4(*source_group,*target_group);
1778 dec4.setMethod("P0");
1779 parafield->getField()->setNature(IntensiveConservation);
1780 if(source_group->containsMyRank())
1782 double *vals=parafield->getField()->getArray()->getPointer();
1784 { vals[0]=7.; vals[1]=8.; }
1786 { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
1787 dec4.attachLocalField(parafield);
1791 const double *res=parafield->getField()->getArray()->getConstPointer();
1794 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1795 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1799 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1800 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
1801 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
1806 dec4.attachLocalField(parafield);
1809 const double *res=parafield->getField()->getArray()->getConstPointer();
1810 CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
1816 delete target_group;
1817 delete source_group;
1819 MPI_Barrier(MPI_COMM_WORLD);
1822 void ParaMEDMEMTest::testInterpKernelDECPartialProcs()
1826 MPI_Comm_size(MPI_COMM_WORLD,&size);
1827 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1831 set<int> procs_source;
1832 set<int> procs_target;
1834 procs_source.insert(0);
1835 procs_target.insert(1);
1837 MEDCoupling::MEDCouplingUMesh *mesh=0;
1838 MEDCoupling::ParaMESH *paramesh=0;
1839 MEDCoupling::ParaFIELD *parafield=0;
1841 MEDCoupling::CommInterface interface;
1843 MPI_Barrier(MPI_COMM_WORLD);
1844 double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
1846 int grpIds[2]={0,1};
1847 MPI_Group grp,group_world;
1848 comm.commGroup(MPI_COMM_WORLD,&group_world);
1849 comm.groupIncl(group_world,2,grpIds,&grp);
1850 MPI_Comm partialComm;
1851 comm.commCreate(MPI_COMM_WORLD,grp,&partialComm);
1853 ProcessorGroup* target_group=0;
1854 ProcessorGroup* source_group=0;
1856 MEDCoupling::InterpKernelDEC *dec=0;
1857 if(rank==0 || rank==1)
1859 target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target,partialComm);
1860 source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source,partialComm);
1861 if(source_group->containsMyRank())
1863 mesh=MEDCouplingUMesh::New();
1864 mesh->setMeshDimension(2);
1865 DataArrayDouble *myCoords=DataArrayDouble::New();
1866 myCoords->alloc(4,2);
1867 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
1868 mesh->setCoords(myCoords);
1869 myCoords->decrRef();
1870 mcIdType targetConn[4]={0,2,3,1};
1871 mesh->allocateCells(1);
1872 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
1873 mesh->finishInsertingCells();
1874 MEDCoupling::ComponentTopology comptopo;
1875 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1876 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1877 parafield->getField()->setNature(IntensiveMaximum);
1878 double *vals=parafield->getField()->getArray()->getPointer();
1880 dec=new MEDCoupling::InterpKernelDEC(*source_group,*target_group);
1881 dec->attachLocalField(parafield);
1888 mesh=MEDCouplingUMesh::New();
1889 mesh->setMeshDimension(2);
1890 DataArrayDouble *myCoords=DataArrayDouble::New();
1891 myCoords->alloc(4,2);
1892 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
1893 mesh->setCoords(myCoords);
1894 myCoords->decrRef();
1895 mcIdType targetConn[6]={0,2,1,2,3,1};
1896 mesh->allocateCells(2);
1897 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
1898 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
1899 mesh->finishInsertingCells();
1900 MEDCoupling::ComponentTopology comptopo;
1901 paramesh=new ParaMESH(mesh,*target_group,"target mesh");
1902 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1903 parafield->getField()->setNature(IntensiveMaximum);
1904 dec=new MEDCoupling::InterpKernelDEC(*source_group,*target_group);
1905 dec->attachLocalField(parafield);
1915 delete target_group;
1916 delete source_group;
1918 if(partialComm != MPI_COMM_NULL)
1919 comm.commFree(&partialComm);
1920 comm.groupFree(&grp);
1921 comm.groupFree(&group_world);
1922 MPI_Barrier(MPI_COMM_WORLD);
1926 * This test reproduces bug of Gauthier on 13/9/2010 concerning 3DSurf meshes.
1927 * 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.
1929 void ParaMEDMEMTest::testInterpKernelDEC3DSurfEmptyBBox()
1933 MPI_Comm_size(MPI_COMM_WORLD,&size);
1934 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
1938 int nproc_source = 1;
1939 set<int> self_procs;
1940 set<int> procs_source;
1941 set<int> procs_target;
1943 for (int i=0; i<nproc_source; i++)
1944 procs_source.insert(i);
1945 for (int i=nproc_source; i<size; i++)
1946 procs_target.insert(i);
1947 self_procs.insert(rank);
1949 MEDCoupling::MEDCouplingUMesh *mesh=0;
1950 MEDCoupling::ParaMESH *paramesh=0;
1951 MEDCoupling::ParaFIELD *parafieldP0=0;
1953 MEDCoupling::CommInterface interface;
1955 ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
1956 ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
1957 ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
1959 MPI_Barrier(MPI_COMM_WORLD);
1960 if(source_group->containsMyRank())
1962 double coords[15]={1.,0.,0., 2.,0.,0., 2.,2.,0., 0.,2.,0., 0.5,0.5,1.};
1963 mcIdType conn[7]={0,1,2,3,0,3,4};
1964 mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
1965 mesh->allocateCells(2);
1966 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1967 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
1968 mesh->finishInsertingCells();
1969 DataArrayDouble *myCoords=DataArrayDouble::New();
1970 myCoords->alloc(5,3);
1971 std::copy(coords,coords+15,myCoords->getPointer());
1972 mesh->setCoords(myCoords);
1973 myCoords->decrRef();
1975 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
1976 MEDCoupling::ComponentTopology comptopo;
1977 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
1978 double *valueP0=parafieldP0->getField()->getArray()->getPointer();
1979 parafieldP0->getField()->setNature(IntensiveMaximum);
1980 valueP0[0]=7.; valueP0[1]=8.;
1984 const char targetMeshName[]="target mesh";
1987 double coords[12]={0.25,0.25,0.5, 0.,0.25,0.5, 0.,0.,0.5, 0.25,0.,0.5};
1988 mcIdType conn[4]={0,1,2,3};
1989 mesh=MEDCouplingUMesh::New("Target mesh Proc1",2);
1990 mesh->allocateCells(1);
1991 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
1992 mesh->finishInsertingCells();
1993 DataArrayDouble *myCoords=DataArrayDouble::New();
1994 myCoords->alloc(4,3);
1995 std::copy(coords,coords+12,myCoords->getPointer());
1996 mesh->setCoords(myCoords);
1997 myCoords->decrRef();
1998 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
2002 double coords[12]={0.,0.25,0.5, 0.,0.,0.5, -1.,0.,0.5, -1.,0.25,0.5};
2003 mcIdType conn[4]={0,1,2,3};
2004 mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
2005 mesh->allocateCells(1);
2006 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
2007 mesh->finishInsertingCells();
2008 DataArrayDouble *myCoords=DataArrayDouble::New();
2009 myCoords->alloc(4,3);
2010 std::copy(coords,coords+12,myCoords->getPointer());
2011 mesh->setCoords(myCoords);
2012 myCoords->decrRef();
2013 paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
2015 MEDCoupling::ComponentTopology comptopo;
2016 parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2017 parafieldP0->getField()->setNature(IntensiveMaximum);
2020 MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
2021 if (source_group->containsMyRank())
2023 dec.setMethod("P0");
2024 dec.attachLocalField(parafieldP0);
2026 // dec.setForcedRenormalization(false);
2029 // const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
2032 // CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
2033 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
2037 // CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
2041 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
2046 dec.setMethod("P0");
2047 dec.attachLocalField(parafieldP0);
2049 // dec.setForcedRenormalization(false);
2051 // const double *res=parafieldP0->getField()->getArray()->getConstPointer();
2054 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
2055 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
2056 // CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
2060 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
2061 // CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
2062 // CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
2071 delete target_group;
2072 delete source_group;
2074 MPI_Barrier(MPI_COMM_WORLD);
2078 * Tests an asynchronous exchange between two codes
2079 * one sends data with dtA as an interval, the max time being tmaxA
2080 * the other one receives with dtB as an interval, the max time being tmaxB
2082 void ParaMEDMEMTest::testAsynchronousInterpKernelDEC_2D(double dtA, double tmaxA,
2083 double dtB, double tmaxB, bool WithPointToPoint, bool Asynchronous,
2084 bool WithInterp, const char *srcMeth, const char *targetMeth)
2086 std::string srcM(srcMeth);
2087 std::string targetM(targetMeth);
2090 MPI_Comm_size(MPI_COMM_WORLD,&size);
2091 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
2093 //the test is meant to run on five processors
2094 if (size !=5) return ;
2096 int nproc_source = 3;
2097 set<int> self_procs;
2098 set<int> procs_source;
2099 set<int> procs_target;
2101 for (int i=0; i<nproc_source; i++)
2102 procs_source.insert(i);
2103 for (int i=nproc_source; i<size; i++)
2104 procs_target.insert(i);
2105 self_procs.insert(rank);
2107 MEDCoupling::CommInterface interface;
2109 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
2110 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
2111 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
2113 //loading the geometry for the source group
2115 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
2117 MEDCoupling::MEDCouplingUMesh* mesh = nullptr;
2118 MEDCoupling::ParaMESH* paramesh = nullptr;
2119 MEDCoupling::ParaFIELD* parafield = nullptr;
2120 ICoCo::MEDDoubleField* icocofield = nullptr;
2122 char * tmp_dir_c = getenv("TMP");
2124 if (tmp_dir_c != NULL)
2125 tmp_dir = string(tmp_dir_c);
2128 string filename_xml1 = "square1_split";
2129 string filename_xml2 = "square2_split";
2130 //string filename_seq_wr = makeTmpFile("");
2131 //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
2133 // To remove tmp files from disk
2134 ParaMEDMEMTest_TmpFilesRemover aRemover;
2136 MPI_Barrier(MPI_COMM_WORLD);
2138 if (source_group->containsMyRank())
2140 string master = filename_xml1;
2142 ostringstream strstream;
2143 strstream <<master<<rank+1<<".med";
2144 string fName = INTERP_TEST::getResourceFile(strstream.str());
2145 ostringstream meshname ;
2146 meshname<< "Mesh_2_"<< rank+1;
2148 mesh=ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
2150 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
2152 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
2153 MEDCoupling::ComponentTopology comptopo;
2156 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2157 parafield->getField()->setNature(IntensiveMaximum);//InvertIntegral);//IntensiveMaximum);
2160 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
2164 nb_local=mesh->getNumberOfCells();
2166 nb_local=mesh->getNumberOfNodes();
2167 // double * value= new double[nb_local];
2168 double *value=parafield->getField()->getArray()->getPointer();
2169 for(int ielem=0; ielem<nb_local;ielem++)
2172 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
2173 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
2175 dec.attachLocalField(icocofield);
2180 //loading the geometry for the target group
2181 if (target_group->containsMyRank())
2183 string master= filename_xml2;
2184 ostringstream strstream;
2185 strstream << master<<(rank-nproc_source+1)<<".med";
2186 string fName = INTERP_TEST::getResourceFile(strstream.str());
2187 ostringstream meshname ;
2188 meshname<< "Mesh_3_"<<rank-nproc_source+1;
2190 mesh = ReadUMeshFromFile(fName.c_str(),meshname.str().c_str(),0);
2192 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
2193 // MEDCoupling::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
2194 MEDCoupling::ComponentTopology comptopo;
2197 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
2198 parafield->getField()->setNature(IntensiveMaximum);//InvertIntegral);//IntensiveMaximum);
2201 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
2205 nb_local=mesh->getNumberOfCells();
2207 nb_local=mesh->getNumberOfNodes();
2209 double *value=parafield->getField()->getArray()->getPointer();
2210 for(int ielem=0; ielem<nb_local;ielem++)
2212 // ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
2213 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
2215 dec.attachLocalField(icocofield);
2219 //attaching a DEC to the source group
2221 if (source_group->containsMyRank())
2223 cout<<"DEC usage"<<endl;
2224 dec.setAsynchronous(Asynchronous);
2226 dec.setTimeInterpolationMethod(LinearTimeInterp);
2228 if ( WithPointToPoint ) {
2229 dec.setAllToAllMethod(PointToPoint);
2232 dec.setAllToAllMethod(Native);
2235 dec.setForcedRenormalization(false);
2236 for (double time=0; time<tmaxA+1e-10; time+=dtA)
2238 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2239 << " dtA " << dtA << " tmaxA " << tmaxA << endl ;
2240 if ( time+dtA < tmaxA+1e-7 ) {
2241 dec.sendData( time , dtA );
2244 dec.sendData( time , 0 );
2246 double* value = parafield->getField()->getArray()->getPointer();
2247 int nb_local=parafield->getField()->getMesh()->getNumberOfCells();
2248 for (int i=0; i<nb_local;i++)
2255 //attaching a DEC to the target group
2256 if (target_group->containsMyRank())
2258 cout<<"DEC usage"<<endl;
2259 dec.setAsynchronous(Asynchronous);
2261 dec.setTimeInterpolationMethod(LinearTimeInterp);
2263 if ( WithPointToPoint ) {
2264 dec.setAllToAllMethod(PointToPoint);
2267 dec.setAllToAllMethod(Native);
2270 dec.setForcedRenormalization(false);
2271 vector<double> times;
2272 for (double time=0; time<tmaxB+1e-10; time+=dtB)
2274 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2275 << " dtB " << dtB << " tmaxB " << tmaxB << endl ;
2276 dec.recvData( time );
2277 double vi = parafield->getVolumeIntegral(0,true);
2278 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
2279 << " VolumeIntegral " << vi
2280 << " time*10000 " << time*10000 << endl ;
2282 CPPUNIT_ASSERT_DOUBLES_EQUAL(vi,time*10000,0.001);
2287 delete source_group;
2288 delete target_group;
2295 cout << "testAsynchronousInterpKernelDEC_2D" << rank << " MPI_Barrier " << endl ;
2297 if (Asynchronous) MPI_Barrier(MPI_COMM_WORLD);
2298 cout << "end of InterpKernelDEC_2D test"<<endl;