1 #ifndef DUNE_AMG_KAMG_HH
2 #define DUNE_AMG_KAMG_HH
28 :
public Preconditioner<typename AMG::Domain,typename AMG::Range>
49 : amg_(amg), coarseSolver_(coarseSolver)
69 bool processFineLevel =
70 amg_.moveToCoarseLevel();
76 coarseSolver_->
apply(x, b, res);
80 amg_.moveToFineLevel(processFineLevel);
118 template<
class M,
class X,
class S,
class PI=SequentialInformation,
165 std::size_t preSmoothingSteps =1, std::size_t postSmoothingSteps = 1,
166 std::size_t maxLevelKrylovSteps = 3 ,
double minDefectReduction =1e-1);
187 std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1,
188 std::size_t maxLevelKrylovSteps=3,
double minDefectReduction=1e-1,
205 std::size_t maxLevelKrylovSteps;
208 double levelDefectReduction;
211 std::vector<typename Amg::ScalarProduct*> scalarproducts;
214 std::vector<KAmgTwoGrid<Amg>*> ksolvers;
217 template<
class M,
class X,
class S,
class P,
class K,
class A>
220 std::size_t gamma, std::size_t preSmoothingSteps,
221 std::size_t postSmoothingSteps,
222 std::size_t ksteps,
double reduction)
223 : amg(matrices, coarseSolver, smootherArgs, gamma, preSmoothingSteps,
224 postSmoothingSteps), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
227 template<
class M,
class X,
class S,
class P,
class K,
class A>
231 std::size_t preSmoothingSteps, std::size_t postSmoothingSteps,
232 std::size_t ksteps,
double reduction,
234 : amg(fineOperator, criterion, smootherArgs, gamma, preSmoothingSteps,
235 postSmoothingSteps, false, pinfo), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
239 template<
class M,
class X,
class S,
class P,
class K,
class A>
243 scalarproducts.reserve(amg.levels());
244 ksolvers.reserve(amg.levels());
246 typename OperatorHierarchy::ParallelMatrixHierarchy::Iterator
247 matrix = amg.matrices_->matrices().coarsest();
249 pinfo = amg.matrices_->parallelInformation().coarsest();
250 bool hasCoarsest=(amg.levels()==amg.maxlevels());
254 if(matrix==amg.matrices_->matrices().finest())
262 std::ostringstream s;
264 if(matrix!=amg.matrices_->matrices().finest())
266 scalarproducts.push_back(Amg::ScalarProductChooser::construct(*pinfo));
267 ks =
new KrylovSolver(*matrix, *(scalarproducts.back()),
268 *(ksolvers.back()), levelDefectReduction,
269 maxLevelKrylovSteps, 0);
273 if(matrix==amg.matrices_->matrices().finest())
279 template<
class M,
class X,
class S,
class P,
class K,
class A>
282 typedef typename std::vector<KAmgTwoGrid<Amg>*>::reverse_iterator KIterator;
283 typedef typename std::vector<ScalarProduct*>::iterator SIterator;
284 for(KIterator kiter = ksolvers.rbegin(); kiter != ksolvers.rend();
287 if((*kiter)->coarseSolver()!=amg.solver_)
288 delete (*kiter)->coarseSolver();
291 for(SIterator siter = scalarproducts.begin(); siter!=scalarproducts.end();
298 template<
class M,
class X,
class S,
class P,
class K,
class A>
301 amg.initIteratorsWithFineLevel();
302 if(ksolvers.size()==0)
306 amg.solver_->apply(x,tb,res);
308 ksolvers.back()->apply(x,b);
311 template<
class M,
class X,
class S,
class P,
class K,
class A>
314 return amg.maxlevels();