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_ILU_H
00031 #define IFPACK_ILU_H
00032
00033 #include "Ifpack_ConfigDefs.h"
00034 #include "Ifpack_Preconditioner.h"
00035 #include "Ifpack_Condest.h"
00036 #include "Ifpack_ScalingType.h"
00037 #include "Ifpack_IlukGraph.h"
00038 #include "Epetra_CompObject.h"
00039 #include "Epetra_MultiVector.h"
00040 #include "Epetra_Vector.h"
00041 #include "Epetra_CrsGraph.h"
00042 #include "Epetra_CrsMatrix.h"
00043 #include "Epetra_BlockMap.h"
00044 #include "Epetra_Map.h"
00045 #include "Epetra_Object.h"
00046 #include "Epetra_Comm.h"
00047 #include "Epetra_RowMatrix.h"
00048 #include "Epetra_Time.h"
00049 #include "Teuchos_RefCountPtr.hpp"
00050
00051 namespace Teuchos {
00052 class ParameterList;
00053 }
00054
00056
00069 class Ifpack_ILU: public Ifpack_Preconditioner {
00070
00071 public:
00072
00074 Ifpack_ILU(Epetra_RowMatrix* A);
00075
00077 ~Ifpack_ILU()
00078 {
00079 Destroy();
00080 }
00081
00082
00083
00084
00086 int Initialize();
00087
00089 bool IsInitialized() const
00090 {
00091 return(IsInitialized_);
00092 }
00093
00095
00103 int Compute();
00104
00106 bool IsComputed() const
00107 {
00108 return(IsComputed_);
00109 }
00110
00112
00113
00114
00115
00116
00117
00118
00119 int SetParameters(Teuchos::ParameterList& parameterlist);
00120
00122
00131 int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
00132
00133
00134
00135
00136 int Apply(const Epetra_MultiVector& X,
00137 Epetra_MultiVector& Y) const
00138 {
00139 return(Multiply(false,X,Y));
00140 }
00141
00142 int Multiply(bool Trans, const Epetra_MultiVector& X,
00143 Epetra_MultiVector& Y) const;
00144
00146
00159 int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00160
00162 double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
00163 const int MaxIters = 1550,
00164 const double Tol = 1e-9,
00165 Epetra_RowMatrix* Matrix_in = 0);
00166
00168 double Condest() const
00169 {
00170 return(Condest_);
00171 }
00172
00173
00174
00175
00177 const Epetra_CrsMatrix & L() const {return(*L_);};
00178
00180 const Epetra_Vector & D() const {return(*D_);};
00181
00183 const Epetra_CrsMatrix & U() const {return(*U_);};
00184
00186 const char* Label() const {return(Label_);}
00187
00189 int SetLabel(const char* Label_in)
00190 {
00191 strcpy(Label_,Label_in);
00192 return(0);
00193 }
00194
00196 double NormInf() const {return(0.0);};
00197
00199 bool HasNormInf() const {return(false);};
00200
00202 bool UseTranspose() const {return(UseTranspose_);};
00203
00205 const Epetra_Map & OperatorDomainMap() const {return(U_->OperatorDomainMap());};
00206
00208 const Epetra_Map & OperatorRangeMap() const{return(L_->OperatorRangeMap());};
00209
00211 const Epetra_Comm & Comm() const{return(Comm_);};
00212
00214 const Epetra_RowMatrix& Matrix() const
00215 {
00216 return(*A_);
00217 }
00218
00220 virtual ostream& Print(ostream& os) const;
00221
00223 virtual int NumInitialize() const
00224 {
00225 return(NumInitialize_);
00226 }
00227
00229 virtual int NumCompute() const
00230 {
00231 return(NumCompute_);
00232 }
00233
00235 virtual int NumApplyInverse() const
00236 {
00237 return(NumApplyInverse_);
00238 }
00239
00241 virtual double InitializeTime() const
00242 {
00243 return(InitializeTime_);
00244 }
00245
00247 virtual double ComputeTime() const
00248 {
00249 return(ComputeTime_);
00250 }
00251
00253 virtual double ApplyInverseTime() const
00254 {
00255 return(ApplyInverseTime_);
00256 }
00257
00259 virtual double InitializeFlops() const
00260 {
00261 return(0.0);
00262 }
00263
00264 virtual double ComputeFlops() const
00265 {
00266 return(ComputeFlops_);
00267 }
00268
00269 virtual double ApplyInverseFlops() const
00270 {
00271 return(ApplyInverseFlops_);
00272 }
00273
00274 private:
00275
00276
00277
00278
00280 Ifpack_ILU(const Ifpack_ILU& RHS) :
00281 Comm_(RHS.Comm()),
00282 Time_(RHS.Comm())
00283 {}
00284
00286 Ifpack_ILU& operator=(const Ifpack_ILU& RHS)
00287 {
00288 return(*this);
00289 }
00290
00292 void Destroy();
00293
00295
00305 int Solve(bool Trans, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00306
00307 int ComputeSetup();
00308 int InitAllValues(const Epetra_RowMatrix & A, int MaxNumEntries);
00309
00311 int LevelOfFill() const {return LevelOfFill_;}
00312
00314 double RelaxValue() const {return RelaxValue_;}
00315
00317 double AbsoluteThreshold() const {return Athresh_;}
00318
00320 double RelativeThreshold() const {return Rthresh_;}
00321
00323 int NumGlobalRows() const {return(Graph().NumGlobalRows());};
00324
00326 int NumGlobalCols() const {return(Graph().NumGlobalCols());};
00327
00329 int NumGlobalNonzeros() const {return(L().NumGlobalNonzeros()+U().NumGlobalNonzeros());};
00330
00332 virtual int NumGlobalBlockDiagonals() const {return(Graph().NumGlobalBlockDiagonals());};
00333
00335 int NumMyRows() const {return(Graph().NumMyRows());};
00336
00338 int NumMyCols() const {return(Graph().NumMyCols());};
00339
00341 int NumMyNonzeros() const {return(L().NumMyNonzeros()+U().NumMyNonzeros());};
00342
00344 virtual int NumMyBlockDiagonals() const {return(Graph().NumMyBlockDiagonals());};
00345
00347 virtual int NumMyDiagonals() const {return(NumMyDiagonals_);};
00348
00350 int IndexBase() const {return(Graph().IndexBase());};
00351
00353 const Ifpack_IlukGraph & Graph() const {return(*Graph_);};
00354
00356 Epetra_RowMatrix& Matrix()
00357 {
00358 return(*A_);
00359 }
00360
00361
00362
00363
00365 Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
00366 Teuchos::RefCountPtr<Ifpack_IlukGraph> Graph_;
00367 Teuchos::RefCountPtr<Epetra_CrsGraph> CrsGraph_;
00368 Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
00369 Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
00370 Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
00371 const Epetra_Map * U_DomainMap_;
00372 const Epetra_Map * L_RangeMap_;
00373 const Epetra_Comm & Comm_;
00375 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
00377 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
00378 Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
00379 Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
00381 Teuchos::RefCountPtr<Epetra_Vector> D_;
00382 bool UseTranspose_;
00383
00384 int NumMyDiagonals_;
00385 bool Allocated_;
00386 bool ValuesInitialized_;
00387 bool Factored_;
00389 double RelaxValue_;
00391 double Athresh_;
00393 double Rthresh_;
00395 double Condest_;
00397 int LevelOfFill_;
00399 bool IsInitialized_;
00401 bool IsComputed_;
00403 char Label_[160];
00405 int NumInitialize_;
00407 int NumCompute_;
00409 mutable int NumApplyInverse_;
00411 double InitializeTime_;
00413 double ComputeTime_;
00415 mutable double ApplyInverseTime_;
00417 double ComputeFlops_;
00419 mutable double ApplyInverseFlops_;
00421 mutable Epetra_Time Time_;
00422
00423 };
00424
00425 #endif