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 "ParaMESH.hxx"
31 #include "ParaFIELD.hxx"
32 #include "ComponentTopology.hxx"
42 using namespace ParaMEDMEM;
43 using namespace ICoCo;
45 void afficheGauthier1( const TrioField& field, const double *vals, int lgth)
47 CPPUNIT_ASSERT_EQUAL(lgth,field._nb_elems);
48 for (int ele=0;ele<field._nb_elems;ele++)
49 CPPUNIT_ASSERT_DOUBLES_EQUAL(vals[ele],field._field[ele],1e-12);
52 void remplit_coordGauthier1(double* coords)
54 double angle,epaisseur;
55 angle=0*45*(asin(1.)/90);
61 coords[1*3+0]=cos(angle);
63 coords[1*3+2]=sin(angle);
66 coords[2*3+0]=-sin(angle);
68 coords[2*3+2]=cos(angle);
71 coords[3*3+d]=coords[1*3+d]+ coords[2*3+d];
76 coords[i*3+d]=coords[(i-4)*3+d];
77 coords[i*3+1]+=epaisseur;
82 void init_quadGauthier1(TrioField& champ_quad,int is_master)
85 champ_quad.setName("champ_quad");
86 champ_quad._space_dim=3;
87 champ_quad._mesh_dim=2;
88 champ_quad._nodes_per_elem=4;
89 champ_quad._itnumber=0;
92 champ_quad._nb_field_components=1;
96 champ_quad._nbnodes=8;
97 champ_quad._nb_elems=2;
99 champ_quad._coords=new double[champ_quad._nbnodes*champ_quad._space_dim];
100 //memcpy(afield._coords,sommets.addr(),champ_quad._nbnodes*champ_quad._space_dim*sizeof(double));
102 remplit_coordGauthier1(champ_quad._coords);
105 champ_quad._connectivity=new int[champ_quad._nb_elems*champ_quad._nodes_per_elem];
106 champ_quad._connectivity[0*champ_quad._nodes_per_elem+0]=0;
107 champ_quad._connectivity[0*champ_quad._nodes_per_elem+1]=1;
108 champ_quad._connectivity[0*champ_quad._nodes_per_elem+2]=3;
109 champ_quad._connectivity[0*champ_quad._nodes_per_elem+3]=2;
110 champ_quad._connectivity[1*champ_quad._nodes_per_elem+0]=4;
111 champ_quad._connectivity[1*champ_quad._nodes_per_elem+1]=5;
112 champ_quad._connectivity[1*champ_quad._nodes_per_elem+2]=7;
113 champ_quad._connectivity[1*champ_quad._nodes_per_elem+3]=6;
118 champ_quad._nbnodes=0;
119 champ_quad._nb_elems=0;
120 champ_quad._coords=new double[champ_quad._nbnodes*champ_quad._space_dim];
123 champ_quad._has_field_ownership=false;
125 //champ_quad._field=new double[champ_quad._nb_elems];
126 // assert(champ_quad._nb_field_components==1);
128 void init_triangleGauthier1(TrioField& champ_triangle,int is_master)
131 champ_triangle.setName("champ_triangle");
132 champ_triangle._space_dim=3;
133 champ_triangle._mesh_dim=2;
134 champ_triangle._nodes_per_elem=3;
135 champ_triangle._itnumber=0;
136 champ_triangle._time1=0;
137 champ_triangle._time2=1;
138 champ_triangle._nb_field_components=1;
142 champ_triangle._nb_elems=4;
143 champ_triangle._nbnodes=8;
145 champ_triangle._coords=new double[champ_triangle._nbnodes*champ_triangle._space_dim];
146 //memcpy(afield._coords,sommets.addr(),champ_triangle._nbnodes*champ_triangle._space_dim*sizeof(double));
147 remplit_coordGauthier1(champ_triangle._coords);
149 champ_triangle._connectivity=new int[champ_triangle._nb_elems*champ_triangle._nodes_per_elem];
150 champ_triangle._connectivity[0*champ_triangle._nodes_per_elem+0]=0;
151 champ_triangle._connectivity[0*champ_triangle._nodes_per_elem+1]=1;
152 champ_triangle._connectivity[0*champ_triangle._nodes_per_elem+2]=2;
153 champ_triangle._connectivity[1*champ_triangle._nodes_per_elem+0]=1;
154 champ_triangle._connectivity[1*champ_triangle._nodes_per_elem+1]=2;
155 champ_triangle._connectivity[1*champ_triangle._nodes_per_elem+2]=3;
157 champ_triangle._connectivity[2*champ_triangle._nodes_per_elem+0]=4;
158 champ_triangle._connectivity[2*champ_triangle._nodes_per_elem+1]=5;
159 champ_triangle._connectivity[2*champ_triangle._nodes_per_elem+2]=7;
160 champ_triangle._connectivity[3*champ_triangle._nodes_per_elem+0]=4;
161 champ_triangle._connectivity[3*champ_triangle._nodes_per_elem+1]=6;
162 champ_triangle._connectivity[3*champ_triangle._nodes_per_elem+2]=7;
166 champ_triangle._nb_elems=0;
167 champ_triangle._nbnodes=0;
168 champ_triangle._coords=new double[champ_triangle._nbnodes*champ_triangle._space_dim];
171 champ_triangle._has_field_ownership=false;
172 // champ_triangle._field=new double[champ_triangle._nb_elems];
173 champ_triangle._field=0;
178 void ParaMEDMEMTest::testGauthier1()
182 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
183 MPI_Comm_size(MPI_COMM_WORLD,&size);
188 set<int> emetteur_ids;
189 set<int> recepteur_ids;
190 emetteur_ids.insert(0);
193 recepteur_ids.insert(1);
195 recepteur_ids.insert(2);
197 emetteur_ids.insert(3);
198 if ((rank==0)||(rank==1))
201 MPIProcessorGroup recepteur_group(comm,recepteur_ids);
202 MPIProcessorGroup emetteur_group(comm,emetteur_ids);
206 if (recepteur_group.containsMyRank())
209 //freopen("recpeteur.out","w",stdout);
210 //freopen("recepteur.err","w",stderr);
216 // freopen("emetteur.out","w",stdout);
217 //freopen("emetteur.err","w",stderr);
219 double expected[8][4]={
223 {40.,1.,1e200,1e200},
227 {20.5,1.,1e200,1e200}
230 int expectedLgth[8]={4,4,2,2,4,4,2,2};
232 for (int send=0;send<2;send++)
233 for (int rec=0;rec<2;rec++)
235 InterpKernelDEC dec_emetteur(emetteur_group, recepteur_group);
236 dec_emetteur.setOrientation(2);
237 TrioField champ_emetteur, champ_recepteur;
240 init_quadGauthier1(champ_emetteur,is_master);
242 init_triangleGauthier1(champ_emetteur,is_master);
244 init_triangleGauthier1(champ_recepteur,is_master);
246 init_quadGauthier1(champ_recepteur,is_master);
250 champ_emetteur._field=new double[champ_emetteur._nb_elems];
251 for (int ele=0;ele<champ_emetteur._nb_elems;ele++)
252 champ_emetteur._field[ele]=1;
254 champ_emetteur._has_field_ownership=true;
258 MPI_Barrier(MPI_COMM_WORLD);
260 //clock_t clock0= clock ();
263 bool init=true; // first time step ??
265 //boucle sur les pas de quads
269 //clock_t clocki= clock ();
270 //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
271 for (int non_unif=0;non_unif<2;non_unif++)
273 // if (champ_recepteur._field)
274 // delete [] champ_recepteur._field;
275 champ_recepteur._field=0;
276 // champ_recepteur._has_field_ownership=false;
284 champ_emetteur._field[0]=40;
286 //bool ok=false; // Is the time interval successfully solved ?
288 // Loop on the time interval tries
293 dec_emetteur.attachLocalField((ICoCo::Field*) &champ_emetteur);
295 dec_emetteur.attachLocalField((ICoCo::Field*) &champ_recepteur);
298 if(init) dec_emetteur.synchronize();
301 if (cas=="emetteur") {
302 // affiche(champ_emetteur);
303 dec_emetteur.sendData();
305 else if (cas=="recepteur")
307 dec_emetteur.recvData();
309 afficheGauthier1(champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
313 MPI_Barrier(MPI_COMM_WORLD);
318 // destruction des champs, des DEC, et des tableaux associés
323 void ParaMEDMEMTest::testGauthier2()
325 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";
327 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";
328 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";
330 const char *save_vit_outs[2]={save_vit_out_1_2,save_vit_out_0_2};
332 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";
334 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";
336 double valuesExpected1[2]={0.,0.};
337 double valuesExpected2[2]={0.95,0.970625};
339 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};
340 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 };
342 double *valuesExpected3[2]={valuesExpected30,valuesExpected31};
345 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
346 MPI_Comm_size(MPI_COMM_WORLD,&size);
351 set<int> entree_chaude_ids;
352 Genepi_ids.insert(0);
353 for (int i=1;i<size;i++)
354 entree_chaude_ids.insert(i);
355 for (int type=0;type<2;type++)
357 MPIProcessorGroup entree_chaude_group(comm,entree_chaude_ids);
358 MPIProcessorGroup Genepi_group(comm,Genepi_ids);
361 InterpKernelDEC dec_vit_in_chaude(entree_chaude_group, Genepi_group);
363 if ( entree_chaude_group.containsMyRank())
365 istringstream save_vit(save_vit_in);
366 vitesse.restore(save_vit);
370 istringstream save_vit(save_vit_out_1_0);
371 vitesse.restore(save_vit);
372 vitesse._has_field_ownership=false;
376 delete [] vitesse._field;
377 // cette ligne est super importante sinon c'est tout faux !!!!!!!
380 // pour tester P1->P0
384 if (vitesse._type==1)
385 dec_vit_in_chaude.setMethod("P1");
389 dec_vit_in_chaude.attachLocalField((ICoCo::Field*) &vitesse);
391 dec_vit_in_chaude.synchronize();
394 // Envois - receptions
395 if (entree_chaude_group.containsMyRank())
397 dec_vit_in_chaude.sendData();
401 dec_vit_in_chaude.recvData();
403 if (entree_chaude_group.containsMyRank() )
407 ostringstream save_vit(save_vit_in_2);
408 vitesse.save(save_vit);
414 double pmin=1e38, pmax=-1e38;
416 for(int i=0;i<vitesse.nb_values()*vitesse._nb_field_components;i++)
418 double p=*(vitesse._field+i);
422 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected1[type],pmin,1e-12);
423 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected2[type],pmax,1e-12);
425 ostringstream save_vit(save_vit_outs[type]);
426 vitesse.save(save_vit);
428 for(int i=0;i<vitesse.nb_values();i++)
430 for(int c=0;c<vitesse._nb_field_components;c++)
432 double p=vitesse._field[i*vitesse._nb_field_components+c];
433 CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected3[type][i*vitesse._nb_field_components+c],p,1e-12);
442 * Non regression test testing copy constructor of InterpKernelDEC.
444 void ParaMEDMEMTest::testGauthier3()
448 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
449 MPI_Comm_size(MPI_COMM_WORLD,&size);
454 set<int> emetteur_ids;
455 set<int> recepteur_ids;
456 emetteur_ids.insert(0);
459 recepteur_ids.insert(1);
461 recepteur_ids.insert(2);
463 emetteur_ids.insert(3);
464 if ((rank==0)||(rank==1))
467 MPIProcessorGroup recepteur_group(comm,recepteur_ids);
468 MPIProcessorGroup emetteur_group(comm,emetteur_ids);
472 if (recepteur_group.containsMyRank())
475 //freopen("recpeteur.out","w",stdout);
476 //freopen("recepteur.err","w",stderr);
482 // freopen("emetteur.out","w",stdout);
483 //freopen("emetteur.err","w",stderr);
485 double expected[8][4]={
489 {40.,1.,1e200,1e200},
493 {20.5,1.,1e200,1e200}
496 int expectedLgth[8]={4,4,2,2,4,4,2,2};
498 for (int send=0;send<2;send++)
499 for (int rec=0;rec<2;rec++)
502 std::vector<InterpKernelDEC> decu(1);
503 decu[0]=InterpKernelDEC(emetteur_group, recepteur_group);
504 InterpKernelDEC& dec_emetteur=decu[0];
506 //InterpKernelDEC dec_emetteur(emetteur_group, recepteur_group);
507 dec_emetteur.setOrientation(2);
508 TrioField champ_emetteur, champ_recepteur;
511 init_quadGauthier1(champ_emetteur,is_master);
513 init_triangleGauthier1(champ_emetteur,is_master);
515 init_triangleGauthier1(champ_recepteur,is_master);
517 init_quadGauthier1(champ_recepteur,is_master);
521 champ_emetteur._field=new double[champ_emetteur._nb_elems];
522 for (int ele=0;ele<champ_emetteur._nb_elems;ele++)
523 champ_emetteur._field[ele]=1;
525 champ_emetteur._has_field_ownership=true;
529 MPI_Barrier(MPI_COMM_WORLD);
531 //clock_t clock0= clock ();
534 bool init=true; // first time step ??
536 //boucle sur les pas de quads
540 //clock_t clocki= clock ();
541 //cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
542 for (int non_unif=0;non_unif<2;non_unif++)
544 // if (champ_recepteur._field)
545 // delete [] champ_recepteur._field;
546 champ_recepteur._field=0;
547 // champ_recepteur._has_field_ownership=false;
555 champ_emetteur._field[0]=40;
557 //bool ok=false; // Is the time interval successfully solved ?
559 // Loop on the time interval tries
564 dec_emetteur.attachLocalField((ICoCo::Field*) &champ_emetteur);
566 dec_emetteur.attachLocalField((ICoCo::Field*) &champ_recepteur);
569 if(init) dec_emetteur.synchronize();
572 if (cas=="emetteur") {
573 // affiche(champ_emetteur);
574 dec_emetteur.sendData();
576 else if (cas=="recepteur")
578 dec_emetteur.recvData();
580 afficheGauthier1(champ_recepteur,expected[num_cas],expectedLgth[num_cas]);
584 MPI_Barrier(MPI_COMM_WORLD);
589 // destruction des champs, des DEC, et des tableaux associés
595 * This test is the parallel version of MEDCouplingBasicsTest.test3D1DOnP1P0_1 test.
597 void ParaMEDMEMTest::testGauthier4()
600 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};
601 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};
602 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};
603 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.};
604 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};
605 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.};
606 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};
610 MPI_Comm_size(MPI_COMM_WORLD,&size);
611 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
615 int nproc_source = 1;
617 set<int> procs_source;
618 set<int> procs_target;
620 for (int i=0; i<nproc_source; i++)
621 procs_source.insert(i);
622 for (int i=nproc_source; i<size; i++)
623 procs_target.insert(i);
624 self_procs.insert(rank);
626 ParaMEDMEM::MEDCouplingUMesh *mesh=0;
627 ParaMEDMEM::ParaMESH *paramesh=0;
628 ParaMEDMEM::ParaFIELD* parafield=0;
630 ParaMEDMEM::CommInterface interface;
632 ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
633 ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
634 ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
636 MPI_Barrier(MPI_COMM_WORLD);
637 if(source_group->containsMyRank())
639 std::ostringstream stream; stream << "sourcemesh2D proc " << rank;
640 mesh=MEDCouplingUMesh::New(stream.str().c_str(),1);
641 mesh->allocateCells();
642 for(int i=0;i<18;i++)
643 mesh->insertNextCell(INTERP_KERNEL::NORM_SEG2,2,sourceConn+2*i);
644 mesh->finishInsertingCells();
645 DataArrayDouble *myCoords=DataArrayDouble::New();
646 myCoords->alloc(19,3);
647 std::copy(sourceCoords,sourceCoords+19*3,myCoords->getPointer());
648 mesh->setCoords(myCoords);
650 paramesh=new ParaMESH(mesh,*source_group,"source mesh");
651 ParaMEDMEM::ComponentTopology comptopo;
652 parafield = new ParaFIELD(ON_NODES,NO_TIME,paramesh,comptopo);
653 double *value=parafield->getField()->getArray()->getPointer();
654 std::copy(sourceVals,sourceVals+19,value);
660 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
661 mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
662 mesh->allocateCells();
664 mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn0+8*i);
665 mesh->finishInsertingCells();
666 DataArrayDouble *myCoords=DataArrayDouble::New();
667 myCoords->alloc(20,3);
668 std::copy(targetCoords0,targetCoords0+20*3,myCoords->getPointer());
669 mesh->setCoords(myCoords);
671 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
672 ParaMEDMEM::ComponentTopology comptopo;
673 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
677 std::ostringstream stream; stream << "targetmesh2D proc " << rank-nproc_source;
678 mesh=MEDCouplingUMesh::New(stream.str().c_str(),3);
679 mesh->allocateCells();
681 mesh->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8,targetConn1+8*i);
682 mesh->finishInsertingCells();
683 DataArrayDouble *myCoords=DataArrayDouble::New();
684 myCoords->alloc(28,3);
685 std::copy(targetCoords1,targetCoords1+28*3,myCoords->getPointer());
686 mesh->setCoords(myCoords);
688 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
689 ParaMEDMEM::ComponentTopology comptopo;
690 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
693 //test 1 - primaire -> secondaire
694 ParaMEDMEM::InterpKernelDEC dec(*source_group,*target_group);
695 dec.setIntersectionType(INTERP_KERNEL::PointLocator);
696 parafield->getField()->setNature(ConservativeVolumic);//very important
697 if (source_group->containsMyRank())
700 dec.attachLocalField(parafield);
702 dec.setForcedRenormalization(false);
708 dec.attachLocalField(parafield);
710 dec.setForcedRenormalization(false);
712 const double *res(parafield->getField()->getArray()->getConstPointer());
715 const double expected0[4]={0.49,7.956666666666667,27.29,0.};
717 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected0[i],res[i],1e-13);
721 const double expected1[6]={59.95666666666667,94.09,0.,125.69,202.89,296.09};
723 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected1[i],res[i],1e-13);
726 MPI_Barrier(MPI_COMM_WORLD);
727 if (source_group->containsMyRank())
730 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.};
731 const double *res(parafield->getField()->getArray()->getConstPointer());
732 for(int i=0;i<19;i++)
733 CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],res[i],1e-13);
746 MPI_Barrier(MPI_COMM_WORLD);