IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
Ifpack_OverlappingPartitioner.cpp
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
43 #include "Ifpack_ConfigDefs.h"
44 #include "Ifpack_Partitioner.h"
45 #include "Ifpack_OverlappingPartitioner.h"
46 #include "Ifpack_Graph.h"
47 
48 #include "Epetra_Comm.h"
49 #include "Epetra_BlockMap.h"
50 #include "Epetra_Map.h"
51 #include "Teuchos_ParameterList.hpp"
52 
53 static const std::string PrintMsg_ = "(Ifpack_OvPartitioner) ";
54 
55 //==============================================================================
58  NumLocalParts_(1),
59  Graph_(Graph),
60  OverlappingLevel_(0),
61  IsComputed_(false),
62  verbose_(false)
63 {
64 }
65 
66 //==============================================================================
68 {
69 }
70 
71 //==============================================================================
72 int Ifpack_OverlappingPartitioner::SetParameters(Teuchos::ParameterList& List)
73 {
74 
75  NumLocalParts_ = List.get("partitioner: local parts", NumLocalParts_);
76  OverlappingLevel_ = List.get("partitioner: overlap", OverlappingLevel_);
77  verbose_ = List.get("partitioner: print level", verbose_);
78 
79  if (NumLocalParts_ < 0)
81  if (NumLocalParts_ == 0)
82  NumLocalParts_ = 1;
83  if (NumLocalParts_ < 0)
84  IFPACK_CHK_ERR(-1);
86  IFPACK_CHK_ERR(-1);
87 
88  if (OverlappingLevel_ < 0)
89  IFPACK_CHK_ERR(-1);
90 
92 
93  return(0);
94 }
95 
96 //==============================================================================
98 {
99  using std::cout;
100  using std::endl;
101 
102  if (NumLocalParts_ < 1)
103  IFPACK_CHK_ERR(-1); // incorrect value
104 
105  if (OverlappingLevel_ < 0)
106  IFPACK_CHK_ERR(-1); // incorrect value
107 
108  // some output
109 
110  if (verbose_ && (Comm().MyPID() == 0)) {
111  cout << PrintMsg_ << "Number of local parts = " << NumLocalParts_ << endl;
112  cout << PrintMsg_ << "Number of global parts = "
113  << NumLocalParts_ * Comm().NumProc() << endl;
114  cout << PrintMsg_ << "Amount of overlap = " << OverlappingLevel_ << endl;
115  }
116 
117  // 1.- allocate memory
118 
119  Partition_.resize(NumMyRows());
120  Parts_.resize(NumLocalParts());
121 
122  // 2.- sanity checks on input graph
123 
124  if (Graph_->Filled() == false)
125  IFPACK_CHK_ERR(-4); // need FillComplete() called
126 
127  if (Graph_->NumGlobalRows64() != Graph_->NumGlobalCols64())
128  IFPACK_CHK_ERR(-3); // can partition square matrices only
129 
130  if (NumLocalParts_ < 1)
131  IFPACK_CHK_ERR(-2); // value not valid
132 
133  // 3.- perform non-overlapping partition
134 
135  IFPACK_CHK_ERR(ComputePartitions());
136 
137  // 4.- compute the partitions with overlapping
138 
139  IFPACK_CHK_ERR(ComputeOverlappingPartitions());
140 
141  // 5.- return to the user
142 
143  IsComputed_ = true;
144 
145  return(0);
146 }
147 
148 // ======================================================================
150 {
151  using std::cerr;
152  using std::endl;
153 
154  // FIXME: the first part of this function should be elsewhere
155  // start defining the subgraphs for no overlap
156 
157  std::vector<int> sizes;
158  sizes.resize(NumLocalParts_);
159 
160  // 1.- compute how many rows are in each subgraph
161  for (int i = 0 ; i < NumLocalParts_ ; ++i)
162  sizes[i] = 0;
163 
164  for (int i = 0 ; i < NumMyRows() ; ++i) {
165  if (Partition_[i] >= NumLocalParts_) {
166  cerr << "ERROR: Partition[" << i << "] = "<< Partition_[i]
167  << ", NumLocalParts = " << NumLocalParts_ << endl;
168  cerr << "(file = " << __FILE__ << ", line = "
169  << __LINE__ << ")" << endl;
170  IFPACK_CHK_ERR(-10);
171  }
172  // no singletons should be here, as the matrix is
173  // supposed to be filtered through Ifpack_SingletonFilter
174  if (Partition_[i] == -1)
175  IFPACK_CHK_ERR(-1);
176  sizes[Partition_[i]]++;
177  }
178 
179  // 2.- allocate space for each subgraph
180  for (int i = 0 ; i < NumLocalParts_ ; ++i)
181  Parts_[i].resize(sizes[i]);
182 
183  // 3.- cycle over all rows and populate the vectors
184 
185  for (int i = 0 ; i < NumLocalParts_ ; ++i)
186  sizes[i] = 0;
187 
188  for (int i = 0 ; i < NumMyRows() ; ++i) {
189  int part = Partition_[i];
190  int count = sizes[part];
191  Parts_[part][count] = i;
192  sizes[part]++;
193  }
194 
195  if (OverlappingLevel_ == 0)
196  return(0);
197 
198  // wider overlap requires further computations
199  for (int level = 1 ; level <= OverlappingLevel_ ; ++level) {
200 
201  std::vector<std::vector<int> > tmp;
202  tmp.resize(NumLocalParts_);
203 
204  // cycle over all rows in the local graph (that is the overlapping
205  // graph). For each row, all columns will belong to the subgraph of
206  // row `i'.
207 
208  int MaxNumEntries_tmp = Graph_->MaxMyNumEntries();
209  std::vector<int> Indices;
210  Indices.resize(MaxNumEntries_tmp);
211 
212  for (int part = 0 ; part < NumLocalParts_ ; ++part) {
213 
214  for (int i = 0; i < (int)Parts_[part].size() ; ++i) {
215 
216  int LRID = Parts_[part][i];
217  int NumIndices;
218  int ierr = Graph_->ExtractMyRowCopy(LRID, MaxNumEntries_tmp,
219  NumIndices, &Indices[0]);
220  IFPACK_CHK_ERR(ierr);
221 
222  for (int j = 0 ; j < NumIndices ; ++j) {
223 
224  // use *local* indices
225  int col = Indices[j];
226  if (col >= NumMyRows())
227  continue;
228 
229  // has this column already been inserted?
230  std::vector<int>::iterator
231  where = find(tmp[part].begin(), tmp[part].end(), col);
232 
233  if (where == tmp[part].end()) {
234  tmp[part].push_back(col);
235  }
236  }
237  }
238  }
239 
240  // now I convert the STL vectors into Epetra_IntSerialDenseVectors.
241  for (int i = 0 ; i < NumLocalParts_ ; ++i) {
242  Parts_[i].resize(tmp[i].size());
243  for (int j = 0 ; j < (int)tmp[i].size() ; ++j)
244  Parts_[i][j] = tmp[i][j];
245  }
246  }
247 
248  return(0);
249 
250 }
251 
252 //============================================================================
254 {
255  return(Graph_->NumMyRows());
256 }
257 
258 //============================================================================
260 {
261  return(Graph_->NumMyNonzeros());
262 }
263 
264 //============================================================================
265 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
267 {
268  return(Graph_->NumGlobalRows());
269 }
270 #endif
271 
272 long long Ifpack_OverlappingPartitioner::NumGlobalRows64() const
273 {
274  return(Graph_->NumGlobalRows64());
275 }
276 //============================================================================
278 {
279  return(Graph_->MaxMyNumEntries());
280 }
281 
282 //============================================================================
284 {
285  return(Graph_->Comm());
286 }
287 
288 // ======================================================================
289 std::ostream& Ifpack_OverlappingPartitioner::Print(std::ostream & os) const
290 {
291  using std::endl;
292 
293  if (Comm().MyPID())
294  return(os);
295 
296  os << "================================================================================" << endl;
297  os << "Ifpack_OverlappingPartitioner" << endl;
298  os << "Number of local rows = " << Graph_->NumMyRows() << endl;
299  os << "Number of global rows = " << Graph_->NumGlobalRows64() << endl;
300  os << "Number of local parts = " << NumLocalParts_ << endl;
301  os << "Overlapping level = " << OverlappingLevel_ << endl;
302  os << "Is computed = " << IsComputed_ << endl;
303  os << "================================================================================" << endl;
304 
305  return(os);
306 }
307 
virtual int SetParameters(Teuchos::ParameterList &List)
Sets all the parameters for the partitioner.
virtual int Compute()
Computes the partitions. Returns 0 if successful.
virtual int ExtractMyRowCopy(int MyRow, int LenOfIndices, int &NumIndices, int *Indices) const =0
Extracts a copy of input local row.
virtual int NumMyRows() const =0
Returns the number of local rows.
virtual int SetPartitionParameters(Teuchos::ParameterList &List)=0
Sets all the parameters for the partitioner.
int MaxNumEntries() const
Returns the max number of local entries in a row.
Ifpack_OverlappingPartitioner(const Ifpack_Graph *Graph)
Constructor.
int NumMyNonzeros() const
Returns the number of local nonzero elements.
int NumMyRows() const
Returns the number of local rows.
std::vector< int > Partition_
Partition_[i] contains the ID of non-overlapping part it belongs to.
virtual int ComputeOverlappingPartitions()
Computes the partitions. Returns 0 if successful.
virtual int NumMyNonzeros() const =0
Returns the number of local nonzero entries.
virtual bool Filled() const =0
Returns true is graph is filled.
bool IsComputed_
If true, the graph has been successfully partitioned.
virtual int NumGlobalRows() const =0
Returns the number of global rows.
bool verbose_
If true, information are reported on cout.
virtual const Epetra_Comm & Comm() const =0
Returns the communicator object of the graph.
const Ifpack_Graph * Graph_
Reference to the graph to be partitioned.
int NumLocalParts_
Number of local subgraphs.
virtual int NumProc() const =0
virtual int MaxMyNumEntries() const =0
Returns the maximun number of entries for row.
virtual int ComputePartitions()=0
Computes the partitions. Returns 0 if successful.
int NumLocalParts() const
Returns the number of computed local partitions.
Ifpack_Graph: a pure virtual class that defines graphs for IFPACK.
Definition: Ifpack_Graph.h:61
virtual std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator&lt;&lt;.
const Epetra_Comm & Comm() const
Returns the communicator object of Graph.
std::vector< std::vector< int > > Parts_
Parts_[i][j] is the ID of the j-th row contained in the (overlapping)
int NumGlobalRows() const
Returns the number of global rows.