1 // Copyright (C) 2007-2014 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "ParaMEDMEMTest.hxx"
21 #include <cppunit/TestAssert.h>
23 #include "CommInterface.hxx"
24 #include "ProcessorGroup.hxx"
25 #include "MPIProcessorGroup.hxx"
27 #include "InterpKernelDEC.hxx"
28 #include "ICoCoTrioField.hxx"
29 #include "MEDCouplingUMesh.hxx"
30 #include "MEDCouplingFieldDouble.hxx"
31 #include "ParaMESH.hxx"
32 #include "ParaFIELD.hxx"
33 #include "ComponentTopology.hxx"
34 #include "BlockTopology.hxx"
44 using namespace ParaMEDMEM;
45 using namespace ICoCo;
47 void afficheGauthier1(const ParaFIELD& field, const double *vals, int lgth)
49 const DataArrayDouble *valsOfField(field.getField()->getArray());
50 CPPUNIT_ASSERT_EQUAL(lgth,valsOfField->getNumberOfTuples());
51 for (int ele=0;ele<valsOfField->getNumberOfTuples();ele++)
52 CPPUNIT_ASSERT_DOUBLES_EQUAL(vals[ele],valsOfField->getIJ(ele,0),1e-12);
55 MEDCouplingUMesh *init_quadGauthier1(int is_master)
57 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_quad",2));
58 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo(DataArrayDouble::New());
61 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};
63 std::copy(dataCoo,dataCoo+24,coo->getPointer());
64 const int conn[8]={0,1,3,2,4,5,7,6};
66 m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
67 m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
78 MEDCouplingUMesh *init_triangleGauthier1(int is_master)
80 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_triangle",2));
81 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo(DataArrayDouble::New());
84 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};
86 std::copy(dataCoo,dataCoo+24,coo->getPointer());
87 const int conn[12]={0,1,2,1,2,3,4,5,7,4,6,7};
90 m->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+3*i);
102 void ParaMEDMEMTest::testGauthier1()
106 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
107 MPI_Comm_size(MPI_COMM_WORLD,&size);
112 set<int> emetteur_ids;
113 set<int> recepteur_ids;
114 emetteur_ids.insert(0);
117 recepteur_ids.insert(1);
119 recepteur_ids.insert(2);
121 emetteur_ids.insert(3);
122 if ((rank==0)||(rank==1))
125 MPIProcessorGroup recepteur_group(comm,recepteur_ids);
126 MPIProcessorGroup emetteur_group(comm,emetteur_ids);
129 if (recepteur_group.containsMyRank())
132 //freopen("recpeteur.out","w",stdout);
133 //freopen("recepteur.err","w",stderr);
138 // freopen("emetteur.out","w",stdout);
139 //freopen("emetteur.err","w",stderr);
141 double expected[8][4]={
145 {40.,1.,1e200,1e200},
149 {20.5,1.,1e200,1e200}
151 int expectedLgth[8]={4,4,2,2,4,4,2,2};
153 for (int send=0;send<2;send++)
154 for (int rec=0;rec<2;rec++)
156 InterpKernelDEC dec_emetteur(emetteur_group, recepteur_group);
157 ParaMEDMEM::ParaFIELD *champ_emetteur(0),*champ_recepteur(0);
158 ParaMEDMEM::ParaMESH *paramesh(0);
159 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh;
160 dec_emetteur.setOrientation(2);
163 mesh=init_quadGauthier1(is_master);
167 mesh=init_triangleGauthier1(is_master);
169 paramesh=new ParaMEDMEM::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"emetteur mesh");
170 ParaMEDMEM::ComponentTopology comptopo;
171 champ_emetteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
172 champ_emetteur->getField()->setNature(ConservativeVolumic);
173 champ_emetteur->setOwnSupport(true);
176 mesh=init_triangleGauthier1(is_master);
180 mesh=init_quadGauthier1(is_master);
182 paramesh=new ParaMEDMEM::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"recepteur mesh");
183 champ_recepteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
184 champ_recepteur->getField()->setNature(ConservativeVolumic);
185 champ_recepteur->setOwnSupport(true);
188 champ_emetteur->getField()->getArray()->fillWithValue(1.);
192 MPI_Barrier(MPI_COMM_WORLD);
194 //clock_t clock0= clock ();
197 bool init=true; // first time step ??
199 //boucle sur les pas de quads
203 //clock_t clocki= clock ();
204 //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
205 for (int non_unif=0;non_unif<2;non_unif++)
211 champ_emetteur->getField()->getArray()->setIJ(0,0,40);
213 //bool ok=false; // Is the time interval successfully solved ?
215 // Loop on the time interval tries
220 dec_emetteur.attachLocalField(champ_emetteur);
222 dec_emetteur.attachLocalField(champ_recepteur);
225 if(init) dec_emetteur.synchronize();
228 if (cas=="emetteur") {
229 // affiche(champ_emetteur);
230 dec_emetteur.sendData();
232 else if (cas=="recepteur")
234 dec_emetteur.recvData();
236 afficheGauthier1(*champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
240 MPI_Barrier(MPI_COMM_WORLD);
246 delete champ_emetteur;
247 delete champ_recepteur;
251 void ParaMEDMEMTest::testGauthier2()
253 const char save_vit_in_2[]="VITESSE_P1_OUT\n1\n2\n3\n63\n3\n80\n0\n 0 1 2\n 3 4 5\n 6 7 8\n 9 10 11\n 12 13 14\n 15 16 17\n 18 19 20\n 21 22 23\n 24 25 26\n 27 28 29\n 30 2 1\n 31 5 4\n 32 8 7\n 33 11 10\n 34 14 13\n 35 17 16\n 36 20 19\n 37 23 22\n 38 26 25\n 39 29 28\n 30 40 2\n 31 41 5\n 32 42 8\n 33 43 11\n 34 44 14\n 35 45 17\n 36 46 20\n 37 47 23\n 38 48 26\n 39 49 29\n 31 2 40\n 32 5 41\n 33 8 42\n 34 11 43\n 35 14 44\n 36 17 45\n 37 20 46\n 38 23 47\n 39 26 48\n 50 29 49\n 3 2 4\n 6 5 7\n 9 8 10\n 12 11 13\n 15 14 16\n 18 17 19\n 21 20 22\n 24 23 25\n 27 26 28\n 51 29 52\n 31 4 2\n 32 7 5\n 33 10 8\n 34 13 11\n 35 16 14\n 36 19 17\n 37 22 20\n 38 25 23\n 39 28 26\n 50 52 29\n 0 2 53\n 3 5 54\n 6 8 55\n 9 11 56\n 12 14 57\n 15 17 58\n 18 20 59\n 21 23 60\n 24 26 61\n 27 29 62\n 3 53 2\n 6 54 5\n 9 55 8\n 12 56 11\n 15 57 14\n 18 58 17\n 21 59 20\n 24 60 23\n 27 61 26\n 51 62 29\n 0 0 0\n 0.5 0 0\n 0.5 0.05 0\n 0 0.1 0\n 0.5 0.1 0\n 0.5 0.15 0\n 0 0.2 0\n 0.5 0.2 0\n 0.5 0.25 0\n 0 0.3 0\n 0.5 0.3 0\n 0.5 0.35 0\n 0 0.4 0\n 0.5 0.4 0\n 0.5 0.45 0\n 0 0.5 0\n 0.5 0.5 0\n 0.5 0.55 0\n 0 0.6 0\n 0.5 0.6 0\n 0.5 0.65 0\n 0 0.7 0\n 0.5 0.7 0\n 0.5 0.75 0\n 0 0.8 0\n 0.5 0.8 0\n 0.5 0.85 0\n 0 0.9 0\n 0.5 0.9 0\n 0.5 0.95 0\n 1 0 0\n 1 0.1 0\n 1 0.2 0\n 1 0.3 0\n 1 0.4 0\n 1 0.5 0\n 1 0.6 0\n 1 0.7 0\n 1 0.8 0\n 1 0.9 0\n 1 0.05 0\n 1 0.15 0\n 1 0.25 0\n 1 0.35 0\n 1 0.45 0\n 1 0.55 0\n 1 0.65 0\n 1 0.75 0\n 1 0.85 0\n 1 0.95 0\n 1 1 0\n 0 1 0\n 0.5 1 0\n 0 0.05 0\n 0 0.15 0\n 0 0.25 0\n 0 0.35 0\n 0 0.45 0\n 0 0.55 0\n 0 0.65 0\n 0 0.75 0\n 0 0.85 0\n 0 0.95 0\n2.9268\n3.1707\n3\n1\n 0 0 0\n 0 0 0\n 0 0 0.05\n 0 0 0.1\n 0 0 0.1\n 0 0 0.15\n 0 0 0.2\n 0 0 0.2\n 0 0 0.25\n 0 0 0.3\n 0 0 0.3\n 0 0 0.35\n 0 0 0.4\n 0 0 0.4\n 0 0 0.45\n 0 0 0.5\n 0 0 0.5\n 0 0 0.55\n 0 0 0.6\n 0 0 0.6\n 0 0 0.65\n 0 0 0.7\n 0 0 0.7\n 0 0 0.75\n 0 0 0.8\n 0 0 0.8\n 0 0 0.85\n 0 0 0.9\n 0 0 0.9\n 0 0 0.95\n 0 0 0\n 0 0 0.1\n 0 0 0.2\n 0 0 0.3\n 0 0 0.4\n 0 0 0.5\n 0 0 0.6\n 0 0 0.7\n 0 0 0.8\n 0 0 0.9\n 0 0 0.05\n 0 0 0.15\n 0 0 0.25\n 0 0 0.35\n 0 0 0.45\n 0 0 0.55\n 0 0 0.65\n 0 0 0.75\n 0 0 0.85\n 0 0 0.95\n 0 0 1\n 0 0 1\n 0 0 1\n 0 0 0.05\n 0 0 0.15\n 0 0 0.25\n 0 0 0.35\n 0 0 0.45\n 0 0 0.55\n 0 0 0.65\n 0 0 0.75\n 0 0 0.85\n 0 0 0.95\n1\n";
255 const char save_vit_out_0_2[]="vitesse_in_chaude\n0\n2\n3\n22\n4\n10\n-1081737852\n 0 1 3 2\n 2 3 5 4\n 4 5 7 6\n 6 7 9 8\n 8 9 11 10\n 10 11 13 12\n 12 13 15 14\n 14 15 17 16\n 16 17 19 18\n 18 19 21 20\n 0 0 0\n 1 0 0\n 0 0.1 0\n 1 0.1 0\n 0 0.2 0\n 1 0.2 0\n 0 0.3 0\n 1 0.3 0\n 0 0.4 0\n 1 0.4 0\n 0 0.5 0\n 1 0.5 0\n 0 0.6 0\n 1 0.6 0\n 0 0.7 0\n 1 0.7 0\n 0 0.8 0\n 1 0.8 0\n 0 0.9 0\n 1 0.9 0\n 0 1 0\n 1 1 0\n2.9268\n3.1707\n3\n1\n 0 0 0.05\n 0 0 0.15\n 0 0 0.25\n 0 0 0.35\n 0 0 0.45\n 0 0 0.55\n 0 0 0.65\n 0 0 0.75\n 0 0 0.85\n 0 0 0.95\n0\n";
256 const char save_vit_out_1_2[]="vitesse_in_chaude\n1\n2\n3\n22\n4\n10\n-1081737852\n 0 1 3 2\n 2 3 5 4\n 4 5 7 6\n 6 7 9 8\n 8 9 11 10\n 10 11 13 12\n 12 13 15 14\n 14 15 17 16\n 16 17 19 18\n 18 19 21 20\n 0 0 0\n 1 0 0\n 0 0.1 0\n 1 0.1 0\n 0 0.2 0\n 1 0.2 0\n 0 0.3 0\n 1 0.3 0\n 0 0.4 0\n 1 0.4 0\n 0 0.5 0\n 1 0.5 0\n 0 0.6 0\n 1 0.6 0\n 0 0.7 0\n 1 0.7 0\n 0 0.8 0\n 1 0.8 0\n 0 0.9 0\n 1 0.9 0\n 0 1 0\n 1 1 0\n2.9268\n3.1707\n3\n1\n 0 0 0.029375\n 0 0 0.029375\n 0 0 0.1\n 0 0 0.1\n 0 0 0.2\n 0 0 0.2\n 0 0 0.3\n 0 0 0.3\n 0 0 0.4\n 0 0 0.4\n 0 0 0.5\n 0 0 0.5\n 0 0 0.6\n 0 0 0.6\n 0 0 0.7\n 0 0 0.7\n 0 0 0.8\n 0 0 0.8\n 0 0 0.9\n 0 0 0.9\n 0 0 0.970625\n 0 0 0.970625\n0\n";
258 const char *save_vit_outs[2]={save_vit_out_1_2,save_vit_out_0_2};
260 const char save_vit_out_1_0[]="vitesse_in_chaude\n1\n2\n3\n22\n4\n10\n-1081737852\n 0 1 3 2\n 2 3 5 4\n 4 5 7 6\n 6 7 9 8\n 8 9 11 10\n 10 11 13 12\n 12 13 15 14\n 14 15 17 16\n 16 17 19 18\n 18 19 21 20\n 0 0 0\n 1 0 0\n 0 0.1 0\n 1 0.1 0\n 0 0.2 0\n 1 0.2 0\n 0 0.3 0\n 1 0.3 0\n 0 0.4 0\n 1 0.4 0\n 0 0.5 0\n 1 0.5 0\n 0 0.6 0\n 1 0.6 0\n 0 0.7 0\n 1 0.7 0\n 0 0.8 0\n 1 0.8 0\n 0 0.9 0\n 1 0.9 0\n 0 1 0\n 1 1 0\n2.9268\n3.1707\n3\n1\n 0 0 0.029375\n 0 0 0.029375\n 0 0 0.1\n 0 0 0.1\n 0 0 0.2\n 0 0 0.2\n 0 0 0.3\n 0 0 0.3\n 0 0 0.4\n 0 0 0.4\n 0 0 0.5\n 0 0 0.5\n 0 0 0.6\n 0 0 0.6\n 0 0 0.7\n 0 0 0.7\n 0 0 0.8\n 0 0 0.8\n 0 0 0.9\n 0 0 0.9\n 0 0 0.970625\n 0 0 0.970625\n0\n";
262 const char save_vit_in[]="VITESSE_P1_OUT\n1\n2\n3\n63\n3\n80\n0\n 0 1 2\n 3 4 5\n 6 7 8\n 9 10 11\n 12 13 14\n 15 16 17\n 18 19 20\n 21 22 23\n 24 25 26\n 27 28 29\n 30 2 1\n 31 5 4\n 32 8 7\n 33 11 10\n 34 14 13\n 35 17 16\n 36 20 19\n 37 23 22\n 38 26 25\n 39 29 28\n 30 40 2\n 31 41 5\n 32 42 8\n 33 43 11\n 34 44 14\n 35 45 17\n 36 46 20\n 37 47 23\n 38 48 26\n 39 49 29\n 31 2 40\n 32 5 41\n 33 8 42\n 34 11 43\n 35 14 44\n 36 17 45\n 37 20 46\n 38 23 47\n 39 26 48\n 50 29 49\n 3 2 4\n 6 5 7\n 9 8 10\n 12 11 13\n 15 14 16\n 18 17 19\n 21 20 22\n 24 23 25\n 27 26 28\n 51 29 52\n 31 4 2\n 32 7 5\n 33 10 8\n 34 13 11\n 35 16 14\n 36 19 17\n 37 22 20\n 38 25 23\n 39 28 26\n 50 52 29\n 0 2 53\n 3 5 54\n 6 8 55\n 9 11 56\n 12 14 57\n 15 17 58\n 18 20 59\n 21 23 60\n 24 26 61\n 27 29 62\n 3 53 2\n 6 54 5\n 9 55 8\n 12 56 11\n 15 57 14\n 18 58 17\n 21 59 20\n 24 60 23\n 27 61 26\n 51 62 29\n 0 0 0\n 0.5 0 0\n 0.5 0.05 0\n 0 0.1 0\n 0.5 0.1 0\n 0.5 0.15 0\n 0 0.2 0\n 0.5 0.2 0\n 0.5 0.25 0\n 0 0.3 0\n 0.5 0.3 0\n 0.5 0.35 0\n 0 0.4 0\n 0.5 0.4 0\n 0.5 0.45 0\n 0 0.5 0\n 0.5 0.5 0\n 0.5 0.55 0\n 0 0.6 0\n 0.5 0.6 0\n 0.5 0.65 0\n 0 0.7 0\n 0.5 0.7 0\n 0.5 0.75 0\n 0 0.8 0\n 0.5 0.8 0\n 0.5 0.85 0\n 0 0.9 0\n 0.5 0.9 0\n 0.5 0.95 0\n 1 0 0\n 1 0.1 0\n 1 0.2 0\n 1 0.3 0\n 1 0.4 0\n 1 0.5 0\n 1 0.6 0\n 1 0.7 0\n 1 0.8 0\n 1 0.9 0\n 1 0.05 0\n 1 0.15 0\n 1 0.25 0\n 1 0.35 0\n 1 0.45 0\n 1 0.55 0\n 1 0.65 0\n 1 0.75 0\n 1 0.85 0\n 1 0.95 0\n 1 1 0\n 0 1 0\n 0.5 1 0\n 0 0.05 0\n 0 0.15 0\n 0 0.25 0\n 0 0.35 0\n 0 0.45 0\n 0 0.55 0\n 0 0.65 0\n 0 0.75 0\n 0 0.85 0\n 0 0.95 0\n2.9268\n3.1707\n3\n1\n 0 0 0\n 0 0 0\n 0 0 0.05\n 0 0 0.1\n 0 0 0.1\n 0 0 0.15\n 0 0 0.2\n 0 0 0.2\n 0 0 0.25\n 0 0 0.3\n 0 0 0.3\n 0 0 0.35\n 0 0 0.4\n 0 0 0.4\n 0 0 0.45\n 0 0 0.5\n 0 0 0.5\n 0 0 0.55\n 0 0 0.6\n 0 0 0.6\n 0 0 0.65\n 0 0 0.7\n 0 0 0.7\n 0 0 0.75\n 0 0 0.8\n 0 0 0.8\n 0 0 0.85\n 0 0 0.9\n 0 0 0.9\n 0 0 0.95\n 0 0 0\n 0 0 0.1\n 0 0 0.2\n 0 0 0.3\n 0 0 0.4\n 0 0 0.5\n 0 0 0.6\n 0 0 0.7\n 0 0 0.8\n 0 0 0.9\n 0 0 0.05\n 0 0 0.15\n 0 0 0.25\n 0 0 0.35\n 0 0 0.45\n 0 0 0.55\n 0 0 0.65\n 0 0 0.75\n 0 0 0.85\n 0 0 0.95\n 0 0 1\n 0 0 1\n 0 0 1\n 0 0 0.05\n 0 0 0.15\n 0 0 0.25\n 0 0 0.35\n 0 0 0.45\n 0 0 0.55\n 0 0 0.65\n 0 0 0.75\n 0 0 0.85\n 0 0 0.95\n1\n";
264 double valuesExpected1[2]={0.,0.};
265 double valuesExpected2[2]={0.95,0.970625};
267 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};
268 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 };
270 double *valuesExpected3[2]={valuesExpected30,valuesExpected31};
273 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
274 MPI_Comm_size(MPI_COMM_WORLD,&size);
279 set<int> entree_chaude_ids;
280 Genepi_ids.insert(0);
281 for (int i=1;i<size;i++)
282 entree_chaude_ids.insert(i);
283 for (int type=0;type<2;type++)
285 MPIProcessorGroup entree_chaude_group(comm,entree_chaude_ids);
286 MPIProcessorGroup Genepi_group(comm,Genepi_ids);
289 InterpKernelDEC dec_vit_in_chaude(entree_chaude_group, Genepi_group);
291 if ( entree_chaude_group.containsMyRank())
293 istringstream save_vit(save_vit_in);
294 vitesse.restore(save_vit);
298 istringstream save_vit(save_vit_out_1_0);
299 vitesse.restore(save_vit);
300 vitesse._has_field_ownership=false;
304 delete [] vitesse._field;
305 // cette ligne est super importante sinon c'est tout faux !!!!!!!
308 // pour tester P1->P0
312 if (vitesse._type==1)
313 dec_vit_in_chaude.setMethod("P1");
317 dec_vit_in_chaude.attachLocalField((ICoCo::Field*) &vitesse);
319 dec_vit_in_chaude.synchronize();
322 // Envois - receptions
323 if (entree_chaude_group.containsMyRank())
325 dec_vit_in_chaude.sendData();
329 dec_vit_in_chaude.recvData();
331 if (entree_chaude_group.containsMyRank() )
335 ostringstream save_vit(save_vit_in_2);
336 vitesse.save(save_vit);
342 double pmin=1e38, pmax=-1e38;
344 for(int i=0;i<vitesse.nb_values()*vitesse._nb_field_components;i++)
346 double p=*(vitesse._field+i);
350 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected1[type],pmin,1e-12);
351 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected2[type],pmax,1e-12);
353 ostringstream save_vit(save_vit_outs[type]);
354 vitesse.save(save_vit);
356 for(int i=0;i<vitesse.nb_values();i++)
358 for(int c=0;c<vitesse._nb_field_components;c++)
360 double p=vitesse._field[i*vitesse._nb_field_components+c];
361 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected3[type][i*vitesse._nb_field_components+c],p,1e-12);
370 * Non regression test testing copy constructor of InterpKernelDEC.
372 void ParaMEDMEMTest::testGauthier3()
376 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
377 MPI_Comm_size(MPI_COMM_WORLD,&size);
382 set<int> emetteur_ids;
383 set<int> recepteur_ids;
384 emetteur_ids.insert(0);
387 recepteur_ids.insert(1);
389 recepteur_ids.insert(2);
391 emetteur_ids.insert(3);
392 if ((rank==0)||(rank==1))
395 MPIProcessorGroup recepteur_group(comm,recepteur_ids);
396 MPIProcessorGroup emetteur_group(comm,emetteur_ids);
399 if (recepteur_group.containsMyRank())
402 //freopen("recpeteur.out","w",stdout);
403 //freopen("recepteur.err","w",stderr);
408 // freopen("emetteur.out","w",stdout);
409 //freopen("emetteur.err","w",stderr);
411 double expected[8][4]={
415 {40.,1.,1e200,1e200},
419 {20.5,1.,1e200,1e200}
421 int expectedLgth[8]={4,4,2,2,4,4,2,2};
423 for (int send=0;send<2;send++)
424 for (int rec=0;rec<2;rec++)
426 std::vector<InterpKernelDEC> decu(1);
427 decu[0]=InterpKernelDEC(emetteur_group,recepteur_group);
428 InterpKernelDEC& dec_emetteur=decu[0];
429 ParaMEDMEM::ParaFIELD *champ_emetteur(0),*champ_recepteur(0);
430 ParaMEDMEM::ParaMESH *paramesh(0);
431 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh;
432 dec_emetteur.setOrientation(2);
435 mesh=init_quadGauthier1(is_master);
439 mesh=init_triangleGauthier1(is_master);
441 paramesh=new ParaMEDMEM::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"emetteur mesh");
442 ParaMEDMEM::ComponentTopology comptopo;
443 champ_emetteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
444 champ_emetteur->getField()->setNature(ConservativeVolumic);
445 champ_emetteur->setOwnSupport(true);
448 mesh=init_triangleGauthier1(is_master);
452 mesh=init_quadGauthier1(is_master);
454 paramesh=new ParaMEDMEM::ParaMESH(mesh,recepteur_group.containsMyRank()?recepteur_group:emetteur_group,"recepteur mesh");
455 champ_recepteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
456 champ_recepteur->getField()->setNature(ConservativeVolumic);
457 champ_recepteur->setOwnSupport(true);
460 champ_emetteur->getField()->getArray()->fillWithValue(1.);
464 MPI_Barrier(MPI_COMM_WORLD);
466 //clock_t clock0= clock ();
469 bool init=true; // first time step ??
471 //boucle sur les pas de quads
475 //clock_t clocki= clock ();
476 //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
477 for (int non_unif=0;non_unif<2;non_unif++)
483 champ_emetteur->getField()->getArray()->setIJ(0,0,40);
485 //bool ok=false; // Is the time interval successfully solved ?
487 // Loop on the time interval tries
492 dec_emetteur.attachLocalField(champ_emetteur);
494 dec_emetteur.attachLocalField(champ_recepteur);
497 if(init) dec_emetteur.synchronize();
500 if (cas=="emetteur") {
501 // affiche(champ_emetteur);
502 dec_emetteur.sendData();
504 else if (cas=="recepteur")
506 dec_emetteur.recvData();
508 afficheGauthier1(*champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
512 MPI_Barrier(MPI_COMM_WORLD);
518 delete champ_emetteur;
519 delete champ_recepteur;
524 * This test is the parallel version of MEDCouplingBasicsTest.test3D1DOnP1P0_1 test.
526 void ParaMEDMEMTest::testGauthier4()
529 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};
530 const int 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};
531 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};
532 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.};
533 const int 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};
534 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.};
535 const int 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};
539 MPI_Comm_size(MPI_COMM_WORLD,&size);
540 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
544 int nproc_source = 1;
546 set<int> procs_source;
547 set<int> procs_target;
549 for (int i=0; i<nproc_source; i++)
550 procs_source.insert(i);
551 for (int i=nproc_source; i<size; i++)
552 procs_target.insert(i);
553 self_procs.insert(rank);
555 ParaMEDMEM::MEDCouplingUMesh *mesh=0;
556 ParaMEDMEM::ParaMESH *paramesh=0;
557 ParaMEDMEM::ParaFIELD* parafield=0;
559 ParaMEDMEM::CommInterface interface;
561 ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
562 ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
563 ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
565 MPI_Barrier(MPI_COMM_WORLD);
566 if(source_group->containsMyRank())
568 std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
569 mesh=MEDCouplingUMesh::New(stream.str().c_str(),1);
570 mesh->allocateCells();
571 for(int i=0;i<18;i++)
572 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,sourceConn+2*i);
573 mesh->finishInsertingCells();
574 DataArrayDouble *myCoords=DataArrayDouble::New();
575 myCoords->alloc(19,3);
576 std::copy(sourceCoords,sourceCoords+19*3,myCoords->getPointer());
577 mesh->setCoords(myCoords);
579 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
580 ParaMEDMEM::ComponentTopology comptopo;
581 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh,comptopo);
582 double *value=parafield->getField()->getArray()->getPointer();
583 std::copy(sourceVals,sourceVals+19,value);
589 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
590 mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
591 mesh->allocateCells();
593 mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn0+8*i);
594 mesh->finishInsertingCells();
595 DataArrayDouble *myCoords=DataArrayDouble::New();
596 myCoords->alloc(20,3);
597 std::copy(targetCoords0,targetCoords0+20*3,myCoords->getPointer());
598 mesh->setCoords(myCoords);
600 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
601 ParaMEDMEM::ComponentTopology comptopo;
602 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
606 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
607 mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
608 mesh->allocateCells();
610 mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn1+8*i);
611 mesh->finishInsertingCells();
612 DataArrayDouble *myCoords=DataArrayDouble::New();
613 myCoords->alloc(28,3);
614 std::copy(targetCoords1,targetCoords1+28*3,myCoords->getPointer());
615 mesh->setCoords(myCoords);
617 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
618 ParaMEDMEM::ComponentTopology comptopo;
619 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
622 //test 1 - primaire -> secondaire
623 ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
624 dec.setIntersectionType(INTERP_KERNEL::PointLocator);
625 parafield->getField()->setNature(ConservativeVolumic);//very important
626 if (source_group->containsMyRank())
629 dec.attachLocalField(parafield);
631 dec.setForcedRenormalization(false);
637 dec.attachLocalField(parafield);
639 dec.setForcedRenormalization(false);
641 const double *res(parafield->getField()->getArray()->getConstPointer());
644 const double expected0[4]={0.49,7.956666666666667,27.29,0.};
646 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected0[i],res[i],1e-13);
650 const double expected1[6]={59.95666666666667,94.09,0.,125.69,202.89,296.09};
652 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected1[i],res[i],1e-13);
655 MPI_Barrier(MPI_COMM_WORLD);
656 if (source_group->containsMyRank())
659 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.};
660 const double *res(parafield->getField()->getArray()->getConstPointer());
661 for(int i=0;i<19;i++)
662 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],res[i],1e-13);
675 MPI_Barrier(MPI_COMM_WORLD);