Salome HOME
Copyright update 2020
[tools/medcoupling.git] / src / ParaMEDMEMTest / ParaMEDMEMTest_Gauthier1.cxx
1 // Copyright (C) 2007-2020  CEA/DEN, EDF R&D
2 //
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
7 //
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
11 // Lesser General Public License for more details.
12 //
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19
20 #include "ParaMEDMEMTest.hxx"
21 #include <cppunit/TestAssert.h>
22
23 #include "CommInterface.hxx"
24 #include "ProcessorGroup.hxx"
25 #include "MPIProcessorGroup.hxx"
26 #include "DEC.hxx"
27 #include "InterpKernelDEC.hxx"
28 #include "MEDCouplingUMesh.hxx"
29 #include "MEDCouplingFieldDouble.hxx"
30 #include "ParaMESH.hxx"
31 #include "ParaFIELD.hxx"
32 #include "ComponentTopology.hxx"
33 #include "BlockTopology.hxx"
34
35 #include <set>
36 #include <time.h>
37 #include <iostream>
38 #include <assert.h>
39 #include <string>
40 #include <math.h>
41
42 using namespace std;
43 using namespace MEDCoupling;
44 using namespace ICoCo;
45
46 void afficheGauthier1(const ParaFIELD& field, const double *vals, int lgth)
47 {
48   const DataArrayDouble *valsOfField(field.getField()->getArray());
49   CPPUNIT_ASSERT_EQUAL(lgth,(int)valsOfField->getNumberOfTuples());
50   for (int ele=0;ele<valsOfField->getNumberOfTuples();ele++)
51     CPPUNIT_ASSERT_DOUBLES_EQUAL(vals[ele],valsOfField->getIJ(ele,0),1e-12);
52 }
53
54 MEDCouplingUMesh *init_quadGauthier1(int is_master)
55 {
56   MCAuto<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_quad",2));
57   MCAuto<DataArrayDouble> coo(DataArrayDouble::New());
58   if(is_master)
59     {
60       const double dataCoo[24]={0,0,0,1,0,0,0,0,1,1,0,1,0,1,0,1,1,0,0,1,1,1,1,1};
61       coo->alloc(8,3);
62       std::copy(dataCoo,dataCoo+24,coo->getPointer());
63       const mcIdType conn[8]={0,1,3,2,4,5,7,6};
64       m->allocateCells(2);
65       m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
66       m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
67     }
68   else
69     {
70       coo->alloc(0,3);
71       m->allocateCells(0);
72     }
73   m->setCoords(coo);
74   return m.retn();
75 }
76
77 MEDCouplingUMesh *init_triangleGauthier1(int is_master)
78 {
79   MCAuto<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_triangle",2));
80   MCAuto<DataArrayDouble> coo(DataArrayDouble::New());
81   if(is_master)
82     {
83       const double dataCoo[24]={0,0,0,1,0,0,0,0,1,1,0,1,0,1,0,1,1,0,0,1,1,1,1,1};
84       coo->alloc(8,3);
85       std::copy(dataCoo,dataCoo+24,coo->getPointer());
86       const mcIdType conn[12]={0,1,2,1,2,3,4,5,7,4,6,7};
87       m->allocateCells(2);
88       for(int i=0;i<4;i++)
89         m->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+3*i);
90     }
91   else
92     {
93       coo->alloc(0,3);
94       m->allocateCells(0);
95     }
96   m->setCoords(coo);
97   return m.retn();
98 }
99
100
101 void ParaMEDMEMTest::testGauthier1()
102 {
103   int num_cas=0;
104   int rank, size;
105   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
106   MPI_Comm_size(MPI_COMM_WORLD,&size);
107   
108   int is_master=0;
109
110   CommInterface comm;
111   set<int> emetteur_ids;
112   set<int> recepteur_ids;
113   emetteur_ids.insert(0);
114   if(size!=4)
115     return;
116   recepteur_ids.insert(1);
117   if (size >2) 
118     recepteur_ids.insert(2);
119   if (size >2) 
120     emetteur_ids.insert(3);
121   if ((rank==0)||(rank==1)) 
122     is_master=1;
123   
124   MPIProcessorGroup recepteur_group(comm,recepteur_ids);
125   MPIProcessorGroup emetteur_group(comm,emetteur_ids);
126
127   string cas;
128   if (recepteur_group.containsMyRank())
129     {
130       cas="recepteur";
131       //freopen("recpeteur.out","w",stdout);
132       //freopen("recepteur.err","w",stderr);
133     }
134   else
135     {
136       cas="emetteur";
137       // freopen("emetteur.out","w",stdout);
138       //freopen("emetteur.err","w",stderr);
139     }
140   double expected[8][4]={
141     {1.,1.,1.,1.},
142     {40., 40., 1., 1.},
143     {1.,1.,1e200,1e200},
144     {40.,1.,1e200,1e200},
145     {1.,1.,1.,1.},
146     {40.,1.,1.,1.},
147     {1.,1.,1e200,1e200},
148     {20.5,1.,1e200,1e200}
149   };
150   int expectedLgth[8]={4,4,2,2,4,4,2,2};
151   
152   for (int send=0;send<2;send++)
153     for (int rec=0;rec<2;rec++)
154       {
155         InterpKernelDEC dec_emetteur(emetteur_group, recepteur_group);
156         MEDCoupling::ParaFIELD *champ_emetteur(0),*champ_recepteur(0);
157         MEDCoupling::ParaMESH *paramesh(0);
158         MCAuto<MEDCouplingUMesh> mesh;
159         dec_emetteur.setOrientation(2);
160         if (send==0)
161           {
162             mesh=init_quadGauthier1(is_master);
163           }
164         else
165           {
166             mesh=init_triangleGauthier1(is_master);
167           }
168         paramesh=new MEDCoupling::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"emetteur mesh");
169         MEDCoupling::ComponentTopology comptopo;
170         champ_emetteur=new MEDCoupling::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
171         champ_emetteur->getField()->setNature(IntensiveMaximum);
172         champ_emetteur->setOwnSupport(true);
173         if (rec==0)
174           {
175             mesh=init_triangleGauthier1(is_master);
176           }
177         else
178           {
179             mesh=init_quadGauthier1(is_master);
180           }
181         paramesh=new MEDCoupling::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"recepteur mesh");
182         champ_recepteur=new MEDCoupling::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
183         champ_recepteur->getField()->setNature(IntensiveMaximum);
184         champ_recepteur->setOwnSupport(true);
185         if (cas=="emetteur") 
186           {
187             champ_emetteur->getField()->getArray()->fillWithValue(1.);
188           }
189   
190   
191         MPI_Barrier(MPI_COMM_WORLD);
192
193         //clock_t clock0= clock ();
194         int compti=0;
195
196         bool init=true; // first time step ??
197         bool stop=false;
198         //boucle sur les pas de quads
199         while (!stop) {
200   
201           compti++;
202           //clock_t clocki= clock ();
203           //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl; 
204           for (int non_unif=0;non_unif<2;non_unif++)
205             {
206               if (cas=="emetteur") 
207                 {
208                   if (non_unif)
209                     if(rank!=3)
210                       champ_emetteur->getField()->getArray()->setIJ(0,0,40);
211                 }
212               //bool ok=false; // Is the time interval successfully solved ?
213     
214               // Loop on the time interval tries
215               if(1) {
216       
217
218                 if (cas=="emetteur")
219                   dec_emetteur.attachLocalField(champ_emetteur);
220                 else
221                   dec_emetteur.attachLocalField(champ_recepteur);
222
223
224                 if(init) dec_emetteur.synchronize();
225                 init=false;
226
227                 if (cas=="emetteur") {
228                   //    affiche(champ_emetteur);
229                   dec_emetteur.sendData();
230                 }
231                 else if (cas=="recepteur")
232                   {
233                     dec_emetteur.recvData();
234                     if (is_master)
235                       afficheGauthier1(*champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
236                   }
237                 else
238                   throw 0;
239                 MPI_Barrier(MPI_COMM_WORLD);
240               }
241               stop=true;
242               num_cas++;
243             }
244         }
245         delete champ_emetteur;
246         delete champ_recepteur;
247       }
248 }
249
250 void ParaMEDMEMTest::testGauthier2()
251 {
252   std::cout << "testGauthier2\n";
253   double valuesExpected1[2]={0.,0.};
254   double valuesExpected2[2]={0.95,0.970625};
255   
256   double valuesExpected30[]={0., 0., 0.05, 0., 0., 0.15, 0., 0., 0.25, 0., 0., 0.35, 0., 0., 0.45, 0., 0., 0.55, 0., 0., 0.65, 0., 0., 0.75, 0., 0., 0.85, 0., 0., 0.95};
257   double valuesExpected31[]={0.,  0.,  0.029375,  0.,  0.,  0.029375,  0.,  0.,  0.1,  0.,  0.,  0.1,  0.,  0.,  0.2,  0.,  0.,  0.2,  0.,  0.,  0.3,  0.,  0.,  0.3,  0.,  0.,  0.4,  0.,  0.,  0.4,  0.,  0.,  0.5,  0.,  0.,  0.5,  0.,  0.,  0.6,  0.,  0.,  0.6,  0.,  0.,  0.7,  0.,  0.,  0.7,  0.,  0.,  0.8,  0.,  0.,  0.8,  0.,  0.,  0.9,  0.,  0.,  0.9,  0.,  0.,  0.970625,  0.,  0.,  0.970625 };
258
259   double *valuesExpected3[2]={valuesExpected30,valuesExpected31};
260
261   int rank, size;
262   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
263   MPI_Comm_size(MPI_COMM_WORLD,&size);
264   if (size <2)
265     return ;
266   CommInterface comm;
267   set<int> Genepi_ids;
268   set<int> entree_chaude_ids;
269   Genepi_ids.insert(0);
270   for (int i=1;i<size;i++)
271     entree_chaude_ids.insert(i);
272   for (int type=0;type<2;type++)
273     {
274       MPIProcessorGroup entree_chaude_group(comm,entree_chaude_ids);
275       MPIProcessorGroup Genepi_group(comm,Genepi_ids);
276
277       MEDCoupling::ParaFIELD *vitesse(0);
278       InterpKernelDEC dec_vit_in_chaude(entree_chaude_group, Genepi_group);
279
280       if ( entree_chaude_group.containsMyRank())
281         {
282           MCAuto<MEDCouplingUMesh> mesh(MEDCouplingUMesh::New("mesh",2));
283           MCAuto<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(63,3);
284           const double cooData[189]={0.,0.,0.,0.5,0.,0.,0.5,0.05,0.,0.,0.1,0.,0.5,0.1,0.,0.5,0.15,0.,0.,0.2,0.,0.5,0.2,0.,0.5,0.25,0.,0.,0.3,0.,0.5,0.3,0.,0.5,0.35,0.,0.,0.4,0.,0.5,0.4,0.,0.5,0.45,0.,0.,0.5,0.,0.5,0.5,0.,0.5,0.55,0.,0.,0.6,0.,0.5,0.6,0.,0.5,0.65,0.,0.,0.7,0.,0.5,0.7,0.,0.5,0.75,0.,0.,0.8,0.,0.5,0.8,0.,0.5,0.85,0.,0.,0.9,0.,0.5,0.9,0.,0.5,0.95,0.,1.,0.,0.,1.,0.1,0.,1.,0.2,0.,1.,0.3,0.,1.,0.4,0.,1.,0.5,0.,1.,0.6,0.,1.,0.7,0.,1.,0.8,0.,1.,0.9,0.,1.,0.05,0.,1.,0.15,0.,1.,0.25,0.,1.,0.35,0.,1.,0.45,0.,1.,0.55,0.,1.,0.65,0.,1.,0.75,0.,1.,0.85,0.,1.,0.95,0.,1.,1.,0.,0.,1.,0.,0.5,1.,0.,0.,0.05,0.,0.,0.15,0.,0.,0.25,0.,0.,0.35,0.,0.,0.45,0.,0.,0.55,0.,0.,0.65,0.,0.,0.75,0.,0.,0.85,0.,0.,0.95,0.};
285           std::copy(cooData,cooData+189,arr->getPointer());
286           mesh->setCoords(arr);
287           mesh->allocateCells(80);
288           const mcIdType conn[240]={0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,2,1,31,5,4,32,8,7,33,11,10,34,14,13,35,17,16,36,20,19,37,23,22,38,26,25,39,29,28,30,40,2,31,41,5,32,42,8,33,43,11,34,44,14,35,45,17,36,46,20,37,47,23,38,48,26,39,49,29,31,2,40,32,5,41,33,8,42,34,11,43,35,14,44,36,17,45,37,20,46,38,23,47,39,26,48,50,29,49,3,2,4,6,5,7,9,8,10,12,11,13,15,14,16,18,17,19,21,20,22,24,23,25,27,26,28,51,29,52,31,4,2,32,7,5,33,10,8,34,13,11,35,16,14,36,19,17,37,22,20,38,25,23,39,28,26,50,52,29,0,2,53,3,5,54,6,8,55,9,11,56,12,14,57,15,17,58,18,20,59,21,23,60,24,26,61,27,29,62,3,53,2,6,54,5,9,55,8,12,56,11,15,57,14,18,58,17,21,59,20,24,60,23,27,61,26,51,62,29};
289           for(int i=0;i<80;i++)
290             mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+3*i);
291           MCAuto<MEDCouplingFieldDouble> f(MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME));
292           const double valsOfField[189]={0.,0.,0.,0.,0.,0.,0.,0.,0.05,0.,0.,0.1,0.,0.,0.1,0.,0.,0.15,0.,0.,0.2,0.,0.,0.2,0.,0.,0.25,0.,0.,0.3,0.,0.,0.3,0.,0.,0.35,0.,0.,0.4,0.,0.,0.4,0.,0.,0.45,0.,0.,0.5,0.,0.,0.5,0.,0.,0.55,0.,0.,0.6,0.,0.,0.6,0.,0.,0.65,0.,0.,0.7,0.,0.,0.7,0.,0.,0.75,0.,0.,0.8,0.,0.,0.8,0.,0.,0.85,0.,0.,0.9,0.,0.,0.9,0.,0.,0.95,0.,0.,0.,0.,0.,0.1,0.,0.,0.2,0.,0.,0.3,0.,0.,0.4,0.,0.,0.5,0.,0.,0.6,0.,0.,0.7,0.,0.,0.8,0.,0.,0.9,0.,0.,0.05,0.,0.,0.15,0.,0.,0.25,0.,0.,0.35,0.,0.,0.45,0.,0.,0.55,0.,0.,0.65,0.,0.,0.75,0.,0.,0.85,0.,0.,0.95,0.,0.,1.,0.,0.,1.,0.,0.,1.,0.,0.,0.05,0.,0.,0.15,0.,0.,0.25,0.,0.,0.35,0.,0.,0.45,0.,0.,0.55,0.,0.,0.65,0.,0.,0.75,0.,0.,0.85,0.,0.,0.95};
293           f->setMesh(mesh); f->setName("VITESSE_P1_OUT");
294           arr=DataArrayDouble::New(); arr->alloc(63,3);
295           std::copy(valsOfField,valsOfField+189,arr->getPointer());
296           f->setArray(arr); f->setNature(IntensiveMaximum);
297           MEDCoupling::ParaMESH *paramesh(new MEDCoupling::ParaMESH(mesh,entree_chaude_group,"emetteur mesh"));
298           vitesse=new MEDCoupling::ParaFIELD(f,paramesh,entree_chaude_group);
299           vitesse->setOwnSupport(true);
300           dec_vit_in_chaude.setMethod("P1");
301         }
302       else
303         {
304           MCAuto<MEDCouplingUMesh> mesh(MEDCouplingUMesh::New("mesh",2));
305           MCAuto<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(22,3);
306           const double cooData[66]={0,0,0,1,0,0,0,0.1,0,1,0.1,0,0,0.2,0,1,0.2,0,0,0.3,0,1,0.3,0,0,0.4,0,1,0.4,0,0,0.5,0,1,0.5,0,0,0.6,0,1,0.6,0,0,0.7,0,1,0.7,0,0,0.8,0,1,0.8,0,0,0.9,0,1,0.9,0,0,1,0,1,1,0};
307           std::copy(cooData,cooData+66,arr->getPointer());
308           mesh->setCoords(arr);
309           mesh->allocateCells(10);
310           const mcIdType conn[40]={0,1,3,2,2,3,5,4,4,5,7,6,6,7,9,8,8,9,11,10,10,11,13,12,12,13,15,14,14,15,17,16,16,17,19,18,18,19,21,20};
311           for(int i=0;i<10;i++)
312             mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4*i);
313           MCAuto<MEDCouplingFieldDouble> f(MEDCouplingFieldDouble::New(type==0?ON_CELLS:ON_NODES,ONE_TIME));
314           f->setMesh(mesh); f->setName("vitesse_in_chaude");
315           arr=DataArrayDouble::New(); arr->alloc(f->getNumberOfTuplesExpected()*3); arr->fillWithZero(); arr->rearrange(3);
316           f->setArray(arr); f->setNature(IntensiveMaximum);
317           MEDCoupling::ParaMESH *paramesh(new MEDCoupling::ParaMESH(mesh,Genepi_group,"recepteur mesh"));
318           vitesse=new MEDCoupling::ParaFIELD(f,paramesh,Genepi_group);
319           vitesse->setOwnSupport(true);
320           dec_vit_in_chaude.setMethod(f->getDiscretization()->getRepr());
321         }
322
323       dec_vit_in_chaude.attachLocalField(vitesse);
324       
325       dec_vit_in_chaude.synchronize();
326   
327   
328       // Envois - receptions
329       if (entree_chaude_group.containsMyRank())
330         {
331           dec_vit_in_chaude.sendData();
332         }
333       else
334         {
335           dec_vit_in_chaude.recvData(); 
336         }
337       if ( !entree_chaude_group.containsMyRank() )
338         {
339           double pmin=1e38, pmax=-1e38;
340           const double *p(vitesse->getField()->getArray()->begin());
341           for(mcIdType i=0;i<vitesse->getField()->getArray()->getNbOfElems();i++,p++)
342             {
343               if (*p<pmin) pmin=*p;
344               if (*p>pmax) pmax=*p;
345             }
346           CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected1[type],pmin,1e-12);
347           CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected2[type],pmax,1e-12);
348       
349           int nbCompo(vitesse->getField()->getNumberOfComponents());
350           p=vitesse->getField()->getArray()->begin();
351           for(int i=0;i<vitesse->getField()->getNumberOfTuples();i++)
352             for(int c=0;c<nbCompo;c++,p++)
353               CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected3[type][i*nbCompo+c],*p,1e-12);
354         }
355       delete vitesse;
356     }
357 }
358
359 void ParaMEDMEMTest::testGauthier3_1()
360 {
361   testGauthier3_GEN(true,4);
362 }
363
364 void ParaMEDMEMTest::testGauthier3_2()
365 {
366   testGauthier3_GEN(false,4);
367 }
368
369 void ParaMEDMEMTest::testGauthier3_3()
370 {
371   testGauthier3_GEN(true,5);
372 }
373
374 void ParaMEDMEMTest::testGauthier3_4()
375 {
376   int size;
377   MPI_Comm_size(MPI_COMM_WORLD,&size);
378
379   if(size!=5)
380     return;
381
382   // Should throw since the two groups (source/target) do not form a partition of
383   // all the procs.
384   CPPUNIT_ASSERT_THROW(testGauthier3_GEN(false,5), INTERP_KERNEL::Exception);
385 }
386
387
388 /*!
389  * Non regression test testing copy constructor of InterpKernelDEC. 
390  */
391 void ParaMEDMEMTest::testGauthier3_GEN(bool withIDs, int nprocs)
392 {
393   int num_cas=0;
394   int rank, size;
395   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
396   MPI_Comm_size(MPI_COMM_WORLD,&size);
397
398   int is_master=0;
399
400   CommInterface comm;
401   set<int> emetteur_ids;
402   set<int> recepteur_ids;
403   emetteur_ids.insert(0);
404   if(size!=nprocs)
405     return;
406   recepteur_ids.insert(1);
407
408   recepteur_ids.insert(size-2);
409
410   emetteur_ids.insert(size-1);
411   if ((rank==0)||(rank==1))
412     is_master=1;
413
414   MPIProcessorGroup recepteur_group(comm,recepteur_ids);
415   MPIProcessorGroup emetteur_group(comm,emetteur_ids);
416
417   string cas;
418   if (recepteur_group.containsMyRank())
419     {
420       cas="recepteur";
421       //freopen("recpeteur.out","w",stdout);
422       //freopen("recepteur.err","w",stderr);
423     }
424   else
425     {
426       if (emetteur_group.containsMyRank())
427         cas="emetteur";
428       else
429         cas="vide";
430       // freopen("emetteur.out","w",stdout);
431       //freopen("emetteur.err","w",stderr);
432     }
433
434   double expected[8][4]={
435     {1.,1.,1.,1.},
436     {40., 40., 1., 1.},
437     {1.,1.,1e200,1e200},
438     {40.,1.,1e200,1e200},
439     {1.,1.,1.,1.},
440     {40.,1.,1.,1.},
441     {1.,1.,1e200,1e200},
442     {20.5,1.,1e200,1e200}
443   };
444   int expectedLgth[8]={4,4,2,2,4,4,2,2};
445
446   for (int send=0;send<2;send++)
447     for (int rec=0;rec<2;rec++)
448       {
449         std::vector<InterpKernelDEC> decu(1);
450         if (withIDs)
451           decu[0] = InterpKernelDEC(emetteur_ids,recepteur_ids);
452         else
453           decu[0] = InterpKernelDEC(emetteur_group,recepteur_group);
454         InterpKernelDEC& dec_emetteur=decu[0];
455         MEDCoupling::ParaFIELD *champ_emetteur(0),*champ_recepteur(0);
456         MEDCoupling::ParaMESH *paramesh(0);
457         MCAuto<MEDCouplingUMesh> mesh;
458         dec_emetteur.setOrientation(2);
459         if (send==0)
460           {
461             mesh=init_quadGauthier1(is_master);
462           }
463         else
464           {
465             mesh=init_triangleGauthier1(is_master);
466           }
467         if (cas!="vide")
468           {
469             paramesh=new MEDCoupling::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"emetteur mesh");
470             MEDCoupling::ComponentTopology comptopo;
471             champ_emetteur=new MEDCoupling::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
472             champ_emetteur->getField()->setNature(IntensiveMaximum);
473             champ_emetteur->setOwnSupport(true);
474             if (rec==0)
475               {
476                 mesh=init_triangleGauthier1(is_master);
477               }
478             else
479               {
480                 mesh=init_quadGauthier1(is_master);
481               }
482             paramesh=new MEDCoupling::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"recepteur mesh");
483             champ_recepteur=new MEDCoupling::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
484             champ_recepteur->getField()->setNature(IntensiveMaximum);
485             champ_recepteur->setOwnSupport(true);
486             if (cas=="emetteur")
487               {
488                 champ_emetteur->getField()->getArray()->fillWithValue(1.);
489               }
490           }
491
492         MPI_Barrier(MPI_COMM_WORLD);
493
494         //clock_t clock0= clock ();
495         int compti=0;
496
497         bool init=true; // first time step ??
498         bool stop=false;
499         //boucle sur les pas de quads
500         while (!stop) {
501
502             compti++;
503             //clock_t clocki= clock ();
504             //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
505             for (int non_unif=0;non_unif<2;non_unif++)
506               {
507                 if (cas=="emetteur")
508                   {
509                     if (non_unif)
510                       if(rank!=3)
511                         champ_emetteur->getField()->getArray()->setIJ(0,0,40);
512                   }
513                 //bool ok=false; // Is the time interval successfully solved ?
514
515                 // Loop on the time interval tries
516                 if(1) {
517
518
519                     if (cas=="emetteur")
520                       dec_emetteur.attachLocalField(champ_emetteur);
521                     else
522                       dec_emetteur.attachLocalField(champ_recepteur);
523
524
525                     if(init) dec_emetteur.synchronize();
526                     init=false;
527
528                     if (cas=="emetteur") {
529                         //    affiche(champ_emetteur);
530                         dec_emetteur.sendData();
531                     }
532                     else if (cas=="recepteur")
533                       {
534                         dec_emetteur.recvData();
535                         if (is_master)
536                           afficheGauthier1(*champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
537                       }
538                     // else
539                     //   throw 0;
540                     MPI_Barrier(MPI_COMM_WORLD);
541                 }
542                 stop=true;
543                 num_cas++;
544               }
545         }
546         delete champ_emetteur;
547         delete champ_recepteur;
548       }
549 }
550
551 /*!
552  * This test is the parallel version of MEDCouplingBasicsTest.test3D1DOnP1P0_1 test.
553  */
554 void ParaMEDMEMTest::testGauthier4()
555 {
556   //
557   const double sourceCoords[19*3]={0.5,0.5,0.1,0.5,0.5,1.2,0.5,0.5,1.6,0.5,0.5,1.8,0.5,0.5,2.43,0.5,0.5,2.55,0.5,0.5,4.1,0.5,0.5,4.4,0.5,0.5,4.9,0.5,0.5,5.1,0.5,0.5,7.6,0.5,0.5,7.7,0.5,0.5,8.2,0.5,0.5,8.4,0.5,0.5,8.6,0.5,0.5,8.8,0.5,0.5,9.2,0.5,0.5,9.6,0.5,0.5,11.5};
558   const mcIdType sourceConn[18*2]={0,1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8,9,9,10,10,11,11,12,12,13,13,14,14,15,15,16,16,17,17,18};
559   const double sourceVals[19]={0.49,2.8899999999999997,7.29,13.69,22.09,32.49,44.89,59.29,75.69,94.09, 114.49,136.89,161.29,187.69,216.09,246.49,278.89,313.29,349.69};
560   const double targetCoords0[20*3]={0.,0.,0.,1.,0.,0.,0.,1.,0.,1.,1.,0.,0.,0.,1.,1.,0.,1.,0.,1.,1.,1.,1.,1.,0.,0.,2.,1.,0.,2.,0.,1.,2.,1.,1.,2.,0.,0.,3.,1.,0.,3.,0.,1.,3.,1.,1.,3.,0.,0.,4.,1.,0.,4.,0.,1.,4.,1.,1.,4.};
561   const mcIdType targetConn0[8*4]={1,0,2,3,5,4,6,7,5,4,6,7,9,8,10,11,9,8,10,11,13,12,14,15,13,12,14,15,17,16,18,19};
562   const double targetCoords1[28*3]={0.,0.,4.,1.,0.,4.,0.,1.,4.,1.,1.,4.,0.,0.,5.,1.,0.,5.,0.,1.,5.,1.,1.,5.,0.,0.,6.,1.,0.,6.,0.,1.,6.,1.,1.,6.,0.,0.,7.,1.,0.,7.,0.,1.,7.,1.,1.,7.,0.,0.,8.,1.,0.,8.,0.,1.,8.,1.,1.,8.,0.,0.,9.,1.,0.,9.,0.,1.,9.,1.,1.,9.,0.,0.,10.,1.,0.,10.,0.,1.,10.,1.,1.,10.};
563   const mcIdType targetConn1[8*6]={1,0,2,3,5,4,6,7,5,4,6,7,9,8,10,11,9,8,10,11,13,12,14,15,13,12,14,15,17,16,18,19,17,16,18,19,21,20,22,23,21,20,22,23,25,24,26,27};
564   //
565   int size;
566   int rank;
567   MPI_Comm_size(MPI_COMM_WORLD,&size);
568   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
569   //
570   if(size!=3)
571     return ;
572   int nproc_source = 1;
573   set<int> self_procs;
574   set<int> procs_source;
575   set<int> procs_target;
576   
577   for (int i=0; i<nproc_source; i++)
578     procs_source.insert(i);
579   for (int i=nproc_source; i<size; i++)
580     procs_target.insert(i);
581   self_procs.insert(rank);
582   //
583   MEDCoupling::MEDCouplingUMesh *mesh=0;
584   MEDCoupling::ParaMESH *paramesh=0;
585   MEDCoupling::ParaFIELD* parafield=0;
586   //
587   MEDCoupling::CommInterface interface;
588   //
589   ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
590   ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
591   ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
592   //
593   MPI_Barrier(MPI_COMM_WORLD);
594   if(source_group->containsMyRank())
595     {
596       std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
597       mesh=MEDCouplingUMesh::New(stream.str().c_str(),1);
598       mesh->allocateCells();
599       for(int i=0;i<18;i++)
600         mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,sourceConn+2*i);
601       mesh->finishInsertingCells();
602       DataArrayDouble *myCoords=DataArrayDouble::New();
603       myCoords->alloc(19,3);
604       std::copy(sourceCoords,sourceCoords+19*3,myCoords->getPointer());
605       mesh->setCoords(myCoords);
606       myCoords->decrRef();
607       paramesh=new ParaMESH(mesh,*source_group,"source mesh");
608       MEDCoupling::ComponentTopology comptopo;
609       parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh,comptopo);
610       double *value=parafield->getField()->getArray()->getPointer();
611       std::copy(sourceVals,sourceVals+19,value);
612     }
613   else
614     {
615       if(rank==1)
616         {
617           std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
618           mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
619           mesh->allocateCells();
620           for(int i=0;i<4;i++)
621             mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn0+8*i);
622           mesh->finishInsertingCells();
623           DataArrayDouble *myCoords=DataArrayDouble::New();
624           myCoords->alloc(20,3);
625           std::copy(targetCoords0,targetCoords0+20*3,myCoords->getPointer());
626           mesh->setCoords(myCoords);
627           myCoords->decrRef();
628           paramesh=new ParaMESH (mesh,*target_group,"target mesh");
629           MEDCoupling::ComponentTopology comptopo;
630           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
631         }
632       else if(rank==2)
633         {
634           std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
635           mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
636           mesh->allocateCells();
637           for(int i=0;i<6;i++)
638             mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn1+8*i);
639           mesh->finishInsertingCells();
640           DataArrayDouble *myCoords=DataArrayDouble::New();
641           myCoords->alloc(28,3);
642           std::copy(targetCoords1,targetCoords1+28*3,myCoords->getPointer());
643           mesh->setCoords(myCoords);
644           myCoords->decrRef();
645           paramesh=new ParaMESH (mesh,*target_group,"target mesh");
646           MEDCoupling::ComponentTopology comptopo;
647           parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
648         }
649     }
650   //test 1 - primaire -> secondaire
651   MEDCoupling::InterpKernelDEC dec(*source_group,*target_group);
652   dec.setIntersectionType(INTERP_KERNEL::PointLocator);
653   parafield->getField()->setNature(IntensiveMaximum);//very important
654   if (source_group->containsMyRank())
655     { 
656       dec.setMethod("P1");
657       dec.attachLocalField(parafield);
658       dec.synchronize();
659       dec.setForcedRenormalization(false);
660       dec.sendData();
661     }
662   else
663     {
664       dec.setMethod("P0");
665       dec.attachLocalField(parafield);
666       dec.synchronize();
667       dec.setForcedRenormalization(false);
668       dec.recvData();
669       const double *res(parafield->getField()->getArray()->getConstPointer());
670       if(rank==1)
671         {
672           const double expected0[4]={0.49,7.956666666666667,27.29,0.};
673           for(int i=0;i<4;i++)
674             CPPUNIT_ASSERT_DOUBLES_EQUAL(expected0[i],res[i],1e-13);
675         }
676       else
677         {
678           const double expected1[6]={59.95666666666667,94.09,0.,125.69,202.89,296.09};
679           for(int i=0;i<6;i++)
680             CPPUNIT_ASSERT_DOUBLES_EQUAL(expected1[i],res[i],1e-13);
681         }
682     }
683   MPI_Barrier(MPI_COMM_WORLD);
684   if (source_group->containsMyRank())
685     {
686       dec.recvData();
687       const double expected2[19]={0.49,7.956666666666667,7.956666666666667,7.956666666666667,27.29,27.29,59.95666666666667,59.95666666666667,59.95666666666667,94.09,125.69,125.69,202.89,202.89,202.89,202.89,296.09,296.09,0.};
688       const double *res(parafield->getField()->getArray()->getConstPointer());
689       for(int i=0;i<19;i++)
690         CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],res[i],1e-13);
691     }
692   else
693     {
694       dec.sendData();
695     }
696   delete parafield;
697   mesh->decrRef();
698   delete paramesh;
699   delete self_group;
700   delete target_group;
701   delete source_group;
702   //
703   MPI_Barrier(MPI_COMM_WORLD);
704 }