00001
00011
00012
00013
00014
00015
00016 #ifndef ML_REFMAXWELL_H
00017 #define ML_REFMAXWELL_H
00018
00019 #define RefMaxwellPreconditioner ML_RMP
00020
00021
00022 #define ENABLE_MS_MATRIX
00023
00024 #if defined(HAVE_ML_EPETRA) && defined(HAVE_ML_TEUCHOS)
00025 #include "ml_common.h"
00026 #include "Epetra_Comm.h"
00027 #include "Epetra_Map.h"
00028 #include "Epetra_Vector.h"
00029 #include "Epetra_CrsMatrix.h"
00030 #include "ml_RefMaxwell_11_Operator.h"
00031 #include "ml_Preconditioner.h"
00032 #include "ml_MultiLevelPreconditioner.h"
00033
00034 #ifdef HAVE_ML_EPETRAEXT
00035 #include "EpetraExt_SolverMap_CrsMatrix.h"
00036 #endif
00037
00038 namespace ML_Epetra
00039 {
00040
00041 int UpdateList(Teuchos::ParameterList &source, Teuchos::ParameterList &dest, bool OverWrite=true);
00043 int SetDefaultsRefMaxwell(Teuchos::ParameterList & inList,bool OverWrite=true);
00044
00045
00046
00059 class RefMaxwellPreconditioner: public virtual ML_Preconditioner
00060 {
00061 public:
00062
00064
00066
00067
00068 RefMaxwellPreconditioner(const Epetra_CrsMatrix& SM_Matrix,
00069 const Epetra_CrsMatrix& D0_Clean_Matrix,
00070 #ifdef ENABLE_MS_MATRIX
00071 const Epetra_CrsMatrix& Ms_Matrix,
00072 #endif
00073 const Epetra_CrsMatrix& M0inv_Matrix,
00074 const Epetra_CrsMatrix& M1_Matrix,
00075
00076 const Teuchos::ParameterList& List,
00077 const bool ComputePrec = true);
00079
00080
00082
00083
00084 ~RefMaxwellPreconditioner();
00086
00087
00089
00091 int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
00092 return(-1);}
00093
00095 int ApplyInverse(const Epetra_MultiVector& B, Epetra_MultiVector& X) const;
00096
00098
00099
00101
00103
00122
00123 int ComputePreconditioner(const bool CheckFiltering = false);
00124
00126 int ReComputePreconditioner(){return(-1);}
00127
00129 void Print(const char *whichHierarchy = "11");
00130
00132 int DestroyPreconditioner();
00133
00135 int SetUseTranspose(bool UseTranspose){return(-1);}
00136
00138 double NormInf() const {return(0.0);};
00139
00141 bool UseTranspose() const {return(false);};
00142
00144 bool HasNormInf() const{return(false);};
00145
00147 const Epetra_Comm& Comm() const{return(*Comm_);};
00148
00150 const Epetra_Map& OperatorDomainMap() const {return(*DomainMap_);};
00151
00153 const Epetra_Map& OperatorRangeMap() const {return(*RangeMap_);};
00155
00156
00157 private:
00158
00160
00162 int SetEdgeSmoother(Teuchos::ParameterList &List_);
00163
00165 int ApplyInverse_Implicit_212(const Epetra_MultiVector& B, Epetra_MultiVector& X) const;
00166
00168 int ApplyInverse_Implicit_Additive(const Epetra_MultiVector& B, Epetra_MultiVector& X) const;
00169
00171 int ApplyInverse_Implicit_121(const Epetra_MultiVector& B, Epetra_MultiVector& X) const;
00172
00173
00175
00177
00179 const Epetra_CrsMatrix * SM_Matrix_;
00181 Epetra_CrsMatrix * D0_Matrix_;
00183 const Epetra_CrsMatrix * D0_Clean_Matrix_;
00184 #ifdef ENABLE_MS_MATRIX
00186 const Epetra_CrsMatrix * Ms_Matrix_;
00187 #endif
00189 const Epetra_CrsMatrix * M0inv_Matrix_;
00191 Epetra_CrsMatrix * M1_Matrix_;
00193 Epetra_CrsMatrix * TMT_Matrix_;
00195 Epetra_CrsMatrix * TMT_Agg_Matrix_;
00196
00197
00198 #ifdef HAVE_ML_EPETRAEXT
00200 EpetraExt::CrsMatrix_SolverMap SM_Matrix_Trans_;
00202 EpetraExt::CrsMatrix_SolverMap D0_Matrix_Trans_;
00204 EpetraExt::CrsMatrix_SolverMap D0_Clean_Matrix_Trans_;
00206 #ifdef ENABLE_MS_MATRIX
00207 EpetraExt::CrsMatrix_SolverMap Ms_Matrix_Trans_;
00208 #endif
00210 EpetraExt::CrsMatrix_SolverMap M0inv_Matrix_Trans_;
00212 EpetraExt::CrsMatrix_SolverMap M1_Matrix_Trans_;
00214 EpetraExt::CrsMatrix_SolverMap TMT_Matrix_Trans_;
00216 EpetraExt::CrsMatrix_SolverMap TMT_Agg_Matrix_Trans_;
00217 #endif
00218
00219
00220
00222 int* BCrows;
00223 int numBCrows;
00224 bool HasOnlyDirichletNodes;
00225
00227 Epetra_Vector* Diagonal_;
00229 ML_RefMaxwell_11_Operator* Operator11_;
00230
00232 ML_Preconditioner * EdgePC;
00234
00235 MultiLevelPreconditioner * NodePC;
00236
00238 MultiLevelPreconditioner *PreEdgeSmoother;
00239 MultiLevelPreconditioner *PostEdgeSmoother;
00240
00242 string mode;
00243
00245 bool aggregate_with_sigma;
00246
00248 bool lump_m1;
00249
00250
00251 bool use_local_nodal_solver;
00252 MultiLevelPreconditioner *LocalNodalSolver;
00253 const Epetra_CrsMatrix* NodesToLocalNodes;
00254 Epetra_CrsMatrix* LocalNodalMatrix;
00255
00257 const Epetra_Map* DomainMap_;
00259 const Epetra_Map* RangeMap_;
00261 const Epetra_Map* NodeMap_;
00262
00264 const Epetra_Comm* Comm_;
00265
00267 bool verbose_;
00268
00270 bool very_verbose_;
00271
00273 bool print_hierarchy;
00275
00276
00278
00279 int NumApplications_;
00281 double ApplicationTime_;
00282 bool FirstApplication_;
00283
00284 double FirstApplicationTime_;
00286 int NumConstructions_;
00288 double ConstructionTime_;
00290 };
00291 }
00292
00293 #endif
00294 #endif