Back to index

salome-med  6.5.0
medpartitioner_para.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 /*
00021   examples of launch
00022   rm ttmp* tttmp*
00023   export verb=11
00024   mpirun -np 2 medpartitioner_para --input-file=blade.med --output-file=ttmp1_ --ndomains=2 --dump-cpu-memory --verbose=$verb
00025   mpirun -np 5 medpartitioner_para --input-file=blade.med --output-file=ttmp1_ --ndomains=2 --dump-cpu-memory --verbose=$verb
00026   mpirun -np 4 medpartitioner_para --input-file=ttmp1_.xml --output-file=tttmp1_ --ndomains=4 --verbose=$verb
00027 
00028   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=tmp_testMesh_20x30x50.med --output-file=ttmp2_ --verbose=$verb
00029   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=ttmp2_.xml --output-file=ttmp3_ --verbose=$verb
00030 
00031   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=tmp_testMeshWithFaces_20x30x50.med --output-file=ttmp2_ --verbose=$verb
00032   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=ttmp2_.xml --output-file=ttmp3_ --verbose=$verb
00033 
00034   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=tmp_testMesh_20x30x50_WithVecFieldOnCells.med --output-file=ttmp2_ --verbose=$verb
00035   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=tmp_testMesh_20x30x50_WithVecFieldOnNodes.med --output-file=ttmp2_ --verbose=$verb
00036   mpirun -np 2 medpartitioner_para --ndomains=2 --input-file=ttmp2_.xml --output-file=ttmp3_ --verbose=$verb
00037 
00038   mpirun -np 4 medpartitioner_para --ndomains=4 --input-file=tmp_testMeshHuge_20x30x50_6.xml --output-file=ttmp3_ --verbose=1
00039 */
00040 
00041 
00042 #include "MEDPARTITIONER_MeshCollection.hxx"
00043 #include "MEDPARTITIONER_ParallelTopology.hxx"
00044 #include "MEDPARTITIONER_ParaDomainSelector.hxx"
00045 #include "MEDPARTITIONER_Utils.hxx"
00046 
00047 #include "MEDLoader.hxx"
00048 
00049 #include <fstream>
00050 #include <iostream>
00051 #include <iomanip>
00052 #include <sstream>
00053 #include <string>
00054 
00055 #ifdef HAVE_MPI2
00056 #include <mpi.h>
00057 #endif
00058 
00059 using namespace std;
00060 using namespace MEDPARTITIONER;
00061 
00062 int main(int argc, char** argv)
00063 {
00064 #if !defined(MED_ENABLE_PARMETIS)
00065   cout << "Sorry, no one split method is available. Please, compile with ParMETIS."<<endl;
00066   return 1;
00067 #else
00068 
00069   // Defining options
00070   // by parsing the command line
00071   
00072   bool split_family=false;
00073   bool empty_groups=false;
00074   bool mesure_memory=false;
00075   bool filter_face=true;
00076 
00077   string input;
00078   string output;
00079   string meshname;
00080   string library="metis";  //default
00081   int ndomains;
00082   int help=0;
00083   int test=0;
00084   MPI_Init(&argc,&argv);
00085   MPI_Comm_size(MPI_COMM_WORLD, &MyGlobals::_World_Size);
00086   MPI_Comm_rank(MPI_COMM_WORLD, &MyGlobals::_Rank);
00087   
00088   //cout<<"proc "<<MyGlobals::_Rank<<" of "<<MyGlobals::_World_Size<<endl; //for debug
00089   //TestVectorOfStringMpi(); //cvw
00090   //TestRandomize();
00091   
00092   //Primitive parsing of command-line options
00093   string desc ("Available options of medpartitioner_para V1.0.3:\n"
00094                "\t--help                   : produces this help message\n"
00095                "\t--verbose                : echoes arguments\n"
00096                "\t--input-file=<string>    : name of the input .med file or .xml master file\n"
00097                "\t--output-file=<string>   : name of the resulting file (without extension)\n"
00098                "\t--ndomains=<number>      : number of subdomains in the output file, default is 1\n"
00099 #if defined(MED_ENABLE_PARMETIS) 
00100 // || defined(MED_ENABLE_PTSCOTCH)
00101 //user can choose! (not yet)
00102                "\t--split-method=<string>  : name of the splitting library (metis/scotch), default is metis\n"
00103 #endif
00104                "\t--creates-boundary-faces : creates boundary faces mesh in the output files\n"
00105                "\t--dump-cpu-memory        : dumps passed CPU time and maximal increase of used memory\n"
00106                //"\t--randomize=<number>     : random seed for other partitionning (only on one proc)\n"
00107                //"\t--atomize                : do the opposite of a good partitionner (only on one proc)\n"
00108                );
00109 
00110   if (argc<=1) help=1;
00111   string value;
00112   for (int i = 1; i < argc; i++)
00113     {
00114       if (strlen(argv[i]) < 3) 
00115         {
00116           if (MyGlobals::_Rank==0) cerr << "bad argument : "<< argv[i] << endl;
00117           MPI_Finalize(); return 1;
00118         }
00119     
00120       if (TestArg(argv[i],"--verbose",value)) 
00121         {
00122           MyGlobals::_Verbose=1;
00123           if (value!="") MyGlobals::_Verbose = atoi(value.c_str());
00124         }
00125       else if (TestArg(argv[i],"--help",value)) help=1;
00126       else if (TestArg(argv[i],"--test",value)) test=1;
00127       else if (TestArg(argv[i],"--input-file",value)) input=value;
00128       else if (TestArg(argv[i],"--output-file",value)) output=value;
00129       else if (TestArg(argv[i],"--split-method",value)) library=value;
00130       else if (TestArg(argv[i],"--ndomains",value)) ndomains=atoi(value.c_str());
00131       else if (TestArg(argv[i],"--randomize",value)) MyGlobals::_Randomize=atoi(value.c_str());
00132       else if (TestArg(argv[i],"--atomize",value)) MyGlobals::_Atomize=atoi(value.c_str());
00133       else if (TestArg(argv[i],"--creates-boundary-faces",value)) MyGlobals::_Creates_Boundary_Faces=1;
00134       else if (TestArg(argv[i],"--dump-cpu-memory",value)) mesure_memory=true;
00135       else 
00136         {
00137           if (MyGlobals::_Rank==0) cerr << "unknown argument : "<< argv[i] << endl;
00138           MPI_Finalize(); return 1;
00139         }
00140     }
00141 
00142 
00143   if (MyGlobals::_Randomize!=0 && MyGlobals::_World_Size!=1)
00144     {
00145       if (MyGlobals::_Rank==0) cerr << "randomize only available on 1 proc (mpirun -np 1)" << endl;
00146       MyGlobals::_Randomize=0;
00147     }
00148 
00149 //no choice
00150 #if defined(MED_ENABLE_PARMETIS) && !defined(MED_ENABLE_SCOTCH)
00151   library = "metis";
00152 #endif
00153 #if !defined(MED_ENABLE_PARMETIS) && defined(MED_ENABLE_SCOTCH)
00154   library = "scotch";
00155 #endif
00156 //user choice
00157 #if defined(MED_ENABLE_PARMETIS) && defined(MED_ENABLE_SCOTCH)
00158   if ((library!="metis") && (library!="scotch"))
00159     {
00160       if (MyGlobals::_Rank==0) cerr << "split-method only available : metis, scotch" << endl;
00161       MPI_Finalize(); return 1;
00162     }
00163 #endif
00164  
00165   if (help==1)
00166     {
00167       if (MyGlobals::_Rank==0) cout<<desc<<"\n";
00168       MPI_Finalize(); return 0;
00169     }
00170   
00171   MyGlobals::_Is0verbose=0;
00172   if (MyGlobals::_Rank==0) MyGlobals::_Is0verbose=MyGlobals::_Verbose;
00173   //MyGlobals::_Is0verbose=((MyGlobals::_Rank==0) && MyGlobals::_Verbose);
00174   if (test==1) //only for debugger
00175     {
00176       if (MyGlobals::_Is0verbose>0) cout<<"tests on "<<MyGlobals::_Atomize<<" "<<ndomains<<endl;
00177       //TestPersistantMpi0To1(MyGlobals::_Atomize, ndomains);
00178       //TestPersistantMpiRing(MyGlobals::_Atomize, ndomains);
00179       TestPersistantMpiRingOnCommSplit(MyGlobals::_Atomize, ndomains);
00180       //MPI_Barrier(MPI_COMM_WORLD);
00181       MPI_Finalize();
00182       return 0;
00183       TestVectorOfStringMpi();
00184       TestMapOfStringIntMpi();
00185       TestMapOfStringVectorOfStringMpi();
00186       TestDataArrayMpi();
00187       MPI_Finalize();
00188       return 1;
00189     }
00190   
00191   if (MyGlobals::_Is0verbose)
00192     {
00193       cout << "medpartitioner_para V1.0 :" << endl;
00194       cout << "  input-file = " << input << endl;
00195       cout << "  output-file = " << output << endl;
00196       cout << "  split-method = " << library << endl;
00197       cout << "  ndomains = " << ndomains << endl;
00198       cout << "  creates_boundary_faces = " << MyGlobals::_Creates_Boundary_Faces << endl;
00199       cout << "  dump-cpu-memory = " << mesure_memory<< endl;
00200       cout << "  verbose = " << MyGlobals::_Verbose << endl;
00201     }
00202   //testing whether it is possible to write a file at the specified location
00203   if (MyGlobals::_Rank==0)
00204     {
00205       string outputtest = output + ".testioms.";
00206       ofstream testfile (outputtest.c_str());
00207       if (testfile.fail())
00208         { 
00209           cerr << "output-file directory does not exist or is in read-only access" << endl;
00210           MPI_Finalize(); return 1;
00211         }
00212       //deletes test file
00213       remove(outputtest.c_str());
00214     }
00215   
00216   // Beginning of the computation
00217 
00218   // Loading the mesh collection
00219   if (MyGlobals::_Is0verbose) cout << "Reading input files "<<endl;
00220   
00221   try
00222     {
00223       MEDPARTITIONER::ParaDomainSelector parallelizer(mesure_memory);
00224       MEDPARTITIONER::MeshCollection collection(input,parallelizer);
00225       MEDPARTITIONER::ParallelTopology* aPT = (MEDPARTITIONER::ParallelTopology*) collection.getTopology();
00226       aPT->setGlobalNumerotationDefault(collection.getParaDomainSelector());
00227       //int nbfiles=MyGlobals::_fileMedNames->size(); //nb domains
00228       //to have unique valid fields names/pointers/descriptions for partitionning
00229       collection.prepareFieldDescriptions();
00230       //int nbfields=collection.getFieldDescriptions().size(); //on all domains
00231       //cout<<ReprVectorOfString(collection.getFieldDescriptions());
00232     
00233       if (MyGlobals::_Is0verbose)
00234         {
00235           cout<<"fileNames :"<<endl
00236               <<ReprVectorOfString(MyGlobals::_File_Names);
00237           cout<<"fieldDescriptions :"<<endl
00238               <<ReprFieldDescriptions(collection.getFieldDescriptions()," ");
00239           cout<<"familyInfo :\n"
00240               <<ReprMapOfStringInt(collection.getFamilyInfo())<<endl;
00241           cout<<"groupInfo :\n"
00242               <<ReprMapOfStringVectorOfString(collection.getGroupInfo())<<endl;
00243         }
00244     
00245       //Creating the graph and partitioning it
00246       if (MyGlobals::_Is0verbose) cout << "Computing partition "<< endl;
00247   
00248       auto_ptr< MEDPARTITIONER::Topology > new_topo;
00249       if (library == "metis")
00250         new_topo.reset( collection.createPartition(ndomains,MEDPARTITIONER::Graph::METIS));
00251       else
00252         new_topo.reset( collection.createPartition(ndomains,MEDPARTITIONER::Graph::SCOTCH));
00253       parallelizer.evaluateMemory();
00254   
00255       //Creating a new mesh collection from the partitioning
00256       if (MyGlobals::_Is0verbose)  cout << "Creating new meshes"<< endl;
00257       MEDPARTITIONER::MeshCollection new_collection(collection,new_topo.get(),split_family,empty_groups);
00258       //cout<<"proc "<<MyGlobals::_Rank<<" : new_collection done"<<endl;
00259       parallelizer.evaluateMemory();
00260   
00261       //if (!xml_output_master) new_collection.setDriverType(MEDPARTITIONER::MedAscii);
00262       if (filter_face) new_collection.filterFaceOnCell();
00263     
00264       //to get infos on all procs
00265     
00266       //see meshName
00267       vector<string> finalInformations;
00268       vector<string> r1,r2;
00269       r1=AllgathervVectorOfString(MyGlobals::_General_Informations);
00270       //if (MyGlobals::_Is0verbose>1000) cout << "generalInformations : \n"<<ReprVectorOfString(r1);
00271       r2=SelectTagsInVectorOfString(r1,"ioldDomain=");
00272       r2=SelectTagsInVectorOfString(r2,"meshName=");
00273       if (r2.size()==(collection.getMesh()).size())
00274         {
00275           for (std::size_t i=0; i<r2.size(); i++)
00276             r2[i]=EraseTagSerialized(r2[i],"ioldDomain=");
00277           r2=DeleteDuplicatesInVectorOfString(r2);
00278           if (r2.size()==1)
00279             {
00280               string finalMesh="finalMeshName="+ExtractFromDescription(r2[0], "meshName=");
00281               finalInformations.push_back(SerializeFromString(finalMesh));
00282             }
00283         }
00284       if (finalInformations.size()==0)
00285         {
00286           if (MyGlobals::_Rank==0)
00287             cerr<<"Problem on final meshName : set at 'Merge'"<<endl;
00288           finalInformations.push_back(SerializeFromString("finalMeshName=Merge"));
00289         }
00290     
00291       //see field info nbComponents & componentInfo (if fields present)
00292       r2=SelectTagsInVectorOfString(r1,"fieldName=");
00293       r2=SelectTagsInVectorOfString(r2,"nbComponents=");
00294       //may be yes? or not?
00295       for (std::size_t i=0; i<r2.size(); i++)
00296         r2[i]=EraseTagSerialized(r2[i],"ioldFieldDouble=");
00297       r2=DeleteDuplicatesInVectorOfString(r2);
00298       for (std::size_t i=0; i<r2.size(); i++)
00299         finalInformations.push_back(r2[i]);
00300     
00301       MyGlobals::_General_Informations=finalInformations;
00302       if (MyGlobals::_Is0verbose) 
00303         cout << "generalInformations : \n"<<ReprVectorOfString(finalInformations);
00304     
00305       //new_collection.setSubdomainBoundaryCreates(creates_boundary_faces);
00306       if (MyGlobals::_Is0verbose) cout << "Writing "<<ndomains<<" output files "<<output<<"xx.med"<<" and "<<output<<".xml"<<endl;
00307       new_collection.write(output);
00308   
00309       if ( mesure_memory )
00310         if ( parallelizer.isOnDifferentHosts() || MyGlobals::_Rank==0 )
00311           {
00312             cout << "Elapsed time = " << parallelizer.getPassedTime()
00313                  << ", max memory usage = " << parallelizer.evaluateMemory() << " KB"
00314                  << endl;
00315           }
00316       if(MyGlobals::_World_Size>1)
00317         MPI_Barrier(MPI_COMM_WORLD);
00318       if (MyGlobals::_Is0verbose>0) cout<<"OK END"<< endl;
00319       MPI_Finalize();
00320       return 0;
00321     }
00322   catch(const char *mess)
00323     {
00324       cerr<<"proc "<<MyGlobals::_Rank<<" : "<<mess<<endl;
00325       fflush(stderr);
00326       MPI_Finalize();
00327       return 1;
00328     }
00329   catch(INTERP_KERNEL::Exception& e)
00330     {
00331       cerr<<"proc "<<MyGlobals::_Rank<<" : INTERP_KERNEL_Exception : "<<e.what()<<endl;
00332       fflush(stderr);
00333       MPI_Finalize();
00334       return 1;
00335     }
00336   catch(std::exception& e)
00337     {
00338       cerr<<"proc "<<MyGlobals::_Rank<<" : std_Exception : "<<e.what()<<endl;
00339       fflush(stderr);
00340       MPI_Finalize();
00341       return 1;
00342     }
00343   catch(...)
00344     {
00345       cerr<<"proc "<<MyGlobals::_Rank<<" : an unknown type exception error was occured"<<endl;
00346       fflush(stderr);
00347       MPI_Finalize();
00348       return 1;
00349     }
00350 #endif
00351 }
00352