IFPACK Development
Loading...
Searching...
No Matches
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//==============================================================================
101SetParameter(const std::string Name, const int Value)
102{
103 if (Name == "reorder: root node")
104 RootNode_ = Value;
105 return(0);
106}
107
108//==============================================================================
110SetParameter(const std::string /* Name */, const double /* Value */)
111{
112 return(0);
113}
114
115//==============================================================================
117SetParameters(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//==============================================================================
237int 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//==============================================================================
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//==============================================================================
294std::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}
int NumVectors() const
Ifpack_Graph_Epetra_RowMatrix: a class to define Ifpack_Graph as a light-weight conversion of Epetra_...
Ifpack_Graph: a pure virtual class that defines graphs for IFPACK.
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<<.
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 SetParameters(Teuchos::ParameterList &List)
Sets all parameters.
virtual int NumMyRows() const
Returns the number of local rows.
virtual int SetParameter(const std::string Name, const int Value)
Sets integer parameters ‘Name’.
virtual bool IsComputed() const
Returns true is the reordering object has been successfully initialized, false otherwise.
Ifpack_RCMReordering & operator=(const Ifpack_RCMReordering &RHS)
Assignment operator.
virtual int InvReorder(const int i) const
Returns the inverse reordered index of row i.
Ifpack_RCMReordering()
Constructor for Ifpack_Graph's.
virtual int Reorder(const int i) const
Returns the reordered index of row i.
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.
virtual int Compute(const Ifpack_Graph &Graph)
Computes all it is necessary to initialize the reordering object.