IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
Ifpack_RCMReordering.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 "Teuchos_ParameterList.hpp"
45 #include "Teuchos_RefCountPtr.hpp"
46 #include "Epetra_MultiVector.h"
47 #include "Ifpack_Graph.h"
48 #include "Epetra_RowMatrix.h"
49 #include "Ifpack_Graph_Epetra_RowMatrix.h"
50 #include "Ifpack_RCMReordering.h"
51 
52 //==============================================================================
55  RootNode_(0),
56  NumMyRows_(0),
57  IsComputed_(false)
58 {
59 }
60 
61 //==============================================================================
64  RootNode_(RHS.RootNode()),
65  NumMyRows_(RHS.NumMyRows()),
66  IsComputed_(RHS.IsComputed())
67 {
68  Reorder_.resize(NumMyRows());
69  InvReorder_.resize(NumMyRows());
70  for (int i = 0 ; i < NumMyRows() ; ++i) {
71  Reorder_[i] = RHS.Reorder(i);
72  InvReorder_[i] = RHS.InvReorder(i);
73  }
74 }
75 
76 //==============================================================================
79 {
80  if (this == &RHS) {
81  return (*this);
82  }
83 
84  NumMyRows_ = RHS.NumMyRows(); // set number of local rows
85  RootNode_ = RHS.RootNode(); // set root node
86  IsComputed_ = RHS.IsComputed();
87  // resize vectors, and copy values from RHS
88  Reorder_.resize(NumMyRows());
89  InvReorder_.resize(NumMyRows());
90  if (IsComputed()) {
91  for (int i = 0 ; i < NumMyRows_ ; ++i) {
92  Reorder_[i] = RHS.Reorder(i);
93  InvReorder_[i] = RHS.InvReorder(i);
94  }
95  }
96  return (*this);
97 }
98 
99 //==============================================================================
101 SetParameter(const std::string Name, const int Value)
102 {
103  if (Name == "reorder: root node")
104  RootNode_ = Value;
105  return(0);
106 }
107 
108 //==============================================================================
110 SetParameter(const std::string /* Name */, const double /* Value */)
111 {
112  return(0);
113 }
114 
115 //==============================================================================
117 SetParameters(Teuchos::ParameterList& List)
118 {
119  RootNode_ = List.get("reorder: root node", RootNode_);
120  return(0);
121 }
122 
123 //==============================================================================
125 {
126  Ifpack_Graph_Epetra_RowMatrix Graph(Teuchos::rcp(&Matrix,false));
127 
128  IFPACK_CHK_ERR(Compute(Graph));
129 
130  return(0);
131 }
132 
133 //==============================================================================
135 {
136  IsComputed_ = false;
137  NumMyRows_ = Graph.NumMyRows();
138 
139  if ((RootNode_ < 0) || (RootNode_ >= NumMyRows_))
140  RootNode_ = 0;
141 
142  Reorder_.resize(NumMyRows_);
143 
144  // the case where one processor holds no chunk of the graph happens...
145  if (!NumMyRows_) {
146  InvReorder_.resize(NumMyRows_);
147  IsComputed_ = true;
148  return(0);
149  }
150 
151  for (int i = 0 ; i < NumMyRows_ ; ++i)
152  Reorder_[i] = -1;
153 
154  std::vector<int> tmp;
155  tmp.push_back(RootNode_);
156 
157  int count = NumMyRows_ - 1;
158  int Length = Graph.MaxMyNumEntries();
159  std::vector<int> Indices(Length);
160 
161  Reorder_[RootNode_] = count;
162  count--;
163 
164  // stop when no nodes have been added in the previous level
165 
166  while (tmp.size()) {
167 
168  std::vector<int> tmp2;
169 
170  // for each node in the previous level, look for non-marked
171  // neighbors.
172  for (int i = 0 ; i < (int)tmp.size() ; ++i) {
173  int NumEntries;
174  IFPACK_CHK_ERR(Graph.ExtractMyRowCopy(tmp[i], Length,
175  NumEntries, &Indices[0]));
176 
177  if (Length > 1)
178  std::sort(Indices.begin(), Indices.begin() + Length);
179 
180  for (int j = 0 ; j < NumEntries ; ++j) {
181  int col = Indices[j];
182  if (col >= NumMyRows_)
183  continue;
184 
185  if (Reorder_[col] == -1) {
186  Reorder_[col] = count;
187  count--;
188  if (col != tmp[i]) {
189  tmp2.push_back(col);
190  }
191  }
192  }
193  }
194 
195  // if no nodes have been found but we still have
196  // rows to walk through, to localize the next -1
197  // and restart.
198  // FIXME: I can replace with STL
199  if ((tmp2.size() == 0) && (count != -1)) {
200  for (int i = 0 ; i < NumMyRows_ ; ++i)
201  if (Reorder_[i] == -1) {
202  tmp2.push_back(i);
203  Reorder_[i] = count--;
204  break;
205  }
206  }
207 
208  // prepare for the next level
209  tmp = tmp2;
210  }
211 
212  // check nothing went wrong
213  for (int i = 0 ; i < NumMyRows_ ; ++i) {
214  if (Reorder_[i] == -1)
215  IFPACK_CHK_ERR(-1);
216  }
217 
218  // build inverse reorder (will be used by ExtractMyRowCopy()
219  InvReorder_.resize(NumMyRows_);
220 
221  for (int i = 0 ; i < NumMyRows_ ; ++i)
222  InvReorder_[i] = -1;
223 
224  for (int i = 0 ; i < NumMyRows_ ; ++i)
225  InvReorder_[Reorder_[i]] = i;
226 
227  for (int i = 0 ; i < NumMyRows_ ; ++i) {
228  if (InvReorder_[i] == -1)
229  IFPACK_CHK_ERR(-1);
230  }
231 
232  IsComputed_ = true;
233  return(0);
234 }
235 
236 //==============================================================================
237 int Ifpack_RCMReordering::Reorder(const int i) const
238 {
239 #ifdef IFPACK_ABC
240  if (!IsComputed())
241  IFPACK_CHK_ERR(-1);
242  if ((i < 0) || (i >= NumMyRows_))
243  IFPACK_CHK_ERR(-1);
244 #endif
245 
246  return(Reorder_[i]);
247 }
248 
249 //==============================================================================
250 int Ifpack_RCMReordering::InvReorder(const int i) const
251 {
252 #ifdef IFPACK_ABC
253  if (!IsComputed())
254  IFPACK_CHK_ERR(-1);
255  if ((i < 0) || (i >= NumMyRows_))
256  IFPACK_CHK_ERR(-1);
257 #endif
258 
259  return(InvReorder_[i]);
260 }
261 //==============================================================================
263  Epetra_MultiVector& X) const
264 {
265  int NumVectors = X.NumVectors();
266 
267  for (int j = 0 ; j < NumVectors ; ++j) {
268  for (int i = 0 ; i < NumMyRows_ ; ++i) {
269  int np = Reorder_[i];
270  X[j][np] = Xorig[j][i];
271  }
272  }
273 
274  return(0);
275 }
276 
277 //==============================================================================
279  Epetra_MultiVector& X) const
280 {
281  int NumVectors = X.NumVectors();
282 
283  for (int j = 0 ; j < NumVectors ; ++j) {
284  for (int i = 0 ; i < NumMyRows_ ; ++i) {
285  int np = Reorder_[i];
286  X[j][i] = Xorig[j][np];
287  }
288  }
289 
290  return(0);
291 }
292 
293 //==============================================================================
294 std::ostream& Ifpack_RCMReordering::Print(std::ostream& os) const
295 {
296  using std::endl;
297 
298  os << "*** Ifpack_RCMReordering" << endl << endl;
299  if (!IsComputed())
300  os << "*** Reordering not yet computed." << endl;
301 
302  os << "*** Number of local rows = " << NumMyRows_ << endl;
303  os << "*** Root node = " << RootNode_ << endl;
304  os << endl;
305  os << "Local Row\tReorder[i]\tInvReorder[i]" << endl;
306  for (int i = 0 ; i < NumMyRows_ ; ++i) {
307  os << '\t' << i << "\t\t" << Reorder_[i] << "\t\t" << InvReorder_[i] << endl;
308  }
309 
310  return(os);
311 }
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 bool IsComputed() const
Returns true is the reordering object has been successfully initialized, false otherwise.
virtual int Pinv(const Epetra_MultiVector &Xorig, Epetra_MultiVector &Xinvreord) const
Applies inverse reordering to multivector X, whose local length equals the number of local rows...
virtual int Compute(const Ifpack_Graph &Graph)
Computes all it is necessary to initialize the reordering object.
virtual int NumMyRows() const
Returns the number of local rows.
virtual int InvReorder(const int i) const
Returns the inverse reordered index of row i.
virtual int RootNode() const
Returns the root node.
Ifpack_RCMReordering: reverse Cuthill-McKee reordering.
virtual std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator&lt;&lt;.
Ifpack_Graph_Epetra_RowMatrix: a class to define Ifpack_Graph as a light-weight conversion of Epetra_...
virtual int MaxMyNumEntries() const =0
Returns the maximun number of entries for row.
Ifpack_Graph: a pure virtual class that defines graphs for IFPACK.
Definition: Ifpack_Graph.h:61
virtual int SetParameters(Teuchos::ParameterList &List)
Sets all parameters.
virtual int Reorder(const int i) const
Returns the reordered index of row i.
virtual int SetParameter(const std::string Name, const int Value)
Sets integer parameters `Name&#39;.
Ifpack_RCMReordering & operator=(const Ifpack_RCMReordering &RHS)
Assignment operator.
virtual int P(const Epetra_MultiVector &Xorig, Epetra_MultiVector &Xreord) const
Applies reordering to multivector X, whose local length equals the number of local rows...
Ifpack_RCMReordering()
Constructor for Ifpack_Graph&#39;s.