Ifpack2 Templated Preconditioning Package  Version 1.0
Ifpack2_OverlappingPartitioner_def.hpp
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack2: Tempated Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2009) 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 #ifndef IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
44 #define IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
45 #include "Ifpack2_ConfigDefs.hpp"
46 #include "Ifpack2_OverlappingPartitioner_decl.hpp"
47 #include "Teuchos_Array.hpp"
48 #include "Teuchos_ArrayRCP.hpp"
49 #include <vector>
50 #include <string>
51 
52 namespace Ifpack2 {
53 
54 template<class GraphType>
57  NumLocalParts_ (1),
58  Graph_ (graph),
59  OverlappingLevel_ (0),
60  IsComputed_ (false),
61  verbose_ (false)
62 {}
63 
64 
65 template<class GraphType>
67 
68 
69 template<class GraphType>
70 int
72 {
73  return NumLocalParts_;
74 }
75 
76 
77 template<class GraphType>
79 {
80  return OverlappingLevel_;
81 }
82 
83 
84 template<class GraphType>
85 typename GraphType::local_ordinal_type
87 operator () (const local_ordinal_type MyRow) const
88 {
90  MyRow < 0 || Teuchos::as<size_t> (MyRow) > Graph_->getNodeNumRows (),
91  std::runtime_error,
92  "Ifpack2::OverlappingPartitioner::operator(): "
93  "Invalid local row index " << MyRow << ".");
94 
95  return Partition_[MyRow];
96 }
97 
98 
99 //==============================================================================
100 template<class GraphType>
101 typename GraphType::local_ordinal_type
103 operator() (const local_ordinal_type i, const local_ordinal_type j) const
104 {
106  i < 0 || i > Teuchos::as<local_ordinal_type> (NumLocalParts_),
107  std::runtime_error,
108  "Ifpack2::OverlappingPartitioner::operator(): "
109  "Invalid local row index i=" << i << ".");
111  j < 0 || j > Teuchos::as<local_ordinal_type> (Parts_[i].size ()),
112  std::runtime_error,
113  "Ifpack2::OverlappingPartitioner::operator(): "
114  "Invalid node index j=" << j << ".");
115  return Parts_[i][j];
116 }
117 
118 //==============================================================================
119 template<class GraphType>
120 size_t
122 numRowsInPart (const local_ordinal_type Part) const
123 {
125  Part < 0 || Part > Teuchos::as<local_ordinal_type> (NumLocalParts_),
126  std::runtime_error,
127  "Ifpack2::OverlappingPartitioner::numRowsInPart: "
128  "Invalid partition index Part=" << Part << ".");
129  return Parts_[Part].size ();
130 }
131 
132 //==============================================================================
133 template<class GraphType>
134 void
136 rowsInPart (const local_ordinal_type Part,
138 {
139  // Let numRowsInPart do the sanity checking...
140  const size_t numRows = numRowsInPart (Part);
141  for (size_t i = 0; i < numRows; ++i) {
142  List[i] = Parts_[Part][i];
143  }
144 }
145 
146 //==============================================================================
147 template<class GraphType>
150 {
151  return Partition_.view (0, Graph_->getNodeNumRows ());
152 }
153 
154 //==============================================================================
155 template<class GraphType>
156 void
159 {
160  NumLocalParts_ = List.get("partitioner: local parts", NumLocalParts_);
161  OverlappingLevel_ = List.get("partitioner: overlap", OverlappingLevel_);
162  verbose_ = List.get("partitioner: print level", verbose_);
163 
164  if (NumLocalParts_ < 0) {
165  NumLocalParts_ = Graph_->getNodeNumRows() / (-NumLocalParts_);
166  }
167  if (NumLocalParts_ == 0) {
168  NumLocalParts_ = 1;
169  }
170 
171  // Sanity checking
173  NumLocalParts_ < 0 ||
174  Teuchos::as<size_t> (NumLocalParts_) > Graph_->getNodeNumRows(),
175  std::runtime_error,
176  "Ifpack2::OverlappingPartitioner::setParameters: "
177  "Invalid NumLocalParts_ = " << NumLocalParts_ << ".");
179  OverlappingLevel_ < 0, std::runtime_error,
180  "Ifpack2::OverlappingPartitioner::setParameters: "
181  "Invalid OverlappingLevel_ = " << OverlappingLevel_ << ".");
182 
183  setPartitionParameters(List);
184 }
185 
186 //==============================================================================
187 template<class GraphType>
189 {
190  using std::cout;
191  using std::endl;
192 
194  NumLocalParts_ < 1 || OverlappingLevel_ < 0,
195  std::runtime_error,
196  "Ifpack2::OverlappingPartitioner::compute: "
197  "Invalid NumLocalParts_ or OverlappingLevel_.");
198 
199  // std::string's constructor has some overhead, so it's better to
200  // use const char[] for local constant strings.
201  const char printMsg[] = "OverlappingPartitioner: ";
202 
203  if (verbose_ && (Graph_->getComm()->getRank() == 0)) {
204  cout << printMsg << "Number of local parts = "
205  << NumLocalParts_ << endl;
206  cout << printMsg << "Approx. Number of global parts = "
207  << NumLocalParts_ * Graph_->getComm ()->getSize () << endl;
208  cout << printMsg << "Amount of overlap = "
209  << OverlappingLevel_ << endl;
210  }
211 
212  // 1.- allocate memory
213  Partition_.resize (Graph_->getNodeNumRows ());
214  //Parts_ is allocated in computeOverlappingPartitions_, where it is used
215 
216  // 2.- sanity checks on input graph
218  ! Graph_->isFillComplete (), std::runtime_error,
219  "Ifpack2::OverlappingPartitioner::compute: "
220  "The input graph must be fill complete.");
221 
223  Graph_->getGlobalNumRows () != Graph_->getGlobalNumCols (),
224  std::runtime_error,
225  "Ifpack2::OverlappingPartitioner::compute: "
226  "The input graph must be (globally) square.");
227 
228  // 3.- perform non-overlapping partition
229  computePartitions ();
230 
231  // 4.- compute the partitions with overlapping
232  computeOverlappingPartitions ();
233 
234  // 5.- mark as computed
235  IsComputed_ = true;
236 }
237 
238 //==============================================================================
239 template<class GraphType>
241 {
242  //If user has explicitly specified parts, then Partition_ size has been set to 0.
243  //In this case, there is no need to compute Parts_.
244  if (Partition_.size() == 0)
245  return;
246 
247  const local_ordinal_type invalid =
249 
250  // Old FIXME from Ifpack: the first part of this function should be elsewhere
251 
252  // start defining the subgraphs for no overlap
253 
254  std::vector<size_t> sizes;
255  sizes.resize (NumLocalParts_);
256 
257  // 1.- compute how many rows are in each subgraph
258  for (int i = 0; i < NumLocalParts_; ++i) {
259  sizes[i] = 0;
260  }
261 
262  for (size_t i = 0; i < Graph_->getNodeNumRows (); ++i) {
264  Partition_[i] >= NumLocalParts_, std::runtime_error,
265  "Ifpack2::OverlappingPartitioner::computeOverlappingPartitions: "
266  "Partition_[i] > NumLocalParts_.");
267  // invalid indicates that this unknown is not in a nonoverlapping
268  // partition
269  if (Partition_[i] != invalid) {
270  sizes[Partition_[i]]++;
271  }
272  }
273 
274  // 2.- allocate space for each subgraph
275  Parts_.resize (NumLocalParts_);
276  for (int i = 0; i < NumLocalParts_; ++i) {
277  Parts_[i].resize (sizes[i]);
278  }
279 
280  // 3.- cycle over all rows and populate the vectors
281  for (int i = 0; i < NumLocalParts_; ++i) {
282  sizes[i] = 0;
283  }
284 
285  for (size_t i = 0; i < Graph_->getNodeNumRows (); ++i) {
286  const local_ordinal_type part = Partition_[i];
287  if (part != invalid) {
288  const size_t count = sizes[part];
289  Parts_[part][count] = i;
290  sizes[part]++;
291  }
292  }
293 
294  // If there is no overlap, we're done, so return
295  if (OverlappingLevel_ == 0) {
296  return;
297  }
298 
299  // wider overlap requires further computations
300  for (int level = 1; level <= OverlappingLevel_; ++level) {
301  std::vector<std::vector<size_t> > tmp;
302  tmp.resize (NumLocalParts_);
303 
304  // cycle over all rows in the local graph (that is the overlapping
305  // graph). For each row, all columns will belong to the subgraph
306  // of row `i'.
307 
308  int MaxNumEntries_tmp = Graph_->getNodeMaxNumRowEntries();
310  Indices.resize (MaxNumEntries_tmp);
311 
312  for (int part = 0; part < NumLocalParts_ ; ++part) {
313  for (size_t i = 0; i < Teuchos::as<size_t> (Parts_[part].size ()); ++i) {
314  const local_ordinal_type LRID = Parts_[part][i];
315 
316  size_t NumIndices;
317  Graph_->getLocalRowCopy (LRID, Indices (), NumIndices);
318 
319  for (size_t j = 0; j < NumIndices; ++j) {
320  // use *local* indices only
321  const local_ordinal_type col = Indices[j];
322  if (Teuchos::as<size_t> (col) >= Graph_->getNodeNumRows ()) {
323  continue;
324  }
325 
326  // has this column already been inserted?
327  std::vector<size_t>::iterator where =
328  std::find (tmp[part].begin (), tmp[part].end (), Teuchos::as<size_t> (col));
329 
330  if (where == tmp[part].end()) {
331  tmp[part].push_back (col);
332  }
333  }
334 
335  // has this column already been inserted?
336  std::vector<size_t>::iterator where =
337  std::find (tmp[part].begin (), tmp[part].end (), Teuchos::as<size_t> (LRID));
338 
339  // This happens here b/c Vanka on Stokes with Stabilized elements will have
340  // a zero pivot entry if this gets pushed back first. So... Last.
341  if (where == tmp[part].end ()) {
342  tmp[part].push_back (LRID);
343  }
344  }
345  }
346 
347  // now I convert the STL vectors into Teuchos Array RCP's
348  //
349  // FIXME (mfh 12 July 2013) You could have started with ArrayRCP
350  // in the first place (which implements push_back and iterators)
351  // and avoided the copy.
352  for (int i = 0; i < NumLocalParts_; ++i) {
353  Parts_[i].resize (tmp[i].size ());
354  for (size_t j = 0; j < tmp[i].size (); ++j) {
355  Parts_[i][j] = tmp[i][j];
356  }
357  }
358  }
359 }
360 
361 //==============================================================================
362 template<class GraphType>
364 {
365  return IsComputed_;
366 }
367 
368 //==============================================================================
369 template<class GraphType>
370 std::ostream&
372 {
373  Teuchos::FancyOStream fos (Teuchos::rcpFromRef (os));
374  fos.setOutputToRootOnly (0);
375  describe (fos);
376  return os;
377 }
378 
379 //==============================================================================
380 template<class GraphType>
382 {
383  std::ostringstream oss;
385  if (isComputed()) {
386  oss << "{status = computed";
387  }
388  else {
389  oss << "{status = is not computed";
390  }
391  oss <<"}";
392  return oss.str();
393 }
394 
395 //==============================================================================
396 template<class GraphType>
398 {
399  using std::endl;
400  if (verbLevel == Teuchos::VERB_NONE) {
401  return;
402  }
403 
404  os << "================================================================================" << endl;
405  os << "Ifpack2::OverlappingPartitioner" << endl;
406  os << "Number of local rows = " << Graph_->getNodeNumRows() << endl;
407  os << "Number of global rows = " << Graph_->getGlobalNumRows() << endl;
408  os << "Number of local parts = " << NumLocalParts_ << endl;
409  os << "Overlapping level = " << OverlappingLevel_ << endl;
410  os << "Is computed = " << IsComputed_ << endl;
411  os << "================================================================================" << endl;
412 }
413 
414 
415 }// namespace Ifpack2
416 
417 #define IFPACK2_OVERLAPPINGPARTITIONER_INSTANT(LO,GO,N) \
418  template class Ifpack2::OverlappingPartitioner<Tpetra::CrsGraph< LO, GO, N > >; \
419  template class Ifpack2::OverlappingPartitioner<Tpetra::RowGraph< LO, GO, N > >;
420 
421 #endif // IFPACK2_OVERLAPPINGPARTITIONER_DEF_HPP
virtual std::ostream & print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator<<.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:371
size_t numRowsInPart(const local_ordinal_type Part) const
the number of rows contained in the given partition.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:122
std::string description() const
Return a simple one-line description of this object.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:381
VERB_NONE
T & get(const std::string &name, T def_value)
OverlappingPartitioner(const Teuchos::RCP< const row_graph_type > &graph)
Constructor.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:56
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual std::string description() const
virtual void computeOverlappingPartitions()
Computes the partitions. Returns 0 if successful.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:240
int overlappingLevel() const
The number of levels of overlap.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:78
virtual ~OverlappingPartitioner()
Destructor.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:66
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
void resize(size_type new_size, const value_type &x=value_type())
virtual Teuchos::ArrayView< const local_ordinal_type > nonOverlappingPartition() const
A view of the local indices of the nonoverlapping partitions of each local row.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:149
virtual void compute()
Computes the partitions. Returns 0 if successful.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:188
void rowsInPart(const local_ordinal_type Part, Teuchos::ArrayRCP< local_ordinal_type > &List) const
Fill List with the local indices of the rows in the (overlapping) partition Part. ...
Definition: Ifpack2_OverlappingPartitioner_def.hpp:136
virtual bool isComputed() const
Returns true if partitions have been computed successfully.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:363
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:397
int numLocalParts() const
Number of computed local partitions.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:71
virtual void setParameters(Teuchos::ParameterList &List)
Set all the parameters for the partitioner.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:158
local_ordinal_type operator()(const local_ordinal_type MyRow) const
Local index of the nonoverlapping partition of the given row.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:87
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72