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_AMDReordering.h"
53#include <trilinos_amd.h>
67 NumMyRows_(RHS.NumMyRows()),
68 IsComputed_(RHS.IsComputed())
73 Reorder_[i] = RHS.Reorder(i);
74 InvReorder_[i] = RHS.InvReorder(i);
86 NumMyRows_ = RHS.NumMyRows();
87 IsComputed_ = RHS.IsComputed();
92 for (
int i = 0 ; i < NumMyRows_ ; ++i) {
93 Reorder_[i] = RHS.Reorder(i);
94 InvReorder_[i] = RHS.InvReorder(i);
126 IFPACK_CHK_ERR(
Compute(Graph));
138 NumMyRows_ = Graph.NumMyRows();
139 int NumNz = Graph.NumMyNonzeros();
145 std::vector<int> ia(NumMyRows_+1,0);
146 std::vector<int> ja(NumNz);
148 for(
int i = 0; i < NumMyRows_; ++i )
150 int * tmpP = &ja[ia[i]];
151 Graph.ExtractMyRowCopy( i, NumNz-ia[i], cnt, tmpP );
152 ia[i+1] = ia[i] + cnt;
156 std::vector<int> iat(NumMyRows_+1);
157 std::vector<int> jat(NumNz);
159 for(
int i = 0; i < NumMyRows_; ++i )
162 for(
int j = ia[i]; j < ia[i+1]; ++j )
164 if( ja[j] < NumMyRows_ )
170 iat[NumMyRows_] = loc;
173 Reorder_.resize(NumMyRows_);
174 std::vector<double> info(TRILINOS_AMD_INFO);
176 trilinos_amd_order( NumMyRows_, &iat[0], &jat[0], &Reorder_[0], NULL, &info[0] );
178 if( info[TRILINOS_AMD_STATUS] == TRILINOS_AMD_INVALID )
179 cout <<
"AMD ORDERING: Invalid!!!!" << endl;
182 InvReorder_.resize(NumMyRows_);
184 for (
int i = 0 ; i < NumMyRows_ ; ++i)
187 for (
int i = 0 ; i < NumMyRows_ ; ++i)
188 InvReorder_[Reorder_[i]] = i;
190 for (
int i = 0 ; i < NumMyRows_ ; ++i) {
191 if (InvReorder_[i] == -1)
205 if ((i < 0) || (i >= NumMyRows_))
218 if ((i < 0) || (i >= NumMyRows_))
222 return(InvReorder_[i]);
230 for (
int j = 0 ; j < NumVectors ; ++j) {
231 for (
int i = 0 ; i < NumMyRows_ ; ++i) {
232 int np = Reorder_[i];
233 X[j][np] = Xorig[j][i];
246 for (
int j = 0 ; j < NumVectors ; ++j) {
247 for (
int i = 0 ; i < NumMyRows_ ; ++i) {
248 int np = Reorder_[i];
249 X[j][i] = Xorig[j][np];
261 os <<
"*** Ifpack_AMDReordering" << endl << endl;
263 os <<
"*** Reordering not yet computed." << endl;
265 os <<
"*** Number of local rows = " << NumMyRows_ << endl;
267 os <<
"Local Row\tReorder[i]\tInvReorder[i]" << endl;
268 for (
int i = 0 ; i < NumMyRows_ ; ++i) {
269 os <<
'\t' << i <<
"\t\t" << Reorder_[i] <<
"\t\t" << InvReorder_[i] << endl;
Ifpack_AMDReordering: approximate minimum degree reordering.
int SetParameters(Teuchos::ParameterList &List)
Sets all parameters.
int Compute(const Ifpack_Graph &Graph)
Computes all it is necessary to initialize the reordering object.
Ifpack_AMDReordering()
Constructor for Ifpack_Graph's.
int InvReorder(const int i) const
Returns the inverse reordered index of row i.
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.
int SetParameter(const std::string Name, const int Value)
Sets integer parameters ‘Name’.
bool IsComputed() const
Returns true is the reordering object has been successfully initialized, false otherwise.
Ifpack_AMDReordering & operator=(const Ifpack_AMDReordering &RHS)
Assignment operator.
std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator<<.
int P(const Epetra_MultiVector &Xorig, Epetra_MultiVector &Xreord) const
Applies reordering to multivector X, whose local length equals the number of local rows.
int Reorder(const int i) const
Returns the reordered index of row i.
int NumMyRows() const
Returns the number of local rows.
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.