5#ifndef DUNE_AMG_KAMG_HH
6#define DUNE_AMG_KAMG_HH
69 *levelContext_->update=0;
70 *levelContext_->rhs = d;
71 *levelContext_->lhs = v;
73 presmooth(*levelContext_, amg_.preSteps_);
74 bool processFineLevel =
75 amg_.moveToCoarseLevel(*levelContext_);
77 if(processFineLevel) {
81 coarseSolver_->apply(x, b, res);
82 *levelContext_->update=x;
85 amg_.moveToFineLevel(*levelContext_, processFineLevel);
88 v=*levelContext_->update;
117 std::shared_ptr<InverseOperator<Domain,Range> > coarseSolver_;
119 std::shared_ptr<typename AMG::LevelContext> levelContext_;
137 template<
class M,
class X,
class S,
class PI=SequentialInformation,
184 std::size_t maxLevelKrylovSteps=3,
double minDefectReduction=1e-1);
202 std::size_t maxLevelKrylovSteps=3,
double minDefectReduction=1e-1,
219 std::size_t maxLevelKrylovSteps;
222 double levelDefectReduction;
225 std::vector<std::shared_ptr<typename Amg::ScalarProduct> > scalarproducts;
228 std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > > ksolvers;
232 template<
class M,
class X,
class S,
class P,
class K,
class A>
235 std::size_t ksteps,
double reduction)
236 : amg(matrices, coarseSolver, smootherArgs, params),
237 maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
241 template<
class M,
class X,
class S,
class P,
class K,
class A>
245 std::size_t ksteps,
double reduction,
247 : amg(fineOperator, criterion, smootherArgs, pinfo),
248 maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
252 template<
class M,
class X,
class S,
class P,
class K,
class A>
256 scalarproducts.reserve(amg.levels());
257 ksolvers.reserve(amg.levels());
259 typename OperatorHierarchy::ParallelMatrixHierarchy::Iterator
260 matrix = amg.matrices_->matrices().coarsest();
262 pinfo = amg.matrices_->parallelInformation().coarsest();
263 bool hasCoarsest=(amg.levels()==amg.maxlevels());
266 if(matrix==amg.matrices_->matrices().finest())
274 std::ostringstream s;
276 if(matrix!=amg.matrices_->matrices().finest())
279 std::shared_ptr<InverseOperator<Domain,Range> > ks =
280 std::shared_ptr<InverseOperator<Domain,Range> >(
new KrylovSolver(*matrix, *(scalarproducts.back()),
281 *(ksolvers.back()), levelDefectReduction,
282 maxLevelKrylovSteps, 0));
286 if(matrix==amg.matrices_->matrices().finest())
292 template<
class M,
class X,
class S,
class P,
class K,
class A>
299 template<
class M,
class X,
class S,
class P,
class K,
class A>
302 if(ksolvers.size()==0)
306 amg.solver_->apply(v,td,res);
309 typedef typename Amg::LevelContext LevelContext;
310 std::shared_ptr<LevelContext> levelContext(
new LevelContext);
311 amg.initIteratorsWithFineLevel(*levelContext);
312 typedef typename std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > >::iterator Iter;
313 for(Iter solver=ksolvers.begin(); solver!=ksolvers.end(); ++solver)
314 (*solver)->setLevelContext(levelContext);
315 ksolvers.back()->apply(v,d);
319 template<
class M,
class X,
class S,
class P,
class K,
class A>
322 return amg.maxlevels();
Define general preconditioner interface.
X Domain
The domain type.
Definition amg.hh:87
void apply(Domain &v, const Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition kamg.hh:300
M Operator
The matrix operator type.
Definition amg.hh:73
KAMG(OperatorHierarchy &matrices, CoarseSolver &coarseSolver, const SmootherArgs &smootherArgs, const Parameters &parms, std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1)
Construct a new amg with a specific coarse solver.
Definition kamg.hh:233
std::size_t maxlevels()
Definition kamg.hh:320
void post(Domain &x)
Clean up.
Definition kamg.hh:293
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition smoother.hh:406
PI ParallelInformation
The type of the parallel information. Either OwnerOverlapCommunication or another type describing the...
Definition amg.hh:80
X Range
The range type.
Definition amg.hh:89
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition smoother.hh:428
void pre(Domain &x, Range &b)
Prepare the preconditioner.
Definition kamg.hh:253
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition amg.hh:194
Definition allocator.hh:11
std::shared_ptr< ScalarProduct< X > > createScalarProduct(const Comm &comm, SolverCategory::Category category)
Definition scalarproducts.hh:242
an algebraic multigrid method using a Krylov-cycle.
Definition kamg.hh:140
AMG< M, X, S, PI, A > Amg
The type of the underlying AMG.
Definition kamg.hh:143
Amg::ParallelInformationHierarchy ParallelInformationHierarchy
The type of the hierarchy of parallel information.
Definition kamg.hh:161
Amg::Operator Operator
the type of the lineatr operator.
Definition kamg.hh:155
Amg::SmootherArgs SmootherArgs
The type of the arguments for construction of the smoothers.
Definition kamg.hh:153
Amg::ParallelInformation ParallelInformation
the type of the parallelinformation to use.
Definition kamg.hh:151
Amg::Domain Domain
the type of the domain.
Definition kamg.hh:157
K KrylovSolver
The type of the Krylov solver for the cycle.
Definition kamg.hh:145
Amg::Range Range
The type of the range.
Definition kamg.hh:159
Amg::OperatorHierarchy OperatorHierarchy
The type of the hierarchy of operators.
Definition kamg.hh:147
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition kamg.hh:166
Amg::CoarseSolver CoarseSolver
The type of the coarse solver.
Definition kamg.hh:149
Amg::ScalarProduct ScalarProduct
The type of the scalar product.
Definition kamg.hh:163
Two grid operator for AMG with Krylov cycle.
Definition kamg.hh:33
void pre(typename AMG::Domain &x, typename AMG::Range &b)
Prepare the preconditioner.
Definition kamg.hh:58
KAmgTwoGrid(AMG &amg, std::shared_ptr< InverseOperator< Domain, Range > > coarseSolver)
Constructor.
Definition kamg.hh:53
~KAmgTwoGrid()
Destructor.
Definition kamg.hh:110
void post(typename AMG::Domain &x)
Clean up.
Definition kamg.hh:62
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition kamg.hh:41
void setLevelContext(std::shared_ptr< typename AMG::LevelContext > p)
Set the level context pointer.
Definition kamg.hh:104
InverseOperator< Domain, Range > * coarseSolver()
Get a pointer to the coarse grid solver.
Definition kamg.hh:95
void apply(typename AMG::Domain &v, const typename AMG::Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition kamg.hh:66
Parallel algebraic multigrid based on agglomeration.
Definition amg.hh:65
LevelIterator< Hierarchy< ParallelInformation, Allocator >, ParallelInformation > Iterator
Definition hierarchy.hh:216
The hierarchies build by the coarsening process.
Definition matrixhierarchy.hh:61
All parameters for AMG.
Definition parameters.hh:416
The default class for the smoother arguments.
Definition smoother.hh:38
Base class for matrix free definition of preconditioners.
Definition preconditioner.hh:33
Base class for scalar product and norm computation.
Definition scalarproducts.hh:52
Statistics about the application of an inverse operator.
Definition solver.hh:50
Abstract base class for all solvers.
Definition solver.hh:101
Category
Definition solvercategory.hh:23
Generalized preconditioned conjugate gradient solver.
Definition solvers.hh:1307