Ifpack_OverlappingPartitioner.h

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 #ifndef IFPACK_OVERLAPPINGPARTITIONER_H
00031 #define IFPACK_OVERLAPPINGPARTITIONER_H
00032 
00033 #include "Ifpack_ConfigDefs.h"
00034 #include "Ifpack_Partitioner.h"
00035 #include "Teuchos_ParameterList.hpp"
00036 class Epetra_Comm;
00037 class Ifpack_Graph;
00038 class Epetra_Map;
00039 class Epetra_BlockMap;
00040 class Epetra_Import;
00041 
00042 /* \brief Ifpack_OverlappingPartitioner: A class to create overlapping
00043     partitions of a local graph.
00044 
00045 Class Ifpack_OverlappingPartitioner enables the extension of 
00046 non-overlapping partitions to an arbitrary value of overlap.
00047 Note that overlap refers to the overlap among \e local parts,
00048 and not the overlap among the processes.
00049 
00050 Supported parameters are:
00051 - \c "partitioner: local parts": the required number of parts;
00052 - \c "partitioner: overlap": the required amount of overlap is set in 
00053   parameter. Default = 0 (integer).
00054 - \c "partitioner: verbose": if \c true, information are reported on
00055   cout. Nothing is reported otherwise.
00056 
00057 This class is a semi-virtual class, that contains the basic utilities
00058 for derived classes Ifpack_LinearPartitioner, Ifpack_GreedyPartitioner,
00059 Ifpack_METISPartitioner, and Ifpack_EquationPartitioner. Graphs in
00060 input to one of these classes are supposed to contain no singletons.
00061 Usually, this means that the graph is derived from an Epetra_RowMatrix,
00062 that has been filtered using Ifpack_SingletonFilter.
00063 
00064 \author Marzio Sala, SNL 9214.
00065 
00066 \date Last update: Oct-04.
00067 */  
00068 class Ifpack_OverlappingPartitioner : public Ifpack_Partitioner {
00069 
00070 public:
00071 
00073   Ifpack_OverlappingPartitioner(const Ifpack_Graph* Graph);
00074 
00076   virtual ~Ifpack_OverlappingPartitioner();
00077 
00079   int NumLocalParts() const 
00080   {
00081     return(NumLocalParts_);
00082   }
00083 
00085   int OverlappingLevel() const 
00086   {
00087     return(OverlappingLevel_);
00088   }
00089 
00091 
00098   int operator() (int MyRow) const
00099   {
00100     if ((MyRow < 0) || (MyRow > NumMyRows()))
00101       IFPACK_CHK_ERR(-1); // input value not valid
00102 
00103     return(Partition_[MyRow]);
00104   }
00105 
00107   int operator() (int i, int j) const
00108   {
00109     if ((i < 0) || (i >= NumLocalParts()))
00110       IFPACK_CHK_ERR(-1);
00111 
00112     if ((j < 0) || (j > (int)Parts_[i].size()))
00113       IFPACK_CHK_ERR(-2);
00114 
00115     return(Parts_[i][j]);
00116   }
00117 
00119   inline int NumRowsInPart(const int Part) const
00120   {
00121     return(Parts_[Part].size());
00122   }
00123     
00124   int RowsInPart(const int Part, int* List) const
00125   {
00126     for (int i = 0 ; i < NumRowsInPart(Part) ; ++i)
00127       List[i] = Parts_[Part][i];
00128 
00129     return(0);
00130   }
00131   
00132   const int* NonOverlappingPartition() const
00133   {
00134     return(&Partition_[0]);
00135   }
00136 
00138 
00143   virtual int SetParameters(Teuchos::ParameterList& List);
00144 
00146 
00150   virtual int SetPartitionParameters(Teuchos::ParameterList& List) = 0;
00151 
00153   virtual int Compute();
00154 
00156   virtual int ComputePartitions() = 0;
00157 
00159   virtual int ComputeOverlappingPartitions();
00160   
00162   bool IsComputed()
00163   {
00164     return(IsComputed_);
00165   }
00166 
00168   virtual ostream& Print(std::ostream& os) const;
00169 
00170 protected:
00171    
00173   int NumMyRows() const;
00175   int NumMyNonzeros() const;
00177   int NumGlobalRows() const;
00179   int MaxNumEntries() const;
00181   const Epetra_Comm& Comm() const;
00183   int NumLocalParts_;
00185   std::vector<int> Partition_; 
00187   // partition i
00188   std::vector<std::vector<int> > Parts_;
00190   const Ifpack_Graph* Graph_;
00192   int OverlappingLevel_;
00194   bool IsComputed_;
00196   bool verbose_;
00197 
00198 }; // class Ifpack_Partitioner
00199 
00200 #endif // IFPACK_OVERLAPPINGPARTITIONER_H

Generated on Tue Jul 13 09:27:13 2010 for IFPACK by  doxygen 1.4.7