00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030 #ifndef IFPACK_IKLU_H
00031 #define IFPACK_IKLU_H
00032
00033 #include "Ifpack_ConfigDefs.h"
00034 #include "Ifpack_CondestType.h"
00035 #include "Ifpack_ScalingType.h"
00036 #include "Ifpack_Preconditioner.h"
00037 #include "Ifpack_IKLU_Utils.h"
00038 #include "Epetra_Vector.h"
00039 #include "Epetra_CrsMatrix.h"
00040 #include "Epetra_Time.h"
00041 #include "Teuchos_RefCountPtr.hpp"
00042
00043 class Epetra_RowMatrix;
00044 class Epetra_SerialComm;
00045 class Epetra_Comm;
00046 class Epetra_Map;
00047 class Epetra_MultiVector;
00048
00049 namespace Teuchos {
00050 class ParameterList;
00051 }
00052
00054
00066 class Ifpack_IKLU: public Ifpack_Preconditioner {
00067
00068 public:
00069
00071 Ifpack_IKLU(const Epetra_RowMatrix* A);
00072
00074 virtual ~Ifpack_IKLU();
00075
00076
00077
00079
00080
00081
00082
00083
00084
00085
00086 int SetParameters(Teuchos::ParameterList& parameterlis);
00087
00089
00095 int Initialize();
00096
00098 bool IsInitialized() const
00099 {
00100 return(IsInitialized_);
00101 }
00102
00104
00112 int Compute();
00113
00115 bool IsComputed() const {return(IsComputed_);};
00116
00117
00118
00120
00128 int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00129
00130 int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00131
00133 double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
00134 const int MaxIters = 1550,
00135 const double Tol = 1e-9,
00136 Epetra_RowMatrix* Matrix_in = 0);
00137
00139 double Condest() const
00140 {
00141 return(Condest_);
00142 }
00143
00145
00154 int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
00155
00157 double NormInf() const {return(0.0);};
00158
00160 bool HasNormInf() const {return(false);};
00161
00163 bool UseTranspose() const {return(UseTranspose_);};
00164
00166 const Epetra_Map & OperatorDomainMap() const {return(A_.OperatorDomainMap());};
00167
00169 const Epetra_Map & OperatorRangeMap() const{return(A_.OperatorRangeMap());};
00170
00172 const Epetra_Comm & Comm() const{return(Comm_);};
00173
00175 const Epetra_RowMatrix& Matrix() const
00176 {
00177 return(A_);
00178 }
00179
00181 const Epetra_CrsMatrix & L() const {return(*L_);};
00182
00184 const Epetra_CrsMatrix & U() const {return(*U_);};
00185
00187 const char* Label() const
00188 {
00189 return(Label_.c_str());
00190 }
00191
00193 int SetLabel(const char* Label_in)
00194 {
00195 Label_ = Label_in;
00196 return(0);
00197 }
00198
00200 virtual ostream& Print(std::ostream& os) const;
00201
00203 virtual int NumInitialize() const
00204 {
00205 return(NumInitialize_);
00206 }
00207
00209 virtual int NumCompute() const
00210 {
00211 return(NumCompute_);
00212 }
00213
00215 virtual int NumApplyInverse() const
00216 {
00217 return(NumApplyInverse_);
00218 }
00219
00221 virtual double InitializeTime() const
00222 {
00223 return(InitializeTime_);
00224 }
00225
00227 virtual double ComputeTime() const
00228 {
00229 return(ComputeTime_);
00230 }
00231
00233 virtual double ApplyInverseTime() const
00234 {
00235 return(ApplyInverseTime_);
00236 }
00237
00239 virtual double InitializeFlops() const
00240 {
00241 return(0.0);
00242 }
00243
00244 virtual double ComputeFlops() const
00245 {
00246 return(ComputeFlops_);
00247 }
00248
00249 virtual double ApplyInverseFlops() const
00250 {
00251 return(ApplyInverseFlops_);
00252 }
00253
00254 inline double LevelOfFill() const {
00255 return(LevelOfFill_);
00256 }
00257
00259 inline double RelaxValue() const {
00260 return(Relax_);
00261 }
00262
00264 inline double AbsoluteThreshold() const
00265 {
00266 return(Athresh_);
00267 }
00268
00270 inline double RelativeThreshold() const
00271 {
00272 return(Rthresh_);
00273 }
00274
00276 inline double DropTolerance() const
00277 {
00278 return(DropTolerance_);
00279 }
00280
00282 int NumGlobalNonzeros() const {
00283
00284 return(L().NumGlobalNonzeros() + U().NumGlobalNonzeros() - L().NumGlobalRows());
00285 }
00286
00288 int NumMyNonzeros() const {
00289 return(L().NumMyNonzeros() + U().NumMyNonzeros());
00290 }
00291
00292 private:
00293
00294
00295
00296
00298 Ifpack_IKLU(const Ifpack_IKLU& RHS) :
00299 A_(RHS.Matrix()),
00300 Comm_(RHS.Comm()),
00301 Time_(Comm())
00302 {};
00303
00305 Ifpack_IKLU& operator=(const Ifpack_IKLU& RHS)
00306 {
00307 return(*this);
00308 }
00309
00311 void Destroy();
00312
00313
00314
00315
00317 const Epetra_RowMatrix& A_;
00319 const Epetra_Comm& Comm_;
00321 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
00323 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
00325 double Condest_;
00327 double Relax_;
00329 double Athresh_;
00331 double Rthresh_;
00333 double LevelOfFill_;
00335 double DropTolerance_;
00337 string Label_;
00339 bool IsInitialized_;
00341 bool IsComputed_;
00343 bool UseTranspose_;
00345 int NumMyRows_;
00347 int NumMyNonzeros_;
00349 int NumInitialize_;
00351 int NumCompute_;
00353 mutable int NumApplyInverse_;
00355 double InitializeTime_;
00357 double ComputeTime_;
00359 mutable double ApplyInverseTime_;
00361 double ComputeFlops_;
00363 mutable double ApplyInverseFlops_;
00365 mutable Epetra_Time Time_;
00367 int GlobalNonzeros_;
00368 Teuchos::RefCountPtr<Epetra_SerialComm> SerialComm_;
00369 Teuchos::RefCountPtr<Epetra_Map> SerialMap_;
00370
00372 csr* csrA_;
00373 css* cssS_;
00375 csrn* csrnN_;
00376
00377 };
00378
00379 #endif