43 #ifndef _IFPACK_CRSRICK_H_
44 #define _IFPACK_CRSRICK_H_
47 #include "Ifpack_IlukGraph.h"
48 #include "Epetra_CompObject.h"
49 #include "Epetra_Operator.h"
50 #include "Epetra_CrsMatrix.h"
51 #include "Epetra_Object.h"
255 int SetParameters(
const Teuchos::ParameterList& parameterlist,
256 bool cerr_warning_if_unused=
false);
324 int Condest(
bool Trans,
double & ConditionNumberEstimate)
const;
379 const char *
Label()
const {
return(Epetra_Object::Label());};
442 void SetFactored(
bool Flag) {Factored_ = Flag;};
443 void SetValuesInitialized(
bool Flag) {ValuesInitialized_ = Flag;};
444 bool Allocated()
const {
return(Allocated_);};
445 int SetAllocated(
bool Flag) {Allocated_ = Flag;
return(0);};
460 bool ValuesInitialized_;
465 mutable double Condest_;
475 std::ostream& operator << (std::ostream& os,
const Ifpack_CrsRick& A);
void SetRelaxValue(double RelaxValue)
Set RIC(k) relaxation parameter.
const Epetra_Map & RangeMap() const
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y.
const Epetra_Vector & D() const
Returns the address of the D factor associated with this factored matrix.
int IndexBase() const
Returns the index base for row and column indices for this graph.
const Epetra_CrsMatrix & U() const
Returns the address of the U factor associated with this factored matrix.
bool ValuesInitialized() const
If values have been initialized, this query returns true, otherwise it returns false.
void SetOverlapMode(Epetra_CombineMode OverlapMode)
Set overlap mode type.
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
int NumMyRows() const
Returns the number of local matrix rows.
int NumGlobalNonzeros() const
Returns the number of nonzero entries in the global graph.
const Ifpack_IlukGraph & Graph() const
Returns the address of the Ifpack_IlukGraph associated with this factored matrix. ...
int InitValues()
Initialize L and U with values from user matrix A.
Epetra_CombineMode GetOverlapMode()
Get overlap mode type.
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 SetParameters(const Teuchos::ParameterList ¶meterlist, bool cerr_warning_if_unused=false)
Set parameters using a Teuchos::ParameterList object.
int NumGlobalRows() const
Returns the number of global matrix rows.
int SetUseTranspose(bool UseTranspose)
If set true, transpose of this operator will be applied.
Ifpack_ScalingType enumerable type.
double GetRelaxValue()
Get RIC(k) relaxation parameter.
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
virtual int NumGlobalDiagonals() const
Returns the number of diagonal entries found in the global input graph.
Ifpack_CrsRick: A class for constructing and using an incomplete Cholesky (IC) factorization of a giv...
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
int NumMyNonzeros() const
Returns the number of nonzero entries in the local graph.
Ifpack_CrsRick(const Epetra_CrsMatrix &A, const Ifpack_IlukGraph &Graph)
Ifpack_CrsRick constuctor with variable number of indices per row.
virtual int NumMyDiagonals() const
Returns the number of diagonal entries found in the local input graph.
bool UseTranspose() const
Returns the current UseTranspose setting.
int NumGlobalCols() const
Returns the number of global matrix columns.
virtual ~Ifpack_CrsRick()
Ifpack_CrsRick Destructor.
void SetAbsoluteThreshold(double Athresh)
Set absolute threshold value.
int NumMyCols() const
Returns the number of local matrix columns.
int Solve(bool Trans, const Epetra_Vector &x, Epetra_Vector &y) const
Returns the result of a Ifpack_CrsRick forward/back solve on a Epetra_Vector x in y...
const Epetra_Map & DomainMap() const
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
int Factor()
Compute IC factor L using the specified graph, diagonal perturbation thresholds and relaxation parame...
friend std::ostream & operator<<(std::ostream &os, const Ifpack_CrsRick &A)
<< operator will work for Ifpack_CrsRick.
const char * Label() const
Returns a character string describing the operator.
int Multiply(bool Trans, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of multiplying U, D and U^T in that order on an Epetra_MultiVector X in Y...
int Condest(bool Trans, double &ConditionNumberEstimate) const
Returns the maximum over all the condition number estimate for each local IC set of factors...
bool Factored() const
If factor is completed, this query returns true, otherwise it returns false.
double GetRelativeThreshold()
Get relative threshold value.
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
void SetRelativeThreshold(double Rthresh)
Set relative threshold value.
double GetAbsoluteThreshold()
Get absolute threshold value.