43 #ifndef IFPACK_CHEBYSHEV_H 44 #define IFPACK_CHEBYSHEV_H 46 #include "Ifpack_ConfigDefs.h" 47 #include "Ifpack_Preconditioner.h" 48 #include "Teuchos_RefCountPtr.hpp" 54 class Epetra_MultiVector;
60 class Epetra_Operator;
61 class Epetra_RowMatrix;
63 #ifdef HAVE_IFPACK_EPETRAEXT 64 class EpetraExt_PointToBlockDiagPermute;
138 UseTranspose_ = UseTranspose_in;
155 virtual inline int Apply(
const Epetra_MultiVector& X, Epetra_MultiVector& Y)
const;
168 virtual int ApplyInverse(
const Epetra_MultiVector& X, Epetra_MultiVector& Y)
const;
179 virtual const char * Label()
const 181 return(Label_.c_str());
187 return(UseTranspose_);
197 virtual const Epetra_Comm &
Comm()
const;
209 return(IsInitialized_);
231 virtual const Epetra_RowMatrix&
Matrix()
const 237 virtual double Condest(
const Ifpack_CondestType CT = Ifpack_Cheap,
238 const int MaxIters = 1550,
239 const double Tol = 1e-9,
240 Epetra_RowMatrix* Matrix_in = 0);
252 virtual std::ostream&
Print(std::ostream & os)
const;
261 return(NumInitialize_);
273 return(NumApplyInverse_);
279 return(InitializeTime_);
285 return(ComputeTime_);
291 return(ApplyInverseTime_);
303 return(ComputeFlops_);
309 return(ApplyInverseFlops_);
316 static int PowerMethod(
const Epetra_Operator& Operator,
317 const Epetra_Vector& InvPointDiagonal,
318 const int MaximumIterations,
319 double& LambdaMax,
const unsigned int * RngSeed=0);
322 static int CG(
const Epetra_Operator& Operator,
323 const Epetra_Vector& InvPointDiagonal,
324 const int MaximumIterations,
325 double& lambda_min,
double& lambda_max,
const unsigned int * RngSeed=0);
327 #ifdef HAVE_IFPACK_EPETRAEXT 330 int CG(
const int MaximumIterations,
331 double& lambda_min,
double& lambda_max,
const unsigned int * RngSeed=0);
334 int PowerMethod(
const int MaximumIterations,
double& lambda_max,
const unsigned int * RngSeed=0);
342 virtual void SetLabel();
364 mutable int NumApplyInverse_;
366 double InitializeTime_;
370 mutable double ApplyInverseTime_;
372 double ComputeFlops_;
374 mutable double ApplyInverseFlops_;
388 bool ComputeCondest_;
402 double MinDiagonalValue_;
411 long long NumGlobalRows_;
413 long long NumGlobalNonzeros_;
415 Teuchos::RefCountPtr<const Epetra_Operator> Operator_;
417 Teuchos::RefCountPtr<const Epetra_RowMatrix> Matrix_;
419 mutable Teuchos::RefCountPtr<Epetra_Vector> InvDiagonal_;
422 #ifdef HAVE_IFPACK_EPETRAEXT 423 Teuchos::ParameterList BlockList_;
425 Teuchos::RefCountPtr<EpetraExt_PointToBlockDiagPermute> InvBlockDiagonal_;
429 bool SolveNormalEquations_;
434 Teuchos::RefCountPtr<Epetra_Time> Time_;
436 bool ZeroStartingSolution_;
443 #endif // IFPACK_CHEBYSHEV_H virtual const Epetra_RowMatrix & Matrix() const
Returns a pointer to the matrix to be preconditioned.
virtual bool UseTranspose() const
Returns the current UseTranspose setting.
virtual int Compute()
Computes the preconditioners.
virtual double ComputeTime() const
Returns the time spent in Compute().
virtual ~Ifpack_Chebyshev()
Destructor.
virtual const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
virtual int SetParameters(Teuchos::ParameterList &List)
Sets all the parameters for the preconditioner.
virtual double GetLambdaMin()
Contains an approximation to the smallest eigenvalue.
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
virtual int NumCompute() const
Returns the number of calls to Compute().
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
Ifpack_Chebyshev: class for preconditioning with Chebyshev polynomials in Ifpack. ...
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
virtual const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the matrix to an Epetra_MultiVector.
virtual double GetLambdaMax()
Returns an approximation to the largest eigenvalue.
virtual bool HasNormInf() const
Returns true if the this object can provide an approximate Inf-norm, false otherwise.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
static int CG(const Epetra_Operator &Operator, const Epetra_Vector &InvPointDiagonal, const int MaximumIterations, double &lambda_min, double &lambda_max, const unsigned int *RngSeed=0)
Uses AztecOO's CG to estimate lambda_min and lambda_max.
virtual double InitializeTime() const
Returns the time spent in Initialize().
virtual double ApplyInverseFlops() const
Returns the number of flops for the application of the preconditioner.
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
virtual bool IsComputed() const
Returns true if the preconditioner has been successfully computed.
virtual int SetUseTranspose(bool UseTranspose_in)
Ifpack_Chebyshev(const Epetra_Operator *Matrix)
Ifpack_Chebyshev constructor with given Epetra_Operator/Epetra_RowMatrix.
virtual const Epetra_Comm & Comm() const
Returns a pointer to the Epetra_Comm communicator associated with this operator.
virtual int Initialize()
Computes all it is necessary to initialize the preconditioner.
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the preconditioner to X, returns the result in Y.
virtual bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized, false otherwise.
virtual double Condest() const
Returns the condition number estimate, or -1.0 if not computed.
virtual std::ostream & Print(std::ostream &os) const
Prints object to an output stream.
static int PowerMethod(const Epetra_Operator &Operator, const Epetra_Vector &InvPointDiagonal, const int MaximumIterations, double &LambdaMax, const unsigned int *RngSeed=0)
Simple power method to compute lambda_max.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
virtual double NormInf() const
Returns the infinity norm of the global matrix (not implemented)