dune-pdelab  2.0.0
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_SEQISTLSOLVERBACKEND_HH
4 #define DUNE_SEQISTLSOLVERBACKEND_HH
5 
6 #include <dune/common/deprecated.hh>
7 #include <dune/common/parallel/mpihelper.hh>
8 
9 #include <dune/istl/owneroverlapcopy.hh>
10 #include <dune/istl/solvercategory.hh>
11 #include <dune/istl/operators.hh>
12 #include <dune/istl/solvers.hh>
13 #include <dune/istl/preconditioners.hh>
14 #include <dune/istl/scalarproducts.hh>
15 #include <dune/istl/paamg/amg.hh>
16 #include <dune/istl/paamg/pinfo.hh>
17 #include <dune/istl/io.hh>
18 #include <dune/istl/superlu.hh>
19 
25 
26 namespace Dune {
27  namespace PDELab {
28 
32 
33  template<typename X, typename Y, typename GOS>
34  class OnTheFlyOperator : public Dune::LinearOperator<X,Y>
35  {
36  public:
37  typedef X domain_type;
38  typedef Y range_type;
39  typedef typename X::field_type field_type;
40 
41  enum {category=Dune::SolverCategory::sequential};
42 
43  OnTheFlyOperator (GOS& gos_)
44  : gos(gos_)
45  {}
46 
47  virtual void apply (const X& x, Y& y) const
48  {
49  y = 0.0;
50  gos.jacobian_apply(x,y);
51  }
52 
53  virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const
54  {
55  Y temp(y);
56  temp = 0.0;
57  gos.jacobian_apply(x,temp);
58  y.axpy(alpha,temp);
59  }
60 
61  private:
62  GOS& gos;
63  };
64 
65  //==============================================================================
66  // Here we add some standard linear solvers conforming to the linear solver
67  // interface required to solve linear and nonlinear problems.
68  //==============================================================================
69 
70  template<template<class,class,class,int> class Preconditioner,
71  template<class> class Solver>
73  : public SequentialNorm, public LinearResultStorage
74  {
75  public:
81  explicit ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
82  : maxiter(maxiter_), verbose(verbose_)
83  {}
84 
85 
86 
94  template<class M, class V, class W>
95  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
96  {
97  Dune::MatrixAdapter<typename M::BaseT,
98  typename V::BaseT,
99  typename W::BaseT> opa(istl::raw(A));
100  Preconditioner<typename M::BaseT,
101  typename V::BaseT,
102  typename W::BaseT,1> prec(istl::raw(A), 3, 1.0);
103  Solver<typename V::BaseT> solver(opa, prec, reduction, maxiter, verbose);
104  Dune::InverseOperatorResult stat;
105  solver.apply(istl::raw(z), istl::raw(r), stat);
106  res.converged = stat.converged;
107  res.iterations = stat.iterations;
108  res.elapsed = stat.elapsed;
109  res.reduction = stat.reduction;
110  res.conv_rate = stat.conv_rate;
111  }
112 
113  private:
114  unsigned maxiter;
115  int verbose;
116  };
117 
118  template<template<typename> class Solver>
120  : public SequentialNorm, public LinearResultStorage
121  {
122  public:
128  explicit ISTLBackend_SEQ_ILU0 (unsigned maxiter_=5000, int verbose_=1)
129  : maxiter(maxiter_), verbose(verbose_)
130  {}
138  template<class M, class V, class W>
139  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
140  {
141  Dune::MatrixAdapter<typename M::BaseT,
142  typename V::BaseT,
143  typename W::BaseT> opa(istl::raw(A));
144  Dune::SeqILU0<typename M::BaseT,
145  typename V::BaseT,
146  typename W::BaseT> ilu0(istl::raw(A), 1.0);
147  Solver<typename V::BaseT> solver(opa, ilu0, reduction, maxiter, verbose);
148  Dune::InverseOperatorResult stat;
149  solver.apply(istl::raw(z), istl::raw(r), stat);
150  res.converged = stat.converged;
151  res.iterations = stat.iterations;
152  res.elapsed = stat.elapsed;
153  res.reduction = stat.reduction;
154  res.conv_rate = stat.conv_rate;
155  }
156  private:
157  unsigned maxiter;
158  int verbose;
159  };
160 
161  template<template<typename> class Solver>
163  : public SequentialNorm, public LinearResultStorage
164  {
165  public:
172  ISTLBackend_SEQ_ILUn (int n, double w, unsigned maxiter_=5000, int verbose_=1)
173  : n_(n), w_(w), maxiter(maxiter_), verbose(verbose_)
174  {}
182  template<class M, class V, class W>
183  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
184  {
185  Dune::MatrixAdapter<typename M::BaseT,
186  typename V::BaseT,
187  typename W::BaseT> opa(istl::raw(A));
188  Dune::SeqILUn<typename M::BaseT,
189  typename V::BaseT,
190  typename W::BaseT> ilun(istl::raw(A), n_, w_);
191  Solver<typename V::BaseT> solver(opa, ilun, reduction, maxiter, verbose);
192  Dune::InverseOperatorResult stat;
193  solver.apply(istl::raw(z), istl::raw(r), stat);
194  res.converged = stat.converged;
195  res.iterations = stat.iterations;
196  res.elapsed = stat.elapsed;
197  res.reduction = stat.reduction;
198  res.conv_rate = stat.conv_rate;
199  }
200  private:
201  int n_;
202  double w_;
203 
204  unsigned maxiter;
205  int verbose;
206  };
207 
210 
215  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>
216  {
217  public:
222  explicit ISTLBackend_SEQ_LOOP_Jac (unsigned maxiter_=5000, int verbose_=1)
223  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>(maxiter_, verbose_)
224  {}
225  };
226 
231  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>
232  {
233  public:
238  explicit ISTLBackend_SEQ_BCGS_Jac (unsigned maxiter_=5000, int verbose_=1)
239  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>(maxiter_, verbose_)
240  {}
241  };
242 
247  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>
248  {
249  public:
255  explicit ISTLBackend_SEQ_BCGS_SSOR (unsigned maxiter_=5000, int verbose_=1)
256  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>(maxiter_, verbose_)
257  {}
258  };
259 
264  : public ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>
265  {
266  public:
272  explicit ISTLBackend_SEQ_BCGS_ILU0 (unsigned maxiter_=5000, int verbose_=1)
273  : ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>(maxiter_, verbose_)
274  {}
275  };
276 
281  : public ISTLBackend_SEQ_ILU0<Dune::CGSolver>
282  {
283  public:
289  explicit ISTLBackend_SEQ_CG_ILU0 (unsigned maxiter_=5000, int verbose_=1)
290  : ISTLBackend_SEQ_ILU0<Dune::CGSolver>(maxiter_, verbose_)
291  {}
292  };
293 
296  : public ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>
297  {
298  public:
307  explicit ISTLBackend_SEQ_BCGS_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
308  : ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>(n_, w_, maxiter_, verbose_)
309  {}
310  };
311 
314  : public ISTLBackend_SEQ_ILUn<Dune::CGSolver>
315  {
316  public:
325  explicit ISTLBackend_SEQ_CG_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
326  : ISTLBackend_SEQ_ILUn<Dune::CGSolver>(n_, w_, maxiter_, verbose_)
327  {}
328  };
329 
334  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>
335  {
336  public:
342  explicit ISTLBackend_SEQ_CG_SSOR (unsigned maxiter_=5000, int verbose_=1)
343  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>(maxiter_, verbose_)
344  {}
345  };
346 
351  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>
352  {
353  public:
359  explicit ISTLBackend_SEQ_MINRES_SSOR (unsigned maxiter_=5000, int verbose_=1)
360  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>(maxiter_, verbose_)
361  {}
362  };
363 
368  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>
369  {
370  public:
375  explicit ISTLBackend_SEQ_CG_Jac (unsigned maxiter_=5000, int verbose_=1)
376  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>(maxiter_, verbose_)
377  {}
378  };
379 
380 #if HAVE_SUPERLU
381 
384  class ISTLBackend_SEQ_SuperLU
385  : public SequentialNorm, public LinearResultStorage
386  {
387  public:
392  explicit ISTLBackend_SEQ_SuperLU (int verbose_=1)
393  : verbose(verbose_)
394  {}
395 
396 
402  ISTLBackend_SEQ_SuperLU (int maxiter, int verbose_)
403  : verbose(verbose_)
404  {}
405 
413  template<class M, class V, class W>
414  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
415  {
416  typedef typename M::Container ISTLM;
417  Dune::SuperLU<ISTLM> solver(istl::raw(A), verbose);
418  Dune::InverseOperatorResult stat;
419  solver.apply(istl::raw(z), istl::raw(r), stat);
420  res.converged = stat.converged;
421  res.iterations = stat.iterations;
422  res.elapsed = stat.elapsed;
423  res.reduction = stat.reduction;
424  res.conv_rate = stat.conv_rate;
425  }
426 
427  private:
428  int verbose;
429  };
430 #endif // HAVE_SUPERLU
431 
434  : public SequentialNorm, public LinearResultStorage
435  {
436  public:
440  {}
441 
449  template<class M, class V, class W>
450  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
451  {
452  Dune::SeqJac<typename M::BaseT,
453  typename V::BaseT,
454  typename W::BaseT> jac(istl::raw(A),1,1.0);
455  jac.pre(z,r);
456  jac.apply(z,r);
457  jac.post(z);
458  res.converged = true;
459  res.iterations = 1;
460  res.elapsed = 0.0;
461  res.reduction = reduction;
462  res.conv_rate = reduction; // pow(reduction,1.0/1)
463  }
464  };
465 
467 
473  {
478  double tprepare;
480  int levels;
482  double tsolve;
484  double tsetup;
489  };
490 
491  template<class GO, template<class,class,class,int> class Preconditioner, template<class> class Solver,
492  bool skipBlocksizeCheck = false>
494  {
495  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
496  typedef typename GO::Traits::Jacobian M;
497  typedef typename M::BaseT MatrixType;
498  typedef typename GO::Traits::Domain V;
499  typedef typename V::BaseT VectorType;
500  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
501  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
502  typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
503  typedef Dune::Amg::AMG<Operator,VectorType,Smoother> AMG;
504  typedef Dune::Amg::Parameters Parameters;
505 
506  public:
507  ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1,
508  bool reuse_=false, bool usesuperlu_=true)
509  : maxiter(maxiter_), params(15,2000), verbose(verbose_),
510  reuse(reuse_), firstapply(true), usesuperlu(usesuperlu_)
511  {
512  params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
513  params.setDebugLevel(verbose_);
514 #if !HAVE_SUPERLU
515  if (usesuperlu == true)
516  {
517  std::cout << "WARNING: You are using AMG without SuperLU!"
518  << " Please consider installing SuperLU,"
519  << " or set the usesuperlu flag to false"
520  << " to suppress this warning." << std::endl;
521  }
522 #endif
523  }
524 
529  void setparams(Parameters params_)
530  {
531  params = params_;
532  }
533 
538  typename V::ElementType norm (const V& v) const
539  {
540  return istl::raw(v).two_norm();
541  }
542 
550  void apply(M& A, V& z, V& r, typename V::ElementType reduction)
551  {
552  Timer watch;
553  MatrixType& mat=istl::raw(A);
554  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
555  Dune::Amg::FirstDiagonal> > Criterion;
556  SmootherArgs smootherArgs;
557  smootherArgs.iterations = 1;
558  smootherArgs.relaxationFactor = 1;
559 
560  Criterion criterion(params);
561  Operator oop(mat);
562  //only construct a new AMG if the matrix changes
563  if (reuse==false || firstapply==true){
564  amg.reset(new AMG(oop, criterion, smootherArgs));
565  firstapply = false;
566  stats.tsetup = watch.elapsed();
567  stats.levels = amg->maxlevels();
568  stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
569  }
570  watch.reset();
571  Dune::InverseOperatorResult stat;
572 
573  Solver<VectorType> solver(oop,*amg,reduction,maxiter,verbose);
574  solver.apply(istl::raw(z),istl::raw(r),stat);
575  stats.tsolve= watch.elapsed();
576  res.converged = stat.converged;
577  res.iterations = stat.iterations;
578  res.elapsed = stat.elapsed;
579  res.reduction = stat.reduction;
580  res.conv_rate = stat.conv_rate;
581  }
582 
583 
589  {
590  return stats;
591  }
592 
593  private:
594  unsigned maxiter;
595  Parameters params;
596  int verbose;
597  bool reuse;
598  bool firstapply;
599  bool usesuperlu;
600  Dune::shared_ptr<AMG> amg;
601  ISTLAMGStatistics stats;
602  };
603 
606 
612  template<class GO>
614  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
615  {
616 
617  public:
626  ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
627  bool reuse_=false, bool usesuperlu_=true)
628  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
629  (maxiter_, verbose_, reuse_, usesuperlu_)
630  {}
631  };
632 
638  template<class GO>
640  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
641  {
642 
643  public:
652  ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
653  bool reuse_=false, bool usesuperlu_=true)
654  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
655  (maxiter_, verbose_, reuse_, usesuperlu_)
656  {}
657  };
658 
664  template<class GO>
666  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
667  {
668 
669  public:
678  ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
679  bool reuse_=false, bool usesuperlu_=true)
680  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
681  (maxiter_, verbose_, reuse_, usesuperlu_)
682  {}
683  };
684 
690  template<class GO>
692  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
693  {
694 
695  public:
704  ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
705  bool reuse_=false, bool usesuperlu_=true)
706  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
707  (maxiter_, verbose_, reuse_, usesuperlu_)
708  {}
709  };
710 
716  template<class GO>
718  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
719  {
720 
721  public:
730  ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
731  bool reuse_=false, bool usesuperlu_=true)
732  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
733  (maxiter_, verbose_, reuse_, usesuperlu_)
734  {}
735  };
736 
739 
740  } // namespace PDELab
741 } // namespace Dune
742 
743 #endif
ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:652
X domain_type
Definition: seqistlsolverbackend.hh:37
Sequential congute gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:313
ISTLBackend_SEQ_MINRES_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:359
Class providing some statistics of the AMG solver.
Definition: seqistlsolverbackend.hh:472
double tsetup
The time needed for building the AMG hierarchy (coarsening).
Definition: seqistlsolverbackend.hh:484
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:183
bool converged
Definition: solver.hh:32
double tsolve
The time spent in solving the system (without building the hierarchy.
Definition: seqistlsolverbackend.hh:482
Definition: seqistlsolverbackend.hh:72
ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:678
Sequential Loop solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:691
ISTLBackend_SEQ_ILUn(int n, double w, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:172
Definition: seqistlsolverbackend.hh:34
Backend for sequential BiCGSTAB solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:230
double elapsed
Definition: solver.hh:34
Backend for sequential BiCGSTAB solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:263
VTKWriter & w
Definition: function.hh:1102
ISTLBackend_SEQ_LOOP_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:222
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:450
int levels
the number of levels in the AMG hierarchy.
Definition: seqistlsolverbackend.hh:480
ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: seqistlsolverbackend.hh:507
Backend for sequential conjugate gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:280
ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:730
Backend for sequential conjugate gradient solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:333
OnTheFlyOperator(GOS &gos_)
Definition: seqistlsolverbackend.hh:43
ISTLBackend_SEQ_BCGS_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:272
ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:626
Sequential Loop solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:717
X::field_type field_type
Definition: seqistlsolverbackend.hh:39
Sequential BiCGStab solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:295
V & raw(V &v)
Returns the raw ISTL object associated with v, or v itself it is already an ISTL object.
Definition: backend/istl/utility.hh:26
void setparams(Parameters params_)
set AMG parameters
Definition: seqistlsolverbackend.hh:529
Sequential BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:639
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:52
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:95
ISTLBackend_SEQ_BCGS_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:238
ISTLBackend_SEQ_CG_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:375
ISTLBackend_SEQ_ExplicitDiagonal()
make a linear solver object
Definition: seqistlsolverbackend.hh:439
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const
Definition: seqistlsolverbackend.hh:53
RFType reduction
Definition: solver.hh:35
Sequential conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:613
ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:81
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: seqistlsolverbackend.hh:588
Definition: solver.hh:42
ISTLBackend_SEQ_BCGS_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:255
Backend for conjugate gradient solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:367
Definition: seqistlsolverbackend.hh:493
bool directCoarseLevelSolver
True if a direct solver was used on the coarset level.
Definition: seqistlsolverbackend.hh:488
ISTLBackend_SEQ_CG_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:325
Backend for sequential loop solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:214
Definition: solver.hh:16
Definition: seqistlsolverbackend.hh:119
Sequential BiCGSTAB solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:665
ISTLBackend_SEQ_BCGS_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:307
Backend using a MINRes solver preconditioned by SSOR.
Definition: seqistlsolverbackend.hh:350
Backend for sequential BiCGSTAB solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:246
Y range_type
Definition: seqistlsolverbackend.hh:38
double tprepare
The needed for computing the parallel information and for adapting the linear system.
Definition: seqistlsolverbackend.hh:478
ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:704
int iterations
The number of iterations performed until convergence was reached.
Definition: seqistlsolverbackend.hh:486
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: seqistlsolverbackend.hh:538
ISTLBackend_SEQ_CG_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:342
RFType conv_rate
Definition: solver.hh:36
ISTLBackend_SEQ_CG_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:289
void apply(M &A, V &z, V &r, typename V::ElementType reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:550
virtual void apply(const X &x, Y &y) const
Definition: seqistlsolverbackend.hh:47
ISTLBackend_SEQ_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:128
Definition: seqistlsolverbackend.hh:41
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:139
Definition: seqistlsolverbackend.hh:162
unsigned int iterations
Definition: solver.hh:33
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: seqistlsolverbackend.hh:433