• Main Page
  • Related Pages
  • Namespaces
  • Classes
  • Files
  • File List
  • File Members

ml_MultiLevelPreconditioner.h

Go to the documentation of this file.
00001 
00009 /* ******************************************************************** */
00010 /* See the file COPYRIGHT for a complete copyright notice, contact      */
00011 /* person and disclaimer.                                               */        
00012 /* ******************************************************************** */
00013 /*#############################################################################
00014 # CVS File Information
00015 #    Current revision: $Revision: 1.57 $
00016 #    Branch:           $Branch$
00017 #    Last modified:    $Date: 2008/08/01 18:14:20 $
00018 #    Modified by:      $Author: jhu $
00019 #############################################################################*/
00020 
00021 #ifndef ML_MULTILEVELPRECONDITIONER_H
00022 #define ML_MULTILEVELPRECONDITIONER_H
00023 
00024 #include "ml_include.h"
00025 
00026 #if defined(HAVE_ML_EPETRA) && defined(HAVE_ML_TEUCHOS) 
00027 // define the following to allow compilation without AztecOO
00028 #ifndef HAVE_ML_AZTECOO
00029 #ifndef AZ_PROC_SIZE
00030 #define AZ_PROC_SIZE 1
00031 #endif
00032 #ifndef AZ_OPTIONS_SIZE
00033 #define AZ_OPTIONS_SIZE 1
00034 #endif
00035 #ifndef AZ_PARAMS_SIZE
00036 #define AZ_PARAMS_SIZE 1
00037 #endif
00038 #ifndef AZ_STATUS_SIZE
00039 #define AZ_STATUS_SIZE 1
00040 #endif
00041 #endif
00042 
00043 class Epetra_Map;
00044 class Epetra_BlockMap;
00045 class Epetra_MultiVector;
00046 class Epetra_Comm;
00047 class Epetra_CrsMatrix;
00048 class Epetra_FECrsMatrix;
00049 class Epetra_VbrMatrix;
00050 
00051 #include "Epetra_SerialDenseMatrix.h"
00052 #include "Epetra_SerialDenseVector.h"
00053 #include "Epetra_SerialDenseSolver.h"
00054 
00055 #define ML_MEM_SIZE      20
00056 #define ML_MEM_INITIAL    0
00057 #define ML_MEM_FINAL      1
00058 #define ML_MEM_SMOOTHER   2
00059 #define ML_MEM_COARSE     3
00060 #define ML_MEM_HIERARCHY  4
00061 #define ML_MEM_PREC_FIRST 5
00062 #define ML_MEM_PREC_OTHER 6
00063 #define ML_MEM_TOT1       7
00064 #define ML_MEM_TOT2       8
00065 #define ML_MEM_INITIAL_MALLOC    10
00066 #define ML_MEM_FINAL_MALLOC      11
00067 #define ML_MEM_SMOOTHER_MALLOC   12
00068 #define ML_MEM_COARSE_MALLOC     13
00069 #define ML_MEM_HIERARCHY_MALLOC  14
00070 #define ML_MEM_PREC_FIRST_MALLOC 15
00071 #define ML_MEM_PREC_OTHER_MALLOC 16
00072 #define ML_MEM_TOT1_MALLOC       17
00073 #define ML_MEM_TOT2_MALLOC       18
00074 
00075 #include "Epetra_Operator.h"
00076 #include "Epetra_RowMatrix.h"
00077 #ifdef HAVE_ML_AZTECOO
00078 #include "Epetra_MultiVector.h"                                                 
00079 #include "Epetra_MsrMatrix.h"
00080 #endif
00081 #include "Teuchos_ParameterList.hpp"
00082 
00083 #ifdef HAVE_ML_EPETRAEXT
00084 #include "EpetraExt_SolverMap_CrsMatrix.h"
00085 #endif
00086 
00087 namespace ML_Epetra
00088 {
00089 
00091 
00115   int SetDefaults(std::string ProblemType, Teuchos::ParameterList & List,
00116       int * options = 0, double * params = 0, const bool OverWrite=true);
00117   
00119   int SetDefaultsDD(Teuchos::ParameterList & List, 
00120                    Teuchos::RCP<std::vector<int> > &options,
00121                        Teuchos::RCP<std::vector<double> > &params,
00122                        bool Overwrite=true);
00123   
00125   int SetDefaultsDD_LU(Teuchos::ParameterList & List, 
00126                    Teuchos::RCP<std::vector<int> > &options,
00127                        Teuchos::RCP<std::vector<double> > &params,
00128                        bool Overwrite=true);
00129   
00131   int SetDefaultsDD_3Levels(Teuchos::ParameterList & List, 
00132                    Teuchos::RCP<std::vector<int> > &options,
00133                        Teuchos::RCP<std::vector<double> > &params,
00134                        bool Overwrite=true);
00135   
00137   int SetDefaultsDD_3Levels_LU(Teuchos::ParameterList & List, 
00138                    Teuchos::RCP<std::vector<int> > &options,
00139                        Teuchos::RCP<std::vector<double> > &params,
00140                        bool Overwrite=true);
00141 
00143   int SetDefaultsMaxwell(Teuchos::ParameterList & List, 
00144                    Teuchos::RCP<std::vector<int> > &options,
00145                        Teuchos::RCP<std::vector<double> > &params,
00146                        bool Overwrite=true);
00147   
00149   int SetDefaultsSA(Teuchos::ParameterList & List, 
00150                    Teuchos::RCP<std::vector<int> > &options,
00151                        Teuchos::RCP<std::vector<double> > &params,
00152                        bool Overwrite=true);
00153   
00155   int SetDefaultsNSSA(Teuchos::ParameterList & List, 
00156                    Teuchos::RCP<std::vector<int> > &options,
00157                        Teuchos::RCP<std::vector<double> > &params,
00158                        bool Overwrite=true);
00159 
00161   int ReadXML(const string &FileName, Teuchos::ParameterList &List,
00162                    const Epetra_Comm &Comm);
00163 
00224 class MultiLevelPreconditioner : public virtual Epetra_Operator {
00225       
00226 public:  
00227 
00229 
00231 
00232   MultiLevelPreconditioner(const Epetra_RowMatrix & RowMatrix,
00233                            const bool ComputePrec = true);
00234 
00236   
00237   MultiLevelPreconditioner(const Epetra_RowMatrix & RowMatrix,
00238          const Teuchos::ParameterList & List,
00239          const bool ComputePrec = true);
00240 
00242   
00243   MultiLevelPreconditioner(ML_Operator* Operator,
00244          const Teuchos::ParameterList& List,
00245          const bool ComputePrec = true);
00246   
00248 
00262   MultiLevelPreconditioner(const Epetra_RowMatrix& EdgeMatrix,
00263          const Epetra_RowMatrix& GradMatrix,
00264          const Epetra_RowMatrix& NodeMatrix,
00265          const Teuchos::ParameterList& List,
00266          const bool ComputePrec = true,
00267                            const bool UseNodeMatrixForSmoother = false);
00268   
00270 
00284   MultiLevelPreconditioner(const Epetra_RowMatrix & CurlCurlMatrix,
00285              const Epetra_RowMatrix & MassMatrix,
00286              const Epetra_RowMatrix & TMatrix,
00287              const Epetra_RowMatrix & NodeMatrix,
00288              const Teuchos::ParameterList & List,
00289              const bool ComputePrec = true);
00290   
00291 #ifdef HAVE_ML_AZTECOO
00293 
00308   MultiLevelPreconditioner(const Epetra_MsrMatrix & EdgeMatrix,
00309                          ML_Operator * GradMatrix,
00310                          AZ_MATRIX * NodeMatrix,
00311                          int       * proc_config,
00312                          const Teuchos::ParameterList & List,
00313                          const bool ComputePrec = true);
00314 #endif
00315 
00317   
00319 
00321   virtual ~MultiLevelPreconditioner() {
00322     if (IsComputePreconditionerOK_) 
00323       DestroyPreconditioner(); 
00324   }
00325 
00327   
00329 
00331   const char* Label() const{return(Label_);};  
00332   
00334   void PrintUnused() const
00335   {
00336     List_.unused(std::cout);
00337   }
00338 
00340   void PrintUnused(std::ostream & os) const
00341   {
00342     List_.unused(os);
00343   }
00344 
00346 
00351   void PrintUnused(const int MyPID) const;
00352 
00354   Teuchos::ParameterList& GetList() 
00355   {
00356     return List_;
00357   }
00358 
00359   // Get a copy of the internally stored output list.
00360   Teuchos::ParameterList GetOutputList() 
00361   {
00362     return OutputList_;
00363   }
00364 
00366   void PrintList();
00367 
00369   int SetParameterList(const Teuchos::ParameterList& List);
00370 
00372   
00374 
00376   int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
00377     return(-1);}
00378 
00380   int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00381 
00383   
00385 
00386 
00388 
00407   int ComputePreconditioner(const bool CheckFiltering = false);
00408 
00410   int ReComputePreconditioner();
00411 
00413   void Print(const char *whichHierarchy = "main");
00414 
00415   int ComputeAdaptivePreconditioner(int TentativeNullSpaceSize,
00416                                     double* TentativeNullSpace);
00417 
00419   int IsPreconditionerComputed()  const
00420   {
00421     return(IsComputePreconditionerOK_);
00422   }
00423 
00424   // following functions are required to derive Epetra_RowMatrix objects.
00425 
00427   int SetOwnership(bool ownership){ ownership_ = ownership; return(-1);};
00428 
00430   int SetUseTranspose(bool UseTranspose){return(-1);}
00431 
00433   double NormInf() const {return(0.0);};
00434 
00436   bool UseTranspose() const {return(false);};
00437   
00439   bool HasNormInf() const{return(false);};
00440   
00442   const Epetra_Comm& Comm() const{return(*Comm_);};
00443   
00445   const Epetra_Map& OperatorDomainMap() const {return(*DomainMap_);};
00446   
00448   const Epetra_Map& OperatorRangeMap() const {return(*RangeMap_);};
00450 
00452   int DestroyPreconditioner();
00453 
00455   const Epetra_RowMatrix& RowMatrix() const
00456   {
00457     return(*RowMatrix_);
00458   }
00459 
00461   const Epetra_BlockMap& Map() const
00462   {
00463     return(RowMatrix_->Map());
00464   }
00465 
00467   int NumGlobalRows() const
00468   {
00469     return(RowMatrix_->NumGlobalRows());
00470   }
00471   
00473   int NumGlobalCols() const
00474   {
00475     return(RowMatrix_->NumGlobalCols());
00476   }
00477   
00479   int NumMyRows() const
00480   {
00481     return(RowMatrix_->NumMyRows());
00482   }
00483   
00485   int NumMyCols() const
00486   {
00487     return(RowMatrix_->NumMyCols());
00488   }
00489 
00491 
00502   int PrintStencil2D(const int nx, const int ny, 
00503          int NodeID = -1,
00504          const int EquationID = 0);
00505 
00507   int AnalyzeHierarchy(const bool AnalyzeMatrices, 
00508                        const int PreCycles, const int PostCycles,
00509                        const int MLCycles);
00510 
00512   int AnalyzeSmoothers(const int NumPreCycles = 1,
00513                        const int NumPostCycles = 1);
00514 
00516   int AnalyzeCoarse();
00517 
00519   int AnalyzeCycle(const int NumCycles = 1);
00520 
00522   int TestSmoothers(Teuchos::ParameterList& InputList,
00523         const bool IsSymmetric = false);
00524 
00526   int TestSmoothers(const bool IsSymmetric = false) {
00527     return(TestSmoothers(List_,IsSymmetric));
00528   }
00529 
00531   const ML* GetML(const int WhichML= -1) const
00532   {
00533     if (WhichML < 0)
00534       return ml_;
00535     else if (WhichML == 0)
00536       return ml_nodes_;
00537     else
00538       return(0);
00539   }
00540 
00541   bool SolvingMaxwell() const
00542   {
00543     return SolvingMaxwell_;
00544   }
00545 
00547   const ML_Aggregate* GetML_Aggregate() const 
00548   {
00549     return agg_;
00550   }
00551 
00553   int Visualize(bool VizAggre, bool VizPreSmoother,
00554     bool VizPostSmoother, bool VizCycle,
00555     int NumApplPreSmoother, int NumApplPostSmoother,
00556     int NumCycleSmoother);
00557 
00559   int VisualizeAggregates();
00560 
00562   int VisualizeSmoothers(int NumPrecCycles = 1,
00563        int NumPostCycles = 1);
00564 
00566   int VisualizeCycle(int NumCycles = 1);
00567 
00571   int CreateLabel();
00572 
00573   void ReportTime();
00574 
00576 
00577 private:
00578 
00580   MultiLevelPreconditioner(const MultiLevelPreconditioner & rhs) 
00581   {};
00582 
00584   MultiLevelPreconditioner & operator = (const MultiLevelPreconditioner & rhs)
00585   {
00586     return *this;
00587   };
00588 
00590 
00591   int Initialize();
00592 
00594   int SetSmoothers();
00595 
00597   int SetCoarse();
00598 
00600   int SetAggregation();
00601 
00603   int SetPreconditioner();
00604 
00606   int SetNullSpace();
00607 
00610   void CheckNullSpace();
00611 
00613   void Apply_BCsToGradient( const Epetra_RowMatrix & EdgeMatrix,
00614                             const Epetra_RowMatrix & T);
00615 
00617 
00625 #ifdef HAVE_ML_EPETRAEXT
00626   Epetra_RowMatrix* ModifyEpetraMatrixColMap( const Epetra_RowMatrix &A,
00627                                    EpetraExt::CrsMatrix_SolverMap &transform,
00628                                    const char* matrixName );
00629 #endif
00630 
00632   int SetSmoothingDamping();
00633 
00635   int SetSmoothingDampingClassic();
00636 
00637 #define OLD_AUX
00638 #ifdef OLD_AUX
00639   int CreateAuxiliaryMatrixCrs(Epetra_FECrsMatrix * & FakeMatrix);
00640 
00641   int CreateAuxiliaryMatrixVbr(Epetra_VbrMatrix * & FakeMatrix);
00642 #endif
00643 
00644   int SetupCoordinates();
00645 
00646   void PrintMem(char *fmt, int size, int, int);
00647 
00648   void PrintMemoryUsage();
00649 
00650   int SetFiltering();
00651 
00652   void RandomAndZero(double *, double *, int);
00653   
00655 
00657   bool CheckPreconditionerKrylov();
00658 
00659   void VectorNorms(double*, int, double*,double*);
00660 
00662 
00664   
00666   ML* ml_;
00668   //  ml_nodes_, one or all of which may be null.
00669   ML_Comm* ml_comm_;
00671   ML_Aggregate* agg_;
00673   char* Label_;
00675   //  of multiple instances of ML_Epetra::MultiLevelPreconditioner.
00676   std::string mlpLabel_;
00677 
00679   const Epetra_RowMatrix* RowMatrix_;
00681   bool IsComputePreconditionerOK_;
00682   
00684   int NumLevels_;
00686   const Epetra_Map* DomainMap_;
00688   const Epetra_Map* RangeMap_;
00690   const Epetra_Comm* Comm_;
00691   bool  ownership_;
00693   int   ProcConfig_[AZ_PROC_SIZE];
00695   Teuchos::RCP<std::vector<int> >    SmootherOptions_;
00697   Teuchos::RCP<std::vector<double> > SmootherParams_;
00699   double SmootherStatus_[AZ_STATUS_SIZE];
00700 
00702   Teuchos::ParameterList List_;
00704   Teuchos::ParameterList OutputList_;      
00705 
00707   int MaxLevels_;
00708 
00710   int CycleApplications_;
00711 
00713   bool ZeroStartingSolution_;
00714 
00716 
00722   std::vector<int> LevelID_;
00723 
00725   double* NullSpaceToFree_;              
00726 
00728   std::string PrintMsg_;
00730   char ErrorMsg_[80];
00732   bool verbose_;
00734   int NumPDEEqns_;
00736   int profileIterations_;
00737 
00739 
00741 
00743   bool SolvingMaxwell_;
00745   const Epetra_RowMatrix* EdgeMatrix_;
00747   const Epetra_RowMatrix* CurlCurlMatrix_;
00749   bool CreatedEdgeMatrix_;
00750   const Epetra_RowMatrix* MassMatrix_;
00752   const Epetra_RowMatrix* NodeMatrix_;
00754   ML_Operator* TtATMatrixML_;
00755   bool UseNodeMatrixForSmoother_;
00756   bool CreatedNodeMatrix_;
00758   ML_Operator* ML_Kn_;
00759   bool CreatedML_Kn_;
00761   const Epetra_RowMatrix* TMatrix_;
00762 #ifdef HAVE_ML_EPETRAEXT
00764   EpetraExt::CrsMatrix_SolverMap RowMatrixColMapTrans_;
00766   EpetraExt::CrsMatrix_SolverMap NodeMatrixColMapTrans_;
00768   EpetraExt::CrsMatrix_SolverMap TMatrixColMapTrans_;
00770   EpetraExt::CrsMatrix_SolverMap CurlCurlMatrixColMapTrans_;
00772   EpetraExt::CrsMatrix_SolverMap MassMatrixColMapTrans_;
00774   EpetraExt::CrsMatrix_SolverMap TtATMatrixColMapTrans_;
00775 #endif
00776   bool CreatedTMatrix_;
00777   ML_Operator* TMatrixML_;
00778   ML_Operator* TMatrixTransposeML_;
00779   ML_Operator** Tmat_array, ** Tmat_trans_array;
00780   ML_Operator** MassMatrix_array; // If curlcurl & mass are separate
00781   ML_Operator** CurlCurlMatrix_array;  // If curlcurl & mass are separate
00783   ML* ml_nodes_;
00784 
00785   void** nodal_args_,** edge_args_;
00786 
00788 
00790 
00791   int NumApplications_;
00793   double ApplicationTime_;
00794   bool FirstApplication_;
00795   //@ CPU time for first application
00796   double FirstApplicationTime_;
00798   int NumConstructions_;
00800   double ConstructionTime_;
00801 
00803   
00804   // other stuff for old ML's compatibility
00805   Epetra_CrsMatrix* RowMatrixAllocated_;
00806 
00807   bool AnalyzeMemory_;
00808   
00809   int memory_[ML_MEM_SIZE];
00810 
00811   // filtering stuff
00812   std::vector<double> flt_NullSpace_;
00813   ML* flt_ml_;
00814   ML_Aggregate* flt_agg_;
00815   
00816   // for reuse of preconditioning
00817   double RateOfConvergence_;
00818 
00819 }; // class MultiLevelPreconditioner
00820  
00821 } // namespace ML_Epetra
00822 
00823 #endif /* defined HAVE_ML_EPETRA and HAVE_ML_TEUCHOS */
00824 
00825 #endif /* define ML_MULTILEVELPRECONDITIONER_H */