• Main Page
  • Related Pages
  • Namespaces
  • Classes
  • Files
  • File List
  • File Members

ml_Ifpack_ML.h

Go to the documentation of this file.
00001 
00013 /* ******************************************************************** */
00014 /* See the file COPYRIGHT for a complete copyright notice, contact      */
00015 /* person and disclaimer.                                               */        
00016 /* ******************************************************************** */
00017 
00018 #ifndef ML_IFPACK_ML_H
00019 #define ML_IFPACK_ML_H
00020 
00021 #include "ml_include.h"
00022 #if defined(HAVE_ML_EPETRA) && defined(HAVE_ML_TEUCHOS) && defined(HAVE_ML_IFPACK)
00023 
00024 #include "ml_epetra.h"
00025 #include "Ifpack_Preconditioner.h"
00026 #include "Ifpack_Amesos.h"
00027 #include "ml_MultiLevelPreconditioner.h"
00028 
00029 namespace ML_Epetra {
00030 
00042 class Ifpack_ML : public Ifpack_Preconditioner {
00043 
00044 public:
00045 
00047   Ifpack_ML(Epetra_RowMatrix* A) :
00048     A_(A),
00049     MLPrec_(0)
00050   {};
00051  
00053   virtual ~Ifpack_ML() 
00054   {
00055     if (MLPrec_)
00056       delete MLPrec_;
00057   }
00058 
00060   virtual int SetParameters(Teuchos::ParameterList& List)
00061   {
00062     string listName = List.get("ML sublist name","ML list");
00063     try{MLList_ = List.sublist(listName,true);}
00064     catch(...) {
00065       if (A_->Comm().MyPID()==0)
00066         cout << "Did not find sublist \"" << listName
00067              << "\" for ML subdomain solver.  Setting \"SA\" defaults." << endl;
00068       SetDefaults("SA",MLList_);
00069     };
00070     return(0);
00071   }
00072 
00074   virtual int Initialize() 
00075   {
00076     return(0);
00077   };
00078 
00080   virtual bool IsInitialized() const
00081   {
00082     return(true);
00083   }
00084 
00086   virtual int Compute()
00087   {
00088     if (MLPrec_)
00089       delete MLPrec_;
00090 
00091     MLPrec_ = new ML_Epetra::MultiLevelPreconditioner(*A_, MLList_);
00092     if (MLPrec_->IsPreconditionerComputed() == false) {
00093       ML_CHK_ERR(-1);
00094     }
00095     else
00096       return(0);
00097   }
00098 
00100   virtual bool IsComputed() const
00101   {
00102     return(MLPrec_->IsPreconditionerComputed());
00103   }
00104 
00106   virtual double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
00107                          const int MaxIters = 1550,
00108                          const double Tol = 1e-9,
00109                          Epetra_RowMatrix* Matrix = 0)
00110   {
00111     return(-1.0);
00112   }
00113 
00115   virtual double Condest() const
00116   {
00117     return(-1.0);
00118   }
00119 
00120 
00122   virtual int ApplyInverse(const Epetra_MultiVector& X,
00123                            Epetra_MultiVector& Y) const
00124   {
00125     ML_RETURN(MLPrec_->ApplyInverse(X, Y));
00126   }
00127 
00129   virtual const Epetra_RowMatrix& Matrix() const
00130   {
00131     return(*A_);
00132   }
00133 
00135   virtual int NumInitialize() const
00136   {
00137     return(-1);
00138   }
00139 
00141   virtual int NumCompute() const
00142   {
00143     return(-1);
00144   }
00145 
00147   virtual int NumApplyInverse() const
00148   {
00149     return(-1);
00150   }
00151 
00153   virtual double InitializeTime() const
00154   {
00155     return(0.0);
00156   }
00157 
00159   virtual double ComputeTime() const
00160   {
00161     return(0.0);
00162   }
00163 
00165   virtual double ApplyInverseTime() const
00166   {
00167     return(0.0);
00168   }
00169 
00171   virtual double InitializeFlops() const
00172   {
00173     return(0.0);
00174   }
00175 
00177   virtual double ComputeFlops() const
00178   {
00179     return(0.0);
00180   }
00181 
00183   virtual double ApplyInverseFlops() const
00184   {
00185     return(0.0);
00186   }
00187 
00189   virtual ostream& Print(std::ostream& os) const
00190   {
00191     return(os);
00192   }
00193 
00195   int SetUseTranspose(bool UseTranspose)
00196   {
00197     if (UseTranspose == true)
00198       ML_CHK_ERR(-1);
00199 
00200     return(0);
00201   }
00202 
00204   int Apply(const Epetra_MultiVector&, Epetra_MultiVector&) const
00205   {
00206     ML_CHK_ERR(-1);
00207   }
00208 
00210   double NormInf() const
00211   {
00212     return(-1.0);
00213   }
00214 
00216   const char* Label() const
00217   {
00218     return("Ifpack_ML");
00219   }
00220 
00222   bool UseTranspose() const
00223   {
00224     ML_CHK_ERR(-1);
00225   }
00226 
00228   bool HasNormInf() const 
00229   {
00230     return(false);
00231   }
00232 
00234   const Epetra_Comm& Comm() const
00235   {
00236     return(A_->Comm());
00237   }
00238 
00240   const Epetra_Map& OperatorDomainMap() const
00241   {
00242     return(A_->OperatorDomainMap());
00243   }
00244 
00246   const Epetra_Map& OperatorRangeMap() const
00247   {
00248     return(A_->OperatorRangeMap());
00249   }
00250 
00251 private:
00253   Epetra_RowMatrix* A_;
00255   ML_Epetra::MultiLevelPreconditioner* MLPrec_;
00257   Teuchos::ParameterList MLList_;
00258 }; // class Ifpack_ML
00259 
00260 } // namespace ML_Epetra
00261 
00262 #endif
00263 #endif