21 #include "Ifpack_PointRelaxation.h"
22 #include "Ifpack_BlockRelaxation.h"
23 #include "Ifpack_DenseContainer.h"
24 #include "Ifpack_AdditiveSchwarz.h"
25 #include "Ifpack_ILU.h"
26 #include "Ifpack_ILUT.h"
27 #include "Ifpack_IC.h"
28 #include "Ifpack_ICT.h"
29 #include "Ifpack_Graph_Epetra_CrsGraph.h"
33 namespace Preconditioners
35 template<
typename Scalar>
36 IfpackPrecond<Scalar>::IfpackPrecond(
const char *cls,
const char *type)
46 template<
typename Scalar>
47 IfpackPrecond<Scalar>::IfpackPrecond(
const char *cls,
const char *type,
int overlap)
55 this->overlap = overlap;
58 template<
typename Scalar>
59 IfpackPrecond<Scalar>::IfpackPrecond(Ifpack_Preconditioner *ipc)
66 template<
typename Scalar>
67 IfpackPrecond<Scalar>::~IfpackPrecond()
69 if(owner)
delete prec;
72 template<
typename Scalar>
73 void IfpackPrecond<Scalar>::set_param(
const char *name,
const char *value)
75 ilist.set(name, value);
78 template<
typename Scalar>
79 void IfpackPrecond<Scalar>::set_param(
const char *name,
int value)
81 ilist.set(name, value);
84 template<
typename Scalar>
85 void IfpackPrecond<Scalar>::set_param(
const char *name,
double value)
87 ilist.set(name, value);
90 template<
typename Scalar>
91 void IfpackPrecond<Scalar>::create(Matrix<Scalar> *m)
93 EpetraMatrix<Scalar> *mt =
static_cast<EpetraMatrix<Scalar> *
>(m);
96 if(strcmp(cls,
"point-relax") == 0)
98 create_point_relax(mat, type);
102 else if(strcmp(cls,
"block-relax") == 0)
104 create_block_relax(mat, type);
107 else if(strcmp(cls,
"add-schwartz") == 0)
109 create_add_schwartz(mat, type, overlap);
115 template<
typename Scalar>
116 void IfpackPrecond<Scalar>::create_point_relax(EpetraMatrix<Scalar> *a,
const char *name)
118 prec =
new Ifpack_PointRelaxation(a->mat);
119 ilist.set(
"relaxation: type", name);
122 template<
typename Scalar>
123 void IfpackPrecond<Scalar>::create_block_relax(EpetraMatrix<Scalar> *a,
const char *name)
125 Teuchos::RCP<const Epetra_CrsGraph> rgraph = Teuchos::rcp(&(a->mat->Graph()));
126 Ifpack_Graph *graph =
new Ifpack_Graph_Epetra_CrsGraph(rgraph);
128 Ifpack_Partitioner *partitioner =
new Ifpack_GreedyPartitioner(graph);
130 Teuchos::ParameterList list;
131 list.set(
"partitioner: local parts", 1000);
132 partitioner->SetParameters(list);
133 partitioner->Compute();
135 prec =
new Ifpack_BlockRelaxation<Ifpack_DenseContainer>(a->mat);
136 ilist.set(
"relaxation: type", name);
141 template<
typename Scalar>
142 void IfpackPrecond<Scalar>::create_add_schwartz(EpetraMatrix<Scalar> *a,
const char *name,
int overlap)
144 if(strcasecmp(name,
"ilu") == 0)
146 prec =
new Ifpack_AdditiveSchwarz<Ifpack_ILU>(a->mat, overlap);
148 else if(strcasecmp(name,
"ilut") == 0)
150 prec =
new Ifpack_AdditiveSchwarz<Ifpack_ILUT>(a->mat, overlap);
152 else if(strcasecmp(name,
"ic") == 0)
154 prec =
new Ifpack_AdditiveSchwarz<Ifpack_IC>(a->mat, overlap);
156 else if(strcasecmp(name,
"ict") == 0)
158 prec =
new Ifpack_AdditiveSchwarz<Ifpack_ICT>(a->mat, overlap);
164 template<
typename Scalar>
165 int IfpackPrecond<Scalar>::initialize()
167 assert(prec != NULL);
168 return prec->Initialize();
171 template<
typename Scalar>
172 void IfpackPrecond<Scalar>::compute()
174 assert(prec != NULL);
178 template<
typename Scalar>
179 void IfpackPrecond<Scalar>::apply_params()
181 prec->SetParameters(ilist);
184 template<
typename Scalar>
185 int IfpackPrecond<Scalar>::ApplyInverse(
const Epetra_MultiVector &r, Epetra_MultiVector &z)
const
187 assert(prec != NULL);
188 return prec->ApplyInverse(r, z);
191 template<
typename Scalar>
192 const Epetra_Comm &IfpackPrecond<Scalar>::Comm()
const
194 return mat->mat->Comm();
197 template<
typename Scalar>
198 const Epetra_Map &IfpackPrecond<Scalar>::OperatorDomainMap()
const
200 return mat->mat->OperatorDomainMap();
203 template<
typename Scalar>
204 const Epetra_Map &IfpackPrecond<Scalar>::OperatorRangeMap()
const
206 return mat->mat->OperatorRangeMap();
208 template class HERMES_API IfpackPrecond<double>;
209 template class HERMES_API IfpackPrecond<std::complex<double> >;