Back to index

salome-med  6.5.0
ParaMEDMEMTest_InterpKernelDEC.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "ParaMEDMEMTest.hxx"
00021 #include <cppunit/TestAssert.h>
00022 
00023 #include "CommInterface.hxx"
00024 #include "ProcessorGroup.hxx"
00025 #include "MPIProcessorGroup.hxx"
00026 #include "Topology.hxx"
00027 #include "DEC.hxx"
00028 #include "MxN_Mapping.hxx"
00029 #include "InterpKernelDEC.hxx"
00030 #include "ParaMESH.hxx"
00031 #include "ParaFIELD.hxx"
00032 #include "ComponentTopology.hxx"
00033 #include "ICoCoMEDField.hxx"
00034 #include "ParaMEDLoader.hxx"
00035 #include "MEDLoader.hxx"
00036 
00037  
00038 #include <string>
00039 #include <iterator>
00040 
00041 // use this define to enable lines, execution of which leads to Segmentation Fault
00042 #define ENABLE_FAULTS
00043 
00044 // use this define to enable CPPUNIT asserts and fails, showing bugs
00045 #define ENABLE_FORCED_FAILURES
00046 
00047 
00048 using namespace std;
00049 using namespace ParaMEDMEM;
00050 
00051 void ParaMEDMEMTest::testInterpKernelDEC_2D()
00052 {
00053   testInterpKernelDEC_2D_("P0","P0");
00054 }
00055 
00056 void ParaMEDMEMTest::testInterpKernelDEC2_2D()
00057 {
00058   testInterpKernelDEC2_2D_("P0","P0");
00059 }
00060 
00061 void ParaMEDMEMTest::testInterpKernelDEC_3D()
00062 {
00063   testInterpKernelDEC_3D_("P0","P0");
00064 }
00065 
00066 void ParaMEDMEMTest::testInterpKernelDEC_2DP0P1()
00067 {
00068   //testInterpKernelDEC_2D_("P0","P1");
00069 }
00070 
00071 void ParaMEDMEMTest::testInterpKernelDEC_1D()
00072 {
00073   int size;
00074   int rank;
00075   MPI_Comm_size(MPI_COMM_WORLD,&size);
00076   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00077   //
00078   if(size!=5)
00079     return ;
00080   int nproc_source = 3;
00081   set<int> self_procs;
00082   set<int> procs_source;
00083   set<int> procs_target;
00084   
00085   for (int i=0; i<nproc_source; i++)
00086     procs_source.insert(i);
00087   for (int i=nproc_source; i<size; i++)
00088     procs_target.insert(i);
00089   self_procs.insert(rank);
00090   //
00091   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
00092   ParaMEDMEM::ParaMESH *paramesh=0;
00093   ParaMEDMEM::ParaFIELD *parafieldP0=0;
00094   //
00095   ParaMEDMEM::CommInterface interface;
00096   //
00097   ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00098   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00099   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00100   //
00101   MPI_Barrier(MPI_COMM_WORLD);
00102   if(source_group->containsMyRank())
00103     {
00104       if(rank==0)
00105         {
00106           double coords[4]={0.3,0.7, 0.9,1.0};
00107           int conn[4]={0,1,2,3};
00108           mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
00109           mesh->allocateCells(2);
00110           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00111           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
00112           mesh->finishInsertingCells();
00113           DataArrayDouble *myCoords=DataArrayDouble::New();
00114           myCoords->alloc(4,1);
00115           std::copy(coords,coords+4,myCoords->getPointer());
00116           mesh->setCoords(myCoords);
00117           myCoords->decrRef();
00118         }
00119       if(rank==1)
00120         {
00121           double coords[2]={0.7,0.9};
00122           int conn[2]={0,1};
00123           mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
00124           mesh->allocateCells(1);
00125           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00126           mesh->finishInsertingCells();
00127           DataArrayDouble *myCoords=DataArrayDouble::New();
00128           myCoords->alloc(2,1);
00129           std::copy(coords,coords+2,myCoords->getPointer());
00130           mesh->setCoords(myCoords);
00131           myCoords->decrRef();
00132         }
00133       if(rank==2)
00134         {
00135           double coords[2]={1.,1.12};
00136           int conn[2]={0,1};
00137           mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
00138           mesh->allocateCells(1);
00139           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00140           mesh->finishInsertingCells();
00141           DataArrayDouble *myCoords=DataArrayDouble::New();
00142           myCoords->alloc(2,1);
00143           std::copy(coords,coords+2,myCoords->getPointer());
00144           mesh->setCoords(myCoords);
00145           myCoords->decrRef();
00146         }
00147       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
00148       ParaMEDMEM::ComponentTopology comptopo;
00149       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00150       double *valueP0=parafieldP0->getField()->getArray()->getPointer();
00151       parafieldP0->getField()->setNature(ConservativeVolumic);
00152       if(rank==0)
00153         {
00154           valueP0[0]=7.; valueP0[1]=8.;
00155         }
00156       if(rank==1)
00157         {
00158           valueP0[0]=9.;
00159         }
00160       if(rank==2)
00161         {
00162           valueP0[0]=10.;
00163         }
00164     }
00165   else
00166     {
00167       const char targetMeshName[]="target mesh";
00168       if(rank==3)
00169         {
00170           double coords[2]={0.5,0.75};
00171           int conn[2]={0,1};
00172           mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
00173           mesh->allocateCells(1);
00174           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00175           mesh->finishInsertingCells();
00176           DataArrayDouble *myCoords=DataArrayDouble::New();
00177           myCoords->alloc(2,1);
00178           std::copy(coords,coords+2,myCoords->getPointer());
00179           mesh->setCoords(myCoords);
00180           myCoords->decrRef();
00181           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
00182         }
00183       if(rank==4)
00184         {
00185           double coords[2]={0.75,1.2};
00186           int conn[2]={0,1};
00187           mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
00188           mesh->allocateCells(1);
00189           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00190           mesh->finishInsertingCells();
00191           DataArrayDouble *myCoords=DataArrayDouble::New();
00192           myCoords->alloc(2,1);
00193           std::copy(coords,coords+2,myCoords->getPointer());
00194           mesh->setCoords(myCoords);
00195           myCoords->decrRef();
00196           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
00197         }
00198       ParaMEDMEM::ComponentTopology comptopo;
00199       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00200       parafieldP0->getField()->setNature(ConservativeVolumic);
00201     }
00202   // test 1
00203   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
00204   if (source_group->containsMyRank())
00205     { 
00206       dec.setMethod("P0");
00207       dec.attachLocalField(parafieldP0);
00208       dec.synchronize();
00209       dec.setForcedRenormalization(false);
00210       dec.sendData();
00211       dec.recvData();
00212       const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
00213       if(rank==0)
00214         {
00215           CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
00216           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
00217         }
00218       if(rank==1)
00219         {
00220           CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
00221         }
00222       if(rank==2)
00223         {
00224           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
00225         }
00226     }
00227   else
00228     {
00229       dec.setMethod("P0");
00230       dec.attachLocalField(parafieldP0);
00231       dec.synchronize();
00232       dec.setForcedRenormalization(false);
00233       dec.recvData();
00234       const double *res=parafieldP0->getField()->getArray()->getConstPointer();
00235       if(rank==3)
00236         {
00237           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
00238           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
00239           CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
00240         }
00241       if(rank==4)
00242         {
00243           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
00244           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
00245           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
00246         }
00247       dec.sendData();
00248     }
00249   //
00250   delete parafieldP0;
00251   mesh->decrRef();
00252   delete paramesh;
00253   delete self_group;
00254   delete target_group;
00255   delete source_group;
00256   //
00257   MPI_Barrier(MPI_COMM_WORLD);
00258 }
00259 
00260 void ParaMEDMEMTest::testInterpKernelDEC_2DCurve()
00261 {
00262   int size;
00263   int rank;
00264   MPI_Comm_size(MPI_COMM_WORLD,&size);
00265   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00266   //
00267   if(size!=5)
00268     return ;
00269   int nproc_source = 3;
00270   set<int> self_procs;
00271   set<int> procs_source;
00272   set<int> procs_target;
00273   
00274   for (int i=0; i<nproc_source; i++)
00275     procs_source.insert(i);
00276   for (int i=nproc_source; i<size; i++)
00277     procs_target.insert(i);
00278   self_procs.insert(rank);
00279   //
00280   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
00281   ParaMEDMEM::ParaMESH *paramesh=0;
00282   ParaMEDMEM::ParaFIELD *parafieldP0=0;
00283   //
00284   ParaMEDMEM::CommInterface interface;
00285   //
00286   ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00287   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00288   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00289   //
00290   MPI_Barrier(MPI_COMM_WORLD);
00291   if(source_group->containsMyRank())
00292     {
00293       if(rank==0)
00294         {
00295           double coords[8]={0.3,0.3,0.7,0.7, 0.9,0.9,1.0,1.0};
00296           int conn[4]={0,1,2,3};
00297           mesh=MEDCouplingUMesh::New("Source mesh Proc0",1);
00298           mesh->allocateCells(2);
00299           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00300           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn+2);
00301           mesh->finishInsertingCells();
00302           DataArrayDouble *myCoords=DataArrayDouble::New();
00303           myCoords->alloc(4,2);
00304           std::copy(coords,coords+8,myCoords->getPointer());
00305           mesh->setCoords(myCoords);
00306           myCoords->decrRef();
00307         }
00308       if(rank==1)
00309         {
00310           double coords[4]={0.7,0.7,0.9,0.9};
00311           int conn[2]={0,1};
00312           mesh=MEDCouplingUMesh::New("Source mesh Proc1",1);
00313           mesh->allocateCells(1);
00314           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00315           mesh->finishInsertingCells();
00316           DataArrayDouble *myCoords=DataArrayDouble::New();
00317           myCoords->alloc(2,2);
00318           std::copy(coords,coords+4,myCoords->getPointer());
00319           mesh->setCoords(myCoords);
00320           myCoords->decrRef();
00321         }
00322       if(rank==2)
00323         {
00324           double coords[4]={1.,1.,1.12,1.12};
00325           int conn[2]={0,1};
00326           mesh=MEDCouplingUMesh::New("Source mesh Proc2",1);
00327           mesh->allocateCells(1);
00328           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00329           mesh->finishInsertingCells();
00330           DataArrayDouble *myCoords=DataArrayDouble::New();
00331           myCoords->alloc(2,2);
00332           std::copy(coords,coords+4,myCoords->getPointer());
00333           mesh->setCoords(myCoords);
00334           myCoords->decrRef();
00335         }
00336       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
00337       ParaMEDMEM::ComponentTopology comptopo;
00338       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00339       double *valueP0=parafieldP0->getField()->getArray()->getPointer();
00340       parafieldP0->getField()->setNature(ConservativeVolumic);
00341       if(rank==0)
00342         {
00343           valueP0[0]=7.; valueP0[1]=8.;
00344         }
00345       if(rank==1)
00346         {
00347           valueP0[0]=9.;
00348         }
00349       if(rank==2)
00350         {
00351           valueP0[0]=10.;
00352         }
00353     }
00354   else
00355     {
00356       const char targetMeshName[]="target mesh";
00357       if(rank==3)
00358         {
00359           double coords[4]={0.5,0.5,0.75,0.75};
00360           int conn[2]={0,1};
00361           mesh=MEDCouplingUMesh::New("Target mesh Proc3",1);
00362           mesh->allocateCells(1);
00363           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00364           mesh->finishInsertingCells();
00365           DataArrayDouble *myCoords=DataArrayDouble::New();
00366           myCoords->alloc(2,2);
00367           std::copy(coords,coords+4,myCoords->getPointer());
00368           mesh->setCoords(myCoords);
00369           myCoords->decrRef();
00370           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
00371         }
00372       if(rank==4)
00373         {
00374           double coords[4]={0.75,0.75,1.2,1.2};
00375           int conn[2]={0,1};
00376           mesh=MEDCouplingUMesh::New("Target mesh Proc4",1);
00377           mesh->allocateCells(1);
00378           mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,conn);
00379           mesh->finishInsertingCells();
00380           DataArrayDouble *myCoords=DataArrayDouble::New();
00381           myCoords->alloc(2,2);
00382           std::copy(coords,coords+4,myCoords->getPointer());
00383           mesh->setCoords(myCoords);
00384           myCoords->decrRef();
00385           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
00386         }
00387       ParaMEDMEM::ComponentTopology comptopo;
00388       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00389       parafieldP0->getField()->setNature(ConservativeVolumic);
00390     }
00391   // test 1
00392   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
00393   if (source_group->containsMyRank())
00394     { 
00395       dec.setMethod("P0");
00396       dec.attachLocalField(parafieldP0);
00397       dec.synchronize();
00398       dec.setForcedRenormalization(false);
00399       dec.sendData();
00400       dec.recvData();
00401       const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
00402       if(rank==0)
00403         {
00404           CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
00405           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
00406         }
00407       if(rank==1)
00408         {
00409           CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
00410         }
00411       if(rank==2)
00412         {
00413           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
00414         }
00415     }
00416   else
00417     {
00418       dec.setMethod("P0");
00419       dec.attachLocalField(parafieldP0);
00420       dec.synchronize();
00421       dec.setForcedRenormalization(false);
00422       dec.recvData();
00423       const double *res=parafieldP0->getField()->getArray()->getConstPointer();
00424       if(rank==3)
00425         {
00426           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
00427           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
00428           CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
00429         }
00430       if(rank==4)
00431         {
00432           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
00433           CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
00434           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
00435         }
00436       dec.sendData();
00437     }
00438   //
00439   delete parafieldP0;
00440   mesh->decrRef();
00441   delete paramesh;
00442   delete self_group;
00443   delete target_group;
00444   delete source_group;
00445   //
00446   MPI_Barrier(MPI_COMM_WORLD);
00447 }
00448 
00449 
00450 /*
00451  * Check methods defined in InterpKernelDEC.hxx
00452  *
00453  InterpKernelDEC();
00454  InterpKernelDEC(ProcessorGroup& local_group, ProcessorGroup& distant_group);
00455  virtual ~InterpKernelDEC();
00456  void synchronize();
00457  void recvData();
00458  void sendData();
00459 */
00460  
00461 void ParaMEDMEMTest::testInterpKernelDEC_2D_(const char *srcMeth, const char *targetMeth)
00462 {
00463   std::string srcM(srcMeth);
00464   std::string targetM(targetMeth);
00465   int size;
00466   int rank;
00467   MPI_Comm_size(MPI_COMM_WORLD,&size);
00468   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00469 
00470   //the test is meant to run on five processors
00471   if (size !=5) return ;
00472    
00473   int nproc_source = 3;
00474   set<int> self_procs;
00475   set<int> procs_source;
00476   set<int> procs_target;
00477   
00478   for (int i=0; i<nproc_source; i++)
00479     procs_source.insert(i);
00480   for (int i=nproc_source; i<size; i++)
00481     procs_target.insert(i);
00482   self_procs.insert(rank);
00483   
00484   ParaMEDMEM::CommInterface interface;
00485     
00486   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00487   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00488   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00489   
00490   //loading the geometry for the source group
00491 
00492   ParaMEDMEM::InterpKernelDEC dec (*source_group,*target_group);
00493 
00494   ParaMEDMEM::MEDCouplingUMesh* mesh;
00495   ParaMEDMEM::ParaMESH* paramesh;
00496   ParaMEDMEM::ParaFIELD* parafield;
00497   ICoCo::Field* icocofield ;
00498   
00499   string filename_xml1              = getResourceFile("square1_split");
00500   string filename_xml2              = getResourceFile("square2_split");
00501   //string filename_seq_wr            = makeTmpFile("");
00502   //string filename_seq_med           = makeTmpFile("myWrField_seq_pointe221.med");
00503   
00504   // To remove tmp files from disk
00505   ParaMEDMEMTest_TmpFilesRemover aRemover;
00506   
00507   MPI_Barrier(MPI_COMM_WORLD);
00508   if (source_group->containsMyRank())
00509     {
00510       string master = filename_xml1;
00511       
00512       ostringstream strstream;
00513       strstream <<master<<rank+1<<".med";
00514       ostringstream meshname ;
00515       meshname<< "Mesh_2_"<< rank+1;
00516       
00517       mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00518       
00519     
00520       paramesh=new ParaMESH (mesh,*source_group,"source mesh");
00521     
00522       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
00523       ParaMEDMEM::ComponentTopology comptopo;
00524       if(srcM=="P0")
00525         {
00526           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00527           parafield->getField()->setNature(ConservativeVolumic);
00528         }
00529       else
00530         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
00531       int nb_local;
00532       if(srcM=="P0")
00533         nb_local=mesh->getNumberOfCells();
00534       else
00535         nb_local=mesh->getNumberOfNodes();
00536       //      double * value= new double[nb_local];
00537       double *value=parafield->getField()->getArray()->getPointer();
00538       for(int ielem=0; ielem<nb_local;ielem++)
00539         value[ielem]=1.0;
00540     
00541       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
00542       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00543       dec.setMethod(srcMeth);
00544       dec.attachLocalField(icocofield);
00545     }
00546   
00547   //loading the geometry for the target group
00548   if (target_group->containsMyRank())
00549     {
00550       string master= filename_xml2;
00551       ostringstream strstream;
00552       strstream << master<<(rank-nproc_source+1)<<".med";
00553       ostringstream meshname ;
00554       meshname<< "Mesh_3_"<<rank-nproc_source+1;
00555       mesh = MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00556       
00557       paramesh=new ParaMESH (mesh,*target_group,"target mesh");
00558       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
00559       ParaMEDMEM::ComponentTopology comptopo;
00560       if(targetM=="P0")
00561         {
00562           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00563           parafield->getField()->setNature(ConservativeVolumic);
00564         }
00565       else
00566         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
00567       int nb_local;
00568       if(targetM=="P0")
00569         nb_local=mesh->getNumberOfCells();
00570       else
00571         nb_local=mesh->getNumberOfNodes();
00572       //      double * value= new double[nb_local];
00573       double *value=parafield->getField()->getArray()->getPointer();
00574       for(int ielem=0; ielem<nb_local;ielem++)
00575         value[ielem]=0.0;
00576       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
00577       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00578       dec.setMethod(targetMeth);
00579       dec.attachLocalField(icocofield);
00580     }
00581     
00582   
00583   //attaching a DEC to the source group 
00584   double field_before_int;
00585   double field_after_int;
00586   
00587   if (source_group->containsMyRank())
00588     { 
00589       field_before_int = parafield->getVolumeIntegral(0,true);
00590       dec.synchronize();
00591       cout<<"DEC usage"<<endl;
00592       dec.setForcedRenormalization(false);
00593 
00594       dec.sendData();
00595       ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
00596       if (source_group->myRank()==0)
00597         aRemover.Register("./sourcesquareb");
00598       ostringstream filename;
00599       filename<<"./sourcesquareb_"<<source_group->myRank()+1;
00600       aRemover.Register(filename.str().c_str());
00601       //MEDLoader::WriteField("./sourcesquareb",parafield->getField());
00602    
00603       dec.recvData();
00604       cout <<"writing"<<endl;
00605       ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
00606       if (source_group->myRank()==0)
00607         aRemover.Register("./sourcesquare");
00608       //MEDLoader::WriteField("./sourcesquare",parafield->getField());
00609       
00610      
00611       filename<<"./sourcesquare_"<<source_group->myRank()+1;
00612       aRemover.Register(filename.str().c_str());
00613       field_after_int = parafield->getVolumeIntegral(0,true);
00614       
00615       
00616       //      MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
00617       //       MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
00618 
00619       CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
00620     
00621     }
00622   
00623   //attaching a DEC to the target group
00624   if (target_group->containsMyRank())
00625     {
00626       dec.synchronize();
00627       dec.setForcedRenormalization(false);
00628 
00629       dec.recvData();
00630       ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
00631       //MEDLoader::WriteField("./targetsquareb",parafield->getField());
00632       if (target_group->myRank()==0)
00633         aRemover.Register("./targetsquareb");
00634       ostringstream filename;
00635       filename<<"./targetsquareb_"<<target_group->myRank()+1;
00636       aRemover.Register(filename.str().c_str());
00637       dec.sendData();
00638       ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
00639       //MEDLoader::WriteField("./targetsquare",parafield->getField());
00640       
00641       if (target_group->myRank()==0)
00642         aRemover.Register("./targetsquareb");
00643       
00644       filename<<"./targetsquareb_"<<target_group->myRank()+1;
00645       aRemover.Register(filename.str().c_str());
00646       //    double field_before_int, field_after_int;
00647       //       MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
00648       //       MPI_Bcast(&field_after_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
00649       
00650       //      CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
00651     
00652     }
00653   
00654   delete source_group;
00655   delete target_group;
00656   delete self_group;
00657   delete parafield;
00658   delete paramesh;
00659   mesh->decrRef();
00660 
00661   delete icocofield;
00662 
00663   MPI_Barrier(MPI_COMM_WORLD);
00664   cout << "end of InterpKernelDEC_2D test"<<endl;
00665 }
00666 
00667 void ParaMEDMEMTest::testInterpKernelDEC2_2D_(const char *srcMeth, const char *targetMeth)
00668 {
00669   std::string srcM(srcMeth);
00670   std::string targetM(targetMeth);
00671   int size;
00672   int rank;
00673   MPI_Comm_size(MPI_COMM_WORLD,&size);
00674   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00675 
00676   //the test is meant to run on five processors
00677   if (size !=5) return ;
00678    
00679   int nproc_source = 3;
00680   set<int> self_procs;
00681   set<int> procs_source;
00682   set<int> procs_target;
00683   
00684   for (int i=0; i<nproc_source; i++)
00685     procs_source.insert(i);
00686   for (int i=nproc_source; i<size; i++)
00687     procs_target.insert(i);
00688   self_procs.insert(rank);
00689   
00690   ParaMEDMEM::CommInterface interface;
00691     
00692   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00693   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00694   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00695   
00696   //loading the geometry for the source group
00697 
00698   ParaMEDMEM::InterpKernelDEC dec (*source_group,*target_group);
00699 
00700   ParaMEDMEM::MEDCouplingUMesh* mesh;
00701   ParaMEDMEM::MEDCouplingFieldDouble* mcfield;
00702   
00703   string filename_xml1              = getResourceFile("square1_split");
00704   string filename_xml2              = getResourceFile("square2_split");
00705   
00706   // To remove tmp files from disk
00707   ParaMEDMEMTest_TmpFilesRemover aRemover;
00708   
00709   MPI_Barrier(MPI_COMM_WORLD);
00710   if (source_group->containsMyRank())
00711     {
00712       string master = filename_xml1;
00713       
00714       ostringstream strstream;
00715       strstream <<master<<rank+1<<".med";
00716       ostringstream meshname ;
00717       meshname<< "Mesh_2_"<< rank+1;
00718       
00719       mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00720       ParaMEDMEM::ComponentTopology comptopo;
00721       if(srcM=="P0")
00722         {
00723           mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
00724           mcfield->setMesh(mesh);
00725           DataArrayDouble *array=DataArrayDouble::New();
00726           array->alloc(mcfield->getNumberOfTuples(),1);
00727           mcfield->setArray(array);
00728           array->decrRef();
00729           mcfield->setNature(ConservativeVolumic);
00730         }
00731       else
00732         {
00733           mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
00734           mcfield->setMesh(mesh);
00735           DataArrayDouble *array=DataArrayDouble::New();
00736           array->alloc(mcfield->getNumberOfTuples(),1);
00737           mcfield->setArray(array);
00738           array->decrRef();
00739         }
00740       int nb_local;
00741       if(srcM=="P0")
00742         nb_local=mesh->getNumberOfCells();
00743       else
00744         nb_local=mesh->getNumberOfNodes();
00745       double *value=mcfield->getArray()->getPointer();
00746       for(int ielem=0; ielem<nb_local;ielem++)
00747         value[ielem]=1.0;
00748       dec.setMethod(srcMeth);
00749       dec.attachLocalField(mcfield);
00750       dec.attachLocalField(mcfield);
00751     }
00752   
00753   //loading the geometry for the target group
00754   if (target_group->containsMyRank())
00755     {
00756       string master= filename_xml2;
00757       ostringstream strstream;
00758       strstream << master<<(rank-nproc_source+1)<<".med";
00759       ostringstream meshname ;
00760       meshname<< "Mesh_3_"<<rank-nproc_source+1;
00761       mesh = MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00762       ParaMEDMEM::ComponentTopology comptopo;
00763       if(targetM=="P0")
00764         {
00765           mcfield = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
00766           mcfield->setMesh(mesh);
00767           DataArrayDouble *array=DataArrayDouble::New();
00768           array->alloc(mcfield->getNumberOfTuples(),1);
00769           mcfield->setArray(array);
00770           array->decrRef();
00771           mcfield->setNature(ConservativeVolumic);
00772         }
00773       else
00774         {
00775           mcfield = MEDCouplingFieldDouble::New(ON_NODES,NO_TIME);
00776           mcfield->setMesh(mesh);
00777           DataArrayDouble *array=DataArrayDouble::New();
00778           array->alloc(mcfield->getNumberOfTuples(),1);
00779           mcfield->setArray(array);
00780           array->decrRef();
00781         }
00782       int nb_local;
00783       if(targetM=="P0")
00784         nb_local=mesh->getNumberOfCells();
00785       else
00786         nb_local=mesh->getNumberOfNodes();
00787       double *value=mcfield->getArray()->getPointer();
00788       for(int ielem=0; ielem<nb_local;ielem++)
00789         value[ielem]=0.0;
00790       dec.setMethod(targetMeth);
00791       dec.attachLocalField(mcfield);
00792       dec.attachLocalField(mcfield);
00793     }
00794     
00795   
00796   //attaching a DEC to the source group 
00797 
00798   if (source_group->containsMyRank())
00799     { 
00800       dec.synchronize();
00801       dec.setForcedRenormalization(false);
00802       dec.sendData();
00803       dec.recvData();
00804     }
00805   
00806   //attaching a DEC to the target group
00807   if (target_group->containsMyRank())
00808     {
00809       dec.synchronize();
00810       dec.setForcedRenormalization(false);
00811       dec.recvData();
00812       dec.sendData();
00813     }
00814   delete source_group;
00815   delete target_group;
00816   delete self_group;
00817   mcfield->decrRef();
00818   mesh->decrRef();
00819 
00820   MPI_Barrier(MPI_COMM_WORLD);
00821   cout << "end of InterpKernelDEC2_2D test"<<endl;
00822 }
00823 
00824 void ParaMEDMEMTest::testInterpKernelDEC_3D_(const char *srcMeth, const char *targetMeth)
00825 {
00826   std::string srcM(srcMeth);
00827   std::string targetM(targetMeth);
00828   int size;
00829   int rank;
00830   MPI_Comm_size(MPI_COMM_WORLD,&size);
00831   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00832 
00833   //the test is meant to run on five processors
00834   if (size !=3) return ;
00835    
00836   int nproc_source = 2;
00837   set<int> self_procs;
00838   set<int> procs_source;
00839   set<int> procs_target;
00840   
00841   for (int i=0; i<nproc_source; i++)
00842     procs_source.insert(i);
00843   for (int i=nproc_source; i<size; i++)
00844     procs_target.insert(i);
00845   self_procs.insert(rank);
00846   
00847   ParaMEDMEM::CommInterface interface;
00848     
00849   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00850   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00851   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00852   
00853   //loading the geometry for the source group
00854 
00855   ParaMEDMEM::InterpKernelDEC dec (*source_group,*target_group);
00856 
00857   ParaMEDMEM::MEDCouplingUMesh* mesh;
00858   ParaMEDMEM::ParaMESH* paramesh;
00859   ParaMEDMEM::ParaFIELD* parafield;
00860   ICoCo::Field* icocofield ;
00861   
00862   string tmp_dir                    = getenv("TMP");
00863   if (tmp_dir == "")
00864     tmp_dir = "/tmp";
00865   string filename_xml1              = getResourceFile("Mesh3D_10_2d");
00866   string filename_xml2              = getResourceFile("Mesh3D_11");
00867   //string filename_seq_wr            = makeTmpFile("");
00868   //string filename_seq_med           = makeTmpFile("myWrField_seq_pointe221.med");
00869   
00870   // To remove tmp files from disk
00871   ParaMEDMEMTest_TmpFilesRemover aRemover;
00872   
00873   MPI_Barrier(MPI_COMM_WORLD);
00874   if (source_group->containsMyRank())
00875     {
00876       string master = filename_xml1;
00877       
00878       ostringstream strstream;
00879       strstream <<master<<rank+1<<".med";
00880       ostringstream meshname ;
00881       meshname<< "Mesh_3_"<< rank+1;
00882       
00883       mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00884       
00885     
00886       paramesh=new ParaMESH (mesh,*source_group,"source mesh");
00887     
00888       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
00889       ParaMEDMEM::ComponentTopology comptopo;
00890       if(srcM=="P0")
00891         {
00892           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00893           parafield->getField()->setNature(ConservativeVolumic);
00894         }
00895       else
00896         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
00897       int nb_local;
00898       if(srcM=="P0")
00899         nb_local=mesh->getNumberOfCells();
00900       else
00901         nb_local=mesh->getNumberOfNodes();
00902       //      double * value= new double[nb_local];
00903       double *value=parafield->getField()->getArray()->getPointer();
00904       for(int ielem=0; ielem<nb_local;ielem++)
00905         value[ielem]=1.0;
00906     
00907       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
00908       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00909       dec.setMethod(srcMeth);
00910       dec.attachLocalField(icocofield);
00911     }
00912   
00913   //loading the geometry for the target group
00914   if (target_group->containsMyRank())
00915     {
00916       string master= filename_xml2;
00917       ostringstream strstream;
00918       strstream << master << ".med";
00919       ostringstream meshname ;
00920       meshname<< "Mesh_6";
00921       mesh = MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00922       
00923       paramesh=new ParaMESH (mesh,*target_group,"target mesh");
00924       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
00925       ParaMEDMEM::ComponentTopology comptopo;
00926       if(targetM=="P0")
00927         {
00928           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00929           parafield->getField()->setNature(ConservativeVolumic);
00930         }
00931       else
00932         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
00933       int nb_local;
00934       if(targetM=="P0")
00935         nb_local=mesh->getNumberOfCells();
00936       else
00937         nb_local=mesh->getNumberOfNodes();
00938       //      double * value= new double[nb_local];
00939       double *value=parafield->getField()->getArray()->getPointer();
00940       for(int ielem=0; ielem<nb_local;ielem++)
00941         value[ielem]=0.0;
00942       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
00943       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00944       dec.setMethod(targetMeth);
00945       dec.attachLocalField(icocofield);
00946     }  
00947   //attaching a DEC to the source group 
00948   double field_before_int;
00949   double field_after_int;
00950   
00951   if (source_group->containsMyRank())
00952     { 
00953       field_before_int = parafield->getVolumeIntegral(0,true);
00954       dec.synchronize();
00955       cout<<"DEC usage"<<endl;
00956       dec.setForcedRenormalization(false);
00957 
00958       dec.sendData();
00959       ParaMEDLoader::WriteParaMesh("./sourcesquareb",paramesh);
00960       if (source_group->myRank()==0)
00961         aRemover.Register("./sourcesquareb");
00962       ostringstream filename;
00963       filename<<"./sourcesquareb_"<<source_group->myRank()+1;
00964       aRemover.Register(filename.str().c_str());
00965       //MEDLoader::WriteField("./sourcesquareb",parafield->getField());
00966    
00967       dec.recvData();
00968       cout <<"writing"<<endl;
00969       ParaMEDLoader::WriteParaMesh("./sourcesquare",paramesh);
00970       if (source_group->myRank()==0)
00971         aRemover.Register("./sourcesquare");
00972       //MEDLoader::WriteField("./sourcesquare",parafield->getField());
00973       
00974      
00975       filename<<"./sourcesquare_"<<source_group->myRank()+1;
00976       aRemover.Register(filename.str().c_str());
00977       field_after_int = parafield->getVolumeIntegral(0,true);
00978 
00979       CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, 1e-6);
00980     
00981     }
00982   
00983   //attaching a DEC to the target group
00984   if (target_group->containsMyRank())
00985     {
00986       dec.synchronize();
00987       dec.setForcedRenormalization(false);
00988 
00989       dec.recvData();
00990       ParaMEDLoader::WriteParaMesh("./targetsquareb",paramesh);
00991       //MEDLoader::WriteField("./targetsquareb",parafield->getField());
00992       if (target_group->myRank()==0)
00993         aRemover.Register("./targetsquareb");
00994       ostringstream filename;
00995       filename<<"./targetsquareb_"<<target_group->myRank()+1;
00996       aRemover.Register(filename.str().c_str());
00997       dec.sendData();
00998       ParaMEDLoader::WriteParaMesh("./targetsquare",paramesh);
00999       //MEDLoader::WriteField("./targetsquare",parafield->getField());
01000       
01001       if (target_group->myRank()==0)
01002         aRemover.Register("./targetsquareb");
01003       
01004       filename<<"./targetsquareb_"<<target_group->myRank()+1;
01005       aRemover.Register(filename.str().c_str());
01006     }
01007   delete source_group;
01008   delete target_group;
01009   delete self_group;
01010   delete parafield;
01011   delete paramesh;
01012   mesh->decrRef();
01013 
01014   delete icocofield;
01015 
01016   MPI_Barrier(MPI_COMM_WORLD);
01017   cout << "end of InterpKernelDEC_3D test"<<endl;
01018 }
01019 
01020 //Synchronous tests without interpolation with native mode (AllToAll(v) from lam/MPI:
01021 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpNativeDEC_2D()
01022 {
01023   testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,false,false,false,"P0","P0");
01024 }
01025 
01026 //Synchronous tests without interpolation :
01027 void ParaMEDMEMTest::testSynchronousEqualInterpKernelWithoutInterpDEC_2D()
01028 {
01029   testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,false,"P0","P0");
01030 }
01031 
01032 //Synchronous tests with interpolation :
01033 void ParaMEDMEMTest::testSynchronousEqualInterpKernelDEC_2D()
01034 {
01035   testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,false,true,"P0","P0");
01036 }
01037 void ParaMEDMEMTest::testSynchronousFasterSourceInterpKernelDEC_2D()
01038 {
01039   testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,false,true,"P0","P0");
01040 }
01041 void ParaMEDMEMTest::testSynchronousSlowerSourceInterpKernelDEC_2D()
01042 {
01043   testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,false,true,"P0","P0");
01044 }
01045 void ParaMEDMEMTest::testSynchronousSlowSourceInterpKernelDEC_2D()
01046 {
01047   testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,false,true,"P0","P0");
01048 }
01049 void ParaMEDMEMTest::testSynchronousFastSourceInterpKernelDEC_2D()
01050 {
01051   testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,false,true,"P0","P0");
01052 }
01053 
01054 //Asynchronous tests with interpolation :
01055 void ParaMEDMEMTest::testAsynchronousEqualInterpKernelDEC_2D()
01056 {
01057   testAsynchronousInterpKernelDEC_2D(0.1,1,0.1,1,true,true,true,"P0","P0");
01058 }
01059 void ParaMEDMEMTest::testAsynchronousFasterSourceInterpKernelDEC_2D()
01060 {
01061   testAsynchronousInterpKernelDEC_2D(0.09,1,0.1,1,true,true,true,"P0","P0");
01062 }
01063 void ParaMEDMEMTest::testAsynchronousSlowerSourceInterpKernelDEC_2D()
01064 {
01065   testAsynchronousInterpKernelDEC_2D(0.11,1,0.1,1,true,true,true,"P0","P0");
01066 }
01067 void ParaMEDMEMTest::testAsynchronousSlowSourceInterpKernelDEC_2D()
01068 {
01069   testAsynchronousInterpKernelDEC_2D(0.11,1,0.01,1,true,true,true,"P0","P0");
01070 }
01071 void ParaMEDMEMTest::testAsynchronousFastSourceInterpKernelDEC_2D()
01072 {
01073   testAsynchronousInterpKernelDEC_2D(0.01,1,0.11,1,true,true,true,"P0","P0");
01074 }
01075 
01076 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P0()
01077 {
01078   //
01079   const double sourceCoordsAll[2][8]={{0.4,0.5,0.4,1.5,1.6,1.5,1.6,0.5},
01080                                       {0.3,-0.5,1.6,-0.5,1.6,-1.5,0.3,-1.5}};
01081   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},
01082                                        {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},
01083                                        {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}};
01084   int conn4All[8]={0,1,2,3,4,5,6,7};
01085   double targetResults[3][2]={{34.,34.},{38.333333333333336,42.666666666666664},{47.,47.}};
01086   double targetResults2[3][2]={{0.28333333333333344,0.56666666666666687},{1.8564102564102569,2.0128205128205132},{1.0846153846153845,0.36153846153846159}};
01087   double targetResults3[3][2]={{3.7777777777777781,7.5555555555555562},{24.511111111111113,26.355555555555558},{14.1,4.7}};
01088   double targetResults4[3][2]={{8.5,17},{8.8461538461538431, 9.8461538461538449},{35.25,11.75}};
01089   //
01090   int size;
01091   int rank;
01092   MPI_Comm_size(MPI_COMM_WORLD,&size);
01093   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01094   //
01095   if(size!=5)
01096     return ;
01097   int nproc_source = 2;
01098   set<int> self_procs;
01099   set<int> procs_source;
01100   set<int> procs_target;
01101   
01102   for (int i=0; i<nproc_source; i++)
01103     procs_source.insert(i);
01104   for (int i=nproc_source; i<size; i++)
01105     procs_target.insert(i);
01106   self_procs.insert(rank);
01107   //
01108   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
01109   ParaMEDMEM::ParaMESH *paramesh=0;
01110   ParaMEDMEM::ParaFIELD* parafield=0;
01111   //
01112   ParaMEDMEM::CommInterface interface;
01113   //
01114   ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
01115   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
01116   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
01117   //
01118   MPI_Barrier(MPI_COMM_WORLD);
01119   if(source_group->containsMyRank())
01120     {
01121       std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
01122       mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
01123       mesh->allocateCells(2);
01124       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
01125       mesh->finishInsertingCells();
01126       DataArrayDouble *myCoords=DataArrayDouble::New();
01127       myCoords->alloc(4,2);
01128       const double *sourceCoords=sourceCoordsAll[rank];
01129       std::copy(sourceCoords,sourceCoords+8,myCoords->getPointer());
01130       mesh->setCoords(myCoords);
01131       myCoords->decrRef();
01132       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
01133       ParaMEDMEM::ComponentTopology comptopo;
01134       parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01135       double *value=parafield->getField()->getArray()->getPointer();
01136       value[0]=34+13*((double)rank);
01137     }
01138   else
01139     {
01140       std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
01141       mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
01142       mesh->allocateCells(2);
01143       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
01144       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All+4);
01145       mesh->finishInsertingCells();
01146       DataArrayDouble *myCoords=DataArrayDouble::New();
01147       myCoords->alloc(8,2);
01148       const double *targetCoords=targetCoordsAll[rank-nproc_source];
01149       std::copy(targetCoords,targetCoords+16,myCoords->getPointer());
01150       mesh->setCoords(myCoords);
01151       myCoords->decrRef();
01152       paramesh=new ParaMESH (mesh,*target_group,"target mesh");
01153       ParaMEDMEM::ComponentTopology comptopo;
01154       parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01155     }
01156   //test 1 - Conservative volumic
01157   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
01158   parafield->getField()->setNature(ConservativeVolumic);
01159   if (source_group->containsMyRank())
01160     { 
01161       dec.setMethod("P0");
01162       dec.attachLocalField(parafield);
01163       dec.synchronize();
01164       dec.setForcedRenormalization(false);
01165       dec.sendData();
01166     }
01167   else
01168     {
01169       dec.setMethod("P0");
01170       dec.attachLocalField(parafield);
01171       dec.synchronize();
01172       dec.setForcedRenormalization(false);
01173       dec.recvData();
01174       const double *res=parafield->getField()->getArray()->getConstPointer();
01175       const double *expected=targetResults[rank-nproc_source];
01176       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
01177       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
01178     }
01179   //test 2 - Integral
01180   ParaMEDMEM::InterpKernelDEC dec2(*source_group,*target_group);
01181   parafield->getField()->setNature(Integral);
01182   if (source_group->containsMyRank())
01183     { 
01184       dec2.setMethod("P0");
01185       dec2.attachLocalField(parafield);
01186       dec2.synchronize();
01187       dec2.setForcedRenormalization(false);
01188       dec2.sendData();
01189     }
01190   else
01191     {
01192       dec2.setMethod("P0");
01193       dec2.attachLocalField(parafield);
01194       dec2.synchronize();
01195       dec2.setForcedRenormalization(false);
01196       dec2.recvData();
01197       const double *res=parafield->getField()->getArray()->getConstPointer();
01198       const double *expected=targetResults2[rank-nproc_source];
01199       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
01200       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
01201     }
01202   //test 3 - Integral with global constraint
01203   ParaMEDMEM::InterpKernelDEC dec3(*source_group,*target_group);
01204   parafield->getField()->setNature(IntegralGlobConstraint);
01205   if (source_group->containsMyRank())
01206     { 
01207       dec3.setMethod("P0");
01208       dec3.attachLocalField(parafield);
01209       dec3.synchronize();
01210       dec3.setForcedRenormalization(false);
01211       dec3.sendData();
01212     }
01213   else
01214     {
01215       dec3.setMethod("P0");
01216       dec3.attachLocalField(parafield);
01217       dec3.synchronize();
01218       dec3.setForcedRenormalization(false);
01219       dec3.recvData();
01220       const double *res=parafield->getField()->getArray()->getConstPointer();
01221       const double *expected=targetResults3[rank-nproc_source];
01222       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
01223       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
01224     }
01225   //test 4 - RevIntegral
01226   ParaMEDMEM::InterpKernelDEC dec4(*source_group,*target_group);
01227   parafield->getField()->setNature(RevIntegral);
01228   if (source_group->containsMyRank())
01229     { 
01230       dec4.setMethod("P0");
01231       dec4.attachLocalField(parafield);
01232       dec4.synchronize();
01233       dec4.setForcedRenormalization(false);
01234       dec4.sendData();
01235     }
01236   else
01237     {
01238       dec4.setMethod("P0");
01239       dec4.attachLocalField(parafield);
01240       dec4.synchronize();
01241       dec4.setForcedRenormalization(false);
01242       dec4.recvData();
01243       const double *res=parafield->getField()->getArray()->getConstPointer();
01244       const double *expected=targetResults4[rank-nproc_source];
01245       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
01246       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
01247     }
01248   //test 5 - Conservative volumic reversed
01249   ParaMEDMEM::InterpKernelDEC dec5(*source_group,*target_group);
01250   parafield->getField()->setNature(ConservativeVolumic);
01251   if (source_group->containsMyRank())
01252     { 
01253       dec5.setMethod("P0");
01254       dec5.attachLocalField(parafield);
01255       dec5.synchronize();
01256       dec5.setForcedRenormalization(false);
01257       dec5.recvData();
01258       const double *res=parafield->getField()->getArray()->getConstPointer();
01259       CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
01260       const double expected[]={37.8518518518519,43.5333333333333};
01261       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
01262     }
01263   else
01264     {
01265       dec5.setMethod("P0");
01266       dec5.attachLocalField(parafield);
01267       dec5.synchronize();
01268       dec5.setForcedRenormalization(false);
01269       double *res=parafield->getField()->getArray()->getPointer();
01270       const double *toSet=targetResults[rank-nproc_source];
01271       res[0]=toSet[0];
01272       res[1]=toSet[1];
01273       dec5.sendData();
01274     }
01275   //test 6 - Integral reversed
01276   ParaMEDMEM::InterpKernelDEC dec6(*source_group,*target_group);
01277   parafield->getField()->setNature(Integral);
01278   if (source_group->containsMyRank())
01279     { 
01280       dec6.setMethod("P0");
01281       dec6.attachLocalField(parafield);
01282       dec6.synchronize();
01283       dec6.setForcedRenormalization(false);
01284       dec6.recvData();
01285       const double *res=parafield->getField()->getArray()->getConstPointer();
01286       CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
01287       const double expected[]={0.794600591715977,1.35631163708087};
01288       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
01289     }
01290   else
01291     {
01292       dec6.setMethod("P0");
01293       dec6.attachLocalField(parafield);
01294       dec6.synchronize();
01295       dec6.setForcedRenormalization(false);
01296       double *res=parafield->getField()->getArray()->getPointer();
01297       const double *toSet=targetResults2[rank-nproc_source];
01298       res[0]=toSet[0];
01299       res[1]=toSet[1];
01300       dec6.sendData();
01301     }
01302   //test 7 - Integral with global constraint reversed
01303   ParaMEDMEM::InterpKernelDEC dec7(*source_group,*target_group);
01304   parafield->getField()->setNature(IntegralGlobConstraint);
01305   if (source_group->containsMyRank())
01306     { 
01307       dec7.setMethod("P0");
01308       dec7.attachLocalField(parafield);
01309       dec7.synchronize();
01310       dec7.setForcedRenormalization(false);
01311       dec7.recvData();
01312       const double *res=parafield->getField()->getArray()->getConstPointer();
01313       CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
01314       const double expected[]={36.4592592592593,44.5407407407407};
01315       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
01316     }
01317   else
01318     {
01319       dec7.setMethod("P0");
01320       dec7.attachLocalField(parafield);
01321       dec7.synchronize();
01322       dec7.setForcedRenormalization(false);
01323       double *res=parafield->getField()->getArray()->getPointer();
01324       const double *toSet=targetResults3[rank-nproc_source];
01325       res[0]=toSet[0];
01326       res[1]=toSet[1];
01327       dec7.sendData();
01328     }
01329   //test 8 - Integral with RevIntegral reversed
01330   ParaMEDMEM::InterpKernelDEC dec8(*source_group,*target_group);
01331   parafield->getField()->setNature(RevIntegral);
01332   if (source_group->containsMyRank())
01333     { 
01334       dec8.setMethod("P0");
01335       dec8.attachLocalField(parafield);
01336       dec8.synchronize();
01337       dec8.setForcedRenormalization(false);
01338       dec8.recvData();
01339       const double *res=parafield->getField()->getArray()->getConstPointer();
01340       CPPUNIT_ASSERT_EQUAL(1,parafield->getField()->getNumberOfTuples());
01341       const double expected[]={0.81314102564102553,1.3428994082840233};
01342       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[rank],res[0],1e-13);
01343     }
01344   else
01345     {
01346       dec8.setMethod("P0");
01347       dec8.attachLocalField(parafield);
01348       dec8.synchronize();
01349       dec8.setForcedRenormalization(false);
01350       double *res=parafield->getField()->getArray()->getPointer();
01351       const double *toSet=targetResults4[rank-nproc_source];
01352       res[0]=toSet[0];
01353       res[1]=toSet[1];
01354       dec8.sendData();
01355     }
01356   //
01357   delete parafield;
01358   mesh->decrRef();
01359   delete paramesh;
01360   delete self_group;
01361   delete target_group;
01362   delete source_group;
01363   //
01364   MPI_Barrier(MPI_COMM_WORLD);
01365 }
01366 
01367 void ParaMEDMEMTest::testInterpKernelDECNonOverlapp_2D_P0P1P1P0()
01368 {
01369   int size;
01370   int rank;
01371   MPI_Comm_size(MPI_COMM_WORLD,&size);
01372   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01373   //
01374   if(size!=5)
01375     return ;
01376   int nproc_source = 2;
01377   set<int> self_procs;
01378   set<int> procs_source;
01379   set<int> procs_target;
01380   
01381   for (int i=0; i<nproc_source; i++)
01382     procs_source.insert(i);
01383   for (int i=nproc_source; i<size; i++)
01384     procs_target.insert(i);
01385   self_procs.insert(rank);
01386   //
01387   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
01388   ParaMEDMEM::ParaMESH *paramesh=0;
01389   ParaMEDMEM::ParaFIELD *parafieldP0=0,*parafieldP1=0;
01390   //
01391   ParaMEDMEM::CommInterface interface;
01392   //
01393   ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
01394   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
01395   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
01396   //
01397   MPI_Barrier(MPI_COMM_WORLD);
01398   if(source_group->containsMyRank())
01399     {
01400       if(rank==0)
01401         {
01402           double coords[6]={-0.3,-0.3, 0.7,0.7, 0.7,-0.3};
01403           int conn[3]={0,1,2};
01404           //int globalNode[3]={1,2,0};
01405           mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
01406           mesh->allocateCells(1);
01407           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
01408           mesh->finishInsertingCells();
01409           DataArrayDouble *myCoords=DataArrayDouble::New();
01410           myCoords->alloc(3,2);
01411           std::copy(coords,coords+6,myCoords->getPointer());
01412           mesh->setCoords(myCoords);
01413           myCoords->decrRef();
01414         }
01415       if(rank==1)
01416         {
01417           double coords[6]={-0.3,-0.3, -0.3,0.7, 0.7,0.7};
01418           int conn[3]={0,1,2};
01419           //int globalNode[3]={1,3,2};
01420           mesh=MEDCouplingUMesh::New("Source mesh Proc1",2);
01421           mesh->allocateCells(1);
01422           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
01423           mesh->finishInsertingCells();
01424           DataArrayDouble *myCoords=DataArrayDouble::New();
01425           myCoords->alloc(3,2);
01426           std::copy(coords,coords+6,myCoords->getPointer());
01427           mesh->setCoords(myCoords);
01428           myCoords->decrRef();
01429         }
01430       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
01431       ParaMEDMEM::ComponentTopology comptopo;
01432       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01433       parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
01434       double *valueP0=parafieldP0->getField()->getArray()->getPointer();
01435       double *valueP1=parafieldP1->getField()->getArray()->getPointer();
01436       parafieldP0->getField()->setNature(ConservativeVolumic);
01437       parafieldP1->getField()->setNature(ConservativeVolumic);
01438       if(rank==0)
01439         {
01440           valueP0[0]=31.;
01441           valueP1[0]=34.; valueP1[1]=77.; valueP1[2]=53.;
01442         }
01443       if(rank==1)
01444         {
01445           valueP0[0]=47.;
01446           valueP1[0]=34.; valueP1[1]=57.; valueP1[2]=77.;
01447         }
01448     }
01449   else
01450     {
01451       const char targetMeshName[]="target mesh";
01452       if(rank==2)
01453         {
01454           double coords[10]={-0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2 };
01455           int conn[7]={0,3,4,1, 1,4,2};
01456           //int globalNode[5]={4,3,0,2,1};
01457           mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
01458           mesh->allocateCells(2);
01459           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
01460           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
01461           mesh->finishInsertingCells();
01462           DataArrayDouble *myCoords=DataArrayDouble::New();
01463           myCoords->alloc(5,2);
01464           std::copy(coords,coords+10,myCoords->getPointer());
01465           mesh->setCoords(myCoords);
01466           myCoords->decrRef();
01467           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
01468           DataArrayInt *da=DataArrayInt::New();
01469           const int globalNumberingP2[5]={0,1,2,3,4};
01470           da->useArray(globalNumberingP2,false,CPP_DEALLOC,5,1);
01471           paramesh->setNodeGlobal(da);
01472           da->decrRef();
01473         }
01474       if(rank==3)
01475         {
01476           double coords[6]={0.2,0.2, 0.7,-0.3, 0.7,0.2};
01477           int conn[3]={0,2,1};
01478           //int globalNode[3]={1,0,5};
01479           mesh=MEDCouplingUMesh::New("Target mesh Proc3",2);
01480           mesh->allocateCells(1);
01481           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn);
01482           mesh->finishInsertingCells();
01483           DataArrayDouble *myCoords=DataArrayDouble::New();
01484           myCoords->alloc(3,2);
01485           std::copy(coords,coords+6,myCoords->getPointer());
01486           mesh->setCoords(myCoords);
01487           myCoords->decrRef();
01488           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
01489           DataArrayInt *da=DataArrayInt::New();
01490           const int globalNumberingP3[3]={4,2,5};
01491           da->useArray(globalNumberingP3,false,CPP_DEALLOC,3,1);
01492           paramesh->setNodeGlobal(da);
01493           da->decrRef();
01494         }
01495       if(rank==4)
01496         {
01497           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};
01498           int conn[8]={0,1,2,3, 3,2,4,5};
01499           //int globalNode[6]={2,6,7,1,8,5};
01500           mesh=MEDCouplingUMesh::New("Target mesh Proc4",2);
01501           mesh->allocateCells(2);
01502           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
01503           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
01504           mesh->finishInsertingCells();
01505           DataArrayDouble *myCoords=DataArrayDouble::New();
01506           myCoords->alloc(6,2);
01507           std::copy(coords,coords+12,myCoords->getPointer());
01508           mesh->setCoords(myCoords);
01509           myCoords->decrRef();
01510           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
01511           DataArrayInt *da=DataArrayInt::New();
01512           const int globalNumberingP4[6]={3,6,7,4,8,5};
01513           da->useArray(globalNumberingP4,false,CPP_DEALLOC,6,1);
01514           paramesh->setNodeGlobal(da);
01515           da->decrRef();
01516         }
01517       ParaMEDMEM::ComponentTopology comptopo;
01518       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01519       parafieldP1 = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
01520       parafieldP0->getField()->setNature(ConservativeVolumic);
01521       parafieldP1->getField()->setNature(ConservativeVolumic);
01522     }
01523   // test 1 - P0 P1
01524   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
01525   if (source_group->containsMyRank())
01526     { 
01527       dec.setMethod("P0");
01528       dec.attachLocalField(parafieldP0);
01529       dec.synchronize();
01530       dec.setForcedRenormalization(false);
01531       dec.sendData();
01532       dec.recvData();
01533       const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
01534       if(rank==0)
01535         {
01536           CPPUNIT_ASSERT_DOUBLES_EQUAL(34.42857143,valueP0[0],1e-7);
01537         }
01538       if(rank==1)
01539         {
01540           CPPUNIT_ASSERT_DOUBLES_EQUAL(44.,valueP0[0],1e-7);
01541         }
01542     }
01543   else
01544     {
01545       dec.setMethod("P1");
01546       dec.attachLocalField(parafieldP1);
01547       dec.synchronize();
01548       dec.setForcedRenormalization(false);
01549       dec.recvData();
01550       const double *res=parafieldP1->getField()->getArray()->getConstPointer();
01551       if(rank==2)
01552         {
01553           const double expectP2[5]={39.0, 31.0, 31.0, 47.0, 39.0};
01554           CPPUNIT_ASSERT_EQUAL(5,parafieldP1->getField()->getNumberOfTuples());
01555           CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
01556           for(int kk=0;kk<5;kk++)
01557             CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP2[kk],res[kk],1e-12);
01558         }
01559       if(rank==3)
01560         {
01561           const double expectP3[3]={39.0, 31.0, 31.0};
01562           CPPUNIT_ASSERT_EQUAL(3,parafieldP1->getField()->getNumberOfTuples());
01563           CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
01564           for(int kk=0;kk<3;kk++)
01565             CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP3[kk],res[kk],1e-12);
01566         }
01567       if(rank==4)
01568         {
01569           const double expectP4[6]={47.0, 47.0, 47.0, 39.0, 39.0, 31.0};
01570           CPPUNIT_ASSERT_EQUAL(6,parafieldP1->getField()->getNumberOfTuples());
01571           CPPUNIT_ASSERT_EQUAL(1,parafieldP1->getField()->getNumberOfComponents());
01572           for(int kk=0;kk<6;kk++)
01573             CPPUNIT_ASSERT_DOUBLES_EQUAL(expectP4[kk],res[kk],1e-12);
01574         }
01575       dec.sendData();
01576     }
01577   //
01578   delete parafieldP0;
01579   delete parafieldP1;
01580   mesh->decrRef();
01581   delete paramesh;
01582   delete self_group;
01583   delete target_group;
01584   delete source_group;
01585   //
01586   MPI_Barrier(MPI_COMM_WORLD);
01587 }
01588 
01589 void ParaMEDMEMTest::testInterpKernelDEC2DM1D_P0P0()
01590 {
01591   int size;
01592   int rank;
01593   MPI_Comm_size(MPI_COMM_WORLD,&size);
01594   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01595   //
01596   if(size!=3)
01597     return ;
01598   int nproc_source=2;
01599   set<int> procs_source;
01600   set<int> procs_target;
01601   //
01602   for (int i=0; i<nproc_source; i++)
01603     procs_source.insert(i);
01604   for (int i=nproc_source;i<size; i++)
01605     procs_target.insert(i);
01606   //
01607   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
01608   ParaMEDMEM::ParaMESH *paramesh=0;
01609   ParaMEDMEM::ParaFIELD *parafield=0;
01610   //
01611   ParaMEDMEM::CommInterface interface;
01612   //
01613   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
01614   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
01615   //
01616   MPI_Barrier(MPI_COMM_WORLD);
01617   if(source_group->containsMyRank())
01618     {
01619       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 };
01620       mesh=MEDCouplingUMesh::New();
01621       mesh->setMeshDimension(2);
01622       DataArrayDouble *myCoords=DataArrayDouble::New();
01623       myCoords->alloc(9,2);
01624       std::copy(targetCoords,targetCoords+18,myCoords->getPointer());
01625       mesh->setCoords(myCoords);
01626       myCoords->decrRef();
01627       if(rank==0)
01628         {
01629           int targetConn[7]={0,3,4,1, 1,4,2};
01630           mesh->allocateCells(2);
01631           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
01632           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+4);
01633           mesh->finishInsertingCells();
01634         }
01635       else
01636         { 
01637           int targetConn[11]={4,5,2, 6,7,4,3, 7,8,5,4};
01638           mesh->allocateCells(3);
01639           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
01640           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+3);
01641           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn+7);
01642           mesh->finishInsertingCells();
01643         }
01644       ParaMEDMEM::ComponentTopology comptopo;
01645       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
01646       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01647       parafield->getField()->setNature(ConservativeVolumic);
01648       double *vals=parafield->getField()->getArray()->getPointer();
01649       if(rank==0)
01650         { vals[0]=7.; vals[1]=8.; }
01651       else
01652         { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
01653     }
01654   else
01655     {
01656       mesh=MEDCouplingUMesh::New("an example of -1 D mesh",-1);
01657       ParaMEDMEM::ComponentTopology comptopo;
01658       paramesh=new ParaMESH(mesh,*target_group,"target mesh");
01659       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01660       parafield->getField()->setNature(ConservativeVolumic);
01661     }
01662   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
01663   if(source_group->containsMyRank())
01664     {
01665       dec.setMethod("P0");
01666       dec.attachLocalField(parafield);
01667       dec.synchronize();
01668       dec.setForcedRenormalization(false);
01669       dec.sendData();
01670       dec.recvData();
01671       const double *res=parafield->getField()->getArray()->getConstPointer();
01672       if(rank==0)
01673         {
01674           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01675           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
01676         }
01677       else
01678         {
01679           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01680           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
01681           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
01682         }
01683     }
01684   else
01685     {
01686       dec.setMethod("P0");
01687       dec.attachLocalField(parafield);
01688       dec.synchronize();
01689       dec.setForcedRenormalization(false);
01690       dec.recvData();
01691       const double *res=parafield->getField()->getArray()->getConstPointer();
01692       CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01693       dec.sendData();
01694     }
01695   ParaMEDMEM::InterpKernelDEC dec2(*source_group,*target_group);
01696   dec2.setMethod("P0");
01697   parafield->getField()->setNature(IntegralGlobConstraint);
01698   if(source_group->containsMyRank())
01699     {
01700       double *vals=parafield->getField()->getArray()->getPointer();
01701       if(rank==0)
01702         { vals[0]=7.; vals[1]=8.; }
01703       else
01704         { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
01705       dec2.attachLocalField(parafield);
01706       dec2.synchronize();
01707       dec2.sendData();
01708       dec2.recvData();
01709       const double *res=parafield->getField()->getArray()->getConstPointer();
01710       if(rank==0)
01711         {
01712           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
01713           CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
01714         }
01715       else
01716         {
01717           CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
01718           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
01719           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
01720         }
01721     }
01722   else
01723     {
01724       dec2.attachLocalField(parafield);
01725       dec2.synchronize();
01726       dec2.recvData();
01727       const double *res=parafield->getField()->getArray()->getConstPointer();
01728       CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
01729       dec2.sendData();
01730     }
01731   //
01732   ParaMEDMEM::InterpKernelDEC dec3(*source_group,*target_group);
01733   dec3.setMethod("P0");
01734   parafield->getField()->setNature(Integral);
01735   if(source_group->containsMyRank())
01736     {
01737       double *vals=parafield->getField()->getArray()->getPointer();
01738       if(rank==0)
01739         { vals[0]=7.; vals[1]=8.; }
01740       else
01741         { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
01742       dec3.attachLocalField(parafield);
01743       dec3.synchronize();
01744       dec3.sendData();
01745       dec3.recvData();
01746       const double *res=parafield->getField()->getArray()->getConstPointer();
01747       if(rank==0)
01748         {
01749           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[0],1e-12);
01750           CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[1],1e-12);
01751         }
01752       else
01753         {
01754           CPPUNIT_ASSERT_DOUBLES_EQUAL(5.625,res[0],1e-12);
01755           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[1],1e-12);
01756           CPPUNIT_ASSERT_DOUBLES_EQUAL(11.25,res[2],1e-12);
01757         }
01758     }
01759   else
01760     {
01761       dec3.attachLocalField(parafield);
01762       dec3.synchronize();
01763       dec3.recvData();
01764       const double *res=parafield->getField()->getArray()->getConstPointer();
01765       CPPUNIT_ASSERT_DOUBLES_EQUAL(45.,res[0],1e-12);
01766       dec3.sendData();
01767     }
01768   //
01769   ParaMEDMEM::InterpKernelDEC dec4(*source_group,*target_group);
01770   dec4.setMethod("P0");
01771   parafield->getField()->setNature(RevIntegral);
01772   if(source_group->containsMyRank())
01773     {
01774       double *vals=parafield->getField()->getArray()->getPointer();
01775       if(rank==0)
01776         { vals[0]=7.; vals[1]=8.; }
01777       else
01778         { vals[0]=9.; vals[1]=10.; vals[2]=11.; }
01779       dec4.attachLocalField(parafield);
01780       dec4.synchronize();
01781       dec4.sendData();
01782       dec4.recvData();
01783       const double *res=parafield->getField()->getArray()->getConstPointer();
01784        if(rank==0)
01785         {
01786           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01787           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
01788         }
01789       else
01790         {
01791           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01792           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[1],1e-12);
01793           CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[2],1e-12);
01794         }
01795     }
01796   else
01797     {
01798       dec4.attachLocalField(parafield);
01799       dec4.synchronize();
01800       dec4.recvData();
01801       const double *res=parafield->getField()->getArray()->getConstPointer();
01802       CPPUNIT_ASSERT_DOUBLES_EQUAL(9.125,res[0],1e-12);
01803       dec4.sendData();
01804     }
01805   delete parafield;
01806   delete paramesh;
01807   mesh->decrRef();
01808   delete target_group;
01809   delete source_group;
01810   //
01811   MPI_Barrier(MPI_COMM_WORLD);
01812 }
01813 
01814 void ParaMEDMEMTest::testInterpKernelDECPartialProcs()
01815 {
01816   int size;
01817   int rank;
01818   MPI_Comm_size(MPI_COMM_WORLD,&size);
01819   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01820   //
01821   if(size!=3)
01822     return ;
01823   set<int> procs_source;
01824   set<int> procs_target;
01825   //
01826   procs_source.insert(0);
01827   procs_target.insert(1);
01828   //
01829   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
01830   ParaMEDMEM::ParaMESH *paramesh=0;
01831   ParaMEDMEM::ParaFIELD *parafield=0;
01832   //
01833   ParaMEDMEM::CommInterface interface;
01834   //
01835   MPI_Barrier(MPI_COMM_WORLD);
01836   double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
01837   CommInterface comm;
01838   int grpIds[2]={0,1};
01839   MPI_Group grp,group_world;
01840   comm.commGroup(MPI_COMM_WORLD,&group_world);
01841   comm.groupIncl(group_world,2,grpIds,&grp);
01842   MPI_Comm partialComm;
01843   comm.commCreate(MPI_COMM_WORLD,grp,&partialComm);
01844   //
01845   ProcessorGroup* target_group=0;
01846   ProcessorGroup* source_group=0;
01847   //
01848   ParaMEDMEM::InterpKernelDEC *dec=0;
01849   if(rank==0 || rank==1)
01850     {
01851       target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target,partialComm);
01852       source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source,partialComm);
01853       if(source_group->containsMyRank())
01854         {    
01855           mesh=MEDCouplingUMesh::New();
01856           mesh->setMeshDimension(2);
01857           DataArrayDouble *myCoords=DataArrayDouble::New();
01858           myCoords->alloc(4,2);
01859           std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
01860           mesh->setCoords(myCoords);
01861           myCoords->decrRef();
01862           int targetConn[4]={0,2,3,1};
01863           mesh->allocateCells(1);
01864           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
01865           mesh->finishInsertingCells();
01866           ParaMEDMEM::ComponentTopology comptopo;
01867           paramesh=new ParaMESH(mesh,*source_group,"source mesh");
01868           parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01869           parafield->getField()->setNature(ConservativeVolumic);
01870           double *vals=parafield->getField()->getArray()->getPointer();
01871           vals[0]=7.;
01872           dec=new ParaMEDMEM::InterpKernelDEC(*source_group,*target_group);
01873           dec->attachLocalField(parafield);
01874           dec->synchronize();
01875           dec->sendData();
01876           dec->recvData();
01877         }
01878       else
01879         {
01880           mesh=MEDCouplingUMesh::New();
01881           mesh->setMeshDimension(2);
01882           DataArrayDouble *myCoords=DataArrayDouble::New();
01883           myCoords->alloc(4,2);
01884           std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
01885           mesh->setCoords(myCoords);
01886           myCoords->decrRef();
01887           int targetConn[6]={0,2,1,2,3,1};
01888           mesh->allocateCells(2);
01889           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
01890           mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
01891           mesh->finishInsertingCells();
01892           ParaMEDMEM::ComponentTopology comptopo;
01893           paramesh=new ParaMESH(mesh,*target_group,"target mesh");
01894           parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01895           parafield->getField()->setNature(ConservativeVolumic);
01896           dec=new ParaMEDMEM::InterpKernelDEC(*source_group,*target_group);
01897           dec->attachLocalField(parafield);
01898           dec->synchronize();
01899           dec->recvData();
01900           dec->sendData();
01901         }
01902     }
01903   delete parafield;
01904   delete paramesh;
01905   if(mesh)
01906     mesh->decrRef();
01907   delete target_group;
01908   delete source_group;
01909   delete dec;
01910   MPI_Barrier(MPI_COMM_WORLD);
01911 }
01912 
01917 void ParaMEDMEMTest::testInterpKernelDEC3DSurfEmptyBBox()
01918 {
01919   int size;
01920   int rank;
01921   MPI_Comm_size(MPI_COMM_WORLD,&size);
01922   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01923   //
01924   if(size!=3)
01925     return ;
01926   int nproc_source = 1;
01927   set<int> self_procs;
01928   set<int> procs_source;
01929   set<int> procs_target;
01930   
01931   for (int i=0; i<nproc_source; i++)
01932     procs_source.insert(i);
01933   for (int i=nproc_source; i<size; i++)
01934     procs_target.insert(i);
01935   self_procs.insert(rank);
01936   //
01937   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
01938   ParaMEDMEM::ParaMESH *paramesh=0;
01939   ParaMEDMEM::ParaFIELD *parafieldP0=0;
01940   //
01941   ParaMEDMEM::CommInterface interface;
01942   //
01943   ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
01944   ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
01945   ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
01946   //
01947   MPI_Barrier(MPI_COMM_WORLD);
01948   if(source_group->containsMyRank())
01949     {
01950       double coords[15]={1.,0.,0., 2.,0.,0., 2.,2.,0., 0.,2.,0., 0.5,0.5,1.};
01951       int conn[7]={0,1,2,3,0,3,4};
01952       mesh=MEDCouplingUMesh::New("Source mesh Proc0",2);
01953       mesh->allocateCells(2);
01954       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
01955       mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+4);
01956       mesh->finishInsertingCells();
01957       DataArrayDouble *myCoords=DataArrayDouble::New();
01958       myCoords->alloc(5,3);
01959       std::copy(coords,coords+15,myCoords->getPointer());
01960       mesh->setCoords(myCoords);
01961       myCoords->decrRef();
01962       //
01963       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
01964       ParaMEDMEM::ComponentTopology comptopo;
01965       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
01966       double *valueP0=parafieldP0->getField()->getArray()->getPointer();
01967       parafieldP0->getField()->setNature(ConservativeVolumic);
01968       valueP0[0]=7.; valueP0[1]=8.;
01969     }
01970   else
01971     {
01972       const char targetMeshName[]="target mesh";
01973       if(rank==1)
01974         {
01975           double coords[12]={0.25,0.25,0.5, 0.,0.25,0.5, 0.,0.,0.5, 0.25,0.,0.5};
01976           int conn[4]={0,1,2,3};
01977           mesh=MEDCouplingUMesh::New("Target mesh Proc1",2);
01978           mesh->allocateCells(1);
01979           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
01980           mesh->finishInsertingCells();
01981           DataArrayDouble *myCoords=DataArrayDouble::New();
01982           myCoords->alloc(4,3);
01983           std::copy(coords,coords+12,myCoords->getPointer());
01984           mesh->setCoords(myCoords);
01985           myCoords->decrRef();
01986           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
01987         }
01988       if(rank==2)
01989         {
01990           double coords[12]={0.,0.25,0.5, 0.,0.,0.5, -1.,0.,0.5, -1.,0.25,0.5};
01991           int conn[4]={0,1,2,3};
01992           mesh=MEDCouplingUMesh::New("Target mesh Proc2",2);
01993           mesh->allocateCells(1);
01994           mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
01995           mesh->finishInsertingCells();
01996           DataArrayDouble *myCoords=DataArrayDouble::New();
01997           myCoords->alloc(4,3);
01998           std::copy(coords,coords+12,myCoords->getPointer());
01999           mesh->setCoords(myCoords);
02000           myCoords->decrRef();
02001           paramesh=new ParaMESH(mesh,*target_group,targetMeshName);
02002         }
02003       ParaMEDMEM::ComponentTopology comptopo;
02004       parafieldP0 = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
02005       parafieldP0->getField()->setNature(ConservativeVolumic);
02006     }
02007   // test 1
02008   ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
02009   if (source_group->containsMyRank())
02010     { 
02011       dec.setMethod("P0");
02012       dec.attachLocalField(parafieldP0);
02013       dec.synchronize();
02014       // dec.setForcedRenormalization(false);
02015       // dec.sendData();
02016       // dec.recvData();
02017       // const double *valueP0=parafieldP0->getField()->getArray()->getPointer();
02018       // if(rank==0)
02019       //   {
02020       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,valueP0[0],1e-7);
02021       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[1],1e-7);
02022       //   }
02023       // if(rank==1)
02024       //   {
02025       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(8.64054054054054,valueP0[0],1e-7);
02026       //   }
02027       // if(rank==2)
02028       //   {
02029       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540544,valueP0[0],1e-7);
02030       //   }
02031     }
02032   else
02033     {
02034       dec.setMethod("P0");
02035       dec.attachLocalField(parafieldP0);
02036       dec.synchronize();
02037       // dec.setForcedRenormalization(false);
02038       // dec.recvData();
02039       // const double *res=parafieldP0->getField()->getArray()->getConstPointer();
02040       // if(rank==3)
02041       //   {
02042       //     CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
02043       //     CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
02044       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(7.4,res[0],1e-12);
02045       //   }
02046       // if(rank==4)
02047       //   {
02048       //     CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfTuples());
02049       //     CPPUNIT_ASSERT_EQUAL(1,parafieldP0->getField()->getNumberOfComponents());
02050       //     CPPUNIT_ASSERT_DOUBLES_EQUAL(9.0540540540540526,res[0],1e-12);
02051       //   }
02052       // dec.sendData();
02053     }
02054   //
02055   delete parafieldP0;
02056   mesh->decrRef();
02057   delete paramesh;
02058   delete self_group;
02059   delete target_group;
02060   delete source_group;
02061   //
02062   MPI_Barrier(MPI_COMM_WORLD);
02063 }
02064 
02070 void ParaMEDMEMTest::testAsynchronousInterpKernelDEC_2D(double dtA, double tmaxA, 
02071                                                         double dtB, double tmaxB, bool WithPointToPoint, bool Asynchronous,
02072                                                         bool WithInterp, const char *srcMeth, const char *targetMeth)
02073 {
02074   std::string srcM(srcMeth);
02075   std::string targetM(targetMeth);
02076   int size;
02077   int rank;
02078   MPI_Comm_size(MPI_COMM_WORLD,&size);
02079   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
02080  
02081   //the test is meant to run on five processors
02082   if (size !=5) return ;
02083    
02084   int nproc_source = 3;
02085   set<int> self_procs;
02086   set<int> procs_source;
02087   set<int> procs_target;
02088   
02089   for (int i=0; i<nproc_source; i++)
02090     procs_source.insert(i);
02091   for (int i=nproc_source; i<size; i++)
02092     procs_target.insert(i);
02093   self_procs.insert(rank);
02094   
02095   ParaMEDMEM::CommInterface interface;
02096     
02097   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
02098   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
02099   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
02100     
02101   //loading the geometry for the source group
02102 
02103   ParaMEDMEM::InterpKernelDEC dec (*source_group,*target_group);
02104   
02105   ParaMEDMEM::MEDCouplingUMesh* mesh;
02106   ParaMEDMEM::ParaMESH* paramesh;
02107   ParaMEDMEM::ParaFIELD* parafield;
02108   
02109   ICoCo::Field* icocofield ;
02110 
02111   string tmp_dir                    = getenv("TMP");
02112   if (tmp_dir == "")
02113     tmp_dir = "/tmp";
02114   string filename_xml1              = getResourceFile("square1_split");
02115   string filename_xml2              = getResourceFile("square2_split"); 
02116   //string filename_seq_wr            = makeTmpFile("");
02117   //string filename_seq_med           = makeTmpFile("myWrField_seq_pointe221.med");
02118   
02119   // To remove tmp files from disk
02120   ParaMEDMEMTest_TmpFilesRemover aRemover;
02121   
02122   MPI_Barrier(MPI_COMM_WORLD);
02123 
02124   if (source_group->containsMyRank())
02125     {
02126       string master = filename_xml1;
02127       
02128       ostringstream strstream;
02129       strstream <<master<<rank+1<<".med";
02130       ostringstream meshname ;
02131       meshname<< "Mesh_2_"<< rank+1;
02132       
02133       mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
02134 
02135       paramesh=new ParaMESH (mesh,*source_group,"source mesh");
02136     
02137       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT( support,*source_group);
02138       ParaMEDMEM::ComponentTopology comptopo;
02139       if(srcM=="P0")
02140         {
02141           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
02142           parafield->getField()->setNature(ConservativeVolumic);//InvertIntegral);//ConservativeVolumic);
02143         }
02144       else
02145         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
02146 
02147       int nb_local;
02148       if(srcM=="P0")
02149         nb_local=mesh->getNumberOfCells();
02150       else
02151         nb_local=mesh->getNumberOfNodes();
02152       //      double * value= new double[nb_local];
02153       double *value=parafield->getField()->getArray()->getPointer();
02154       for(int ielem=0; ielem<nb_local;ielem++)
02155         value[ielem]=0.0;
02156     
02157       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
02158       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
02159      
02160       dec.attachLocalField(icocofield);
02161 
02162 
02163     }
02164   
02165   //loading the geometry for the target group
02166   if (target_group->containsMyRank())
02167     {
02168       string master= filename_xml2;
02169       ostringstream strstream;
02170       strstream << master<<(rank-nproc_source+1)<<".med";
02171       ostringstream meshname ;
02172       meshname<< "Mesh_3_"<<rank-nproc_source+1;
02173       
02174       mesh = MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
02175 
02176       paramesh=new ParaMESH (mesh,*target_group,"target mesh");
02177       //      ParaMEDMEM::ParaSUPPORT* parasupport=new UnstructuredParaSUPPORT(support,*target_group);
02178       ParaMEDMEM::ComponentTopology comptopo;
02179       if(targetM=="P0")
02180         {
02181           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
02182           parafield->getField()->setNature(ConservativeVolumic);//InvertIntegral);//ConservativeVolumic);
02183         }
02184       else
02185         parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh, comptopo);
02186       
02187       int nb_local;
02188       if(targetM=="P0")
02189         nb_local=mesh->getNumberOfCells();
02190       else
02191         nb_local=mesh->getNumberOfNodes();
02192                         
02193       double *value=parafield->getField()->getArray()->getPointer();
02194       for(int ielem=0; ielem<nb_local;ielem++)
02195         value[ielem]=0.0;
02196       //      ICoCo::Field* icocofield=new ICoCo::MEDField(paramesh,parafield);
02197       icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
02198       
02199       dec.attachLocalField(icocofield);
02200     }
02201     
02202   
02203   //attaching a DEC to the source group 
02204   
02205   if (source_group->containsMyRank())
02206     { 
02207       cout<<"DEC usage"<<endl;
02208       dec.setAsynchronous(Asynchronous);
02209       if ( WithInterp ) {
02210         dec.setTimeInterpolationMethod(LinearTimeInterp);
02211       }
02212       if ( WithPointToPoint ) {
02213         dec.setAllToAllMethod(PointToPoint);
02214       }
02215       else {
02216         dec.setAllToAllMethod(Native);
02217       }
02218       dec.synchronize();
02219       dec.setForcedRenormalization(false);
02220       for (double time=0; time<tmaxA+1e-10; time+=dtA)
02221         {
02222           cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
02223                << " dtA " << dtA << " tmaxA " << tmaxA << endl ;
02224           if ( time+dtA < tmaxA+1e-7 ) {
02225             dec.sendData( time , dtA );
02226           }
02227           else {
02228             dec.sendData( time , 0 );
02229           }
02230           double* value = parafield->getField()->getArray()->getPointer();
02231           int nb_local=parafield->getField()->getMesh()->getNumberOfCells();
02232           for (int i=0; i<nb_local;i++)
02233             value[i]= time+dtA;
02234 
02235        
02236         }
02237     }
02238   
02239   //attaching a DEC to the target group
02240   if (target_group->containsMyRank())
02241     {
02242       cout<<"DEC usage"<<endl;
02243       dec.setAsynchronous(Asynchronous);
02244       if ( WithInterp ) {
02245         dec.setTimeInterpolationMethod(LinearTimeInterp);
02246       }
02247       if ( WithPointToPoint ) {
02248         dec.setAllToAllMethod(PointToPoint);
02249       }
02250       else {
02251         dec.setAllToAllMethod(Native);
02252       }
02253       dec.synchronize();
02254       dec.setForcedRenormalization(false);
02255       vector<double> times;
02256       for (double time=0; time<tmaxB+1e-10; time+=dtB)
02257         {
02258           cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
02259                << " dtB " << dtB << " tmaxB " << tmaxB << endl ;
02260           dec.recvData( time );
02261           double vi = parafield->getVolumeIntegral(0,true);
02262           cout << "testAsynchronousInterpKernelDEC_2D" << rank << " time " << time
02263                << " VolumeIntegral " << vi
02264                << " time*10000 " << time*10000 << endl ;
02265           
02266           CPPUNIT_ASSERT_DOUBLES_EQUAL(vi,time*10000,0.001);
02267         }
02268       
02269     }
02270   
02271   delete source_group;
02272   delete target_group;
02273   delete self_group;
02274   delete parafield ;
02275   delete paramesh ;
02276   mesh->decrRef() ;
02277   delete icocofield ;
02278 
02279   cout << "testAsynchronousInterpKernelDEC_2D" << rank << " MPI_Barrier " << endl ;
02280  
02281   if (Asynchronous) MPI_Barrier(MPI_COMM_WORLD);
02282   cout << "end of InterpKernelDEC_2D test"<<endl;
02283 }