Ifpack_RCMReordering.cpp

Go to the documentation of this file.
00001 #include "Ifpack_ConfigDefs.h"
00002 #if defined(HAVE_IFPACK_TEUCHOS)
00003 #include "Teuchos_ParameterList.hpp"
00004 #include "Epetra_MultiVector.h"
00005 #include "Ifpack_Graph.h"
00006 #include "Ifpack_Graph_Epetra_RowMatrix.h"
00007 #include "Ifpack_RCMReordering.h"
00008 
00009 //==============================================================================
00010 Ifpack_RCMReordering::
00011 Ifpack_RCMReordering() :
00012   RootNode_(0),
00013   NumMyRows_(0),
00014   IsComputed_(false)
00015 {
00016 }
00017 
00018 //==============================================================================
00019 Ifpack_RCMReordering::
00020 Ifpack_RCMReordering(const Ifpack_RCMReordering& RHS) :
00021   RootNode_(RHS.RootNode()),
00022   NumMyRows_(RHS.NumMyRows()),
00023   IsComputed_(RHS.IsComputed())
00024 {
00025   Reorder_.resize(NumMyRows());
00026   InvReorder_.resize(NumMyRows());
00027   for (int i = 0 ; i < NumMyRows() ; ++i) {
00028     Reorder_[i] = RHS.Reorder(i);
00029     InvReorder_[i] = RHS.InvReorder(i);
00030   }
00031 }
00032 
00033 //==============================================================================
00034 Ifpack_RCMReordering::~Ifpack_RCMReordering()
00035 {}
00036   
00037 //==============================================================================
00038 Ifpack_RCMReordering& Ifpack_RCMReordering::
00039 operator=(const Ifpack_RCMReordering& RHS)
00040 {
00041   if (this == &RHS) {
00042     return (*this);
00043   }
00044 
00045   NumMyRows_ = RHS.NumMyRows(); // set number of local rows
00046   RootNode_ = RHS.RootNode(); // set root node
00047   IsComputed_ = RHS.IsComputed();
00048   // resize vectors, and copy values from RHS
00049   Reorder_.resize(NumMyRows()); 
00050   InvReorder_.resize(NumMyRows());
00051   if (IsComputed()) {
00052     for (int i = 0 ; i < NumMyRows_ ; ++i) {
00053       Reorder_[i] = RHS.Reorder(i);
00054       InvReorder_[i] = RHS.InvReorder(i);
00055     }
00056   }
00057   return (*this);
00058 }
00059 
00060 //==============================================================================
00061 int Ifpack_RCMReordering::
00062 SetParameter(const string Name, const int Value)
00063 {
00064   if (Name == "reorder: root node")
00065     RootNode_ = Value;
00066   return(0);
00067 }
00068 
00069 //==============================================================================
00070 int Ifpack_RCMReordering::
00071 SetParameter(const string Name, const double Value)
00072 {
00073   return(0);
00074 }
00075 
00076 //==============================================================================
00077 #ifdef HAVE_IFPACK_TEUCHOS  
00078 int Ifpack_RCMReordering::
00079 SetParameters(Teuchos::ParameterList& List)
00080 {
00081   RootNode_ = List.get("reorder: root node", RootNode_);
00082   return(0);
00083 }
00084 #endif
00085 
00086 //==============================================================================
00087 int Ifpack_RCMReordering::Compute(const Epetra_RowMatrix& Matrix)
00088 {
00089   Ifpack_Graph_Epetra_RowMatrix Graph(&Matrix);
00090 
00091   IFPACK_CHK_ERR(Compute(Graph));
00092 
00093   return(0);
00094 }
00095 
00096 //==============================================================================
00097 int Ifpack_RCMReordering::Compute(const Ifpack_Graph& Graph)
00098 {
00099   IsComputed_ = false;
00100   NumMyRows_ = Graph.NumMyRows();
00101   
00102   if (NumMyRows_ == 0)
00103     IFPACK_CHK_ERR(-1); // strange graph this one
00104   
00105   if ((RootNode_ < 0) || (RootNode_ >= NumMyRows_))
00106     RootNode_ = 0;
00107     
00108   Reorder_.resize(NumMyRows_);
00109 
00110   for (int i = 0 ; i < NumMyRows_ ; ++i)
00111     Reorder_[i] = -1;
00112 
00113   vector<int> tmp;
00114   tmp.push_back(RootNode_);
00115 
00116   int count = NumMyRows_ - 1;
00117   int Length = Graph.MaxMyNumEntries();
00118   vector<int> Indices(Length);
00119   
00120   Reorder_[RootNode_] = count;
00121   count--;
00122 
00123   // stop when no nodes have been added in the previous level
00124 
00125   while (tmp.size()) {
00126 
00127     vector<int> tmp2;
00128 
00129     // for each node in the previous level, look for non-marked
00130     // neighbors. 
00131     for (int i = 0 ; i < (int)tmp.size() ; ++i) {
00132       int NumEntries;
00133       IFPACK_CHK_ERR(Graph.ExtractMyRowCopy(tmp[i], Length,
00134                NumEntries, &Indices[0]));
00135 
00136       if (Length > 1)
00137   sort(Indices.begin(), Indices.begin() + Length);
00138 
00139       for (int j = 0 ; j < NumEntries ; ++j) {
00140   int col = Indices[j];
00141   if (col >= NumMyRows_) 
00142     continue;
00143 
00144   if (Reorder_[col] == -1) {
00145     Reorder_[col] = count;
00146     count--;
00147     if (col != tmp[i]) {
00148       tmp2.push_back(col);
00149     }
00150   }
00151       }
00152     }
00153 
00154     // if no nodes have been found but we still have
00155     // rows to walk through, to localize the next -1 
00156     // and restart.
00157     // FIXME: I can replace with STL
00158     if ((tmp2.size() == 0) && (count != -1)) {
00159       for (int i = 0 ; i < NumMyRows_ ; ++i)
00160   if (Reorder_[i] == -1) {
00161     tmp2.push_back(i);
00162     Reorder_[i] = count--;
00163     break;
00164   }
00165     }
00166 
00167     // prepare for the next level
00168     tmp = tmp2;
00169   }
00170 
00171   // check nothing went wrong
00172   for (int i = 0 ; i < NumMyRows_ ; ++i) {
00173     if (Reorder_[i] == -1)
00174       IFPACK_CHK_ERR(-1);
00175   }
00176   
00177   // build inverse reorder (will be used by ExtractMyRowCopy() 
00178   InvReorder_.resize(NumMyRows_);
00179 
00180   for (int i = 0 ; i < NumMyRows_ ; ++i)
00181     InvReorder_[i] = -1;
00182 
00183   for (int i = 0 ; i < NumMyRows_ ; ++i)
00184     InvReorder_[Reorder_[i]] = i;
00185 
00186   for (int i = 0 ; i < NumMyRows_ ; ++i) {
00187     if (InvReorder_[i] == -1)
00188       IFPACK_CHK_ERR(-1);
00189   }
00190 
00191   IsComputed_ = true;
00192   return(0);
00193 }
00194 
00195 //==============================================================================
00196 int Ifpack_RCMReordering::Reorder(const int i) const
00197 {
00198 #ifdef IFPACK_ABC
00199   if (!IsComputed())
00200     IFPACK_CHK_ERR(-1);
00201   if ((i < 0) || (i >= NumMyRows_))
00202     IFPACK_CHK_ERR(-1);
00203 #endif
00204 
00205   return(Reorder_[i]);
00206 }
00207 
00208 //==============================================================================
00209 int Ifpack_RCMReordering::InvReorder(const int i) const
00210 {
00211 #ifdef IFPACK_ABC
00212   if (!IsComputed())
00213     IFPACK_CHK_ERR(-1);
00214   if ((i < 0) || (i >= NumMyRows_))
00215     IFPACK_CHK_ERR(-1);
00216 #endif
00217 
00218   return(InvReorder_[i]);
00219 }
00220 //==============================================================================
00221 int Ifpack_RCMReordering::P(const Epetra_MultiVector& Xorig,
00222           Epetra_MultiVector& X) const
00223 {  
00224   int NumVectors = X.NumVectors();
00225 
00226   for (int j = 0 ; j < NumVectors ; ++j) {
00227     for (int i = 0 ; i < NumMyRows_ ; ++i) {
00228       int np = Reorder_[i];
00229       X[j][np] = Xorig[j][i];
00230     }
00231   }
00232 
00233   return(0);
00234 }
00235 
00236 //==============================================================================
00237 int Ifpack_RCMReordering::Pinv(const Epetra_MultiVector& Xorig,
00238              Epetra_MultiVector& X) const
00239 {
00240   int NumVectors = X.NumVectors();
00241 
00242   for (int j = 0 ; j < NumVectors ; ++j) {
00243     for (int i = 0 ; i < NumMyRows_ ; ++i) {
00244       int np = Reorder_[i];
00245       X[j][i] = Xorig[j][np];
00246     }
00247   }
00248 
00249   return(0);
00250 }
00251 
00252 //==============================================================================
00253 ostream& Ifpack_RCMReordering::Print(std::ostream& os) const
00254 {
00255   os << "*** Ifpack_RCMReordering" << endl << endl;
00256   if (!IsComputed())
00257     os << "*** Reordering not yet computed." << endl;
00258   
00259   os << "*** Number of local rows = " << NumMyRows_ << endl;
00260   os << "*** Root node = " << RootNode_ << endl;
00261   os << endl;
00262   os << "Local Row\tReorder[i]\tInvReorder[i]" << endl;
00263   for (int i = 0 ; i < NumMyRows_ ; ++i) {
00264     os << '\t' << i << "\t\t" << Reorder_[i] << "\t\t" << InvReorder_[i] << endl;
00265   }
00266    
00267   return(os);
00268 }
00269 #endif // HAVE_IFPACK_TEUCHOS

Generated on Thu Sep 18 12:37:22 2008 for Ifpack Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1