Ifpack_OverlappingPartitioner.cpp

Go to the documentation of this file.
00001 #include "Ifpack_ConfigDefs.h"
00002 #ifdef HAVE_IFPACK_TEUCHOS
00003 #include "Ifpack_Partitioner.h"
00004 #include "Ifpack_OverlappingPartitioner.h"
00005 #include "Ifpack_Graph.h"
00006 
00007 #include "Epetra_Comm.h"
00008 #include "Epetra_BlockMap.h"
00009 #include "Epetra_Map.h"
00010 #include "Teuchos_ParameterList.hpp"
00011 
00012 static const string PrintMsg_ = "(Ifpack_OvPartitioner) ";
00013 
00014 //==============================================================================
00015 Ifpack_OverlappingPartitioner::
00016 Ifpack_OverlappingPartitioner(const Ifpack_Graph* Graph) :
00017   NumLocalParts_(1),
00018   Graph_(Graph),
00019   OverlappingLevel_(0),
00020   IsComputed_(false),
00021   verbose_(false)
00022 {
00023 }
00024 
00025 //==============================================================================
00026 Ifpack_OverlappingPartitioner::~Ifpack_OverlappingPartitioner()
00027 {
00028 }
00029 
00030 //==============================================================================
00031 int Ifpack_OverlappingPartitioner::SetParameters(Teuchos::ParameterList& List)
00032 {
00033 
00034   NumLocalParts_ = List.get("partitioner: local parts", NumLocalParts_);
00035   OverlappingLevel_ = List.get("partitioner: overlap", OverlappingLevel_);
00036   verbose_ = List.get("partitioner: print level", verbose_);
00037 
00038   if (NumLocalParts_ < 0)
00039     NumLocalParts_ = Graph_->NumMyRows() / (-NumLocalParts_);
00040   if (NumLocalParts_ == 0)
00041     NumLocalParts_ = 1;
00042   if (NumLocalParts_ < 0)
00043     IFPACK_CHK_ERR(-1);
00044   if (NumLocalParts_ > Graph_->NumMyRows())
00045     IFPACK_CHK_ERR(-1);
00046 
00047   if (OverlappingLevel_ < 0)
00048     IFPACK_CHK_ERR(-1);
00049 
00050   SetPartitionParameters(List);
00051 
00052   return(0);
00053 }
00054 
00055 //==============================================================================
00056 int Ifpack_OverlappingPartitioner::Compute()
00057 {
00058 
00059   if (NumLocalParts_ < 1)
00060     IFPACK_CHK_ERR(-1); // incorrect value
00061 
00062   if (OverlappingLevel_ < 0)
00063     IFPACK_CHK_ERR(-1); // incorrect value
00064 
00065   // some output
00066 
00067   if (verbose_ && (Comm().MyPID() == 0)) {
00068     cout << PrintMsg_ << "Number of local parts  = " << NumLocalParts_ << endl;
00069     cout << PrintMsg_ << "Number of global parts = " 
00070          << NumLocalParts_ * Comm().NumProc() << endl;
00071     cout << PrintMsg_ << "Amount of overlap      = " << OverlappingLevel_ << endl;
00072   }
00073 
00074   // 1.- allocate memory 
00075 
00076   Partition_.resize(NumMyRows());
00077   Parts_.resize(NumLocalParts());
00078 
00079   // 2.- sanity checks on input graph
00080  
00081   if (Graph_->Filled() == false)
00082     IFPACK_CHK_ERR(-4); // need FillComplete() called
00083 
00084   if (Graph_->NumGlobalRows() != Graph_->NumGlobalCols())
00085     IFPACK_CHK_ERR(-3); // can partition square matrices only
00086 
00087   if (NumLocalParts_ < 1)
00088     IFPACK_CHK_ERR(-2); // value not valid
00089  
00090   // 3.- perform non-overlapping partition
00091  
00092   IFPACK_CHK_ERR(ComputePartitions());
00093 
00094   // 4.- compute the partitions with overlapping
00095   
00096   IFPACK_CHK_ERR(ComputeOverlappingPartitions());
00097 
00098   // 5.- return to the user
00099  
00100   IsComputed_ = true;
00101 
00102   return(0);
00103 }
00104 
00105 // ======================================================================
00106 int Ifpack_OverlappingPartitioner::ComputeOverlappingPartitions()
00107 {
00108 
00109   // FIXME: the first part of this function should be elsewhere
00110   // start defining the subgraphs for no overlap
00111 
00112   vector<int> sizes;
00113   sizes.resize(NumLocalParts_);
00114 
00115   // 1.- compute how many rows are in each subgraph
00116   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00117     sizes[i] = 0;
00118 
00119   for (int i = 0 ; i < NumMyRows() ; ++i) {
00120     if (Partition_[i] >= NumLocalParts_) {
00121       cerr << "ERROR: Partition[" << i << "] = "<< Partition_[i] 
00122      << ", NumLocalParts = " << NumLocalParts_ << endl;
00123       cerr << "(file = " << __FILE__ << ", line = "
00124            << __LINE__ << ")" << endl;
00125       IFPACK_CHK_ERR(-10);
00126     }
00127     // no singletons should be here, as the matrix is
00128     // supposed to be filtered through Ifpack_SingletonFilter
00129     if (Partition_[i] == -1)
00130       IFPACK_CHK_ERR(-1);
00131     sizes[Partition_[i]]++;
00132   }
00133 
00134   // 2.- allocate space for each subgraph
00135   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00136     Parts_[i].resize(sizes[i]);
00137 
00138   // 3.- cycle over all rows and populate the vectors
00139 
00140   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00141     sizes[i] = 0;
00142 
00143   for (int i = 0 ; i < NumMyRows() ; ++i) {
00144     int part = Partition_[i];
00145     int count = sizes[part];
00146     Parts_[part][count] = i;
00147     sizes[part]++;
00148   }
00149 
00150   if (OverlappingLevel_ == 0)
00151     return(0);
00152 
00153   // wider overlap requires further computations
00154   for (int level = 1 ; level <= OverlappingLevel_ ; ++level) {
00155 
00156     vector<vector<int> > tmp;
00157     tmp.resize(NumLocalParts_);
00158 
00159     // cycle over all rows in the local graph (that is the overlapping
00160     // graph). For each row, all columns will belong to the subgraph of
00161     // row `i'.
00162 
00163     int MaxNumEntries = Graph_->MaxMyNumEntries();
00164     vector<int> Indices;
00165     Indices.resize(MaxNumEntries);
00166 
00167     for (int part = 0 ; part < NumLocalParts_ ; ++part) {
00168 
00169       for (int i = 0; i < (int)Parts_[part].size() ; ++i) {  
00170 
00171   int LRID = Parts_[part][i];
00172   int NumIndices;
00173   int ierr = Graph_->ExtractMyRowCopy(LRID, MaxNumEntries, 
00174                                             NumIndices, &Indices[0]);
00175   IFPACK_CHK_ERR(ierr);
00176 
00177   for (int j = 0 ; j < NumIndices ; ++j) {
00178 
00179     // use *local* indices
00180     int col = Indices[j];
00181           if (col >= NumMyRows())
00182             continue;
00183 
00184     // has this column already been inserted?
00185     vector<int>::iterator
00186       where = find(tmp[part].begin(), tmp[part].end(), col);
00187 
00188     if (where == tmp[part].end()) {
00189       tmp[part].push_back(col);
00190     }
00191   }
00192       }
00193     }
00194 
00195     // now I convert the STL vectors into Epetra_IntSerialDenseVectors.
00196     for (int i = 0 ; i < NumLocalParts_ ; ++i) {
00197       Parts_[i].resize(tmp[i].size());
00198       for (int j = 0 ; j < (int)tmp[i].size() ; ++j)
00199   Parts_[i][j] = tmp[i][j];
00200     }
00201   }
00202 
00203   return(0);
00204 
00205 }
00206 
00207 //============================================================================
00208 int Ifpack_OverlappingPartitioner::NumMyRows() const
00209 {
00210   return(Graph_->NumMyRows());
00211 }
00212 
00213 //============================================================================
00214 int Ifpack_OverlappingPartitioner::NumMyNonzeros() const
00215 {
00216   return(Graph_->NumMyNonzeros());
00217 }
00218 
00219 //============================================================================
00220 int Ifpack_OverlappingPartitioner::NumGlobalRows() const
00221 {
00222   return(Graph_->NumGlobalRows());
00223 }
00224 
00225 //============================================================================
00226 int Ifpack_OverlappingPartitioner::MaxNumEntries() const
00227 {
00228   return(Graph_->MaxMyNumEntries());
00229 }
00230 
00231 //============================================================================
00232 const Epetra_Comm& Ifpack_OverlappingPartitioner::Comm() const
00233 {
00234   return(Graph_->Comm());
00235 }
00236 
00237 // ======================================================================
00238 ostream& Ifpack_OverlappingPartitioner::Print(ostream & os) const
00239 {
00240 
00241   if (Comm().MyPID()) 
00242     return(os);
00243 
00244   os << "================================================================================" << endl;
00245   os << "Ifpack_OverlappingPartitioner" << endl;
00246   os << "Number of local rows  = " << Graph_->NumMyRows() << endl;
00247   os << "Number of global rows = " << Graph_->NumGlobalRows() << endl;
00248   os << "Number of local parts = " << NumLocalParts_ << endl;
00249   os << "Overlapping level     = " << OverlappingLevel_ << endl;
00250   os << "Is computed           = " << IsComputed_ << endl;
00251   os << "================================================================================" << endl;
00252 
00253   return(os);
00254 }
00255 
00256 #endif // HAVE_IFPACK_TEUCHOS

Generated on Thu Sep 18 12:37:22 2008 for Ifpack Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1