Ifpack_ILU.h

Go to the documentation of this file.
00001 /*@HEADER
00002 // ***********************************************************************
00003 // 
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2002) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 //@HEADER
00028 */
00029 
00030 #ifndef IFPACK_ILU_H
00031 #define IFPACK_ILU_H
00032 
00033 #include "Ifpack_ConfigDefs.h"
00034 #ifdef HAVE_IFPACK_TEUCHOS
00035 #include "Ifpack_Preconditioner.h"
00036 #include "Ifpack_Condest.h"
00037 #include "Ifpack_ScalingType.h"
00038 #include "Ifpack_IlukGraph.h"
00039 #include "Epetra_CompObject.h"
00040 #include "Epetra_MultiVector.h"
00041 #include "Epetra_Vector.h"
00042 #include "Epetra_CrsGraph.h"
00043 #include "Epetra_CrsMatrix.h"
00044 #include "Epetra_BlockMap.h"
00045 #include "Epetra_Map.h"
00046 #include "Epetra_Object.h"
00047 #include "Epetra_Comm.h"
00048 #include "Epetra_RowMatrix.h"
00049 #include "Epetra_Time.h"
00050 
00051 #ifdef HAVE_IFPACK_TEUCHOS
00052 namespace Teuchos {
00053   class ParameterList;
00054   template<class T>
00055   class RefCountPtr;
00056 }
00057 #endif
00058 
00060 
00073 class Ifpack_ILU: public Ifpack_Preconditioner {
00074       
00075 public:
00076   // @{ Constructors and destructors.
00078   Ifpack_ILU(Epetra_RowMatrix* A);
00079   
00081   ~Ifpack_ILU()
00082   {
00083     Destroy();
00084   }
00085 
00086   // @}
00087   // @{ Construction methods
00088   
00090   int Initialize();
00091   
00093   bool IsInitialized() const
00094   {
00095     return(IsInitialized_);
00096   }
00097 
00099 
00107   int Compute();
00108 
00110   bool IsComputed() const 
00111   {
00112     return(IsComputed_);
00113   }
00114 
00116   /* This method is only available if the Teuchos package is enabled.
00117      This method recognizes four parameter names: relax_value,
00118      absolute_threshold, relative_threshold and overlap_mode. These names are
00119      case insensitive, and in each case except overlap_mode, the ParameterEntry
00120      must have type double. For overlap_mode, the ParameterEntry must have
00121      type Epetra_CombineMode.
00122   */
00123   int SetParameters(Teuchos::ParameterList& parameterlist);
00124 
00126 
00135   int SetUseTranspose(bool UseTranspose) {UseTranspose_ = UseTranspose; return(0);};
00136   // @}
00137 
00138   // @{ Mathematical functions.
00139   // Applies the matrix to X, returns the result in Y.
00140   int Apply(const Epetra_MultiVector& X, 
00141          Epetra_MultiVector& Y) const
00142   {
00143     return(Multiply(false,X,Y));
00144   }
00145 
00146   int Multiply(bool Trans, const Epetra_MultiVector& X, 
00147          Epetra_MultiVector& Y) const;
00148 
00150 
00163   int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00164 
00166   double Condest(const Ifpack_CondestType CT = Ifpack_Cheap, 
00167                  const int MaxIters = 1550,
00168                  const double Tol = 1e-9,
00169      Epetra_RowMatrix* Matrix = 0);
00170 
00172   double Condest() const
00173   {
00174     return(Condest_);
00175   }
00176 
00177   // @}
00178   // @{ Query methods
00179   
00181   const Epetra_CrsMatrix & L() const {return(*L_);};
00182     
00184   const Epetra_Vector & D() const {return(*D_);};
00185     
00187   const Epetra_CrsMatrix & U() const {return(*U_);};
00188 
00190   const char* Label() const {return(Label_);}
00191 
00193   int SetLabel(const char* Label)
00194   {
00195     strcpy(Label_,Label);
00196     return(0);
00197   }
00198   
00200   double NormInf() const {return(0.0);};
00201 
00203   bool HasNormInf() const {return(false);};
00204 
00206   bool UseTranspose() const {return(UseTranspose_);};
00207 
00209   const Epetra_Map & OperatorDomainMap() const {return(U_->OperatorDomainMap());};
00210 
00212   const Epetra_Map & OperatorRangeMap() const{return(L_->OperatorRangeMap());};
00213 
00215   const Epetra_Comm & Comm() const{return(Comm_);};
00216 
00218   const Epetra_RowMatrix& Matrix() const
00219   { 
00220     return(*A_);
00221   }
00222 
00224   virtual ostream& Print(ostream& os) const;
00225 
00227   virtual int NumInitialize() const
00228   {
00229     return(NumInitialize_);
00230   }
00231 
00233   virtual int NumCompute() const
00234   {
00235     return(NumCompute_);
00236   }
00237 
00239   virtual int NumApplyInverse() const
00240   {
00241     return(NumApplyInverse_);
00242   }
00243 
00245   virtual double InitializeTime() const
00246   {
00247     return(InitializeTime_);
00248   }
00249 
00251   virtual double ComputeTime() const
00252   {
00253     return(ComputeTime_);
00254   }
00255 
00257   virtual double ApplyInverseTime() const
00258   {
00259     return(ApplyInverseTime_);
00260   }
00261 
00263   virtual double InitializeFlops() const
00264   {
00265     return(0.0);
00266   }
00267 
00268   virtual double ComputeFlops() const
00269   {
00270     return(ComputeFlops_);
00271   }
00272 
00273   virtual double ApplyInverseFlops() const
00274   {
00275     return(ApplyInverseFlops_);
00276   }
00277 
00278 private:
00279 
00280   // @}
00281   // @{ Private methods
00282 
00284   Ifpack_ILU(const Ifpack_ILU& RHS) :
00285     Comm_(RHS.Comm()),
00286     Time_(RHS.Comm())
00287   {}
00288 
00290   Ifpack_ILU& operator=(const Ifpack_ILU& RHS)
00291   {
00292     return(*this);
00293   }
00294 
00296   void Destroy();
00297 
00299 
00309   int Solve(bool Trans, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00310 
00311   int ComputeSetup();
00312   int InitAllValues(const Epetra_RowMatrix & A, int MaxNumEntries);
00313 
00315   int LevelOfFill() const {return LevelOfFill_;}
00316 
00318   double RelaxValue() const {return RelaxValue_;}
00319 
00321   double AbsoluteThreshold() const {return Athresh_;}
00322 
00324   double RelativeThreshold() const {return Rthresh_;}
00325 
00327   int NumGlobalRows() const {return(Graph().NumGlobalRows());};
00328   
00330   int NumGlobalCols() const {return(Graph().NumGlobalCols());};
00331   
00333   int NumGlobalNonzeros() const {return(L().NumGlobalNonzeros()+U().NumGlobalNonzeros());};
00334   
00336   virtual int NumGlobalBlockDiagonals() const {return(Graph().NumGlobalBlockDiagonals());};
00337   
00339   int NumMyRows() const {return(Graph().NumMyRows());};
00340   
00342   int NumMyCols() const {return(Graph().NumMyCols());};
00343   
00345   int NumMyNonzeros() const {return(L().NumMyNonzeros()+U().NumMyNonzeros());};
00346   
00348   virtual int NumMyBlockDiagonals() const {return(Graph().NumMyBlockDiagonals());};
00349   
00351   virtual int NumMyDiagonals() const {return(NumMyDiagonals_);};
00352   
00354   int IndexBase() const {return(Graph().IndexBase());};
00355   
00357   const Ifpack_IlukGraph & Graph() const {return(*Graph_);};
00358   
00360   Epetra_RowMatrix& Matrix()
00361   {
00362     return(*A_);
00363   }
00364 
00365   // @}
00366   // @{ Internal data
00367   
00369   Epetra_RowMatrix* A_;
00370   Teuchos::RefCountPtr<Ifpack_IlukGraph> Graph_;
00371   Teuchos::RefCountPtr<Epetra_CrsGraph> CrsGraph_;
00372   Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
00373   Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
00374   Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
00375   const Epetra_Map * U_DomainMap_;
00376   const Epetra_Map * L_RangeMap_;
00377   const Epetra_Comm & Comm_;
00379   Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
00381   Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
00382   Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
00383   Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
00385   Teuchos::RefCountPtr<Epetra_Vector> D_;
00386   bool UseTranspose_;
00387 
00388   int NumMyDiagonals_;
00389   bool Allocated_;
00390   bool ValuesInitialized_;
00391   bool Factored_;
00393   double RelaxValue_;
00395   double Athresh_;
00397   double Rthresh_;
00399   double Condest_;
00401   int LevelOfFill_;
00403   bool IsInitialized_;
00405   bool IsComputed_;
00407   char Label_[160];
00409   int NumInitialize_;
00411   int NumCompute_;
00413   mutable int NumApplyInverse_;
00415   double InitializeTime_;
00417   double ComputeTime_;
00419   mutable double ApplyInverseTime_;
00421   double ComputeFlops_;
00423   mutable double ApplyInverseFlops_;
00425   mutable Epetra_Time Time_;
00426 
00427 };
00428 
00429 #endif // HAVE_IFPACK_TEUCHOS
00430 #endif /* IFPACK_ILU_H */

Generated on Thu Sep 18 12:37:21 2008 for Ifpack Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1