46 #include "Ifpack_ConfigDefs.h"
48 #ifdef HAVE_IFPACK_SUPERLU
49 #include "Ifpack_Preconditioner.h"
50 #include "Ifpack_Condest.h"
52 #include "Ifpack_IlukGraph.h"
53 #include "Epetra_CompObject.h"
54 #include "Epetra_MultiVector.h"
55 #include "Epetra_Vector.h"
56 #include "Epetra_CrsGraph.h"
57 #include "Epetra_CrsMatrix.h"
58 #include "Epetra_BlockMap.h"
59 #include "Epetra_Map.h"
60 #include "Epetra_Object.h"
61 #include "Epetra_Comm.h"
62 #include "Epetra_RowMatrix.h"
63 #include "Epetra_Time.h"
64 #include "Teuchos_RefCountPtr.hpp"
71 #include "slu_ddefs.h"
107 return(IsInitialized_);
139 int SetUseTranspose(
bool UseTranspose_in) {UseTranspose_ = UseTranspose_in;
return(0);};
148 return(Multiply(
false,X,Y));
171 double Condest(
const Ifpack_CondestType CT = Ifpack_Cheap,
172 const int MaxIters = 1550,
173 const double Tol = 1e-9,
188 const char*
Label()
const {
return(Label_);}
193 strcpy(Label_,Label_in);
222 virtual std::ostream&
Print(std::ostream& os)
const;
227 return(NumInitialize_);
239 return(NumApplyInverse_);
245 return(InitializeTime_);
251 return(ComputeTime_);
257 return(ApplyInverseTime_);
312 double DropTol()
const {
return DropTol_;}
315 double FillTol()
const{
return FillTol_;}
318 double FillFactor()
const{
return FillFactor_;}
321 int DropRule()
const{
return DropRule_;}
323 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
324 int NumGlobalRows()
const {
return(Graph().NumGlobalRows());};
328 int NumGlobalCols()
const {
return(Graph().NumGlobalCols());};
331 int NumGlobalNonzeros()
const {
return(Graph().NumGlobalNonzeros());};
334 virtual int NumGlobalBlockDiagonals()
const {
return(Graph().NumGlobalBlockDiagonals());};
338 long long NumGlobalRows64()
const {
return(Graph().NumGlobalRows64());};
341 long long NumGlobalCols64()
const {
return(Graph().NumGlobalCols64());};
344 long long NumGlobalNonzeros64()
const {
return(Graph().NumGlobalNonzeros64());};
347 virtual long long NumGlobalBlockDiagonals64()
const {
return(Graph().NumGlobalBlockDiagonals64());};
350 int NumMyRows()
const {
return(Graph().NumMyRows());};
353 int NumMyCols()
const {
return(Graph().NumMyCols());};
356 int NumMyNonzeros()
const {
return(Graph().NumMyNonzeros());};
359 virtual int NumMyBlockDiagonals()
const {
return(Graph().NumMyBlockDiagonals());};
362 virtual int NumMyDiagonals()
const {
return(NumMyDiagonals_);};
364 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
365 int IndexBase()
const {
return(Graph().IndexBase());};
368 long long IndexBase64()
const {
return(Graph().IndexBase64());};
385 Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
386 Teuchos::RefCountPtr<Epetra_CrsGraph> Graph_;
387 Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
388 Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
389 Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
392 Teuchos::RefCountPtr<Epetra_CrsMatrix> Aover_;
397 bool ValuesInitialized_;
422 mutable int NumApplyInverse_;
424 double InitializeTime_;
428 mutable double ApplyInverseTime_;
432 #ifdef HAVE_IFPACK_SUPERLU5_API
433 mutable GlobalLU_t lu_;
435 mutable SuperLUStat_t stat_;
438 mutable superlu_options_t options_;
440 mutable SuperMatrix SA_,SAc_,SL_,SU_,SY_;
442 int *etree_,*perm_r_,*perm_c_;
446 template<
typename int_type>
int SetLabel(const char *Label_in)
Sets label for this object.
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
const char * Label() const
Returns a character string describing the operator.
Ifpack_SILU(Epetra_RowMatrix *A)
Constructor.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
bool UseTranspose() const
Returns the current UseTranspose setting.
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
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...
A wrapper to SuperLU 4.0's supernodal ILUT w/ partial pivoting.
virtual double InitializeTime() const
Returns the time spent in Initialize().
~Ifpack_SILU()
Destructor.
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
virtual std::ostream & Print(std::ostream &os) const
Prints on stream basic information about this object.
Ifpack_ScalingType enumerable type.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
int Compute()
Compute ILU factors L and U using the specified graph, diagonal perturbation thresholds and relaxatio...
int Initialize()
Initialize the preconditioner, does not touch matrix values.
int SetParameters(Teuchos::ParameterList ¶meterlist)
Set parameters using a Teuchos::ParameterList object.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
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.
virtual double ComputeTime() const
Returns the time spent in Compute().
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
virtual int NumCompute() const
Returns the number of calls to Compute().
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
virtual int NumApplyInverse() const
Returns the number of calls to 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.
double Condest() const
Returns the computed estimated condition number, or -1.0 if not computed.