Ifpack2 Templated Preconditioning Package Version 1.0
Ifpack2_OverlappingPartitioner_def.hpp
00001 /*@HEADER
00002 // ***********************************************************************
00003 //
00004 //       Ifpack2: Tempated Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2009) 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 #ifndef IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
00031 #define IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
00032 #include "Ifpack2_ConfigDefs.hpp"
00033 #include "Ifpack2_OverlappingPartitioner_decl.hpp"
00034 #include "Teuchos_Array.hpp"
00035 #include "Teuchos_ArrayRCP.hpp"
00036 #include <vector>
00037 #include <string>
00038 
00039 namespace Ifpack2 {
00040 //==============================================================================
00041 template<class GraphType>
00042 OverlappingPartitioner<GraphType>::OverlappingPartitioner(const Teuchos::RCP<const Tpetra::RowGraph<LocalOrdinal,GlobalOrdinal,Node> >& Graph) :
00043   NumLocalParts_(1),
00044   Graph_(Graph),
00045   OverlappingLevel_(0),
00046   IsComputed_(false),
00047   verbose_(false)
00048 {
00049 }
00050 
00051 //==============================================================================
00052 template<class GraphType>
00053 OverlappingPartitioner<GraphType>::~OverlappingPartitioner()
00054 {
00055 }
00056 
00057 //==============================================================================
00058 template<class GraphType>
00059 typename GraphType::local_ordinal_type OverlappingPartitioner<GraphType>::numLocalParts() const 
00060 {
00061   return( NumLocalParts_);
00062 }
00063 
00064 //==============================================================================
00065 template<class GraphType>
00066 size_t OverlappingPartitioner<GraphType>::overlappingLevel() const
00067 {
00068   return(OverlappingLevel_);
00069 }
00070 
00071 //==============================================================================
00072 template<class GraphType>
00073 typename GraphType::local_ordinal_type OverlappingPartitioner<GraphType>::operator() (LocalOrdinal MyRow) const
00074 {
00075   TEUCHOS_TEST_FOR_EXCEPTION(((MyRow < 0) || ((size_t)MyRow > Graph_->getNodeNumRows())), std::runtime_error, "Ifpack2::OverlappingPartitioner::operator() invalid row.");
00076    
00077   return(Partition_[MyRow]);
00078 }
00079 
00080 
00081 //==============================================================================
00082 template<class GraphType>
00083 typename GraphType::local_ordinal_type OverlappingPartitioner<GraphType>::operator() (LocalOrdinal i, LocalOrdinal j) const
00084 {
00085   TEUCHOS_TEST_FOR_EXCEPTION( ((i < 0) || ((i > NumLocalParts_) || (j < 0) || (j > Parts_[i].size()))),
00086             std::runtime_error, "Ifpack2::OverlappingPartitioner::operator() invalid row or node.");
00087   return(Parts_[i][j]);
00088 }
00089 
00090 //==============================================================================
00091 template<class GraphType>
00092 size_t OverlappingPartitioner<GraphType>::numRowsInPart(const LocalOrdinal Part) const
00093 {
00094   TEUCHOS_TEST_FOR_EXCEPTION( ((Part < 0) || (Part > NumLocalParts_)),
00095             std::runtime_error, "Ifpack2::OverlappingPartitioner::numRowsInPart() invalid partition.");
00096   return(Parts_[Part].size());
00097 }
00098 
00099 //==============================================================================
00100 template<class GraphType>
00101 void OverlappingPartitioner<GraphType>::rowsInPart(const LocalOrdinal Part,  Teuchos::ArrayRCP<LocalOrdinal> & List) const
00102 {
00103   // Let numRowsInPart do the sanity checking...
00104   size_t numrows= numRowsInPart(Part); 
00105   for (size_t i = 0 ; i < numrows; i++)
00106     List[i] = Parts_[Part][i];
00107 }
00108 
00109 //==============================================================================
00110 template<class GraphType>
00111 Teuchos::ArrayView<const typename GraphType::local_ordinal_type>  OverlappingPartitioner<GraphType>::nonOverlappingPartition() const
00112 {
00113   return(Partition_.view(0, Graph_->getNodeNumRows()));
00114 }
00115 
00116 //==============================================================================
00117 template<class GraphType>
00118 void OverlappingPartitioner<GraphType>::setParameters(Teuchos::ParameterList& List)
00119 {
00120   NumLocalParts_    = List.get("partitioner: local parts", NumLocalParts_);
00121   OverlappingLevel_ = List.get("partitioner: overlap", OverlappingLevel_);
00122   verbose_          = List.get("partitioner: print level", verbose_);
00123 
00124   if (NumLocalParts_ < 0)
00125     NumLocalParts_ = Graph_->getNodeNumRows() / (-NumLocalParts_);
00126   if (NumLocalParts_ == 0)
00127     NumLocalParts_ = 1;
00128   
00129   // Sanity checking
00130   TEUCHOS_TEST_FOR_EXCEPTION( ((NumLocalParts_ < 0) || (size_t)NumLocalParts_ > Graph_->getNodeNumRows()),
00131             std::runtime_error, "Ifpack2::OverlappingPartitioner::setParameters() invalid NumLocalParts_");
00132   TEUCHOS_TEST_FOR_EXCEPTION( (OverlappingLevel_ < 0),
00133             std::runtime_error, "Ifpack2::OverlappingPartitioner::setParameters() invalid OverlappingLevel_");
00134 
00135   setPartitionParameters(List);
00136 }
00137 
00138 //==============================================================================
00139 template<class GraphType>
00140 void OverlappingPartitioner<GraphType>::compute()
00141 {
00142   TEUCHOS_TEST_FOR_EXCEPTION( ((NumLocalParts_ < 1) ||  (OverlappingLevel_ < 0)),
00143             std::runtime_error, "Ifpack2::OverlappingPartitioner::compute() invalid NumLocalParts_ or OverlappingLevel_");
00144 
00145   std::string PrintMsg_("OverlappingPartitioner: ");
00146 
00147   // some output
00148   if (verbose_ && (Graph_->getComm()->getRank() == 0)) {
00149     std::cout << PrintMsg_ << "Number of local parts          = " << NumLocalParts_ << std::endl;
00150     std::cout << PrintMsg_ << "Approx. Number of global parts = " 
00151         << NumLocalParts_ * Graph_->getComm()->getSize() << std::endl;
00152     std::cout << PrintMsg_ << "Amount of overlap              = " << OverlappingLevel_ << std::endl;
00153   }
00154 
00155   // 1.- allocate memory 
00156   Partition_.resize(Graph_->getNodeNumRows());
00157   Parts_.resize(NumLocalParts_);
00158 
00159   // 2.- sanity checks on input graph
00160   TEUCHOS_TEST_FOR_EXCEPTION( (!Graph_->isFillComplete() || Graph_->getGlobalNumRows() != Graph_->getGlobalNumCols()),
00161             std::runtime_error, "Ifpack2::OverlappingPartitioner::compute() input graph error");
00162  
00163   // 3.- perform non-overlapping partition 
00164   computePartitions();
00165 
00166   // 4.- compute the partitions with overlapping
00167   computeOverlappingPartitions();
00168 
00169   // 5.- mark as computed
00170   IsComputed_ = true;
00171 }
00172 
00173 //==============================================================================
00174 template<class GraphType>
00175 void OverlappingPartitioner<GraphType>::computeOverlappingPartitions()
00176 {
00177   using std::vector;
00178 
00179   // Old FIXME from Ifpack: the first part of this function should be elsewhere
00180  
00181   // start defining the subgraphs for no overlap
00182 
00183   vector<size_t> sizes;
00184   sizes.resize(NumLocalParts_);
00185 
00186   // 1.- compute how many rows are in each subgraph
00187   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00188     sizes[i] = 0;
00189 
00190   for (size_t i = 0 ; i < Graph_->getNodeNumRows() ; ++i) {
00191     TEUCHOS_TEST_FOR_EXCEPTION( (Partition_[i] >= NumLocalParts_),
00192         std::runtime_error, "Ifpack2::OverlappingPartitioner::computeOverlappingPartitions() Partition_[i] > NumLocalParts_");
00193 
00194     // no singletons should be here, as the matrix is
00195     // supposed to be filtered through Ifpack2_SingletonFilter
00196     TEUCHOS_TEST_FOR_EXCEPTION((Partition_[i] == -1),
00197              std::runtime_error, "Ifpack2::OverlappingPartitioner::computeOverlappingPartitions() singletons are forbidden");
00198     sizes[Partition_[i]]++;
00199   }
00200 
00201   // 2.- allocate space for each subgraph
00202   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00203     Parts_[i].resize(sizes[i]);
00204 
00205   // 3.- cycle over all rows and populate the vectors
00206   for (int i = 0 ; i < NumLocalParts_ ; ++i)
00207     sizes[i] = 0;
00208 
00209   for (size_t i = 0 ; i < Graph_->getNodeNumRows() ; ++i) {
00210     LocalOrdinal part = Partition_[i];
00211     size_t count = sizes[part];
00212     Parts_[part][count] = i;
00213     sizes[part]++;
00214   }
00215 
00216   // If there is no overlap, we're done, so return
00217   if (OverlappingLevel_ == 0) return;
00218 
00219   // wider overlap requires further computations
00220   for (size_t level = 1 ; level <= OverlappingLevel_ ; ++level) {
00221     vector<vector<size_t> > tmp;
00222     tmp.resize(NumLocalParts_);
00223 
00224     // cycle over all rows in the local graph (that is the overlapping
00225     // graph). For each row, all columns will belong to the subgraph of
00226     // row `i'.
00227 
00228     int MaxNumEntries_tmp = Graph_->getNodeMaxNumRowEntries();
00229     Teuchos::Array<LocalOrdinal> Indices;
00230     Indices.resize(MaxNumEntries_tmp);
00231 
00232     for (int part = 0 ; part < NumLocalParts_ ; ++part) {
00233       for (size_t i = 0; i < (size_t)Parts_[part].size() ; ++i) {  
00234   LocalOrdinal LRID = Parts_[part][i];
00235   size_t NumIndices;
00236   Graph_->getLocalRowCopy(LRID,Indices(),NumIndices);
00237 
00238   for (size_t j = 0 ; j < NumIndices ; ++j) {
00239     // use *local* indices only
00240     LocalOrdinal col = Indices[j];
00241           if ((size_t)col >= Graph_->getNodeNumRows())
00242             continue;
00243 
00244     // has this column already been inserted?
00245     vector<size_t>::iterator
00246       where = std::find(tmp[part].begin(), tmp[part].end(), (size_t) col);
00247 
00248     if (where == tmp[part].end()) {
00249       tmp[part].push_back(col);
00250     }
00251   }
00252       }
00253     }
00254 
00255     // now I convert the STL vectors into Teuchos Array RCP's
00256     for (int i = 0 ; i < NumLocalParts_ ; ++i) {
00257       Parts_[i].resize(tmp[i].size());
00258       for (size_t j = 0 ; j < tmp[i].size() ; ++j)
00259   Parts_[i][j] = tmp[i][j];
00260     }
00261   }
00262 }
00263 
00264 //==============================================================================
00265 template<class GraphType>
00266 bool OverlappingPartitioner<GraphType>::isComputed() const
00267 {
00268   return(IsComputed_);
00269 }
00270 
00271 //==============================================================================
00272 template<class GraphType>
00273 std::ostream& OverlappingPartitioner<GraphType>::print(std::ostream& os) const
00274 {
00275   Teuchos::FancyOStream fos(Teuchos::rcp(&os,false));
00276   fos.setOutputToRootOnly(0);
00277   describe(fos);
00278   return(os);
00279 }
00280 
00281 //==============================================================================
00282 template<class GraphType>
00283 std::string OverlappingPartitioner<GraphType>::description() const
00284 {
00285   std::ostringstream oss;
00286   oss << Teuchos::Describable::description();
00287   if (isComputed()) {
00288     oss << "{status = computed";
00289   }
00290   else {
00291     oss << "{status = is not computed";
00292   }
00293   oss <<"}";
00294   return oss.str();
00295 }
00296 
00297 //==============================================================================
00298 template<class GraphType>
00299 void  OverlappingPartitioner<GraphType>::describe(Teuchos::FancyOStream &os, const Teuchos::EVerbosityLevel verbLevel) const
00300 {
00301   using std::endl;
00302   if(verbLevel==Teuchos::VERB_NONE) return;
00303 
00304   os << "================================================================================" << endl;
00305   os << "Ifpack2::OverlappingPartitioner" << endl;
00306   os << "Number of local rows  = " << Graph_->getNodeNumRows() << endl;
00307   os << "Number of global rows = " << Graph_->getGlobalNumRows() << endl;
00308   os << "Number of local parts = " << NumLocalParts_ << endl;
00309   os << "Overlapping level     = " << OverlappingLevel_ << endl;
00310   os << "Is computed           = " << IsComputed_ << endl;
00311   os << "================================================================================" << endl;
00312 }
00313 
00314 
00315 }// namespace Ifpack2
00316 
00317 #endif // IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends