Back to index

salome-med  6.5.0
MEDPARTITIONER_ParMetisGraph.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_MetisGraph.hxx"
00021 #include "MEDPARTITIONER_ParaDomainSelector.hxx"
00022 #include "MEDPARTITIONER_Utils.hxx"
00023 
00024 #include "InterpKernelException.hxx"
00025 
00026 #include <iostream>
00027 
00028 #ifdef MED_ENABLE_PARMETIS
00029 #include <mpi.h>
00030 #include "parmetis.h"
00031 #endif
00032 
00033 using namespace MEDPARTITIONER;
00034 
00035 METISGraph::METISGraph():Graph()
00036 {
00037 }
00038 
00039 METISGraph::METISGraph(MEDPARTITIONER::SkyLineArray* graph, int* edgeweight)
00040   :Graph(graph,edgeweight)
00041 {
00042 }
00043 
00044 METISGraph::~METISGraph()
00045 {
00046 }
00047 
00048 void METISGraph::partGraph(int ndomain,
00049                            const std::string& options_string,
00050                            ParaDomainSelector *parallelizer)
00051 {
00052   using std::vector;
00053   vector<int> ran,vx,va; //for randomize
00054   
00055   if (MyGlobals::_Verbose>10)
00056     std::cout << "proc " << MyGlobals::_Rank << " : METISGraph::partGraph" << std::endl;
00057   
00058   // number of graph vertices
00059   int n=_graph->getNumberOf();
00060   //graph
00061   int * xadj=const_cast<int*>(_graph->getIndex());
00062   int * adjncy=const_cast<int*>(_graph->getValue());
00063   //constraints
00064   int * vwgt=_cell_weight;
00065   int * adjwgt=_edge_weight;
00066   int wgtflag=(_edge_weight!=0)?1:0+(_cell_weight!=0)?2:0;
00067   //base 0 or 1
00068   int base=0;
00069   //ndomain
00070   int nparts=ndomain;
00071   //options
00072   /*
00073     (0=default_option,option,random_seed) see defs.h
00074     #define PMV3_OPTION_DBGLVL 1
00075     #define PMV3_OPTION_SEED 2
00076     #define PMV3_OPTION_IPART 3
00077     #define PMV3_OPTION_PSR 3
00078     seems no changes int options[4]={1,0,33,0}; //test for a random seed of 33
00079   */
00080   int options[4]={0,0,0,0};
00081   // output parameters
00082   int edgecut;
00083 #if !defined(MED_ENABLE_PARMETIS)
00084   throw INTERP_KERNEL::Exception("METISGraph::partGraph : PARMETIS is not available. Check your products, please.");
00085 #else
00086   int* partition=new int[n];
00087   
00088   if (MyGlobals::_Verbose>10) 
00089     std::cout << "proc " << MyGlobals::_Rank << " : METISGraph::partGraph ParMETIS_PartKway new" << std::endl;
00090   int * vtxdist=parallelizer->getProcVtxdist();
00091   MPI_Comm comm=MPI_COMM_WORLD;
00092   ParMETIS_PartKway(vtxdist, xadj, adjncy, vwgt, 
00093                                     adjwgt, &wgtflag, &base, &nparts, options, 
00094                                     &edgecut, partition, &comm );
00095 
00096 
00097   /*doc from parmetis.h
00098     void __cdecl ParMETIS_PartKway(
00099     idxtype *vtxdist, idxtype *xadj, idxtype *adjncy, idxtype *vwgt, 
00100     idxtype *adjwgt, int *wgtflag, int *numflag, int *nparts, int *options, 
00101     int *edgecut, idxtype *part, MPI_Comm *comm);
00102 
00103     void __cdecl ParMETIS_V3_PartKway(
00104     idxtype *vtxdist, idxtype *xadj, idxtype *adjncy, idxtype *vwgt, 
00105     idxtype *adjwgt, int *wgtflag, int *numflag, int *ncon, int *nparts, 
00106     float *tpwgts, float *ubvec, int *options, int *edgecut, idxtype *part, 
00107     MPI_Comm *comm);
00108   */
00109 
00110   vector<int> index(n+1);
00111   vector<int> value(n);
00112   index[0]=0;
00113   if (ran.size()>0 && MyGlobals::_Atomize==0) //there is randomize
00114     {
00115       if (MyGlobals::_Is0verbose>100)
00116         std::cout << "randomize" << std::endl;
00117       for (int i=0; i<n; i++)
00118         {
00119           index[i+1]=index[i]+1;
00120           value[ran[i]]=partition[i];
00121         }
00122     }
00123   else
00124     {
00125       for (int i=0; i<n; i++)
00126         {
00127           index[i+1]=index[i]+1;
00128           value[i]=partition[i];
00129         }
00130     }
00131   delete [] partition;
00132 
00133   //creating a skylinearray with no copy of the index and partition array
00134   //the fifth argument true specifies that only the pointers are passed 
00135   //to the object
00136   
00137   _partition = new MEDPARTITIONER::SkyLineArray(index,value);
00138 #endif
00139 }
00140