dune-istl  2.4
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 
9 namespace Dune
10 {
11  namespace Amg
12  {
13 
28  template<class AMG>
29  class KAmgTwoGrid
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  {
57  DUNE_UNUSED_PARAMETER(x); DUNE_UNUSED_PARAMETER(b);
58  }
59 
61  void post(typename AMG::Domain& x)
62  {
63  DUNE_UNUSED_PARAMETER(x);
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 
135  template<class M, class X, class S, class PI=SequentialInformation,
136  class K=GeneralizedPCGSolver<X>, class A=std::allocator<X> >
137  class KAMG : public Preconditioner<X,X>
138  {
139  public:
143  typedef K KrylovSolver;
147  typedef typename Amg::CoarseSolver CoarseSolver;
151  typedef typename Amg::SmootherArgs SmootherArgs;
153  typedef typename Amg::Operator Operator;
155  typedef typename Amg::Domain Domain;
157  typedef typename Amg::Range Range;
161  typedef typename Amg::ScalarProduct ScalarProduct;
162 
163  enum {
166  };
180  KAMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver,
181  const SmootherArgs& smootherArgs, std::size_t gamma,
182  std::size_t preSmoothingSteps =1, std::size_t postSmoothingSteps = 1,
183  std::size_t maxLevelKrylovSteps = 3 , double minDefectReduction =1e-1);
184 
201  template<class C>
202  KAMG(const Operator& fineOperator, const C& criterion,
203  const SmootherArgs& smootherArgs, std::size_t gamma=1,
204  std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1,
205  std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1,
206  const ParallelInformation& pinfo=ParallelInformation());
207 
209  void pre(Domain& x, Range& b);
211  void post(Domain& x);
213  void apply(Domain& v, const Range& d);
214 
215  std::size_t maxlevels();
216 
217  private:
219  Amg amg;
220 
222  std::size_t maxLevelKrylovSteps;
223 
225  double levelDefectReduction;
226 
228  std::vector<std::shared_ptr<typename Amg::ScalarProduct> > scalarproducts;
229 
231  std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > > ksolvers;
232  };
233 
234  template<class M, class X, class S, class P, class K, class A>
235  KAMG<M,X,S,P,K,A>::KAMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver,
236  const SmootherArgs& smootherArgs,
237  std::size_t gamma, std::size_t preSmoothingSteps,
238  std::size_t postSmoothingSteps,
239  std::size_t ksteps, double reduction)
240  : amg(matrices, coarseSolver, smootherArgs, gamma, preSmoothingSteps,
241  postSmoothingSteps), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
242  {}
243 
244  template<class M, class X, class S, class P, class K, class A>
245  template<class C>
246  KAMG<M,X,S,P,K,A>::KAMG(const Operator& fineOperator, const C& criterion,
247  const SmootherArgs& smootherArgs, std::size_t gamma,
248  std::size_t preSmoothingSteps, std::size_t postSmoothingSteps,
249  std::size_t ksteps, double reduction,
250  const ParallelInformation& pinfo)
251  : amg(fineOperator, criterion, smootherArgs, gamma, preSmoothingSteps,
252  postSmoothingSteps, false, pinfo), maxLevelKrylovSteps(ksteps), levelDefectReduction(reduction)
253  {}
254 
255 
256  template<class M, class X, class S, class P, class K, class A>
257  void KAMG<M,X,S,P,K,A>::pre(Domain& x, Range& b)
258  {
259  amg.pre(x,b);
260  scalarproducts.reserve(amg.levels());
261  ksolvers.reserve(amg.levels());
262 
263  typename OperatorHierarchy::ParallelMatrixHierarchy::Iterator
264  matrix = amg.matrices_->matrices().coarsest();
266  pinfo = amg.matrices_->parallelInformation().coarsest();
267  bool hasCoarsest=(amg.levels()==amg.maxlevels());
268 
269  if(hasCoarsest) {
270  if(matrix==amg.matrices_->matrices().finest())
271  return;
272  --matrix;
273  --pinfo;
274  ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, amg.solver_)));
275  }else
276  ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, std::shared_ptr<InverseOperator<Domain,Range> >())));
277 
278  std::ostringstream s;
279 
280  if(matrix!=amg.matrices_->matrices().finest())
281  while(true) {
282  scalarproducts.push_back(std::shared_ptr<typename Amg::ScalarProduct>(Amg::ScalarProductChooser::construct(*pinfo)));
283  std::shared_ptr<InverseOperator<Domain,Range> > ks =
284  std::shared_ptr<InverseOperator<Domain,Range> >(new KrylovSolver(*matrix, *(scalarproducts.back()),
285  *(ksolvers.back()), levelDefectReduction,
286  maxLevelKrylovSteps, 0));
287  ksolvers.push_back(std::shared_ptr<KAmgTwoGrid<Amg> >(new KAmgTwoGrid<Amg>(amg, ks)));
288  --matrix;
289  --pinfo;
290  if(matrix==amg.matrices_->matrices().finest())
291  break;
292  }
293  }
294 
295 
296  template<class M, class X, class S, class P, class K, class A>
297  void KAMG<M,X,S,P,K,A>::post(Domain& x)
298  {
299  amg.post(x);
300 
301  }
302 
303  template<class M, class X, class S, class P, class K, class A>
304  void KAMG<M,X,S,P,K,A>::apply(Domain& v, const Range& d)
305  {
306  if(ksolvers.size()==0)
307  {
308  Range td=d;
310  amg.solver_->apply(v,td,res);
311  }else
312  {
313  typedef typename Amg::LevelContext LevelContext;
314  std::shared_ptr<LevelContext> levelContext(new LevelContext);
315  amg.initIteratorsWithFineLevel(*levelContext);
316  typedef typename std::vector<std::shared_ptr<KAmgTwoGrid<Amg> > >::iterator Iter;
317  for(Iter solver=ksolvers.begin(); solver!=ksolvers.end(); ++solver)
318  (*solver)->setLevelContext(levelContext);
319  ksolvers.back()->apply(v,d);
320  }
321  }
322 
323  template<class M, class X, class S, class P, class K, class A>
325  {
326  return amg.maxlevels();
327  }
328 
330  } // Amg
331 } // Dune
332 
333 #endif
~KAmgTwoGrid()
Destructor.
Definition: kamg.hh:111
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition: smoother.hh:408
Amg::OperatorHierarchy OperatorHierarchy
The type of the hierarchy of operators.
Definition: kamg.hh:145
Definition: basearray.hh:19
AMG< M, X, S, PI, A > Amg
The type of the underlying AMG.
Definition: kamg.hh:141
Amg::ScalarProduct ScalarProduct
The type of the scalar product.
Definition: kamg.hh:161
void apply(Domain &v, const Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition: kamg.hh:304
Two grid operator for AMG with Krylov cycle.
Definition: amg.hh:44
Amg::ParallelInformationHierarchy ParallelInformationHierarchy
The type of the hierarchy of parallel information.
Definition: kamg.hh:159
Amg::Domain Domain
the type of the domain.
Definition: kamg.hh:155
X Range
The range type.
Definition: amg.hh:80
M Operator
The matrix operator type.
Definition: amg.hh:64
void pre(typename AMG::Domain &x, typename AMG::Range &b)
Prepare the preconditioner.
Definition: kamg.hh:55
void post(typename AMG::Domain &x)
Clean up.
Definition: kamg.hh:61
Matrix & A
Definition: matrixmatrix.hh:216
K KrylovSolver
The type of the Krylov solver for the cycle.
Definition: kamg.hh:143
X Domain
The domain type.
Definition: amg.hh:78
The solver category.
Definition: amg.hh:95
PI ParallelInformation
The type of the parallel information. Either OwnerOverlapCommunication or another type describing the...
Definition: amg.hh:71
KAmgTwoGrid(AMG &amg, std::shared_ptr< InverseOperator< Domain, Range > > coarseSolver)
Constructor.
Definition: kamg.hh:50
void pre(Domain &x, Range &b)
Prepare the preconditioner.
Definition: kamg.hh:257
SmootherTraits< Smoother >::Arguments SmootherArgs
The argument type for the construction of the smoother.
Definition: amg.hh:91
The solver category.
Definition: kamg.hh:165
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
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
InverseOperator< Domain, Range > * coarseSolver()
Get a pointer to the coarse grid solver.
Definition: kamg.hh:96
Abstract base class for all solvers.
Definition: solver.hh:79
The AMG preconditioner.
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:235
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition: smoother.hh:430
Amg::Operator Operator
the type of the lineatr operator.
Definition: kamg.hh:153
Amg::Range Range
The type of the range.
Definition: kamg.hh:157
Amg::CoarseSolver CoarseSolver
The type of the coarse solver.
Definition: kamg.hh:147
std::size_t maxlevels()
Definition: kamg.hh:324
Amg::SmootherArgs SmootherArgs
The type of the arguments for construction of the smoothers.
Definition: kamg.hh:151
Parallel algebraic multigrid based on agglomeration.
Definition: amg.hh:55
void setLevelContext(std::shared_ptr< typename AMG::LevelContext > p)
Set the level context pointer.
Definition: kamg.hh:105
void post(Domain &x)
Clean up.
Definition: kamg.hh:297
Generalized preconditioned conjugate gradient solver.
Definition: solvers.hh:1445
The solver category.
Definition: kamg.hh:40
Statistics about the application of an inverse operator.
Definition: solver.hh:31
Define general preconditioner interface.
The hierarchies build by the coarsening process.
Definition: hierarchy.hh:316
Amg::ParallelInformation ParallelInformation
the type of the parallelinformation to use.
Definition: kamg.hh:149
LevelIterator< Hierarchy< ParallelInformation, Allocator >, ParallelInformation > Iterator
Type of the mutable iterator.
Definition: hierarchy.hh:257