dune-pdelab  2.5-dev
seqistlsolverbackend.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=8 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
4 #define DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
5 
6 // this is here for backwards compatibility and deprecation warnings, remove after 2.5.0
7 #include "ensureistlinclude.hh"
8 
9 #include <dune/common/deprecated.hh>
10 #include <dune/common/parallel/mpihelper.hh>
11 
12 #include <dune/istl/owneroverlapcopy.hh>
13 #include <dune/istl/solvercategory.hh>
14 #include <dune/istl/operators.hh>
15 #include <dune/istl/solvers.hh>
16 #include <dune/istl/preconditioners.hh>
17 #include <dune/istl/scalarproducts.hh>
18 #include <dune/istl/paamg/amg.hh>
19 #include <dune/istl/paamg/pinfo.hh>
20 #include <dune/istl/io.hh>
21 #include <dune/istl/superlu.hh>
22 #include <dune/istl/umfpack.hh>
23 
29 
30 namespace Dune {
31  namespace PDELab {
32 
36 
37  template<typename X, typename Y, typename GOS>
38  class OnTheFlyOperator : public Dune::LinearOperator<X,Y>
39  {
40  public:
41  typedef X domain_type;
42  typedef Y range_type;
43  typedef typename X::field_type field_type;
44 
45  enum {category=Dune::SolverCategory::sequential};
46 
47  OnTheFlyOperator (GOS& gos_)
48  : gos(gos_)
49  {}
50 
51  virtual void apply (const X& x, Y& y) const
52  {
53  y = 0.0;
54  gos.jacobian_apply(x,y);
55  }
56 
57  virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const
58  {
59  Y temp(y);
60  temp = 0.0;
61  gos.jacobian_apply(x,temp);
62  y.axpy(alpha,temp);
63  }
64 
65  private:
66  GOS& gos;
67  };
68 
69  //==============================================================================
70  // Here we add some standard linear solvers conforming to the linear solver
71  // interface required to solve linear and nonlinear problems.
72  //==============================================================================
73 
74  template<template<class,class,class,int> class Preconditioner,
75  template<class> class Solver>
77  : public SequentialNorm, public LinearResultStorage
78  {
79  public:
85  explicit ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
86  : maxiter(maxiter_), verbose(verbose_)
87  {}
88 
89 
90 
98  template<class M, class V, class W>
99  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
100  {
101  using Backend::Native;
102  using Backend::native;
103 
104  Dune::MatrixAdapter<Native<M>,
105  Native<V>,
106  Native<W>> opa(native(A));
107  Preconditioner<Native<M>,
108  Native<V>,
109  Native<W>,
110  1> prec(native(A), 3, 1.0);
111  Solver<Native<V>> solver(opa, prec, reduction, maxiter, verbose);
112  Dune::InverseOperatorResult stat;
113  solver.apply(native(z), native(r), stat);
114  res.converged = stat.converged;
115  res.iterations = stat.iterations;
116  res.elapsed = stat.elapsed;
117  res.reduction = stat.reduction;
118  res.conv_rate = stat.conv_rate;
119  }
120 
121  private:
122  unsigned maxiter;
123  int verbose;
124  };
125 
126  template<template<typename> class Solver>
128  : public SequentialNorm, public LinearResultStorage
129  {
130  public:
136  explicit ISTLBackend_SEQ_ILU0 (unsigned maxiter_=5000, int verbose_=1)
137  : maxiter(maxiter_), verbose(verbose_)
138  {}
146  template<class M, class V, class W>
147  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
148  {
149  using Backend::Native;
150  using Backend::native;
151  Dune::MatrixAdapter<Native<M>,
152  Native<V>,
153  Native<W>> opa(native(A));
154  Dune::SeqILU0<Native<M>,
155  Native<V>,
156  Native<W>
157  > ilu0(native(A), 1.0);
158  Solver<Native<V>> solver(opa, ilu0, reduction, maxiter, verbose);
159  Dune::InverseOperatorResult stat;
160  solver.apply(native(z), native(r), stat);
161  res.converged = stat.converged;
162  res.iterations = stat.iterations;
163  res.elapsed = stat.elapsed;
164  res.reduction = stat.reduction;
165  res.conv_rate = stat.conv_rate;
166  }
167  private:
168  unsigned maxiter;
169  int verbose;
170  };
171 
172  template<template<typename> class Solver>
174  : public SequentialNorm, public LinearResultStorage
175  {
176  public:
183  ISTLBackend_SEQ_ILUn (int n, double w, unsigned maxiter_=5000, int verbose_=1)
184  : n_(n), w_(w), maxiter(maxiter_), verbose(verbose_)
185  {}
193  template<class M, class V, class W>
194  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
195  {
196  using Backend::Native;
197  using Backend::native;
198  Dune::MatrixAdapter<Native<M>,
199  Native<V>,
200  Native<W>
201  > opa(native(A));
202  Dune::SeqILUn<Native<M>,
203  Native<V>,
204  Native<W>
205  > ilun(native(A), n_, w_);
206  Solver<Native<V>> solver(opa, ilun, reduction, maxiter, verbose);
207  Dune::InverseOperatorResult stat;
208  solver.apply(native(z), native(r), stat);
209  res.converged = stat.converged;
210  res.iterations = stat.iterations;
211  res.elapsed = stat.elapsed;
212  res.reduction = stat.reduction;
213  res.conv_rate = stat.conv_rate;
214  }
215  private:
216  int n_;
217  double w_;
218 
219  unsigned maxiter;
220  int verbose;
221  };
222 
225 
230  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>
231  {
232  public:
237  explicit ISTLBackend_SEQ_LOOP_Jac (unsigned maxiter_=5000, int verbose_=1)
238  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>(maxiter_, verbose_)
239  {}
240  };
241 
246  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>
247  {
248  public:
253  explicit ISTLBackend_SEQ_BCGS_Jac (unsigned maxiter_=5000, int verbose_=1)
254  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>(maxiter_, verbose_)
255  {}
256  };
257 
262  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>
263  {
264  public:
270  explicit ISTLBackend_SEQ_BCGS_SSOR (unsigned maxiter_=5000, int verbose_=1)
271  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>(maxiter_, verbose_)
272  {}
273  };
274 
279  : public ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>
280  {
281  public:
287  explicit ISTLBackend_SEQ_BCGS_ILU0 (unsigned maxiter_=5000, int verbose_=1)
288  : ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>(maxiter_, verbose_)
289  {}
290  };
291 
296  : public ISTLBackend_SEQ_ILU0<Dune::CGSolver>
297  {
298  public:
304  explicit ISTLBackend_SEQ_CG_ILU0 (unsigned maxiter_=5000, int verbose_=1)
305  : ISTLBackend_SEQ_ILU0<Dune::CGSolver>(maxiter_, verbose_)
306  {}
307  };
308 
311  : public ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>
312  {
313  public:
322  explicit ISTLBackend_SEQ_BCGS_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
323  : ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>(n_, w_, maxiter_, verbose_)
324  {}
325  };
326 
329  : public ISTLBackend_SEQ_ILUn<Dune::CGSolver>
330  {
331  public:
340  explicit ISTLBackend_SEQ_CG_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
341  : ISTLBackend_SEQ_ILUn<Dune::CGSolver>(n_, w_, maxiter_, verbose_)
342  {}
343  };
344 
349  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>
350  {
351  public:
357  explicit ISTLBackend_SEQ_CG_SSOR (unsigned maxiter_=5000, int verbose_=1)
358  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>(maxiter_, verbose_)
359  {}
360  };
361 
366  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>
367  {
368  public:
374  explicit ISTLBackend_SEQ_MINRES_SSOR (unsigned maxiter_=5000, int verbose_=1)
375  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>(maxiter_, verbose_)
376  {}
377  };
378 
383  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>
384  {
385  public:
390  explicit ISTLBackend_SEQ_CG_Jac (unsigned maxiter_=5000, int verbose_=1)
391  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>(maxiter_, verbose_)
392  {}
393  };
394 
395 #if HAVE_SUPERLU || DOXYGEN
396 
400  : public SequentialNorm, public LinearResultStorage
401  {
402  public:
407  explicit ISTLBackend_SEQ_SuperLU (int verbose_=1)
408  : verbose(verbose_)
409  {}
410 
411 
417  ISTLBackend_SEQ_SuperLU (int maxiter, int verbose_)
418  : verbose(verbose_)
419  {}
420 
428  template<class M, class V, class W>
429  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
430  {
431  using Backend::Native;
432  using Backend::native;
433  using ISTLM = Native<M>;
434  Dune::SuperLU<ISTLM> solver(native(A), verbose);
435  Dune::InverseOperatorResult stat;
436  solver.apply(native(z), native(r), stat);
437  res.converged = stat.converged;
438  res.iterations = stat.iterations;
439  res.elapsed = stat.elapsed;
440  res.reduction = stat.reduction;
441  res.conv_rate = stat.conv_rate;
442  }
443 
444  private:
445  int verbose;
446  };
447 #endif // HAVE_SUPERLU || DOXYGEN
448 
449 #if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
450 
454  : public SequentialNorm, public LinearResultStorage
455  {
456  public:
461  explicit ISTLBackend_SEQ_UMFPack (int verbose_=1)
462  : verbose(verbose_)
463  {}
464 
465 
471  ISTLBackend_SEQ_UMFPack (int maxiter, int verbose_)
472  : verbose(verbose_)
473  {}
474 
482  template<class M, class V, class W>
483  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
484  {
485  using Backend::native;
486  using ISTLM = Backend::Native<M>;
487  Dune::UMFPack<ISTLM> solver(native(A), verbose);
488  Dune::InverseOperatorResult stat;
489  solver.apply(native(z), native(r), stat);
490  res.converged = stat.converged;
491  res.iterations = stat.iterations;
492  res.elapsed = stat.elapsed;
493  res.reduction = stat.reduction;
494  res.conv_rate = stat.conv_rate;
495  }
496 
497  private:
498  int verbose;
499  };
500 #endif // HAVE_SUITESPARSE_UMFPACK || DOXYGEN
501 
504  : public SequentialNorm, public LinearResultStorage
505  {
506  public:
510  {}
511 
519  template<class M, class V, class W>
520  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
521  {
522  using Backend::Native;
523  Dune::SeqJac<Native<M>,
524  Native<V>,
525  Native<W>
526  > jac(Backend::native(A),1,1.0);
527  jac.pre(z,r);
528  jac.apply(z,r);
529  jac.post(z);
530  res.converged = true;
531  res.iterations = 1;
532  res.elapsed = 0.0;
533  res.reduction = reduction;
534  res.conv_rate = reduction; // pow(reduction,1.0/1)
535  }
536  };
537 
539 
545  {
550  double tprepare;
552  int levels;
554  double tsolve;
556  double tsetup;
561  };
562 
563  template<class GO, template<class,class,class,int> class Preconditioner, template<class> class Solver,
564  bool skipBlocksizeCheck = false>
566  {
567  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
568  typedef typename GO::Traits::Jacobian M;
569  typedef Backend::Native<M> MatrixType;
570  typedef typename GO::Traits::Domain V;
571  typedef Backend::Native<V> VectorType;
572  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
573  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
574  typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
575  typedef Dune::Amg::AMG<Operator,VectorType,Smoother> AMG;
576  typedef Dune::Amg::Parameters Parameters;
577 
578  public:
579  ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1,
580  bool reuse_=false, bool usesuperlu_=true)
581  : maxiter(maxiter_), params(15,2000), verbose(verbose_),
582  reuse(reuse_), firstapply(true), usesuperlu(usesuperlu_)
583  {
584  params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
585  params.setDebugLevel(verbose_);
586 #if !HAVE_SUPERLU
587  if (usesuperlu == true)
588  {
589  std::cout << "WARNING: You are using AMG without SuperLU!"
590  << " Please consider installing SuperLU,"
591  << " or set the usesuperlu flag to false"
592  << " to suppress this warning." << std::endl;
593  }
594 #endif
595  }
596 
601  void setparams(Parameters params_)
602  {
603  params = params_;
604  }
605 
607  void setReuse(bool reuse_)
608  {
609  reuse = reuse_;
610  }
611 
613  bool getReuse() const
614  {
615  return reuse;
616  }
617 
622  typename V::ElementType norm (const V& v) const
623  {
624  return Backend::native(v).two_norm();
625  }
626 
634  void apply(M& A, V& z, V& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
635  {
636  Timer watch;
637  MatrixType& mat = Backend::native(A);
638  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
639  Dune::Amg::FirstDiagonal> > Criterion;
640  SmootherArgs smootherArgs;
641  smootherArgs.iterations = 1;
642  smootherArgs.relaxationFactor = 1;
643 
644  Criterion criterion(params);
645  //only construct a new AMG if the matrix changes
646  if (reuse==false || firstapply==true){
647  oop.reset(new Operator(mat));
648  amg.reset(new AMG(*oop, criterion, smootherArgs));
649  firstapply = false;
650  stats.tsetup = watch.elapsed();
651  stats.levels = amg->maxlevels();
652  stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
653  }
654  watch.reset();
655  Dune::InverseOperatorResult stat;
656 
657  Solver<VectorType> solver(*oop,*amg,reduction,maxiter,verbose);
658  solver.apply(Backend::native(z),Backend::native(r),stat);
659  stats.tsolve= watch.elapsed();
660  res.converged = stat.converged;
661  res.iterations = stat.iterations;
662  res.elapsed = stat.elapsed;
663  res.reduction = stat.reduction;
664  res.conv_rate = stat.conv_rate;
665  }
666 
667 
673  {
674  return stats;
675  }
676 
677  private:
678  unsigned maxiter;
679  Parameters params;
680  int verbose;
681  bool reuse;
682  bool firstapply;
683  bool usesuperlu;
684  std::shared_ptr<Operator> oop;
685  std::shared_ptr<AMG> amg;
686  ISTLAMGStatistics stats;
687  };
688 
691 
697  template<class GO>
699  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
700  {
701 
702  public:
711  ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
712  bool reuse_=false, bool usesuperlu_=true)
713  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
714  (maxiter_, verbose_, reuse_, usesuperlu_)
715  {}
716  };
717 
723  template<class GO>
725  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
726  {
727 
728  public:
737  ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
738  bool reuse_=false, bool usesuperlu_=true)
739  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
740  (maxiter_, verbose_, reuse_, usesuperlu_)
741  {}
742  };
743 
749  template<class GO>
751  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
752  {
753 
754  public:
763  ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
764  bool reuse_=false, bool usesuperlu_=true)
765  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
766  (maxiter_, verbose_, reuse_, usesuperlu_)
767  {}
768  };
769 
775  template<class GO>
777  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
778  {
779 
780  public:
789  ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
790  bool reuse_=false, bool usesuperlu_=true)
791  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
792  (maxiter_, verbose_, reuse_, usesuperlu_)
793  {}
794  };
795 
801  template<class GO>
803  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
804  {
805 
806  public:
815  ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
816  bool reuse_=false, bool usesuperlu_=true)
817  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
818  (maxiter_, verbose_, reuse_, usesuperlu_)
819  {}
820  };
821 
828  : public SequentialNorm, public LinearResultStorage
829  {
830  public :
831 
838  explicit ISTLBackend_SEQ_GMRES_ILU0(int restart_ = 200, int maxiter_ = 5000, int verbose_ = 1)
839  : restart(restart_), maxiter(maxiter_), verbose(verbose_)
840  {}
841 
849  template<class M, class V, class W>
850  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType>::real_type reduction)
851  {
852  using Backend::Native;
853  using Backend::native;
854  Dune::MatrixAdapter<
855  Native<M>,
856  Native<V>,
857  Native<W>
858  > opa(native(A));
859  Dune::SeqILU0<
860  Native<M>,
861  Native<V>,
862  Native<W>
863  > ilu0(native(A), 1.0);
864  Dune::RestartedGMResSolver<Native<V>> solver(opa,ilu0,reduction,restart,maxiter,verbose);
865  Dune::InverseOperatorResult stat;
866  solver.apply(native(z), native(r), stat);
867  res.converged = stat.converged;
868  res.iterations = stat.iterations;
869  res.elapsed = stat.elapsed;
870  res.reduction = stat.reduction;
871  res.conv_rate = stat.conv_rate;
872  }
873 
874  private :
875  int restart, maxiter, verbose;
876  };
877 
880 
881  } // namespace PDELab
882 } // namespace Dune
883 
884 #endif // DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: seqistlsolverbackend.hh:503
ISTLBackend_SEQ_CG_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:357
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:483
Definition: seqistlsolverbackend.hh:127
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > &>::type native(T &t)
Definition: backend/interface.hh:192
Linear solver backend for Restarted GMRes preconditioned with ILU(0)
Definition: seqistlsolverbackend.hh:827
Sequential BiCGSTAB solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:750
ISTLBackend_SEQ_BCGS_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:270
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: seqistlsolverbackend.hh:622
typename native_type< T >::type Native
Alias of the native container type associated with T or T itself if it is not a backend wrapper...
Definition: backend/interface.hh:176
Definition: seqistlsolverbackend.hh:76
ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: seqistlsolverbackend.hh:579
ISTLBackend_SEQ_CG_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:340
ISTLBackend_SEQ_CG_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:304
Y range_type
Definition: seqistlsolverbackend.hh:42
virtual void apply(const X &x, Y &y) const
Definition: seqistlsolverbackend.hh:51
ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:737
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:429
Definition: seqistlsolverbackend.hh:173
OnTheFlyOperator(GOS &gos_)
Definition: seqistlsolverbackend.hh:47
ISTLBackend_SEQ_LOOP_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:237
ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:815
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:194
ISTLBackend_SEQ_MINRES_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:374
Backend for sequential conjugate gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:295
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: seqistlsolverbackend.hh:672
ISTLBackend_SEQ_ILUn(int n, double w, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:183
ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:763
double tsolve
The time spent in solving the system (without building the hierarchy.
Definition: seqistlsolverbackend.hh:554
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:634
double tsetup
The time needed for building the AMG hierarchy (coarsening).
Definition: seqistlsolverbackend.hh:556
void setReuse(bool reuse_)
Set whether the AMG should be reused again during call to apply().
Definition: seqistlsolverbackend.hh:607
Backend for conjugate gradient solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:382
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:850
Backend for sequential BiCGSTAB solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:245
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:520
Definition: seqistlsolverbackend.hh:565
ISTLBackend_SEQ_SuperLU(int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:407
int levels
the number of levels in the AMG hierarchy.
Definition: seqistlsolverbackend.hh:552
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:147
ISTLBackend_SEQ_UMFPack(int maxiter, int verbose_)
make a linear solver object
Definition: seqistlsolverbackend.hh:471
ISTLBackend_SEQ_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:136
Sequential congute gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:328
ISTLBackend_SEQ_ExplicitDiagonal()
make a linear solver object
Definition: seqistlsolverbackend.hh:509
bool getReuse() const
Return whether the AMG is reused during call to apply()
Definition: seqistlsolverbackend.hh:613
Sequential Loop solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:802
Backend for sequential BiCGSTAB solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:278
bool directCoarseLevelSolver
True if a direct solver was used on the coarset level.
Definition: seqistlsolverbackend.hh:560
Definition: seqistlsolverbackend.hh:45
Backend for sequential loop solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:229
double tprepare
The needed for computing the parallel information and for adapting the linear system.
Definition: seqistlsolverbackend.hh:550
ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:789
ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:711
Solver backend using UMFPack as a direct solver.
Definition: seqistlsolverbackend.hh:453
ISTLBackend_SEQ_SuperLU(int maxiter, int verbose_)
make a linear solver object
Definition: seqistlsolverbackend.hh:417
Backend using a MINRes solver preconditioned by SSOR.
Definition: seqistlsolverbackend.hh:365
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:99
ISTLBackend_SEQ_BCGS_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:253
Sequential BiCGStab solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:310
X::field_type field_type
Definition: seqistlsolverbackend.hh:43
X domain_type
Definition: seqistlsolverbackend.hh:41
int iterations
The number of iterations performed until convergence was reached.
Definition: seqistlsolverbackend.hh:558
Backend for sequential BiCGSTAB solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:261
ISTLBackend_SEQ_BCGS_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:287
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const
Definition: seqistlsolverbackend.hh:57
ISTLBackend_SEQ_CG_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:390
ISTLBackend_SEQ_UMFPack(int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:461
VTKWriter & w
Definition: function.hh:829
Sequential BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:724
ISTLBackend_SEQ_GMRES_ILU0(int restart_=200, int maxiter_=5000, int verbose_=1)
make linear solver object
Definition: seqistlsolverbackend.hh:838
Definition: solver.hh:16
Class providing some statistics of the AMG solver.
Definition: seqistlsolverbackend.hh:544
ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:85
Solver backend using SuperLU as a direct solver.
Definition: seqistlsolverbackend.hh:399
Backend for sequential conjugate gradient solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:348
ISTLBackend_SEQ_BCGS_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:322
Definition: seqistlsolverbackend.hh:38
void setparams(Parameters params_)
set AMG parameters
Definition: seqistlsolverbackend.hh:601
Sequential Loop solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:776
Definition: solver.hh:42
Sequential conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:698