1 // Copyright (C) 2007-2022 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 "MEDLoader.hxx"
37 #include "ICoCoMEDDoubleField.hxx"
38 #include "MEDCouplingUMesh.hxx"
44 // use this define to enable lines, execution of which leads to Segmentation Fault
47 // use this define to enable CPPUNIT asserts and fails, showing bugs
48 #define ENABLE_FORCED_FAILURES
52 #define CLK_TCK sysconf(_SC_CLK_TCK);
56 using namespace MEDCoupling;
58 void testInterpKernelDEC_2D(const string& filename1, const string& meshname1,
59 const string& filename2, const string& meshname2,
60 int nproc_source, double epsilon, bool tri, bool all);
61 void get_time( float *telps, float *tuser, float *tsys, float *tcpu );
63 int main(int argc, char *argv[])
65 string filename1, filename2;
66 string meshname1, meshname2;
67 int nproc_source=1, rank;
73 MPI_Init(&argc,&argv);
75 for(int i=1;i<argc;i++){
76 if( strcmp(argv[i],"-f1") == 0 ){
77 filename1 = argv[++i];
80 else if( strcmp(argv[i],"-f2") == 0 ){
81 filename2 = argv[++i];
84 else if( strcmp(argv[i],"-m1") == 0 ){
85 meshname1 = argv[++i];
88 else if( strcmp(argv[i],"-m2") == 0 ){
89 meshname2 = argv[++i];
92 else if( strcmp(argv[i],"-ns") == 0 ){
93 nproc_source = atoi(argv[++i]);
95 else if( strcmp(argv[i],"-eps") == 0 ){
96 epsilon = atof(argv[++i]);
98 else if( strcmp(argv[i],"-tri") == 0 ){
101 else if( strcmp(argv[i],"-all") == 0 ){
107 cout << "usage test_perf -f1 filename1 -m1 meshname1 -f2 filename2 -m2 meshname2 (-ns nproc_source -eps epsilon -tri -all)" << endl;
111 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
112 testInterpKernelDEC_2D(filename1,meshname1,filename2,meshname2,nproc_source,epsilon,tri,all);
117 void testInterpKernelDEC_2D(const string& filename_xml1, const string& meshname1,
118 const string& filename_xml2, const string& meshname2,
119 int nproc_source, double epsilon, bool tri, bool all)
121 float tcpu, tcpu_u, tcpu_s, telps;
124 MPI_Comm_size(MPI_COMM_WORLD,&size);
125 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
128 set<int> procs_source;
129 set<int> procs_target;
131 for (int i=0; i<nproc_source; i++)
132 procs_source.insert(i);
133 for (int i=nproc_source; i<size; i++)
134 procs_target.insert(i);
135 self_procs.insert(rank);
137 MEDCoupling::CommInterface interface;
139 MEDCoupling::ProcessorGroup* self_group = new MEDCoupling::MPIProcessorGroup(interface,self_procs);
140 MEDCoupling::ProcessorGroup* target_group = new MEDCoupling::MPIProcessorGroup(interface,procs_target);
141 MEDCoupling::ProcessorGroup* source_group = new MEDCoupling::MPIProcessorGroup(interface,procs_source);
143 //loading the geometry for the source group
145 MEDCoupling::InterpKernelDEC dec (*source_group,*target_group);
147 dec.setIntersectionType(INTERP_KERNEL::Triangulation);
149 dec.setIntersectionType(INTERP_KERNEL::Convex);
151 MEDCoupling::MEDCouplingUMesh* mesh = nullptr;
152 MEDCoupling::ParaMESH* paramesh = nullptr;
153 MEDCoupling::ParaFIELD* parafield = nullptr;
154 ICoCo::MEDDoubleField* icocofield = nullptr;
156 // To remove tmp files from disk
157 ParaMEDMEMTest_TmpFilesRemover aRemover;
159 MPI_Barrier(MPI_COMM_WORLD);
160 if (source_group->containsMyRank()){
161 string master = filename_xml1;
163 ostringstream strstream;
164 if( nproc_source == 1 )
165 strstream <<master<<".med";
167 strstream <<master<<rank+1<<".med";
169 ostringstream meshname ;
170 if( nproc_source == 1 )
171 meshname<< meshname1;
173 meshname<< meshname1<<"_"<< rank+1;
175 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
176 mesh=ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
177 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
179 cout << "IO : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
182 paramesh=new ParaMESH (mesh,*source_group,"source mesh");
184 MEDCoupling::ComponentTopology comptopo;
185 parafield = new ParaFIELD(ON_CELLS, NO_TIME, paramesh, comptopo);
187 int nb_local=mesh->getNumberOfCells();
188 double *value=parafield->getField()->getArray()->getPointer();
189 for(int ielem=0; ielem<nb_local;ielem++)
192 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
194 dec.attachLocalField(icocofield);
197 //loading the geometry for the target group
198 if (target_group->containsMyRank()){
199 string master= filename_xml2;
200 ostringstream strstream;
201 if( (size-nproc_source) == 1 )
202 strstream << master<<".med";
204 strstream << master<<(rank-nproc_source+1)<<".med";
205 ostringstream meshname ;
206 if( (size-nproc_source) == 1 )
207 meshname<< meshname2;
209 meshname<< meshname2<<"_"<<rank-nproc_source+1;
211 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
212 mesh = ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
213 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
216 paramesh=new ParaMESH (mesh,*target_group,"target mesh");
217 MEDCoupling::ComponentTopology comptopo;
218 parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
220 int nb_local=mesh->getNumberOfCells();
221 double *value=parafield->getField()->getArray()->getPointer();
222 for(int ielem=0; ielem<nb_local;ielem++)
224 icocofield=new ICoCo::MEDDoubleField(parafield->getField());
226 dec.attachLocalField(icocofield);
230 //attaching a DEC to the source group
231 double field_before_int;
232 double field_after_int;
234 if (source_group->containsMyRank()){
235 field_before_int = parafield->getVolumeIntegral(0,true);
236 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
238 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
240 cout << "SYNCHRONIZE : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
241 cout<<"DEC usage"<<endl;
242 dec.setForcedRenormalization(false);
244 dec.setAllToAllMethod(PointToPoint);
246 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
249 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
251 cout << "SEND DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
254 field_after_int = parafield->getVolumeIntegral(0,true);
255 // CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, epsilon);
259 //attaching a DEC to the target group
260 if (target_group->containsMyRank()){
261 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
263 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
264 dec.setForcedRenormalization(false);
266 dec.setAllToAllMethod(PointToPoint);
268 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
270 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
274 get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
276 cout << "RECV DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
286 MPI_Barrier(MPI_COMM_WORLD);
287 cout << "end of InterpKernelDEC_2D test"<<endl;
290 void get_time( float *telps, float *tuser, float *tsys, float *tcpu )
293 /* Variables declaration */
294 static time_t zsec = 0;
295 static long zusec = 0;
298 static clock_t zclock = 0;
300 static clock_t zuser = 0;
301 static clock_t zsys = 0;
308 MPI_Barrier(MPI_COMM_WORLD);
310 /* Elapsed time reading */
312 gettimeofday(&tp,&tzp);
315 *telps = (float)(nsec-zsec) + (float)(nusec-zusec)/(float)CLOCKS_PER_SEC;
320 /* User and system CPU time reading */
323 nuser = local.tms_utime;
324 nsys = local.tms_stime;
325 *tuser = (float)(nuser-zuser) / (float)CLK_TCK;
326 *tsys = (float)(nsys-zsys) / (float)CLK_TCK;
331 /* CPU time reading */
334 *tcpu = (float)(nclock-zclock) / (float)CLOCKS_PER_SEC;