00001
00011
00012
00013
00014
00015
00016 #ifndef ML_EDGE_MATRIX_FREE_PRECONDITIONER_H
00017 #define ML_EDGE_MATRIX_FREE_PRECONDITIONER_H
00018 #if defined(HAVE_ML_EPETRA) && defined(HAVE_ML_TEUCHOS)
00019 #include "Epetra_Comm.h"
00020 #include "Epetra_Map.h"
00021 #include "Epetra_Operator.h"
00022 #include "Epetra_Vector.h"
00023 #include "Epetra_IntVector.h"
00024 #include "Epetra_CrsMatrix.h"
00025 #include "Epetra_RowMatrix.h"
00026 #include "ml_Preconditioner.h"
00027 #include "Epetra_Operator_With_MatMat.h"
00028
00029 #ifdef HAVE_ML_IFPACK
00030 #include "Ifpack_Chebyshev.h"
00031 #endif
00032
00033 #include "ml_include.h"
00034 #include "ml_MultiLevelPreconditioner.h"
00035 #ifdef HAVE_ML_EPETRAEXT
00036 #include "EpetraExt_SolverMap_CrsMatrix.h"
00037 #endif
00038
00039 #ifdef HAVE_ML_IFPACK
00040 class Ifpack_Chebyshev;
00041 #endif
00042
00043 namespace ML_Epetra
00044 {
00045
00048 class EdgeMatrixFreePreconditioner: public virtual ML_Preconditioner
00049 {
00050 public:
00052
00053 EdgeMatrixFreePreconditioner(const Epetra_Operator_With_MatMat & Operator, const Epetra_Vector& Diagonal,
00054 const Epetra_CrsMatrix & D0_Matrix,const Epetra_CrsMatrix & D0_Clean_Matrix,
00055 const Epetra_CrsMatrix &TMT_Matrix,
00056 const int* BCedges, const int numBCedges,
00057 const Teuchos::ParameterList &List,const bool ComputePrec = true);
00059
00060
00062
00063 ~EdgeMatrixFreePreconditioner();
00065
00066
00068
00070
00089
00090 int ComputePreconditioner(const bool CheckFiltering = false);
00091
00093 int ReComputePreconditioner(){return(-1);}
00094
00096 int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
00097 return(-1);}
00098
00100 int ApplyInverse(const Epetra_MultiVector& B, Epetra_MultiVector& X) const;
00101
00103 void Print(const char *whichHierarchy = "main");
00104
00106 int DestroyPreconditioner();
00107
00109 int SetUseTranspose(bool UseTranspose){return(-1);}
00110
00112 double NormInf() const {return(0.0);};
00113
00115 bool UseTranspose() const {return(false);};
00116
00118 bool HasNormInf() const{return(false);};
00119
00121 const Epetra_Comm& Comm() const{return(*Comm_);};
00122
00124 const Epetra_Map& OperatorDomainMap() const {return(*EdgeDomainMap_);};
00125
00127 const Epetra_Map& OperatorRangeMap() const {return(*EdgeRangeMap_);};
00128
00130
00131 private:
00133
00135 int SetupSmoother();
00136
00138 Epetra_MultiVector * BuildNullspace();
00139
00141 int BuildProlongator(const Epetra_MultiVector & nullspace);
00142
00144 int FormCoarseMatrix();
00145
00146
00147
00149
00150
00152
00153 int dim;
00154
00156 ML_Comm* ml_comm_;
00157
00159 const Epetra_Operator_With_MatMat * Operator_;
00160
00162 const Epetra_CrsMatrix * D0_Matrix_;
00163
00165 const Epetra_CrsMatrix * D0_Clean_Matrix_;
00166
00168 const Epetra_CrsMatrix * TMT_Matrix_;
00169
00171 const int* BCedges_;
00172 const int numBCedges_;
00173
00175 Epetra_CrsMatrix * Prolongator_;
00176 #ifdef HAVE_ML_EPETRAEXT
00177 EpetraExt::CrsMatrix_SolverMap ProlongatorColMapTrans_;
00178 #endif
00179
00181 Epetra_Vector * InvDiagonal_;
00182
00184 Epetra_CrsMatrix * CoarseMatrix;
00185 ML_Operator * CoarseMat_ML;
00186
00188 MultiLevelPreconditioner * CoarsePC;
00189
00190 #ifdef HAVE_ML_IFPACK
00192 Ifpack_Chebyshev* Smoother_;
00193 #endif
00194
00196 const Epetra_Map* EdgeDomainMap_;
00198 const Epetra_Map* EdgeRangeMap_;
00200 const Epetra_Comm* Comm_;
00201
00202
00204 const Epetra_Map* NodeDomainMap_;
00206 const Epetra_Map* NodeRangeMap_;
00208 Epetra_Map* CoarseMap_;
00209
00211 int num_cycles;
00212 int MaxLevels;
00213 bool verbose_;
00214 bool very_verbose_;
00215 bool print_hierarchy;
00216
00217
00219 };
00220
00221 }
00222
00223
00224 #endif
00225
00226 #endif