Back to index

salome-med  6.5.0
MEDPARTITIONER_UtilsPara.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 "MEDPARTITIONER_Utils.hxx"
00021 
00022 #include "MEDLoader.hxx"
00023 #include "MEDLoaderBase.hxx"
00024 #include "MEDFileUtilities.hxx"
00025 #include "CellModel.hxx"
00026 #include "MEDCouplingUMesh.hxx"
00027 #include "MEDCouplingFieldDouble.hxx"
00028 #include "InterpKernelException.hxx"
00029 #include "MEDCouplingAutoRefCountObjectPtr.hxx"
00030 #include "InterpKernelAutoPtr.hxx"
00031 
00032 #include <fstream>
00033 #include <iostream>
00034 #include <iomanip>
00035 #include <sstream>
00036 #include <string>
00037 
00038 #ifdef HAVE_MPI2
00039 #include <mpi.h>
00040 #endif
00041 
00042 using namespace MEDPARTITIONER;
00043 
00048 std::vector<std::string> MEDPARTITIONER::SendAndReceiveVectorOfString(const std::vector<std::string>& vec, const int source, const int target)
00049 {
00050   int rank=MyGlobals::_Rank;
00051 
00052   MPI_Status status;
00053   int tag = 111001;
00054   if (rank == source)
00055     {
00056       std::string str=SerializeFromVectorOfString(vec);
00057       int size=str.length();
00058       MPI_Send( &size, 1, MPI_INT, target, tag, MPI_COMM_WORLD );
00059       MPI_Send( (void*)str.data(), str.length(), MPI_CHAR, target, tag+100, MPI_COMM_WORLD );
00060     }
00061   
00062   int recSize=0;
00063   if (rank == target)
00064     {
00065       MPI_Recv(&recSize, 1, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00066       std::string recData(recSize,'x');
00067       MPI_Recv((void*)recData.data(), recSize, MPI_CHAR, source, tag+100, MPI_COMM_WORLD, &status);
00068       return DeserializeToVectorOfString(recData); //not empty one for target proc
00069     }
00070   std::vector<std::string> res;
00071   return res; //empty one for other proc
00072 }
00073 
00077 std::vector<std::string> MEDPARTITIONER::AllgathervVectorOfString(const std::vector<std::string>& vec)
00078 {
00079   if (MyGlobals::_World_Size==1) //nothing to do
00080     return vec;
00081 
00082   int world_size=MyGlobals::_World_Size;
00083   std::string str=SerializeFromVectorOfString(vec);
00084   
00085   std::vector<int> indexes(world_size);
00086   int size=str.length();
00087   MPI_Allgather(&size, 1, MPI_INT, 
00088                 &indexes[0], 1, MPI_INT, MPI_COMM_WORLD);
00089   
00090   //calcul of displacement
00091   std::vector< int > disp(1,0);
00092   for (int i=0; i<world_size; i++) disp.push_back( disp.back() + indexes[i] );
00093   
00094   std::string recData(disp.back(),'x');
00095   MPI_Allgatherv((void*)str.data(), str.length(), MPI_CHAR,
00096                  (void*)recData.data(), &indexes[0], &disp[0], MPI_CHAR,
00097                  MPI_COMM_WORLD);
00098   
00099   //really extraordinary verbose for debug
00100   std::vector<std::string> deserial=DeserializeToVectorOfString(recData);
00101   if (MyGlobals::_Verbose>1000) 
00102     {
00103       std::cout << "proc "<<MyGlobals::_Rank<<" : receive '" << recData << "'" << std::endl;
00104       std::cout << "deserialize is : a vector of size " << deserial.size() << std::endl;
00105       std::cout << ReprVectorOfString(deserial) << std::endl;
00106     }
00107   return deserial;
00108 }
00109 
00115 void MEDPARTITIONER::SendDoubleVec(const std::vector<double>& vec, const int target)
00116 {
00117   int tag = 111002;
00118   int size=vec.size();
00119   if (MyGlobals::_Verbose>1000) 
00120     std::cout << "proc " << MyGlobals::_Rank << " : --> SendDoubleVec " << size << std::endl;
00121 #ifdef HAVE_MPI2
00122   MPI_Send(&size, 1, MPI_INT, target, tag, MPI_COMM_WORLD);
00123   MPI_Send(const_cast<double*>(&vec[0]), size, MPI_DOUBLE, target, tag+100, MPI_COMM_WORLD);
00124 #endif
00125 }
00126 
00133 std::vector<double>* MEDPARTITIONER::RecvDoubleVec(const int source)
00134 {
00135   int tag = 111002;
00136   int size;
00137 #ifdef HAVE_MPI2
00138   MPI_Status status;  
00139   MPI_Recv(&size, 1, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00140   if (MyGlobals::_Verbose>1000) 
00141     std::cout << "proc " << MyGlobals::_Rank << " : <-- RecvDoubleVec " << size << std::endl;
00142   std::vector<double>* vec=new std::vector<double>;
00143   vec->resize(size);
00144   MPI_Recv(&vec[0], size, MPI_DOUBLE, source, tag+100, MPI_COMM_WORLD, &status);
00145 #endif
00146   return vec;
00147 }
00148 
00149 void MEDPARTITIONER::RecvDoubleVec(std::vector<double>& vec, const int source)
00150 {
00151   int tag = 111002;
00152   int size;
00153 #ifdef HAVE_MPI2
00154   MPI_Status status;  
00155   MPI_Recv(&size, 1, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00156   if (MyGlobals::_Verbose>1000)
00157     std::cout<< "proc " << MyGlobals::_Rank << " : <-- RecvDoubleVec " << size << std::endl;;
00158   vec.resize(size);
00159   MPI_Recv(&vec[0], size, MPI_DOUBLE, source, tag+100, MPI_COMM_WORLD, &status);
00160 #endif
00161 }
00167 void MEDPARTITIONER::SendIntVec(const std::vector<int>& vec, const int target)
00168 {
00169   int tag = 111003;
00170   int size=vec.size();
00171   if (MyGlobals::_Verbose>1000)
00172     std::cout << "proc " << MyGlobals::_Rank << " : --> SendIntVec " << size << std::endl;
00173 #ifdef HAVE_MPI2
00174   MPI_Send(&size, 1, MPI_INT, target, tag, MPI_COMM_WORLD);
00175   MPI_Send(const_cast<int*>(&vec[0]), size,MPI_INT, target, tag+100, MPI_COMM_WORLD);
00176 #endif
00177 }
00178 
00184 std::vector<int> *MEDPARTITIONER::RecvIntVec(const int source)
00185 {
00186   int tag = 111003;
00187   int size;
00188 #ifdef HAVE_MPI2
00189   MPI_Status status;  
00190   MPI_Recv(&size, 1, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00191   if (MyGlobals::_Verbose>1000)
00192     std::cout << "proc " << MyGlobals::_Rank << " : <-- RecvIntVec " << size << std::endl;
00193   std::vector<int> *vec=new std::vector<int>;
00194   vec->resize(size);
00195   MPI_Recv(&vec[0], size, MPI_INT, source, tag+100, MPI_COMM_WORLD, &status);
00196 #endif
00197   return vec;
00198 }
00199 
00200 void MEDPARTITIONER::RecvIntVec(std::vector<int>& vec, const int source)
00201 {
00202   int tag = 111003;
00203   int size;
00204 #ifdef HAVE_MPI2
00205   MPI_Status status;  
00206   MPI_Recv(&size, 1, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00207   if (MyGlobals::_Verbose>1000)
00208     std::cout << "proc " << MyGlobals::_Rank << " : <-- RecvIntVec " << size << std::endl;
00209   vec.resize(size);
00210   MPI_Recv(&vec[0], size, MPI_INT, source, tag+100, MPI_COMM_WORLD,&status);
00211 #endif
00212 }
00213 
00220 void MEDPARTITIONER::SendDataArrayInt(const ParaMEDMEM::DataArrayInt *da, const int target)
00221 {
00222   if (da==0)
00223     throw INTERP_KERNEL::Exception("Problem send DataArrayInt* NULL");
00224   int tag = 111004;
00225   int size[3];
00226   size[0]=da->getNbOfElems();
00227   size[1]=da->getNumberOfTuples();
00228   size[2]=da->getNumberOfComponents();
00229   if (MyGlobals::_Verbose>1000) 
00230     std::cout << "proc " << MyGlobals::_Rank << " : --> SendDataArrayInt " << size[0] << std::endl;
00231 #ifdef HAVE_MPI2
00232   MPI_Send(&size, 3, MPI_INT, target, tag, MPI_COMM_WORLD);
00233   const int *p=da->getConstPointer();
00234   MPI_Send(const_cast<int*>(&p[0]), size[0] ,MPI_INT, target, tag+100, MPI_COMM_WORLD);
00235 #endif
00236 }
00237 
00243 ParaMEDMEM::DataArrayInt *MEDPARTITIONER::RecvDataArrayInt(const int source)
00244 {
00245   int tag = 111004;
00246   int size[3];
00247 #ifdef HAVE_MPI2
00248   MPI_Status status;
00249   MPI_Recv(size, 3, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00250   if (MyGlobals::_Verbose>1000)
00251     std::cout << "proc " << MyGlobals::_Rank << " : <-- RecvDataArrayInt " << size[0] << std::endl;
00252   if (size[0]!=(size[1]*size[2]))
00253     throw INTERP_KERNEL::Exception("Problem in RecvDataArrayInt incoherent sizes");
00254   ParaMEDMEM::DataArrayInt* da=ParaMEDMEM::DataArrayInt::New();
00255   da->alloc(size[1],size[2]);
00256   int *p=da->getPointer();
00257   MPI_Recv(const_cast<int*>(&p[0]), size[0], MPI_INT, source, tag+100, MPI_COMM_WORLD, &status);
00258 #endif
00259   return da;
00260 }
00261 
00268 void MEDPARTITIONER::SendDataArrayDouble(const ParaMEDMEM::DataArrayDouble *da, const int target)
00269 {
00270   if (da==0)
00271     throw INTERP_KERNEL::Exception("Problem send DataArrayDouble* NULL");
00272   int tag = 111005;
00273   int size[3];
00274   size[0]=da->getNbOfElems();
00275   size[1]=da->getNumberOfTuples();
00276   size[2]=da->getNumberOfComponents();
00277   if (MyGlobals::_Verbose>1000) 
00278     std::cout << "proc " << MyGlobals::_Rank << " : --> SendDataArrayDouble " << size[0] << std::endl;
00279 #ifdef HAVE_MPI2
00280   MPI_Send(&size, 3, MPI_INT, target, tag, MPI_COMM_WORLD);
00281   const double *p=da->getConstPointer();
00282   MPI_Send(const_cast<double*>(&p[0]), size[0] ,MPI_DOUBLE, target, tag+100, MPI_COMM_WORLD);
00283 #endif
00284 }
00285 
00291 ParaMEDMEM::DataArrayDouble* MEDPARTITIONER::RecvDataArrayDouble(const int source)
00292 {
00293   int tag = 111005;
00294   int size[3];
00295 #ifdef HAVE_MPI2
00296   MPI_Status status;
00297   MPI_Recv(size, 3, MPI_INT, source, tag, MPI_COMM_WORLD, &status);
00298   if (MyGlobals::_Verbose>1000)
00299     std::cout << "proc " << MyGlobals::_Rank << " : <-- RecvDataArrayDouble " << size[0] << std::endl;
00300   if (size[0]!=(size[1]*size[2]))
00301     throw INTERP_KERNEL::Exception("Problem in RecvDataArrayDouble incoherent sizes");
00302   ParaMEDMEM::DataArrayDouble* da=ParaMEDMEM::DataArrayDouble::New();
00303   da->alloc(size[1],size[2]);
00304   double *p=da->getPointer();
00305   MPI_Recv(const_cast<double*>(&p[0]), size[0], MPI_DOUBLE, source, tag+100, MPI_COMM_WORLD, &status);
00306 #endif
00307   return da;
00308 }
00309 
00310 void MEDPARTITIONER::TestVectorOfStringMpi()
00311 {
00312   int rank=MyGlobals::_Rank;
00313   int world_size=MyGlobals::_World_Size;
00314   std::vector<std::string> myVector;
00315   std::ostringstream oss;
00316   oss << "hello from " << std::setw(5) << rank << " " << std::string(rank+1,'n') <<
00317     " next is an empty one";
00318   myVector.push_back(oss.str());
00319   myVector.push_back("");
00320   myVector.push_back("next is an singleton");
00321   myVector.push_back("1");
00322   
00323   if (rank==0)
00324     {
00325       std::string s0=SerializeFromVectorOfString(myVector);
00326       std::vector<std::string> res=DeserializeToVectorOfString(s0);
00327       if (res.size()!=myVector.size()) 
00328         throw INTERP_KERNEL::Exception("Problem in (de)serialise VectorOfString incoherent sizes");
00329       for (std::size_t i=0; i<myVector.size(); i++)
00330         if (res[i]!=myVector[i])
00331           throw INTERP_KERNEL::Exception("Problem in (de)serialise VectorOfString incoherent elements");
00332     }
00333 
00334   for (int i=0; i<world_size; i++)
00335     {
00336       for (int j=0; j<world_size; j++)
00337         {
00338           std::vector<std::string> res=SendAndReceiveVectorOfString(myVector, i, j);
00339           if ((rank==j) && MyGlobals::_Verbose>20)
00340             std::cout << "proc " << rank << " : receive \n" << ReprVectorOfString(res) << std::endl;
00341           if (rank==j)
00342             {
00343               if (res.size()!=myVector.size()) 
00344                 throw INTERP_KERNEL::Exception("Problem in SendAndReceiveVectorOfString incoherent sizes");
00345               for (std::size_t ii=1; ii<myVector.size(); ii++) //first is different
00346                 if (res[i]!=myVector[ii])
00347                   throw INTERP_KERNEL::Exception("Problem in SendAndReceiveVectorOfString incoherent elements");
00348             }
00349           else 
00350             {
00351               if (res.size()!=0) 
00352                 throw INTERP_KERNEL::Exception("Problem in SendAndReceiveVectorOfString size have to be 0");
00353             }
00354         }
00355     }
00356   std::vector<std::string> res=AllgathervVectorOfString(myVector);
00357   //sometimes for test
00358   res=AllgathervVectorOfString(myVector);
00359   res=AllgathervVectorOfString(myVector);
00360   if (rank==0 && MyGlobals::_Verbose>20)
00361     std::cout << "proc " << rank << " : receive \n" << ReprVectorOfString(res) << std::endl;
00362   if (res.size()!=myVector.size()*world_size) 
00363     throw INTERP_KERNEL::Exception("Problem in AllgathervVectorOfString incoherent sizes");
00364   int jj=-1;
00365   for (int j=0; j<world_size; j++)
00366     {
00367       for (int i=0; i<(int)myVector.size(); i++)
00368         {
00369           jj=jj+1;
00370           if (i==0)
00371             continue; //first is different
00372           if (res[jj]!=myVector[i])
00373             throw INTERP_KERNEL::Exception("Problem in AllgathervVectorOfString incoherent elements");
00374         }
00375     }
00376   if (MyGlobals::_Verbose)
00377     std::cout << "proc " << rank << " : OK TestVectorOfStringMpi END" << std::endl;
00378 }
00379 
00380 void MEDPARTITIONER::TestMapOfStringIntMpi()
00381 {
00382   int rank=MyGlobals::_Rank;
00383   std::map<std::string,int> myMap;
00384   myMap["one"]=1;
00385   myMap["two"]=22;  //a bug
00386   myMap["three"]=3;
00387   myMap["two"]=2; //last speaking override
00388   
00389   if (rank==0)
00390     {
00391       std::vector<std::string> v2=VectorizeFromMapOfStringInt(myMap);
00392       std::map<std::string,int> m3=DevectorizeToMapOfStringInt(v2);
00393       if (ReprMapOfStringInt(m3)!=ReprMapOfStringInt(myMap))
00394         throw INTERP_KERNEL::Exception("Problem in (de)vectorize MapOfStringInt");
00395     }
00396     
00397   std::vector<std::string> v2=AllgathervVectorOfString(VectorizeFromMapOfStringInt(myMap));
00398   if (rank==0 && MyGlobals::_Verbose>20)
00399     {
00400       std::cout << "v2 is : a vector of size " << v2.size() << std::endl;
00401       std::cout << ReprVectorOfString(v2) << std::endl;
00402       std::map<std::string,int> m2=DevectorizeToMapOfStringInt(v2);
00403       std::cout << "m2 is : a map of size " << m2.size() << std::endl;
00404       std::cout << ReprMapOfStringInt(m2) << std::endl;
00405     }
00406   if (MyGlobals::_Verbose)
00407     std::cout << "proc " << rank << " : OK TestMapOfStringIntMpi END" << std::endl;
00408 }
00409 
00410 void MEDPARTITIONER::TestMapOfStringVectorOfStringMpi()
00411 {
00412   int rank=MyGlobals::_Rank;
00413   std::vector<std::string> myVector;
00414   std::ostringstream oss;
00415   oss << "hello from " << std::setw(5) << MyGlobals::_Rank << " " << std::string(rank+1,'n') << " next is an empty one";
00416   myVector.push_back(oss.str());
00417   myVector.push_back("");
00418   myVector.push_back("next is an singleton");
00419   myVector.push_back("1");
00420   
00421   if (rank==0)
00422     {
00423       std::map< std::string,std::vector<std::string> > m2;
00424       m2["first key"]=myVector;
00425       m2["second key"]=myVector;
00426       std::vector<std::string> v2=VectorizeFromMapOfStringVectorOfString(m2);
00427       std::map< std::string,std::vector<std::string> > m3=DevectorizeToMapOfStringVectorOfString(v2);
00428       if (rank==0 && MyGlobals::_Verbose>20)
00429         std::cout << "m2 is : a MapOfStringVectorOfString of size " << m2.size() << std::endl;
00430       std::cout << ReprMapOfStringVectorOfString(m2) << std::endl;
00431       std::cout << "v2 is : a vector of size " << v2.size() << std::endl;
00432       std::cout << ReprVectorOfString(v2) << std::endl;
00433       std::cout << "m3 is : a map of size "<<m3.size() << std::endl;
00434       std::cout << ReprMapOfStringVectorOfString(m3) << std::endl;
00435       if (ReprMapOfStringVectorOfString(m3)!=ReprMapOfStringVectorOfString(m2))
00436         throw INTERP_KERNEL::Exception("Problem in (de)vectorize MapOfStringVectorOfString");
00437     }
00438     
00439   std::map< std::string,std::vector<std::string> > m4;
00440   m4["1rst key"]=myVector;
00441   m4["2snd key"]=myVector;
00442   std::vector<std::string> v4=AllgathervVectorOfString(VectorizeFromMapOfStringVectorOfString(m4));
00443   if (rank==0 && MyGlobals::_Verbose>20)
00444     {
00445       std::map< std::string,std::vector<std::string> > m5=DevectorizeToMapOfStringVectorOfString(v4);
00446       std::map< std::string,std::vector<std::string> > m6=DeleteDuplicatesInMapOfStringVectorOfString(m5);
00447       std::cout<< "m5 is : a map of size "<<m5.size() << std::endl;
00448       std::cout<< ReprMapOfStringVectorOfString(m5) << std::endl;
00449       std::cout<< "m6 is : a map from m5 with deleteDuplicates of size " << m6.size() << std::endl;
00450       std::cout<< ReprMapOfStringVectorOfString(m6) << std::endl;
00451     }
00452   if (MyGlobals::_Verbose)
00453     std::cout<<"proc " << rank << " : OK TestMapOfStringVectorOfStringMpi END" << std::endl;
00454 }
00455 
00456 void MEDPARTITIONER::TestDataArrayMpi()
00457 {
00458   int rank=MyGlobals::_Rank;
00459   //int
00460   {
00461     ParaMEDMEM::DataArrayInt* send=ParaMEDMEM::DataArrayInt::New();
00462     ParaMEDMEM::DataArrayInt* recv=0;
00463     int nbOfTuples=5;
00464     int numberOfComponents=3;
00465     send->alloc(nbOfTuples,numberOfComponents);
00466     std::vector<int> vals;
00467     for (int j=0; j<nbOfTuples; j++)
00468       for (int i=0; i<numberOfComponents; i++) vals.push_back((j+1)*10+i+1);
00469     std::copy(vals.begin(),vals.end(),send->getPointer());
00470     if (rank==0)
00471       SendDataArrayInt(send, 1);
00472     if (rank==1)
00473       recv=RecvDataArrayInt(0);
00474     if (rank==1 && MyGlobals::_Verbose>20)
00475       {
00476         std::cout << send->repr() << std::endl;
00477         std::cout << recv->repr() << std::endl;
00478       }
00479     if (rank==1)
00480       {
00481         if (send->repr()!=recv->repr())
00482           throw INTERP_KERNEL::Exception("Problem in send&recv DataArrayInt");
00483       }
00484     send->decrRef();
00485     if (rank==1)
00486       recv->decrRef();
00487   }
00488   //double
00489   {
00490     ParaMEDMEM::DataArrayDouble* send=ParaMEDMEM::DataArrayDouble::New();
00491     ParaMEDMEM::DataArrayDouble* recv=0;
00492     int nbOfTuples=5;
00493     int numberOfComponents=3;
00494     send->alloc(nbOfTuples,numberOfComponents);
00495     std::vector<double> vals;
00496     for (int j=0; j<nbOfTuples; j++)
00497       for (int i=0; i<numberOfComponents; i++) vals.push_back(double(j+1)+double(i+1)/10);
00498     std::copy(vals.begin(),vals.end(),send->getPointer());
00499     if (rank==0) SendDataArrayDouble(send, 1);
00500     if (rank==1) recv=RecvDataArrayDouble(0);
00501     if (rank==1 && MyGlobals::_Verbose>20)
00502       {
00503         std::cout << send->repr() << std::endl;
00504         std::cout << recv->repr() << std::endl;
00505       }
00506     if (rank==1)
00507       {
00508         if (send->repr()!=recv->repr())
00509           throw INTERP_KERNEL::Exception("Problem in send&recv DataArrayDouble");
00510       }
00511     send->decrRef();
00512     if (rank==1) recv->decrRef();
00513   }
00514   
00515   if (MyGlobals::_Verbose)
00516     std::cout << "proc " << rank << " : OK TestDataArrayMpi END" << std::endl;
00517 }
00518 
00519 void MEDPARTITIONER::TestPersistantMpi0To1(int taille, int nb)
00520 {
00521   double temps_debut=MPI_Wtime();
00522   int rank=MyGlobals::_Rank;
00523   std::vector<int> x, y;
00524   int tag=111111;
00525   MPI_Request requete0, requete1;
00526   MPI_Status statut;
00527   int ok=0;
00528   std::string res;
00529   if (rank==0)
00530     {
00531       x.resize(taille);
00532       MPI_Ssend_init(&x[0], taille, MPI_INT, 1, tag, MPI_COMM_WORLD , &requete0);
00533       for(int k=0; k<nb; k++)
00534         {
00535           for (int i=0; i<taille; ++i) x[i]=k;
00536           //Envoi d’un gros message --> cela peut prendre du temps
00537           MPI_Start(&requete0);
00538           //Traitement sequentiel independant de "x"
00539           //...
00540           MPI_Wait(&requete0, &statut);
00541           //Traitement sequentiel impliquant une modification de "x" en memoire
00542           //x=...
00543         }
00544       MPI_Request_free(&requete0);
00545     }
00546   else if (rank == 1)
00547     {
00548       y.resize(taille);
00549       MPI_Recv_init(&y[0], taille,  MPI_INT, 0, tag, MPI_COMM_WORLD , &requete1);
00550       for(int k=0; k<nb; k++)
00551         {
00552           //Pre-traitement sequentiel
00553           //...
00554           for (int i=0; i<taille; ++i) y[i]=(-1);
00555           //Reception du gros message --> cela peut prendre du temps
00556           MPI_Start(&requete1);
00557           //Traitement sequentiel independant de "y"
00558           //...
00559           MPI_Wait(&requete1, &statut);
00560           //Traitement sequentiel dependant de "y"
00561           //...=f(y)
00562           int nbb=0;
00563           for (int i=0; i<taille; ++i)
00564             if (y[i]==k)
00565               nbb++;
00566           if (nbb==taille)
00567             ok++;
00568           if (MyGlobals::_Verbose>9)
00569             {
00570               res="0K";
00571               if (nbb!=taille)
00572                 res="KO";
00573               std::cout << res << k << " ";
00574             }
00575         }
00576       res="0K";
00577       if (ok!=nb)
00578         res="BAD";
00579       if (MyGlobals::_Verbose>1) 
00580         std::cout << "result " << res << " time(sec) " << MPI_Wtime()-temps_debut << std::endl;
00581       MPI_Request_free(&requete1);
00582     }
00583   //end_time=(MPI_WTIME()-start_time);
00584 }
00585 
00586 void MEDPARTITIONER::TestPersistantMpiRing(int taille, int nb)
00587 {
00588   double temps_debut=MPI_Wtime();
00589   int befo, next, rank, wsize, tagbefo, tagnext;
00590   rank=MyGlobals::_Rank;
00591   wsize=MyGlobals::_World_Size;
00592   befo=rank-1; if (befo<0) befo=wsize-1;
00593   next=rank+1; if (next>=wsize) next=0;
00594   std::vector<int> x, y;
00595   tagbefo=111111+befo;
00596   tagnext=111111+rank;
00597   MPI_Request requete0, requete1;
00598   MPI_Status statut1, statut2;
00599   int ok=0;
00600   std::string res;
00601   //cout<<"ini|"<<rank<<'|'<<befo<<'|'<<next<<' ';
00602   {
00603     x.resize(taille);
00604     y.resize(taille);
00605     MPI_Ssend_init(&x[0], taille, MPI_INT, next, tagnext, MPI_COMM_WORLD , &requete0);
00606     MPI_Recv_init(&y[0], taille,  MPI_INT, befo, tagbefo, MPI_COMM_WORLD , &requete1);
00607     for(int k=0; k<nb; k++)
00608       {
00609         for (int i=0; i<taille; ++i) x[i]=k+rank;
00610         //Envoi d’un gros message --> cela peut prendre du temps
00611         MPI_Start(&requete0);
00612         //Reception du gros message --> cela peut prendre du temps
00613         for (int i=0; i<taille; ++i) y[i]=(-1);
00614         MPI_Start(&requete1);
00615         //Traitement sequentiel independant de "x"
00616         //...
00617         //Traitement sequentiel independant de "y"
00618         //...
00619         MPI_Wait(&requete1, &statut1);
00620         //Traitement sequentiel dependant de "y"
00621         //...=f(y)
00622         int nbb=0;
00623         for (int i=0; i<taille; ++i)
00624           if (y[i]==k+befo)
00625             nbb++;
00626         if (nbb==taille)
00627           ok++;
00628         if (MyGlobals::_Verbose>9)
00629           {
00630             res="0K"+IntToStr(rank);
00631             if (nbb!=taille)
00632               res="KO"+IntToStr(rank);
00633             std::cout << res << k << " ";
00634           }
00635         MPI_Wait(&requete0, &statut2);
00636         //Traitement sequentiel impliquant une modification de "x" en memoire
00637         //x=...
00638       }
00639     res="0K"; if (ok!=nb) res="MAUVAIS";
00640     temps_debut=MPI_Wtime()-temps_debut;
00641     MPI_Request_free(&requete1);
00642     MPI_Request_free(&requete0);
00643   }
00644   //end_time=(MPI_WTIME()-start_time);
00645   if (MyGlobals::_Verbose>1) 
00646     std::cout << "result on proc " << rank << " " << res << " time(sec) " << temps_debut << std::endl;
00647 }
00648 
00649 void MEDPARTITIONER::TestPersistantMpiRingOnCommSplit(int size, int nb)
00650 {
00651   double temps_debut=MPI_Wtime();
00652   int rank=MyGlobals::_Rank;
00653   MPI_Comm newcomm;
00654   int color=1;
00655   int rankMax=4;
00656   if (rank>=rankMax)
00657     color=MPI_UNDEFINED;
00658   //MPI_Comm_dup (MPI_COMM_WORLD, &newcomm) ;
00659   MPI_Comm_split(MPI_COMM_WORLD, color, rank, &newcomm);
00660   
00661   int befo, next, wsize, tagbefo, tagnext;
00662   wsize=rankMax;
00663   if (wsize>MyGlobals::_World_Size)
00664     wsize=MyGlobals::_World_Size;
00665   befo=rank-1;
00666   if (befo<0)
00667     befo=wsize-1;
00668   next=rank+1;
00669   if (next>=wsize)
00670     next=0;
00671   std::vector<int> x, y;
00672   tagbefo=111111+befo;
00673   tagnext=111111+rank;
00674   MPI_Request requete0, requete1;
00675   MPI_Status statut1, statut2;
00676   int ok=0;
00677   std::string res;
00678   
00679   if (color==1)
00680     {
00681       x.resize(size);
00682       y.resize(size);
00683       MPI_Ssend_init(&x[0], size, MPI_INT, next, tagnext, newcomm , &requete0);
00684       MPI_Recv_init(&y[0], size,  MPI_INT, befo, tagbefo, newcomm , &requete1);
00685       for(int k=0; k<nb; k++)
00686         {
00687           for (int i=0; i<size; ++i)
00688             x[i]=k+rank;
00689           //Send of big message --> time consuming
00690           MPI_Start(&requete0);
00691           //Reception of big message --> time consuming
00692           for (int i=0; i<size; ++i)
00693             y[i]=-1;
00694           MPI_Start(&requete1);
00695           //Traitement sequentiel independant de "x"
00696           //...
00697           //Traitement sequentiel independant de "y"
00698           //...
00699           //cout<<"dsr|"<<rank<<' ';
00700           MPI_Wait(&requete1, &statut1);
00701           //Traitement sequentiel dependant de "y"
00702           //...=f(y)
00703           int nbb=0;
00704           for (int i=0; i<size; ++i)
00705             if (y[i]==k+befo)
00706               nbb++;
00707           if (nbb==size)
00708             ok++;
00709           if (MyGlobals::_Verbose>9)
00710             {
00711               res="0K"+IntToStr(rank);
00712               if (nbb!=size)
00713                 res="KO"+IntToStr(rank);
00714               std::cout << res << k << " ";
00715             }
00716           MPI_Wait(&requete0, &statut2);
00717           //Traitement sequentiel impliquant une modification de "x" en memoire
00718           //x=...
00719         }
00720       res="0K";
00721       if (ok!=nb)
00722         res="MAUVAIS";
00723       temps_debut=MPI_Wtime()-temps_debut;
00724       MPI_Request_free(&requete1);
00725       MPI_Request_free(&requete0);
00726     }
00727   //MPI_Barrier(MPI_COMM_WORLD);
00728   if (color==1)
00729     MPI_Comm_free(&newcomm);
00730   if (MyGlobals::_Verbose>1) 
00731     std::cout << "resultat proc " << rank <<" " << res << " time(sec) " << temps_debut << std::endl;
00732 }