dune-istl  2.4.1
fastamg.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_ISTL_FASTAMG_HH
4 #define DUNE_ISTL_FASTAMG_HH
5 
6 #include <memory>
7 #include <dune/common/exceptions.hh>
8 #include <dune/common/typetraits.hh>
9 #include <dune/common/unused.hh>
13 #include <dune/istl/solvers.hh>
15 #include <dune/istl/superlu.hh>
16 #include <dune/istl/umfpack.hh>
17 #include <dune/istl/solvertype.hh>
18 #include <dune/istl/io.hh>
20 
21 #include "fastamgsmoother.hh"
22 
31 namespace Dune
32 {
33  namespace Amg
34  {
57  template<class M, class X, class PI=SequentialInformation, class A=std::allocator<X> >
58  class FastAMG : public Preconditioner<X,X>
59  {
60  public:
62  typedef M Operator;
69  typedef PI ParallelInformation;
74 
76  typedef X Domain;
78  typedef X Range;
81 
82  enum {
85  };
86 
94  FastAMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver,
95  const Parameters& parms,
96  bool symmetric=true);
97 
109  template<class C>
110  FastAMG(const Operator& fineOperator, const C& criterion,
111  const Parameters& parms=Parameters(),
112  bool symmetric=true,
113  const ParallelInformation& pinfo=ParallelInformation());
114 
118  FastAMG(const FastAMG& amg);
119 
120  ~FastAMG();
121 
123  void pre(Domain& x, Range& b);
124 
126  void apply(Domain& v, const Range& d);
127 
129  void post(Domain& x);
130 
135  template<class A1>
136  void getCoarsestAggregateNumbers(std::vector<std::size_t,A1>& cont);
137 
138  std::size_t levels();
139 
140  std::size_t maxlevels();
141 
151  {
152  matrices_->recalculateGalerkin(NegateSet<typename PI::OwnerSet>());
153  }
154 
159  bool usesDirectCoarseLevelSolver() const;
160 
161  private:
168  template<class C>
169  void createHierarchies(C& criterion, Operator& matrix,
170  const PI& pinfo);
171 
178  struct LevelContext
179  {
191  typename OperatorHierarchy::RedistributeInfoList::const_iterator redist;
195  typename OperatorHierarchy::AggregatesMapList::const_iterator aggregates;
211  std::size_t level;
212  };
213 
215  void mgc(LevelContext& levelContext, Domain& x, const Range& b);
216 
223  void presmooth(LevelContext& levelContext, Domain& x, const Range& b);
224 
231  void postsmooth(LevelContext& levelContext, Domain& x, const Range& b);
232 
239  void moveToFineLevel(LevelContext& levelContext, bool processedFineLevel,
240  Domain& fineX);
241 
246  bool moveToCoarseLevel(LevelContext& levelContext);
247 
252  void initIteratorsWithFineLevel(LevelContext& levelContext);
253 
255  std::shared_ptr<OperatorHierarchy> matrices_;
257  std::shared_ptr<CoarseSolver> solver_;
259  Hierarchy<Range,A>* rhs_;
261  Hierarchy<Domain,A>* lhs_;
263  Hierarchy<Domain,A>* residual_;
264 
268  typedef typename ScalarProductChooserType::ScalarProduct ScalarProduct;
269  typedef std::shared_ptr<ScalarProduct> ScalarProductPointer;
271  ScalarProductPointer scalarProduct_;
273  std::size_t gamma_;
275  std::size_t preSteps_;
277  std::size_t postSteps_;
278  std::size_t level;
279  bool buildHierarchy_;
280  bool symmetric;
281  bool coarsesolverconverged;
283  typedef std::shared_ptr<Smoother> SmootherPointer;
284  SmootherPointer coarseSmoother_;
286  std::size_t verbosity_;
287  };
288 
289  template<class M, class X, class PI, class A>
291  : matrices_(amg.matrices_), solver_(amg.solver_),
292  rhs_(), lhs_(), residual_(), scalarProduct_(amg.scalarProduct_),
293  gamma_(amg.gamma_), preSteps_(amg.preSteps_), postSteps_(amg.postSteps_),
294  symmetric(amg.symmetric), coarsesolverconverged(amg.coarsesolverconverged),
295  coarseSmoother_(amg.coarseSmoother_), verbosity_(amg.verbosity_)
296  {
297  if(amg.rhs_)
298  rhs_=new Hierarchy<Range,A>(*amg.rhs_);
299  if(amg.lhs_)
300  lhs_=new Hierarchy<Domain,A>(*amg.lhs_);
301  if(amg.residual_)
302  residual_=new Hierarchy<Domain,A>(*amg.residual_);
303  }
304 
305  template<class M, class X, class PI, class A>
306  FastAMG<M,X,PI,A>::FastAMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver,
307  const Parameters& parms, bool symmetric_)
308  : matrices_(&matrices), solver_(&coarseSolver),
309  rhs_(), lhs_(), residual_(), scalarProduct_(),
310  gamma_(parms.getGamma()), preSteps_(parms.getNoPreSmoothSteps()),
311  postSteps_(parms.getNoPostSmoothSteps()), buildHierarchy_(false),
312  symmetric(symmetric_), coarsesolverconverged(true),
313  coarseSmoother_(), verbosity_(parms.debugLevel())
314  {
315  if(preSteps_>1||postSteps_>1)
316  {
317  std::cerr<<"WARNING only one step of smoothing is supported!"<<std::endl;
318  preSteps_=postSteps_=0;
319  }
320  assert(matrices_->isBuilt());
321  static_assert(is_same<PI,SequentialInformation>::value,
322  "Currently only sequential runs are supported");
323  }
324  template<class M, class X, class PI, class A>
325  template<class C>
326  FastAMG<M,X,PI,A>::FastAMG(const Operator& matrix,
327  const C& criterion,
328  const Parameters& parms,
329  bool symmetric_,
330  const PI& pinfo)
331  : solver_(), rhs_(), lhs_(), residual_(), scalarProduct_(), gamma_(parms.getGamma()),
332  preSteps_(parms.getNoPreSmoothSteps()), postSteps_(parms.getNoPostSmoothSteps()),
333  buildHierarchy_(true),
334  symmetric(symmetric_), coarsesolverconverged(true),
335  coarseSmoother_(), verbosity_(criterion.debugLevel())
336  {
337  if(preSteps_>1||postSteps_>1)
338  {
339  std::cerr<<"WARNING only one step of smoothing is supported!"<<std::endl;
340  preSteps_=postSteps_=1;
341  }
342  static_assert(is_same<PI,SequentialInformation>::value,
343  "Currently only sequential runs are supported");
344  // TODO: reestablish compile time checks.
345  //static_assert(static_cast<int>(PI::category)==static_cast<int>(S::category),
346  // "Matrix and Solver must match in terms of category!");
347  createHierarchies(criterion, const_cast<Operator&>(matrix), pinfo);
348  }
349 
350  template<class M, class X, class PI, class A>
352  {
353  if(buildHierarchy_) {
354  if(solver_)
355  solver_.reset();
356  if(coarseSmoother_)
357  coarseSmoother_.reset();
358  }
359  if(lhs_)
360  delete lhs_;
361  lhs_=nullptr;
362  if(residual_)
363  delete residual_;
364  residual_=nullptr;
365  if(rhs_)
366  delete rhs_;
367  rhs_=nullptr;
368  }
369 
370  template<class M, class X, class PI, class A>
371  template<class C>
372  void FastAMG<M,X,PI,A>::createHierarchies(C& criterion, Operator& matrix,
373  const PI& pinfo)
374  {
375  Timer watch;
376  matrices_.reset(new OperatorHierarchy(matrix, pinfo));
377 
378  matrices_->template build<NegateSet<typename PI::OwnerSet> >(criterion);
379 
380  if(verbosity_>0 && matrices_->parallelInformation().finest()->communicator().rank()==0)
381  std::cout<<"Building Hierarchy of "<<matrices_->maxlevels()<<" levels took "<<watch.elapsed()<<" seconds."<<std::endl;
382 
383  if(buildHierarchy_ && matrices_->levels()==matrices_->maxlevels()) {
384  // We have the carsest level. Create the coarse Solver
385  typedef typename SmootherTraits<Smoother>::Arguments SmootherArgs;
386  SmootherArgs sargs;
387  sargs.iterations = 1;
388 
390  cargs.setArgs(sargs);
391  if(matrices_->redistributeInformation().back().isSetup()) {
392  // Solve on the redistributed partitioning
393  cargs.setMatrix(matrices_->matrices().coarsest().getRedistributed().getmat());
394  cargs.setComm(matrices_->parallelInformation().coarsest().getRedistributed());
395  }else{
396  cargs.setMatrix(matrices_->matrices().coarsest()->getmat());
397  cargs.setComm(*matrices_->parallelInformation().coarsest());
398  }
399 
400  coarseSmoother_.reset(ConstructionTraits<Smoother>::construct(cargs));
401  scalarProduct_.reset(ScalarProductChooserType::construct(cargs.getComm()));
402 
403 #if HAVE_SUPERLU|| HAVE_UMFPACK
404 #if HAVE_UMFPACK
405 #define DIRECTSOLVER UMFPack
406 #else
407 #define DIRECTSOLVER SuperLU
408 #endif
409  // Use superlu if we are purely sequential or with only one processor on the coarsest level.
410  if(is_same<ParallelInformation,SequentialInformation>::value // sequential mode
411  || matrices_->parallelInformation().coarsest()->communicator().size()==1 //parallel mode and only one processor
412  || (matrices_->parallelInformation().coarsest().isRedistributed()
413  && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()==1
414  && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()>0)) { // redistribute and 1 proc
415  if(verbosity_>0 && matrices_->parallelInformation().coarsest()->communicator().rank()==0)
416  std::cout<<"Using superlu"<<std::endl;
417  if(matrices_->parallelInformation().coarsest().isRedistributed())
418  {
419  if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
420  // We are still participating on this level
421  solver_.reset(new DIRECTSOLVER<typename M::matrix_type>(matrices_->matrices().coarsest().getRedistributed().getmat(), false, false));
422  else
423  solver_.reset();
424  }else
425  solver_.reset(new DIRECTSOLVER<typename M::matrix_type>(matrices_->matrices().coarsest()->getmat(), false, false));
426  }else
427 #undef DIRECTSOLVER
428 #endif
429  {
430  if(matrices_->parallelInformation().coarsest().isRedistributed())
431  {
432  if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
433  // We are still participating on this level
434  solver_.reset(new BiCGSTABSolver<X>(const_cast<M&>(matrices_->matrices().coarsest().getRedistributed()),
435  *scalarProduct_,
436  *coarseSmoother_, 1E-2, 1000, 0));
437  else
438  solver_.reset();
439  }else
440  solver_.reset(new BiCGSTABSolver<X>(const_cast<M&>(*matrices_->matrices().coarsest()),
441  *scalarProduct_,
442  *coarseSmoother_, 1E-2, 1000, 0));
443  }
444  }
445 
446  if(verbosity_>0 && matrices_->parallelInformation().finest()->communicator().rank()==0)
447  std::cout<<"Building Hierarchy of "<<matrices_->maxlevels()<<" levels took "<<watch.elapsed()<<" seconds."<<std::endl;
448  }
449 
450 
451  template<class M, class X, class PI, class A>
452  void FastAMG<M,X,PI,A>::pre(Domain& x, Range& b)
453  {
454  Timer watch, watch1;
455  // Detect Matrix rows where all offdiagonal entries are
456  // zero and set x such that A_dd*x_d=b_d
457  // Thus users can be more careless when setting up their linear
458  // systems.
459  typedef typename M::matrix_type Matrix;
460  typedef typename Matrix::ConstRowIterator RowIter;
461  typedef typename Matrix::ConstColIterator ColIter;
462  typedef typename Matrix::block_type Block;
463  Block zero;
464  zero=typename Matrix::field_type();
465 
466  const Matrix& mat=matrices_->matrices().finest()->getmat();
467  for(RowIter row=mat.begin(); row!=mat.end(); ++row) {
468  bool isDirichlet = true;
469  bool hasDiagonal = false;
470  ColIter diag;
471  for(ColIter col=row->begin(); col!=row->end(); ++col) {
472  if(row.index()==col.index()) {
473  diag = col;
474  hasDiagonal = false;
475  }else{
476  if(*col!=zero)
477  isDirichlet = false;
478  }
479  }
480  if(isDirichlet && hasDiagonal)
481  diag->solve(x[row.index()], b[row.index()]);
482  }
483  std::cout<<" Preprocessing Dirichlet took "<<watch1.elapsed()<<std::endl;
484  watch1.reset();
485  // No smoother to make x consistent! Do it by hand
486  matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
487  Range* copy = new Range(b);
488  if(rhs_)
489  delete rhs_;
490  rhs_ = new Hierarchy<Range,A>(copy);
491  Domain* dcopy = new Domain(x);
492  if(lhs_)
493  delete lhs_;
494  lhs_ = new Hierarchy<Domain,A>(dcopy);
495  dcopy = new Domain(x);
496  residual_ = new Hierarchy<Domain,A>(dcopy);
497  matrices_->coarsenVector(*rhs_);
498  matrices_->coarsenVector(*lhs_);
499  matrices_->coarsenVector(*residual_);
500 
501  // The preconditioner might change x and b. So we have to
502  // copy the changes to the original vectors.
503  x = *lhs_->finest();
504  b = *rhs_->finest();
505  }
506  template<class M, class X, class PI, class A>
508  {
509  return matrices_->levels();
510  }
511  template<class M, class X, class PI, class A>
513  {
514  return matrices_->maxlevels();
515  }
516 
518  template<class M, class X, class PI, class A>
519  void FastAMG<M,X,PI,A>::apply(Domain& v, const Range& d)
520  {
521  LevelContext levelContext;
522  // Init all iterators for the current level
523  initIteratorsWithFineLevel(levelContext);
524 
525  assert(v.two_norm()==0);
526 
527  level=0;
528  if(matrices_->maxlevels()==1){
529  // The coarse solver might modify the d!
530  Range b(d);
531  mgc(levelContext, v, b);
532  }else
533  mgc(levelContext, v, d);
534  if(postSteps_==0||matrices_->maxlevels()==1)
535  levelContext.pinfo->copyOwnerToAll(v, v);
536  }
537 
538  template<class M, class X, class PI, class A>
539  void FastAMG<M,X,PI,A>::initIteratorsWithFineLevel(LevelContext& levelContext)
540  {
541  levelContext.matrix = matrices_->matrices().finest();
542  levelContext.pinfo = matrices_->parallelInformation().finest();
543  levelContext.redist =
544  matrices_->redistributeInformation().begin();
545  levelContext.aggregates = matrices_->aggregatesMaps().begin();
546  levelContext.lhs = lhs_->finest();
547  levelContext.residual = residual_->finest();
548  levelContext.rhs = rhs_->finest();
549  levelContext.level=0;
550  }
551 
552  template<class M, class X, class PI, class A>
553  bool FastAMG<M,X,PI,A>
554  ::moveToCoarseLevel(LevelContext& levelContext)
555  {
556  bool processNextLevel=true;
557 
558  if(levelContext.redist->isSetup()) {
559  throw "bla";
560  levelContext.redist->redistribute(static_cast<const Range&>(*levelContext.residual),
561  levelContext.residual.getRedistributed());
562  processNextLevel = levelContext.residual.getRedistributed().size()>0;
563  if(processNextLevel) {
564  //restrict defect to coarse level right hand side.
565  ++levelContext.pinfo;
567  ::restrictVector(*(*levelContext.aggregates), *levelContext.rhs,
568  static_cast<const Range&>(levelContext.residual.getRedistributed()),
569  *levelContext.pinfo);
570  }
571  }else{
572  //restrict defect to coarse level right hand side.
573  ++levelContext.rhs;
574  ++levelContext.pinfo;
576  ::restrictVector(*(*levelContext.aggregates), *levelContext.rhs,
577  static_cast<const Range&>(*levelContext.residual), *levelContext.pinfo);
578  }
579 
580  if(processNextLevel) {
581  // prepare coarse system
582  ++levelContext.residual;
583  ++levelContext.lhs;
584  ++levelContext.matrix;
585  ++levelContext.level;
586  ++levelContext.redist;
587 
588  if(levelContext.matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()) {
589  // next level is not the globally coarsest one
590  ++levelContext.aggregates;
591  }
592  // prepare the lhs on the next level
593  *levelContext.lhs=0;
594  *levelContext.residual=0;
595  }
596  return processNextLevel;
597  }
598 
599  template<class M, class X, class PI, class A>
600  void FastAMG<M,X,PI,A>
601  ::moveToFineLevel(LevelContext& levelContext, bool processNextLevel, Domain& x)
602  {
603  if(processNextLevel) {
604  if(levelContext.matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()) {
605  // previous level is not the globally coarsest one
606  --levelContext.aggregates;
607  }
608  --levelContext.redist;
609  --levelContext.level;
610  //prolongate and add the correction (update is in coarse left hand side)
611  --levelContext.matrix;
612  --levelContext.residual;
613 
614  }
615 
616  typename Hierarchy<Domain,A>::Iterator coarseLhs = levelContext.lhs--;
617  if(levelContext.redist->isSetup()) {
618 
619  // Need to redistribute during prolongate
621  ::prolongateVector(*(*levelContext.aggregates), *coarseLhs, x,
622  levelContext.lhs.getRedistributed(),
623  matrices_->getProlongationDampingFactor(),
624  *levelContext.pinfo, *levelContext.redist);
625  }else{
627  ::prolongateVector(*(*levelContext.aggregates), *coarseLhs, x,
628  matrices_->getProlongationDampingFactor(), *levelContext.pinfo);
629 
630  // printvector(std::cout, *lhs, "prolongated coarse grid correction", "lhs", 10, 10, 10);
631  }
632 
633 
634  if(processNextLevel) {
635  --levelContext.rhs;
636  }
637 
638  }
639 
640 
641  template<class M, class X, class PI, class A>
642  void FastAMG<M,X,PI,A>
643  ::presmooth(LevelContext& levelContext, Domain& x, const Range& b)
644  {
646  x,
647  *levelContext.residual,
648  b);
649  }
650 
651  template<class M, class X, class PI, class A>
652  void FastAMG<M,X,PI,A>
653  ::postsmooth(LevelContext& levelContext, Domain& x, const Range& b)
654  {
656  ::apply(levelContext.matrix->getmat(), x, *levelContext.residual, b);
657  }
658 
659 
660  template<class M, class X, class PI, class A>
662  {
664  }
665 
666  template<class M, class X, class PI, class A>
667  void FastAMG<M,X,PI,A>::mgc(LevelContext& levelContext, Domain& v, const Range& b){
668 
669  if(levelContext.matrix == matrices_->matrices().coarsest() && levels()==maxlevels()) {
670  // Solve directly
672  res.converged=true; // If we do not compute this flag will not get updated
673  if(levelContext.redist->isSetup()) {
674  levelContext.redist->redistribute(b, levelContext.rhs.getRedistributed());
675  if(levelContext.rhs.getRedistributed().size()>0) {
676  // We are still participating in the computation
677  levelContext.pinfo.getRedistributed().copyOwnerToAll(levelContext.rhs.getRedistributed(),
678  levelContext.rhs.getRedistributed());
679  solver_->apply(levelContext.lhs.getRedistributed(), levelContext.rhs.getRedistributed(), res);
680  }
681  levelContext.redist->redistributeBackward(v, levelContext.lhs.getRedistributed());
682  levelContext.pinfo->copyOwnerToAll(v, v);
683  }else{
684  levelContext.pinfo->copyOwnerToAll(b, b);
685  solver_->apply(v, const_cast<Range&>(b), res);
686  }
687 
688  // printvector(std::cout, *lhs, "coarse level update", "u", 10, 10, 10);
689  // printvector(std::cout, *rhs, "coarse level rhs", "rhs", 10, 10, 10);
690  if (!res.converged)
691  coarsesolverconverged = false;
692  }else{
693  // presmoothing
694  presmooth(levelContext, v, b);
695  // printvector(std::cout, *lhs, "update", "u", 10, 10, 10);
696  // printvector(std::cout, *residual, "post presmooth residual", "r", 10);
697 #ifndef DUNE_AMG_NO_COARSEGRIDCORRECTION
698  bool processNextLevel = moveToCoarseLevel(levelContext);
699 
700  if(processNextLevel) {
701  // next level
702  for(std::size_t i=0; i<gamma_; i++)
703  mgc(levelContext, *levelContext.lhs, *levelContext.rhs);
704  }
705 
706  moveToFineLevel(levelContext, processNextLevel, v);
707 #else
708  *lhs=0;
709 #endif
710 
711  if(levelContext.matrix == matrices_->matrices().finest()) {
712  coarsesolverconverged = matrices_->parallelInformation().finest()->communicator().prod(coarsesolverconverged);
713  if(!coarsesolverconverged)
714  DUNE_THROW(MathError, "Coarse solver did not converge");
715  }
716 
717  // printvector(std::cout, *lhs, "update corrected", "u", 10, 10, 10);
718  // postsmoothing
719  postsmooth(levelContext, v, b);
720  // printvector(std::cout, *lhs, "update postsmoothed", "u", 10, 10, 10);
721 
722  }
723  }
724 
725 
727  template<class M, class X, class PI, class A>
728  void FastAMG<M,X,PI,A>::post(Domain& x)
729  {
730  DUNE_UNUSED_PARAMETER(x);
731  delete lhs_;
732  lhs_=nullptr;
733  delete rhs_;
734  rhs_=nullptr;
735  delete residual_;
736  residual_=nullptr;
737  }
738 
739  template<class M, class X, class PI, class A>
740  template<class A1>
741  void FastAMG<M,X,PI,A>::getCoarsestAggregateNumbers(std::vector<std::size_t,A1>& cont)
742  {
743  matrices_->getCoarsestAggregatesOnFinest(cont);
744  }
745 
746  } // end namespace Amg
747 } // end namespace Dune
748 
749 #endif
OperatorHierarchy::ParallelInformationHierarchy ParallelInformationHierarchy
The parallal data distribution hierarchy type.
Definition: fastamg.hh:73
M Operator
The matrix operator type.
Definition: fastamg.hh:62
OperatorHierarchy::AggregatesMapList::const_iterator aggregates
The iterator over the aggregates maps.
Definition: fastamg.hh:195
A fast (sequential) algebraic multigrid based on agglomeration that saves memory bandwidth.
Definition: fastamg.hh:58
void pre(Domain &x, Range &b)
Prepare the preconditioner.
Definition: fastamg.hh:452
static void apply(const M &A, X &x, Y &d, const Y &b)
Definition: fastamgsmoother.hh:53
Matrix & mat
Definition: matrixmatrix.hh:343
FastAMG(const OperatorHierarchy &matrices, CoarseSolver &coarseSolver, const Parameters &parms, bool symmetric=true)
Construct a new amg with a specific coarse solver.
Definition: fastamg.hh:306
InverseOperator< X, X > CoarseSolver
the type of the coarse solver.
Definition: fastamg.hh:80
Prolongation and restriction for amg.
std::size_t levels()
Definition: fastamg.hh:507
Define general preconditioner interface.
void recalculateHierarchy()
Recalculate the matrix hierarchy.
Definition: fastamg.hh:150
row_type::const_iterator ConstColIterator
Const iterator for the entries of each row.
Definition: matrix.hh:53
Provides a classes representing the hierarchies in AMG.
Hierarchy< Domain, A >::Iterator lhs
The iterator over the left hand side.
Definition: fastamg.hh:199
Definition: example.cc:34
static void restrictVector(const AggregatesMap< Vertex > &aggregates, Vector &coarse, const Vector &fine, T &comm)
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
void post(Domain &x)
Clean up.
Definition: fastamg.hh:728
Definition: solvertype.hh:13
PI ParallelInformation
The type of the parallel information. Either OwnerOverlapCommunication or another type describing the...
Definition: fastamg.hh:69
bool usesDirectCoarseLevelSolver() const
Check whether the coarse solver used is a direct solver.
Definition: fastamg.hh:661
Traits class for getting the attribute class of a smoother.
Definition: smoother.hh:63
A hierarchy of coantainers (e.g. matrices or vectors)
Definition: hierarchy.hh:66
Hierarchy< Domain, A >::Iterator residual
The iterator over the residuals.
Definition: fastamg.hh:203
The solver category.
Definition: fastamg.hh:84
T block_type
Export the type representing the components.
Definition: matrix.hh:32
A generic dynamic dense matrix.
Definition: matrix.hh:24
void apply(Domain &v, const Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition: fastamg.hh:519
void getCoarsestAggregateNumbers(std::vector< std::size_t, A1 > &cont)
Get the aggregate number of each unknown on the coarsest level.
Definition: fastamg.hh:741
LevelIterator< const Hierarchy< MatrixOperator, Allocator >, const MatrixOperator > ConstIterator
Type of the const iterator.
Definition: hierarchy.hh:260
Iterator finest()
Get an iterator positioned at the finest level.
Definition: hierarchy.hh:1379
Templates characterizing the type of a solver.
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition: smoother.hh:430
Hierarchy< Range, A >::Iterator rhs
The iterator over the right hand sided.
Definition: fastamg.hh:207
All parameters for AMG.
Definition: parameters.hh:390
Classes for using UMFPack with ISTL matrices.
Define base class for scalar product and norm.
std::size_t level
The level index.
Definition: fastamg.hh:211
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition: solvers.hh:540
LevelIterator< Hierarchy< ParallelInformation, Allocator >, ParallelInformation > Iterator
Type of the mutable iterator.
Definition: hierarchy.hh:257
~FastAMG()
Definition: fastamg.hh:351
MatrixHierarchy< M, ParallelInformation, A > OperatorHierarchy
The operator hierarchy type.
Definition: fastamg.hh:71
Statistics about the application of an inverse operator.
Definition: solver.hh:31
std::size_t maxlevels()
Definition: fastamg.hh:512
Classes for using SuperLU with ISTL matrices.
Implementations of the inverse operator interface.
Sequential SSOR preconditioner.
Definition: preconditioners.hh:127
ParallelInformationHierarchy::Iterator pinfo
The iterator over the parallel information.
Definition: fastamg.hh:187
RowIterator end()
Get iterator to one beyond last row.
Definition: matrix.hh:85
VariableBlockVector< T, A >::ConstIterator ConstRowIterator
Const iterator over the matrix rows.
Definition: matrix.hh:50
Traits class for generically constructing non default constructable types.
Definition: novlpschwarz.hh:326
Col col
Definition: matrixmatrix.hh:347
Classes for the generic construction and application of the smoothers.
const void * Arguments
A type holding all the arguments needed to call the constructor.
Definition: construction.hh:44
Category for sequential solvers.
Definition: solvercategory.hh:21
static void apply(const M &A, X &x, Y &d, const Y &b)
Definition: fastamgsmoother.hh:17
OperatorHierarchy::RedistributeInfoList::const_iterator redist
The iterator over the redistribution information.
Definition: fastamg.hh:191
bool converged
True if convergence criterion has been met.
Definition: solver.hh:56
The hierarchies build by the coarsening process.
Definition: hierarchy.hh:316
X Domain
The domain type.
Definition: fastamg.hh:76
T::field_type field_type
Export the type representing the underlying field.
Definition: matrix.hh:29
Some generic functions for pretty printing vectors and matrices.
static void prolongateVector(const AggregatesMap< Vertex > &aggregates, Vector &coarse, Vector &fine, Vector &fineRedist, T1 damp, R &redistributor=R())
X Range
The range type.
Definition: fastamg.hh:78
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition: smoother.hh:408
Choose the approriate scalar product for a solver category.
Definition: scalarproducts.hh:76
Definition: basearray.hh:19
OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator matrix
The iterator over the matrices.
Definition: fastamg.hh:183
RowIterator begin()
Get iterator to first row.
Definition: matrix.hh:79