1 // Copyright (C) 2007-2015 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 "CommInterface.hxx"
22 #include "ProcessorGroup.hxx"
23 #include "MPIProcessorGroup.hxx"
24 #include "ComponentTopology.hxx"
25 #include "ParaMESH.hxx"
26 #include "ParaFIELD.hxx"
27 #include "InterpKernelDEC.hxx"
29 #include "MEDCouplingUMesh.hxx"
38 using namespace ParaMEDMEM;
39 using namespace ICoCo;
41 typedef enum {sync_and,sync_or} synctype;
42 void synchronize_bool(bool& stop, synctype s)
45 int my_stop_temp = stop?1:0;
47 MPI_Allreduce(&my_stop_temp,&my_stop,1,MPI_INTEGER,MPI_MIN,MPI_COMM_WORLD);
49 MPI_Allreduce(&my_stop_temp,&my_stop,1,MPI_INTEGER,MPI_MAX,MPI_COMM_WORLD);
53 void synchronize_dt(double& dt)
56 MPI_Allreduce(&dttemp,&dt,1,MPI_DOUBLE,MPI_MIN,MPI_COMM_WORLD);
60 void affiche(const ParaFIELD& field)
62 cout <<field.getField()->getName()<<endl;
63 const double *vals(field.getField()->getArray()->begin());
64 for(int ele=0;ele<field.getField()->getNumberOfTuples();ele++)
65 cout << ele <<": "<< vals[ele] << endl;
68 MEDCouplingUMesh *init_quad()
70 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_quad",2));
71 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo(DataArrayDouble::New());
72 const double dataCoo[24]={0.,0.,0.,1.,0.,0.,0.,0.,1.,1.,0.,1.,0.,1e-05,0.,1.,1e-05,0.,0.,1e-05,1.,1.,1e-05,1.};
74 std::copy(dataCoo,dataCoo+24,coo->getPointer());
75 const int conn[8]={0,1,3,2,4,5,7,6};
77 m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn);
78 m->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn+4);
83 MEDCouplingUMesh *init_triangle()
85 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m(MEDCouplingUMesh::New("champ_triangle",2));
86 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo(DataArrayDouble::New());
87 const double dataCoo[24]={0.,0.,0.,1.,0.,0.,0.,0.,1.,1.,0.,1.,0.,1e-05,0.,1.,1e-05,0.,0.,1e-05,1.,1.,1e-05,1.};
89 std::copy(dataCoo,dataCoo+24,coo->getPointer());
90 const int conn[12]={0,1,2,1,2,3,4,5,7,4,6,7};
93 m->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn+3*i);
98 void ParaMEDMEMTest::testICoco1()
102 MPI_Comm_size(MPI_COMM_WORLD,&size);
103 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
105 //the test is meant to run on 2 processors
106 if (size !=2) return ;
109 set<int> emetteur_ids;
110 set<int> recepteur_ids;
111 emetteur_ids.insert(0);
112 recepteur_ids.insert(1);
114 MPIProcessorGroup recepteur_group(comm,recepteur_ids);
115 MPIProcessorGroup emetteur_group(comm,emetteur_ids);
118 if (recepteur_group.containsMyRank())
123 InterpKernelDEC dec_emetteur(emetteur_group,recepteur_group);
124 dec_emetteur.setOrientation(2);
125 ParaMEDMEM::ParaFIELD *champ_emetteur(0),*champ_recepteur(0);
126 ParaMEDMEM::ParaMESH *paramesh(0);
129 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingUMesh> mesh_emetteur(init_triangle());
130 paramesh=new ParaMEDMEM::ParaMESH(mesh_emetteur,emetteur_group,"emetteur mesh");
131 ParaMEDMEM::ComponentTopology comptopo;
132 champ_emetteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
133 champ_emetteur->getField()->setNature(ConservativeVolumic);
134 champ_emetteur->setOwnSupport(true);
135 champ_emetteur->getField()->getArray()->fillWithValue(1.);
139 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingUMesh> mesh_recepteur(init_quad());
140 paramesh=new ParaMEDMEM::ParaMESH(mesh_recepteur,recepteur_group,"recepteur mesh");
141 ParaMEDMEM::ComponentTopology comptopo;
142 champ_recepteur=new ParaMEDMEM::ParaFIELD(ON_CELLS,ONE_TIME,paramesh,comptopo);
143 champ_recepteur->getField()->setNature(ConservativeVolumic);
144 champ_recepteur->setOwnSupport(true);
148 MPI_Barrier(MPI_COMM_WORLD);
150 clock_t clock0(clock());
153 bool init(true),stop(false);
154 //boucle sur les pas de quads
158 clock_t clocki= clock ();
159 cout << compti << " CLOCK " << (clocki-clock0)*1.e-6 << endl;
160 for (int non_unif=0;non_unif<2;non_unif++)
164 champ_emetteur->getField()->getArray()->setIJ(0,0,40.);
165 //bool ok=false; // Is the time interval successfully solved ?
167 // Loop on the time interval tries
169 dec_emetteur.attachLocalField(champ_emetteur);
171 dec_emetteur.attachLocalField(champ_recepteur);
174 dec_emetteur.synchronize();
179 dec_emetteur.sendData();
180 affiche(*champ_emetteur);
182 else if (cas=="recepteur")
184 dec_emetteur.recvData();
185 affiche(*champ_recepteur);
192 delete champ_recepteur;
193 delete champ_emetteur;