Back to index

salome-med  6.5.0
test_perf.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include <time.h>
00021 #include <sys/times.h>
00022 #include <sys/time.h>
00023 #include "ParaMEDMEMTest.hxx"
00024 #include <cppunit/TestAssert.h>
00025 
00026 #include "CommInterface.hxx"
00027 #include "ProcessorGroup.hxx"
00028 #include "MPIProcessorGroup.hxx"
00029 #include "Topology.hxx"
00030 #include "DEC.hxx"
00031 #include "MxN_Mapping.hxx"
00032 #include "InterpKernelDEC.hxx"
00033 #include "ParaMESH.hxx"
00034 #include "ParaFIELD.hxx"
00035 #include "ComponentTopology.hxx"
00036 #include "ICoCoMEDField.hxx"
00037 #include "MEDLoader.hxx"
00038  
00039 #include <string>
00040 #include <cstring>
00041 
00042 // use this define to enable lines, execution of which leads to Segmentation Fault
00043 #define ENABLE_FAULTS
00044 
00045 // use this define to enable CPPUNIT asserts and fails, showing bugs
00046 #define ENABLE_FORCED_FAILURES
00047 
00048 #ifndef CLK_TCK 
00049 #include <unistd.h>
00050 #define CLK_TCK sysconf(_SC_CLK_TCK);
00051 #endif 
00052 
00053 using namespace std;
00054 using namespace ParaMEDMEM;
00055  
00056 void testInterpKernelDEC_2D(const string& filename1, const string& meshname1,
00057                             const string& filename2, const string& meshname2,
00058                             int nproc_source, double epsilon, bool tri, bool all);
00059 void get_time( float *telps, float *tuser, float *tsys, float *tcpu );
00060 
00061 int main(int argc, char *argv[])
00062 {
00063   string filename1, filename2;
00064   string meshname1, meshname2;
00065   int nproc_source=1, rank;
00066   double epsilon=1.e-6;
00067   int count=0;
00068   bool tri=false;
00069   bool all=false;
00070 
00071   MPI_Init(&argc,&argv);
00072 
00073   for(int i=1;i<argc;i++){
00074     if( strcmp(argv[i],"-f1") == 0 ){
00075       filename1 = argv[++i];
00076       count++;
00077     }
00078     else if( strcmp(argv[i],"-f2") == 0 ){
00079       filename2 = argv[++i];
00080       count++;
00081     }
00082     else if( strcmp(argv[i],"-m1") == 0 ){
00083       meshname1 = argv[++i];
00084       count++;
00085     }
00086     else if( strcmp(argv[i],"-m2") == 0 ){
00087       meshname2 = argv[++i];
00088       count++;
00089     }
00090     else if( strcmp(argv[i],"-ns") == 0 ){
00091       nproc_source = atoi(argv[++i]);
00092     }
00093     else if( strcmp(argv[i],"-eps") == 0 ){
00094       epsilon = atof(argv[++i]);
00095     }
00096     else if( strcmp(argv[i],"-tri") == 0 ){
00097       tri = true;
00098     }
00099     else if( strcmp(argv[i],"-all") == 0 ){
00100       all = true;
00101     }
00102   }
00103 
00104   if( count != 4 ){
00105     cout << "usage test_perf -f1 filename1 -m1 meshname1 -f2 filename2 -m2 meshname2 (-ns nproc_source -eps epsilon -tri -all)" << endl;
00106     exit(0);
00107   }
00108 
00109   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00110   testInterpKernelDEC_2D(filename1,meshname1,filename2,meshname2,nproc_source,epsilon,tri,all);
00111 
00112   MPI_Finalize();
00113 }
00114 
00115 void testInterpKernelDEC_2D(const string& filename_xml1, const string& meshname1,
00116                             const string& filename_xml2, const string& meshname2,
00117                             int nproc_source, double epsilon, bool tri, bool all)
00118 {
00119   float tcpu, tcpu_u, tcpu_s, telps;
00120   int size;
00121   int rank;
00122   MPI_Comm_size(MPI_COMM_WORLD,&size);
00123   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00124  
00125   set<int> self_procs;
00126   set<int> procs_source;
00127   set<int> procs_target;
00128   
00129   for (int i=0; i<nproc_source; i++)
00130     procs_source.insert(i);
00131   for (int i=nproc_source; i<size; i++)
00132     procs_target.insert(i);
00133   self_procs.insert(rank);
00134   
00135   ParaMEDMEM::CommInterface interface;
00136     
00137   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00138   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00139   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00140   
00141   //loading the geometry for the source group
00142 
00143   ParaMEDMEM::InterpKernelDEC dec (*source_group,*target_group);
00144   if(tri)
00145     dec.setIntersectionType(INTERP_KERNEL::Triangulation);
00146   else
00147     dec.setIntersectionType(INTERP_KERNEL::Convex);
00148 
00149   ParaMEDMEM::MEDCouplingUMesh* mesh;
00150   ParaMEDMEM::ParaMESH* paramesh;
00151   ParaMEDMEM::ParaFIELD* parafield;
00152   ICoCo::Field* icocofield ;
00153   
00154   // To remove tmp files from disk
00155   ParaMEDMEMTest_TmpFilesRemover aRemover;
00156   
00157   MPI_Barrier(MPI_COMM_WORLD);
00158   if (source_group->containsMyRank()){
00159     string master = filename_xml1;
00160       
00161     ostringstream strstream;
00162     if( nproc_source == 1 )
00163       strstream <<master<<".med";
00164     else
00165       strstream <<master<<rank+1<<".med";
00166 
00167     ostringstream meshname ;
00168     if( nproc_source == 1 )
00169       meshname<< meshname1;
00170     else
00171       meshname<< meshname1<<"_"<< rank+1;
00172       
00173     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00174     mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00175     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00176     if( rank == 0 )
00177       cout << "IO : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
00178     mesh->incrRef();
00179     
00180     paramesh=new ParaMESH (mesh,*source_group,"source mesh");
00181     
00182     ParaMEDMEM::ComponentTopology comptopo;
00183     parafield = new ParaFIELD(ON_CELLS, NO_TIME, paramesh, comptopo);
00184 
00185     int nb_local=mesh->getNumberOfCells();
00186     double *value=parafield->getField()->getArray()->getPointer();
00187     for(int ielem=0; ielem<nb_local;ielem++)
00188       value[ielem]=1.0;
00189     
00190     icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00191      
00192     dec.attachLocalField(icocofield);
00193   }
00194   
00195   //loading the geometry for the target group
00196   if (target_group->containsMyRank()){
00197     string master= filename_xml2;
00198     ostringstream strstream;
00199     if( (size-nproc_source) == 1 )
00200       strstream << master<<".med";
00201     else
00202       strstream << master<<(rank-nproc_source+1)<<".med";
00203     ostringstream meshname ;
00204     if( (size-nproc_source) == 1 )
00205       meshname<< meshname2;
00206     else
00207       meshname<< meshname2<<"_"<<rank-nproc_source+1;
00208       
00209     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00210     mesh = MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00211     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00212     mesh->incrRef();
00213 
00214     paramesh=new ParaMESH (mesh,*target_group,"target mesh");
00215     ParaMEDMEM::ComponentTopology comptopo;
00216     parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00217 
00218     int nb_local=mesh->getNumberOfCells();
00219     double *value=parafield->getField()->getArray()->getPointer();
00220     for(int ielem=0; ielem<nb_local;ielem++)
00221       value[ielem]=0.0;
00222     icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00223       
00224     dec.attachLocalField(icocofield);
00225   }
00226     
00227   
00228   //attaching a DEC to the source group 
00229   double field_before_int;
00230   double field_after_int;
00231   
00232   if (source_group->containsMyRank()){ 
00233     field_before_int = parafield->getVolumeIntegral(0,true);
00234     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00235     dec.synchronize();
00236     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00237     if( rank == 0 )
00238       cout << "SYNCHRONIZE : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
00239     cout<<"DEC usage"<<endl;
00240     dec.setForcedRenormalization(false);
00241     if(all)
00242       dec.setAllToAllMethod(PointToPoint);
00243 
00244     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00245     dec.sendData();
00246     
00247     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00248     if( rank == 0 )
00249       cout << "SEND DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
00250     dec.recvData();
00251      
00252     field_after_int = parafield->getVolumeIntegral(0,true);
00253 //    CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, epsilon);
00254       
00255   }
00256   
00257   //attaching a DEC to the target group
00258   if (target_group->containsMyRank()){
00259     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00260     dec.synchronize();
00261     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00262     dec.setForcedRenormalization(false);
00263     if(all)
00264       dec.setAllToAllMethod(PointToPoint);
00265 
00266     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00267     dec.recvData();
00268     get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00269     dec.sendData();
00270   }
00271   
00272   get_time( &telps, &tcpu_u, &tcpu_s, &tcpu );
00273   if( rank == 0 )
00274     cout << "RECV DATA : Telapse = " << telps << " TuserCPU = " << tcpu_u << " TsysCPU = " << tcpu_s << " TCPU = " << tcpu << endl;
00275 
00276   delete source_group;
00277   delete target_group;
00278   delete self_group;
00279   delete paramesh;
00280   delete parafield;
00281   mesh->decrRef() ;
00282   delete icocofield;
00283 
00284   MPI_Barrier(MPI_COMM_WORLD);
00285   cout << "end of InterpKernelDEC_2D test"<<endl;
00286 }
00287 
00288 void get_time( float *telps, float *tuser, float *tsys, float *tcpu )
00289 {
00290 
00291   /* Variables declaration */
00292   static time_t zsec = 0;
00293   static long zusec = 0;
00294   time_t nsec;
00295   long nusec;
00296   static clock_t zclock = 0;
00297   clock_t nclock;
00298   static clock_t zuser = 0;
00299   static clock_t zsys = 0;
00300   clock_t nuser, nsys;
00301 
00302   struct timeval tp;
00303   struct timezone tzp;
00304   struct tms local;
00305 
00306   MPI_Barrier(MPI_COMM_WORLD);
00307 
00308   /* Elapsed time reading */
00309 
00310   gettimeofday(&tp,&tzp);
00311   nsec = tp.tv_sec;
00312   nusec = tp.tv_usec;
00313   *telps = (float)(nsec-zsec) + (float)(nusec-zusec)/(float)CLOCKS_PER_SEC;
00314   
00315   zsec = nsec;
00316   zusec = nusec;
00317 
00318   /* User and system CPU time reading */
00319 
00320   times(&local);
00321   nuser = local.tms_utime;
00322   nsys = local.tms_stime;
00323   *tuser = (float)(nuser-zuser) / (float)CLK_TCK;
00324   *tsys = (float)(nsys-zsys) / (float)CLK_TCK;
00325 
00326   zuser = nuser;
00327   zsys = nsys;
00328 
00329   /* CPU time reading */
00330 
00331   nclock = clock();
00332   *tcpu = (float)(nclock-zclock) / (float)CLOCKS_PER_SEC;
00333   zclock = nclock;
00334 
00335 }
00336 
00337