46#include "Ifpack_ConfigDefs.h"
47#include "Ifpack_Preconditioner.h"
48#include "Ifpack_Condest.h"
50#include "Ifpack_IlukGraph.h"
51#include "Epetra_CompObject.h"
52#include "Epetra_MultiVector.h"
53#include "Epetra_Vector.h"
54#include "Epetra_CrsGraph.h"
55#include "Epetra_CrsMatrix.h"
56#include "Epetra_BlockMap.h"
57#include "Epetra_Map.h"
58#include "Epetra_Object.h"
59#include "Epetra_Comm.h"
60#include "Epetra_RowMatrix.h"
61#include "Epetra_Time.h"
62#include "Teuchos_RefCountPtr.hpp"
105 return(IsInitialized_);
144 int SetUseTranspose(
bool UseTranspose_in) {UseTranspose_ = UseTranspose_in;
return(0);};
152 return(Multiply(
false,X,Y));
175 double Condest(
const Ifpack_CondestType CT = Ifpack_Cheap,
176 const int MaxIters = 1550,
177 const double Tol = 1e-9,
199 const char*
Label()
const {
return(Label_);}
204 strcpy(Label_,Label_in);
233 virtual std::ostream&
Print(std::ostream& os)
const;
238 return(NumInitialize_);
250 return(NumApplyInverse_);
256 return(InitializeTime_);
262 return(ComputeTime_);
268 return(ApplyInverseTime_);
279 return(ComputeFlops_);
284 return(ApplyInverseFlops_);
324 int LevelOfFill()
const {
return LevelOfFill_;}
327 double RelaxValue()
const {
return RelaxValue_;}
330 double AbsoluteThreshold()
const {
return Athresh_;}
333 double RelativeThreshold()
const {
return Rthresh_;}
335#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
337 int NumGlobalRows()
const {
return(Graph().NumGlobalRows());};
340 int NumGlobalCols()
const {
return(Graph().NumGlobalCols());};
343 int NumGlobalNonzeros()
const {
return(
L().NumGlobalNonzeros()+
U().NumGlobalNonzeros());};
346 virtual int NumGlobalBlockDiagonals()
const {
return(Graph().NumGlobalBlockDiagonals());};
349 long long NumGlobalRows64()
const {
return(Graph().NumGlobalRows64());};
351 long long NumGlobalCols64()
const {
return(Graph().NumGlobalCols64());};
353 long long NumGlobalNonzeros64()
const {
return(
L().NumGlobalNonzeros64()+
U().NumGlobalNonzeros64());};
355 virtual long long NumGlobalBlockDiagonals64()
const {
return(Graph().NumGlobalBlockDiagonals64());};
358 int NumMyRows()
const {
return(Graph().NumMyRows());};
361 int NumMyCols()
const {
return(Graph().NumMyCols());};
364 int NumMyNonzeros()
const {
return(
L().NumMyNonzeros()+
U().NumMyNonzeros());};
367 virtual int NumMyBlockDiagonals()
const {
return(Graph().NumMyBlockDiagonals());};
370 virtual int NumMyDiagonals()
const {
return(NumMyDiagonals_);};
373#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
374 int IndexBase()
const {
return(Graph().IndexBase());};
376 long long IndexBase64()
const {
return(Graph().IndexBase64());};
391 Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
392 Teuchos::RefCountPtr<Ifpack_IlukGraph> Graph_;
393 Teuchos::RefCountPtr<Epetra_CrsGraph> CrsGraph_;
394 Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
395 Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
396 Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
401 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
403 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
404 Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
405 Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
407 Teuchos::RefCountPtr<Epetra_Vector> D_;
412 bool ValuesInitialized_;
435 mutable int NumApplyInverse_;
437 double InitializeTime_;
441 mutable double ApplyInverseTime_;
443 double ComputeFlops_;
445 mutable double ApplyInverseFlops_;
Ifpack_ScalingType enumerable type.
Ifpack_ILU: A class for constructing and using an incomplete lower/upper (ILU) factorization of a giv...
virtual double InitializeTime() const
Returns the time spent in Initialize().
const Epetra_CrsMatrix & U() const
Returns the address of the L factor associated with this factored matrix.
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
int SetLabel(const char *Label_in)
Sets label for this object.
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
virtual double ComputeTime() const
Returns the time spent in Compute().
int SetParameters(Teuchos::ParameterList ¶meterlist)
Set parameters using a Teuchos::ParameterList object.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
int Compute()
Compute ILU factors L and U using the specified graph, diagonal perturbation thresholds and relaxatio...
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
virtual std::ostream & Print(std::ostream &os) const
Prints on stream basic information about this object.
const Epetra_CrsMatrix & L() const
Returns the address of the L factor associated with this factored matrix.
const Epetra_Vector & D() const
Returns the address of the D factor associated with this factored matrix.
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
bool UseTranspose() const
Returns the current UseTranspose setting.
double Condest() const
Returns the computed estimated condition number, or -1.0 if not computed.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied.
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y.
int Initialize()
Initialize the preconditioner, does not touch matrix values.
virtual int NumCompute() const
Returns the number of calls to Compute().
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
Ifpack_ILU(Epetra_RowMatrix *A)
Constructor.
const char * Label() const
Returns a character string describing the operator.
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.