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_);
145 int SetUseTranspose(
bool UseTranspose_in) {UseTranspose_ = UseTranspose_in;
return(0);};
153 return(Multiply(
false,X,Y));
176 double Condest(
const Ifpack_CondestType CT = Ifpack_Cheap,
177 const int MaxIters = 1550,
178 const double Tol = 1e-9,
200 const char*
Label()
const {
return(Label_);}
205 strcpy(Label_,Label_in);
234 virtual std::ostream&
Print(std::ostream& os)
const;
239 return(NumInitialize_);
251 return(NumApplyInverse_);
257 return(InitializeTime_);
263 return(ComputeTime_);
269 return(ApplyInverseTime_);
280 return(ComputeFlops_);
285 return(ApplyInverseFlops_);
325 int LevelOfFill()
const {
return LevelOfFill_;}
328 double RelaxValue()
const {
return RelaxValue_;}
331 double AbsoluteThreshold()
const {
return Athresh_;}
334 double RelativeThreshold()
const {
return Rthresh_;}
336 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 337 int NumGlobalRows()
const {
return(Graph().NumGlobalRows());};
341 int NumGlobalCols()
const {
return(Graph().NumGlobalCols());};
344 int NumGlobalNonzeros()
const {
return(
L().NumGlobalNonzeros()+
U().NumGlobalNonzeros());};
347 virtual int NumGlobalBlockDiagonals()
const {
return(Graph().NumGlobalBlockDiagonals());};
350 long long NumGlobalRows64()
const {
return(Graph().NumGlobalRows64());};
352 long long NumGlobalCols64()
const {
return(Graph().NumGlobalCols64());};
354 long long NumGlobalNonzeros64()
const {
return(
L().NumGlobalNonzeros64()+
U().NumGlobalNonzeros64());};
356 virtual long long NumGlobalBlockDiagonals64()
const {
return(Graph().NumGlobalBlockDiagonals64());};
359 int NumMyRows()
const {
return(Graph().NumMyRows());};
362 int NumMyCols()
const {
return(Graph().NumMyCols());};
365 int NumMyNonzeros()
const {
return(
L().NumMyNonzeros()+
U().NumMyNonzeros());};
368 virtual int NumMyBlockDiagonals()
const {
return(Graph().NumMyBlockDiagonals());};
371 virtual int NumMyDiagonals()
const {
return(NumMyDiagonals_);};
374 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 375 int IndexBase()
const {
return(Graph().IndexBase());};
377 long long IndexBase64()
const {
return(Graph().IndexBase64());};
392 Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
393 Teuchos::RefCountPtr<Ifpack_IlukGraph> Graph_;
394 Teuchos::RefCountPtr<Epetra_CrsGraph> CrsGraph_;
395 Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
396 Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
397 Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
402 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
404 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
405 Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
406 Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
408 Teuchos::RefCountPtr<Epetra_Vector> D_;
413 bool ValuesInitialized_;
436 mutable int NumApplyInverse_;
438 double InitializeTime_;
442 mutable double ApplyInverseTime_;
444 double ComputeFlops_;
446 mutable double ApplyInverseFlops_;
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
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.
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
const Epetra_CrsMatrix & L() const
Returns the address of the L factor associated with this factored matrix.
const Epetra_CrsMatrix & U() const
Returns the address of the L factor associated with this factored matrix.
virtual double InitializeTime() const
Returns the time spent in Initialize().
const char * Label() const
Returns a character string describing the operator.
Ifpack_ILU: A class for constructing and using an incomplete lower/upper (ILU) factorization of a giv...
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
int SetLabel(const char *Label_in)
Sets label for this object.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
int Compute()
Compute ILU factors L and U using the specified graph, diagonal perturbation thresholds and relaxatio...
Ifpack_ScalingType enumerable type.
int SetParameters(Teuchos::ParameterList ¶meterlist)
Set parameters using a Teuchos::ParameterList object.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
virtual double ComputeTime() const
Returns the time spent in Compute().
int Initialize()
Initialize the preconditioner, does not touch matrix values.
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
Ifpack_ILU(Epetra_RowMatrix *A)
Constructor.
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
virtual std::ostream & Print(std::ostream &os) const
Prints on stream basic information about this object.
bool UseTranspose() const
Returns the current UseTranspose setting.
const Epetra_Vector & D() const
Returns the address of the D factor associated with this factored matrix.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
double Condest() const
Returns the computed estimated condition number, or -1.0 if not computed.
virtual int NumCompute() const
Returns the number of calls to Compute().
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...