Dune Core Modules (2.4.1)

kamg.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_AMG_KAMG_HH
4#define DUNE_AMG_KAMG_HH
5
7#include "amg.hh"
8
9namespace Dune
10{
11 namespace Amg
12 {
13
28 template<class AMG>
30 : public Preconditioner<typename AMG::Domain,typename AMG::Range>
31 {
33 typedef typename AMG::Domain Domain;
35 typedef typename AMG::Range Range;
36 public:
37
38 enum {
41 };
42
51 : amg_(amg), coarseSolver_(coarseSolver)
52 {}
53
55 void pre(typename AMG::Domain& x, typename AMG::Range& b)
56 {
58 }
59
61 void post(typename AMG::Domain& x)
62 {
64 }
65
67 void apply(typename AMG::Domain& v, const typename AMG::Range& d)
68 {
69 // Copy data
70 *levelContext_->update=0;
71 *levelContext_->rhs = d;
72 *levelContext_->lhs = v;
73
74 presmooth(*levelContext_, amg_.preSteps_);
75 bool processFineLevel =
76 amg_.moveToCoarseLevel(*levelContext_);
77
78 if(processFineLevel) {
79 typename AMG::Range b=*levelContext_->rhs;
80 typename AMG::Domain x=*levelContext_->update;
82 coarseSolver_->apply(x, b, res);
83 *levelContext_->update=x;
84 }
85
86 amg_.moveToFineLevel(*levelContext_, processFineLevel);
87
88 postsmooth(*levelContext_, amg_.postSteps_);
89 v=*levelContext_->update;
90 }
91
97 {
98 return coarseSolver_;
99 }
100
105 void setLevelContext(std::shared_ptr<typename AMG::LevelContext> p)
106 {
107 levelContext_=p;
108 }
109
112 {}
113
114 private:
116 AMG& amg_;
118 std::shared_ptr<InverseOperator<Domain,Range> > coarseSolver_;
120 std::shared_ptr<typename AMG::LevelContext> levelContext_;
121 };
122
123
124
138 template<class M, class X, class S, class PI=SequentialInformation,
139 class K=GeneralizedPCGSolver<X>, class A=std::allocator<X> >
140 class KAMG : public Preconditioner<X,X>
141 {
142 public:
146 typedef K KrylovSolver;
156 typedef typename Amg::Operator Operator;
158 typedef typename Amg::Domain Domain;
160 typedef typename Amg::Range Range;
164 typedef typename Amg::ScalarProduct ScalarProduct;
165
166 enum {
169 };
183 KAMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver,
184 const SmootherArgs& smootherArgs, std::size_t gamma,
185 std::size_t preSmoothingSteps =1, std::size_t postSmoothingSteps = 1,
186 std::size_t maxLevelKrylovSteps = 3 , double minDefectReduction =1e-1);
187
204 template<class C>
205 KAMG(const Operator& fineOperator, const C& criterion,
206 const SmootherArgs& smootherArgs, std::size_t gamma=1,
207 std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1,
208 std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1,
210
212 void pre(Domain& x, Range& b);
214 void post(Domain& x);
216 void apply(Domain& v, const Range& d);
217
218 std::size_t maxlevels();
219
220 private:
222 Amg amg;
223
225 std::size_t maxLevelKrylovSteps;
226
228 double levelDefectReduction;
229
231 std::vector<std::shared_ptr<typename Amg::ScalarProduct> > scalarproducts;
232
234 std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > > ksolvers;
235 };
236
237 template<class M, class X, class S, class P, class K, class A>
239 const SmootherArgs& smootherArgs,
240 std::size_t gamma, std::size_t preSmoothingSteps,
241 std::size_t postSmoothingSteps,
242 std::size_t ksteps, double reduction)
243 : amg(matrices, coarseSolver, smootherArgs, gamma, preSmoothingSteps,
244 postSmoothingSteps), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
245 {}
246
247 template<class M, class X, class S, class P, class K, class A>
248 template<class C>
249 KAMG<M,X,S,P,K,A>::KAMG(const Operator& fineOperator, const C& criterion,
250 const SmootherArgs& smootherArgs, std::size_t gamma,
251 std::size_t preSmoothingSteps, std::size_t postSmoothingSteps,
252 std::size_t ksteps, double reduction,
253 const ParallelInformation& pinfo)
254 : amg(fineOperator, criterion, smootherArgs, gamma, preSmoothingSteps,
255 postSmoothingSteps, false, pinfo), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
256 {}
257
258
259 template<class M, class X, class S, class P, class K, class A>
261 {
262 amg.pre(x,b);
263 scalarproducts.reserve(amg.levels());
264 ksolvers.reserve(amg.levels());
265
266 typename OperatorHierarchy::ParallelMatrixHierarchy::Iterator
267 matrix = amg.matrices_->matrices().coarsest();
269 pinfo = amg.matrices_->parallelInformation().coarsest();
270 bool hasCoarsest=(amg.levels()==amg.maxlevels());
271
272 if(hasCoarsest) {
273 if(matrix==amg.matrices_->matrices().finest())
274 return;
275 --matrix;
276 --pinfo;
277 ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, amg.solver_)));
278 }else
279 ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, std::shared_ptr<InverseOperator<Domain,Range> >())));
280
281 std::ostringstream s;
282
283 if(matrix!=amg.matrices_->matrices().finest())
284 while(true) {
285 scalarproducts.push_back(std::shared_ptr<typename Amg::ScalarProduct>(Amg::ScalarProductChooser::construct(*pinfo)));
286 std::shared_ptr<InverseOperator<Domain,Range> > ks =
287 std::shared_ptr<InverseOperator<Domain,Range> >(new KrylovSolver(*matrix, *(scalarproducts.back()),
288 *(ksolvers.back()), levelDefectReduction,
289 maxLevelKrylovSteps, 0));
290 ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, ks)));
291 --matrix;
292 --pinfo;
293 if(matrix==amg.matrices_->matrices().finest())
294 break;
295 }
296 }
297
298
299 template<class M, class X, class S, class P, class K, class A>
301 {
302 amg.post(x);
303
304 }
305
306 template<class M, class X, class S, class P, class K, class A>
308 {
309 if(ksolvers.size()==0)
310 {
311 Range td=d;
313 amg.solver_->apply(v,td,res);
314 }else
315 {
316 typedef typename Amg::LevelContext LevelContext;
317 std::shared_ptr<LevelContext> levelContext(new LevelContext);
318 amg.initIteratorsWithFineLevel(*levelContext);
319 typedef typename std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > >::iterator Iter;
320 for(Iter solver=ksolvers.begin(); solver!=ksolvers.end(); ++solver)
321 (*solver)->setLevelContext(levelContext);
322 ksolvers.back()->apply(v,d);
323 }
324 }
325
326 template<class M, class X, class S, class P, class K, class A>
327 std::size_t KAMG<M,X,S,P,K,A>::maxlevels()
328 {
329 return amg.maxlevels();
330 }
331
333 } // Amg
334} // Dune
335
336#endif
The AMG preconditioner.
Parallel algebraic multigrid based on agglomeration.
Definition: amg.hh:56
LevelIterator< Hierarchy< ParallelInformation, Allocator >, ParallelInformation > Iterator
Type of the mutable iterator.
Definition: hierarchy.hh:257
an algebraic multigrid method using a Krylov-cycle.
Definition: kamg.hh:141
Amg::Domain Domain
the type of the domain.
Definition: kamg.hh:158
Amg::SmootherArgs SmootherArgs
The type of the arguments for construction of the smoothers.
Definition: kamg.hh:154
Amg::ParallelInformation ParallelInformation
the type of the parallelinformation to use.
Definition: kamg.hh:152
Amg::CoarseSolver CoarseSolver
The type of the coarse solver.
Definition: kamg.hh:150
Amg::OperatorHierarchy OperatorHierarchy
The type of the hierarchy of operators.
Definition: kamg.hh:148
Amg::Range Range
The type of the range.
Definition: kamg.hh:160
Amg::ScalarProduct ScalarProduct
The type of the scalar product.
Definition: kamg.hh:164
AMG< M, X, S, PI, A > Amg
The type of the underlying AMG.
Definition: kamg.hh:144
Amg::Operator Operator
the type of the lineatr operator.
Definition: kamg.hh:156
Amg::ParallelInformationHierarchy ParallelInformationHierarchy
The type of the hierarchy of parallel information.
Definition: kamg.hh:162
K KrylovSolver
The type of the Krylov solver for the cycle.
Definition: kamg.hh:146
@ category
The solver category.
Definition: kamg.hh:168
Two grid operator for AMG with Krylov cycle.
Definition: kamg.hh:31
void pre(typename AMG::Domain &x, typename AMG::Range &b)
Prepare the preconditioner.
Definition: kamg.hh:55
KAmgTwoGrid(AMG &amg, std::shared_ptr< InverseOperator< Domain, Range > > coarseSolver)
Constructor.
Definition: kamg.hh:50
~KAmgTwoGrid()
Destructor.
Definition: kamg.hh:111
void post(typename AMG::Domain &x)
Clean up.
Definition: kamg.hh:61
void setLevelContext(std::shared_ptr< typename AMG::LevelContext > p)
Set the level context pointer.
Definition: kamg.hh:105
@ category
The solver category.
Definition: kamg.hh:40
InverseOperator< Domain, Range > * coarseSolver()
Get a pointer to the coarse grid solver.
Definition: kamg.hh:96
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:67
The hierarchies build by the coarsening process.
Definition: hierarchy.hh:317
Generalized preconditioned conjugate gradient solver.
Definition: solvers.hh:1456
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
void apply(Domain &v, const Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition: kamg.hh:307
X Domain
The domain type.
Definition: amg.hh:78
SmootherTraits< Smoother >::Arguments SmootherArgs
The argument type for the construction of the smoother.
Definition: amg.hh:91
M Operator
The matrix operator type.
Definition: amg.hh:64
void post(Domain &x)
Clean up.
Definition: kamg.hh:300
X Range
The range type.
Definition: amg.hh:80
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition: smoother.hh:408
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition: smoother.hh:430
KAMG(const OperatorHierarchy &matrices, CoarseSolver &coarseSolver, const SmootherArgs &smootherArgs, std::size_t gamma, std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1, std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1)
Construct a new amg with a specific coarse solver.
Definition: kamg.hh:238
void pre(Domain &x, Range &b)
Prepare the preconditioner.
Definition: kamg.hh:260
PI ParallelInformation
The type of the parallel information. Either OwnerOverlapCommunication or another type describing the...
Definition: amg.hh:71
@ category
The solver category.
Definition: amg.hh:95
Dune namespace.
Definition: alignment.hh:10
Define general preconditioner interface.
Statistics about the application of an inverse operator.
Definition: solver.hh:32
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentional unused function parameters with.
Definition: unused.hh:18
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Dec 22, 23:30, 2024)