Ifpack_OverlappingPartitioner.cpp

00001 /*@HEADER
00002 // ***********************************************************************
00003 //
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2002) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00025 //
00026 // ***********************************************************************
00027 //@HEADER
00028 */
00029 
00030 #include "Ifpack_ConfigDefs.h"
00031 #include "Ifpack_Partitioner.h"
00032 #include "Ifpack_OverlappingPartitioner.h"
00033 #include "Ifpack_Graph.h"
00034 
00035 #include "Epetra_Comm.h"
00036 #include "Epetra_BlockMap.h"
00037 #include "Epetra_Map.h"
00038 #include "Teuchos_ParameterList.hpp"
00039 
00040 static const string PrintMsg_ = "(Ifpack_OvPartitioner) ";
00041 
00042 //==============================================================================
00043 Ifpack_OverlappingPartitioner::
00044 Ifpack_OverlappingPartitioner(const Ifpack_Graph* Graph) :
00045   NumLocalParts_(1),
00046   Graph_(Graph),
00047   OverlappingLevel_(0),
00048   IsComputed_(false),
00049   verbose_(false)
00050 {
00051 }
00052 
00053 //==============================================================================
00054 Ifpack_OverlappingPartitioner::~Ifpack_OverlappingPartitioner()
00055 {
00056 }
00057 
00058 //==============================================================================
00059 int Ifpack_OverlappingPartitioner::SetParameters(Teuchos::ParameterList& List)
00060 {
00061 
00062   NumLocalParts_ = List.get("partitioner: local parts", NumLocalParts_);
00063   OverlappingLevel_ = List.get("partitioner: overlap", OverlappingLevel_);
00064   verbose_ = List.get("partitioner: print level", verbose_);
00065 
00066   if (NumLocalParts_ < 0)
00067     NumLocalParts_ = Graph_->NumMyRows() / (-NumLocalParts_);
00068   if (NumLocalParts_ == 0)
00069     NumLocalParts_ = 1;
00070   if (NumLocalParts_ < 0)
00071     IFPACK_CHK_ERR(-1);
00072   if (NumLocalParts_ > Graph_->NumMyRows())
00073     IFPACK_CHK_ERR(-1);
00074 
00075   if (OverlappingLevel_ < 0)
00076     IFPACK_CHK_ERR(-1);
00077 
00078   SetPartitionParameters(List);
00079 
00080   return(0);
00081 }
00082 
00083 //==============================================================================
00084 int Ifpack_OverlappingPartitioner::Compute()
00085 {
00086 
00087   if (NumLocalParts_ < 1)
00088     IFPACK_CHK_ERR(-1); // incorrect value
00089 
00090   if (OverlappingLevel_ < 0)
00091     IFPACK_CHK_ERR(-1); // incorrect value
00092 
00093   // some output
00094 
00095   if (verbose_ && (Comm().MyPID() == 0)) {
00096     cout << PrintMsg_ << "Number of local parts  = " << NumLocalParts_ << endl;
00097     cout << PrintMsg_ << "Number of global parts = " 
00098          << NumLocalParts_ * Comm().NumProc() << endl;
00099     cout << PrintMsg_ << "Amount of overlap      = " << OverlappingLevel_ << endl;
00100   }
00101 
00102   // 1.- allocate memory 
00103 
00104   Partition_.resize(NumMyRows());
00105   Parts_.resize(NumLocalParts());
00106 
00107   // 2.- sanity checks on input graph
00108  
00109   if (Graph_->Filled() == false)
00110     IFPACK_CHK_ERR(-4); // need FillComplete() called
00111 
00112   if (Graph_->NumGlobalRows() != Graph_->NumGlobalCols())
00113     IFPACK_CHK_ERR(-3); // can partition square matrices only
00114 
00115   if (NumLocalParts_ < 1)
00116     IFPACK_CHK_ERR(-2); // value not valid
00117  
00118   // 3.- perform non-overlapping partition
00119  
00120   IFPACK_CHK_ERR(ComputePartitions());
00121 
00122   // 4.- compute the partitions with overlapping
00123   
00124   IFPACK_CHK_ERR(ComputeOverlappingPartitions());
00125 
00126   // 5.- return to the user
00127  
00128   IsComputed_ = true;
00129 
00130   return(0);
00131 }
00132 
00133 // ======================================================================
00134 int Ifpack_OverlappingPartitioner::ComputeOverlappingPartitions()
00135 {
00136 
00137   // FIXME: the first part of this function should be elsewhere
00138   // start defining the subgraphs for no overlap
00139 
00140   vector<int> sizes;
00141   sizes.resize(NumLocalParts_);
00142 
00143   // 1.- compute how many rows are in each subgraph
00144   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00145     sizes[i] = 0;
00146 
00147   for (int i = 0 ; i < NumMyRows() ; ++i) {
00148     if (Partition_[i] >= NumLocalParts_) {
00149       cerr << "ERROR: Partition[" << i << "] = "<< Partition_[i] 
00150        << ", NumLocalParts = " << NumLocalParts_ << endl;
00151       cerr << "(file = " << __FILE__ << ", line = "
00152            << __LINE__ << ")" << endl;
00153       IFPACK_CHK_ERR(-10);
00154     }
00155     // no singletons should be here, as the matrix is
00156     // supposed to be filtered through Ifpack_SingletonFilter
00157     if (Partition_[i] == -1)
00158       IFPACK_CHK_ERR(-1);
00159     sizes[Partition_[i]]++;
00160   }
00161 
00162   // 2.- allocate space for each subgraph
00163   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00164     Parts_[i].resize(sizes[i]);
00165 
00166   // 3.- cycle over all rows and populate the vectors
00167 
00168   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00169     sizes[i] = 0;
00170 
00171   for (int i = 0 ; i < NumMyRows() ; ++i) {
00172     int part = Partition_[i];
00173     int count = sizes[part];
00174     Parts_[part][count] = i;
00175     sizes[part]++;
00176   }
00177 
00178   if (OverlappingLevel_ == 0)
00179     return(0);
00180 
00181   // wider overlap requires further computations
00182   for (int level = 1 ; level <= OverlappingLevel_ ; ++level) {
00183 
00184     vector<vector<int> > tmp;
00185     tmp.resize(NumLocalParts_);
00186 
00187     // cycle over all rows in the local graph (that is the overlapping
00188     // graph). For each row, all columns will belong to the subgraph of
00189     // row `i'.
00190 
00191     int MaxNumEntries_tmp = Graph_->MaxMyNumEntries();
00192     vector<int> Indices;
00193     Indices.resize(MaxNumEntries_tmp);
00194 
00195     for (int part = 0 ; part < NumLocalParts_ ; ++part) {
00196 
00197       for (int i = 0; i < (int)Parts_[part].size() ; ++i) {  
00198 
00199     int LRID = Parts_[part][i];
00200     int NumIndices;
00201     int ierr = Graph_->ExtractMyRowCopy(LRID, MaxNumEntries_tmp, 
00202                                             NumIndices, &Indices[0]);
00203     IFPACK_CHK_ERR(ierr);
00204 
00205     for (int j = 0 ; j < NumIndices ; ++j) {
00206 
00207       // use *local* indices
00208       int col = Indices[j];
00209           if (col >= NumMyRows())
00210             continue;
00211 
00212       // has this column already been inserted?
00213       vector<int>::iterator
00214         where = find(tmp[part].begin(), tmp[part].end(), col);
00215 
00216       if (where == tmp[part].end()) {
00217         tmp[part].push_back(col);
00218       }
00219     }
00220       }
00221     }
00222 
00223     // now I convert the STL vectors into Epetra_IntSerialDenseVectors.
00224     for (int i = 0 ; i < NumLocalParts_ ; ++i) {
00225       Parts_[i].resize(tmp[i].size());
00226       for (int j = 0 ; j < (int)tmp[i].size() ; ++j)
00227     Parts_[i][j] = tmp[i][j];
00228     }
00229   }
00230 
00231   return(0);
00232 
00233 }
00234 
00235 //============================================================================
00236 int Ifpack_OverlappingPartitioner::NumMyRows() const
00237 {
00238   return(Graph_->NumMyRows());
00239 }
00240 
00241 //============================================================================
00242 int Ifpack_OverlappingPartitioner::NumMyNonzeros() const
00243 {
00244   return(Graph_->NumMyNonzeros());
00245 }
00246 
00247 //============================================================================
00248 int Ifpack_OverlappingPartitioner::NumGlobalRows() const
00249 {
00250   return(Graph_->NumGlobalRows());
00251 }
00252 
00253 //============================================================================
00254 int Ifpack_OverlappingPartitioner::MaxNumEntries() const
00255 {
00256   return(Graph_->MaxMyNumEntries());
00257 }
00258 
00259 //============================================================================
00260 const Epetra_Comm& Ifpack_OverlappingPartitioner::Comm() const
00261 {
00262   return(Graph_->Comm());
00263 }
00264 
00265 // ======================================================================
00266 ostream& Ifpack_OverlappingPartitioner::Print(ostream & os) const
00267 {
00268 
00269   if (Comm().MyPID()) 
00270     return(os);
00271 
00272   os << "================================================================================" << endl;
00273   os << "Ifpack_OverlappingPartitioner" << endl;
00274   os << "Number of local rows  = " << Graph_->NumMyRows() << endl;
00275   os << "Number of global rows = " << Graph_->NumGlobalRows() << endl;
00276   os << "Number of local parts = " << NumLocalParts_ << endl;
00277   os << "Overlapping level     = " << OverlappingLevel_ << endl;
00278   os << "Is computed           = " << IsComputed_ << endl;
00279   os << "================================================================================" << endl;
00280 
00281   return(os);
00282 }
00283 

Generated on Wed May 12 21:30:18 2010 for IFPACK by  doxygen 1.4.7