1 // Copyright (C) 2007-2016 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
21 #include <sys/times.h>
23 #include "ParaMEDMEMTest.hxx"
24 #include <cppunit/TestAssert.h>
26 #include "CommInterface.hxx"
27 #include "ProcessorGroup.hxx"
28 #include "MPIProcessorGroup.hxx"
29 #include "Topology.hxx"
31 #include "MxN_Mapping.hxx"
32 #include "InterpKernelDEC.hxx"
33 #include "ParaMESH.hxx"
34 #include "ParaFIELD.hxx"
35 #include "ComponentTopology.hxx"
36 #include "ICoCoMEDField.hxx"
37 #include "MEDLoader.hxx"
42 // use this define to enable lines, execution of which leads to Segmentation Fault
45 // use this define to enable CPPUNIT asserts and fails, showing bugs
46 #define ENABLE_FORCED_FAILURES
50 #define CLK_TCK sysconf(_SC_CLK_TCK);
54 using namespace MEDCoupling;
56 void testInterpKernelDEC_2D(const string& filename1, const string& meshname1,
57 const string& filename2, const string& meshname2,
58 int nproc_source, double epsilon, bool tri, bool all);
59 void get_time( float *telps, float *tuser, float *tsys, float *tcpu );
61 int main(int argc, char *argv[])
63 string filename1, filename2;
64 string meshname1, meshname2;
65 int nproc_source=1, rank;
71 MPI_Init(&argc,&argv);
73 for(int i=1;i<argc;i++){
74 if( strcmp(argv[i],"-f1") == 0 ){
75 filename1 = argv[++i];
78 else if( strcmp(argv[i],"-f2") == 0 ){
79 filename2 = argv[++i];
82 else if( strcmp(argv[i],"-m1") == 0 ){
83 meshname1 = argv[++i];
86 else if( strcmp(argv[i],"-m2") == 0 ){
87 meshname2 = argv[++i];
90 else if( strcmp(argv[i],"-ns") == 0 ){
91 nproc_source = atoi(argv[++i]);
93 else if( strcmp(argv[i],"-eps") == 0 ){
94 epsilon = atof(argv[++i]);
96 else if( strcmp(argv[i],"-tri") == 0 ){
99 else if( strcmp(argv[i],"-all") == 0 ){
105 cout << "usage test_perf -f1 filename1 -m1 meshname1 -f2 filename2 -m2 meshname2 (-ns nproc_source -eps epsilon -tri -all)" << endl;
109 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
110 testInterpKernelDEC_2D(filename1,meshname1,filename2,meshname2,nproc_source,epsilon,tri,all);
115 void testInterpKernelDEC_2D(const string& filename_xml1, const string& meshname1,
116 const string& filename_xml2, const string& meshname2,
117 int nproc_source, double epsilon, bool tri, bool all)
119 float tcpu, tcpu_u, tcpu_s, telps;
122 MPI_Comm_size(MPI_COMM_WORLD,&size);
123 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
126 set<int> procs_source;
127 set<int> procs_target;
129 for (int i=0; i<nproc_source; i++)
130 procs_source.insert(i);
131 for (int i=nproc_source; i<size; i++)
132 procs_target.insert(i);
133 self_procs.insert(rank);
135 MEDCoupling::CommInterface interface;
137 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
138 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
139 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
141 //loading the geometry for the source group
143 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
145 dec.setIntersectionType(INTERP_KERNEL::Triangulation);
147 dec.setIntersectionType(INTERP_KERNEL::Convex);
149 MEDCoupling::MEDCouplingUMesh* mesh;
150 MEDCoupling::ParaMESH* paramesh;
151 MEDCoupling::ParaFIELD* parafield;
152 ICoCo::MEDField* icocofield ;
154 // To remove tmp files from disk
155 ParaMEDMEMTest_TmpFilesRemover aRemover;
157 MPI_Barrier(MPI_COMM_WORLD);
158 if (source_group->containsMyRank()){
159 string master = filename_xml1;
161 ostringstream strstream;
162 if( nproc_source == 1 )
163 strstream <<master<<".med";
165 strstream <<master<<rank+1<<".med";
167 ostringstream meshname ;
168 if( nproc_source == 1 )
169 meshname<< meshname1;
171 meshname<< meshname1<<"_"<< rank+1;
173 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
174 mesh=ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
175 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
177 cout << "IO : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
180 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
182 MEDCoupling::ComponentTopology comptopo;
183 parafield = new ParaFIELD(ON_CELLS, NO_TIME, paramesh, comptopo);
185 int nb_local=mesh->getNumberOfCells();
186 double *value=parafield->getField()->getArray()->getPointer();
187 for(int ielem=0; ielem<nb_local;ielem++)
190 icocofield=new ICoCo::MEDField(parafield->getField());
192 dec.attachLocalField(icocofield);
195 //loading the geometry for the target group
196 if (target_group->containsMyRank()){
197 string master= filename_xml2;
198 ostringstream strstream;
199 if( (size-nproc_source) == 1 )
200 strstream << master<<".med";
202 strstream << master<<(rank-nproc_source+1)<<".med";
203 ostringstream meshname ;
204 if( (size-nproc_source) == 1 )
205 meshname<< meshname2;
207 meshname<< meshname2<<"_"<<rank-nproc_source+1;
209 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
210 mesh = ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
211 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
214 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
215 MEDCoupling::ComponentTopology comptopo;
216 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
218 int nb_local=mesh->getNumberOfCells();
219 double *value=parafield->getField()->getArray()->getPointer();
220 for(int ielem=0; ielem<nb_local;ielem++)
222 icocofield=new ICoCo::MEDField(parafield->getField());
224 dec.attachLocalField(icocofield);
228 //attaching a DEC to the source group
229 double field_before_int;
230 double field_after_int;
232 if (source_group->containsMyRank()){
233 field_before_int = parafield->getVolumeIntegral(0,true);
234 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
236 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
238 cout << "SYNCHRONIZE : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
239 cout<<"DEC usage"<<endl;
240 dec.setForcedRenormalization(false);
242 dec.setAllToAllMethod(PointToPoint);
244 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
247 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
249 cout << "SEND DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
252 field_after_int = parafield->getVolumeIntegral(0,true);
253 // CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, epsilon);
257 //attaching a DEC to the target group
258 if (target_group->containsMyRank()){
259 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
261 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
262 dec.setForcedRenormalization(false);
264 dec.setAllToAllMethod(PointToPoint);
266 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
268 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
272 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
274 cout << "RECV DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
284 MPI_Barrier(MPI_COMM_WORLD);
285 cout << "end of InterpKernelDEC_2D test"<<endl;
288 void get_time( float *telps, float *tuser, float *tsys, float *tcpu )
291 /* Variables declaration */
292 static time_t zsec = 0;
293 static long zusec = 0;
296 static clock_t zclock = 0;
298 static clock_t zuser = 0;
299 static clock_t zsys = 0;
306 MPI_Barrier(MPI_COMM_WORLD);
308 /* Elapsed time reading */
310 gettimeofday(&tp,&tzp);
313 *telps = (float)(nsec-zsec) + (float)(nusec-zusec)/(float)CLOCKS_PER_SEC;
318 /* User and system CPU time reading */
321 nuser = local.tms_utime;
322 nsys = local.tms_stime;
323 *tuser = (float)(nuser-zuser) / (float)CLK_TCK;
324 *tsys = (float)(nsys-zsys) / (float)CLK_TCK;
329 /* CPU time reading */
332 *tcpu = (float)(nclock-zclock) / (float)CLOCKS_PER_SEC;