dune-pdelab  2.4-dev
istl/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 #include <dune/istl/umfpack.hh>
20 
26 
27 namespace Dune {
28  namespace PDELab {
29 
33 
34  template<typename X, typename Y, typename GOS>
35  class OnTheFlyOperator : public Dune::LinearOperator<X,Y>
36  {
37  public:
38  typedef X domain_type;
39  typedef Y range_type;
40  typedef typename X::field_type field_type;
41 
42  enum {category=Dune::SolverCategory::sequential};
43 
44  OnTheFlyOperator (GOS& gos_)
45  : gos(gos_)
46  {}
47 
48  virtual void apply (const X& x, Y& y) const
49  {
50  y = 0.0;
51  gos.jacobian_apply(x,y);
52  }
53 
54  virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const
55  {
56  Y temp(y);
57  temp = 0.0;
58  gos.jacobian_apply(x,temp);
59  y.axpy(alpha,temp);
60  }
61 
62  private:
63  GOS& gos;
64  };
65 
66  //==============================================================================
67  // Here we add some standard linear solvers conforming to the linear solver
68  // interface required to solve linear and nonlinear problems.
69  //==============================================================================
70 
71  template<template<class,class,class,int> class Preconditioner,
72  template<class> class Solver>
74  : public SequentialNorm, public LinearResultStorage
75  {
76  public:
82  explicit ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
83  : maxiter(maxiter_), verbose(verbose_)
84  {}
85 
86 
87 
95  template<class M, class V, class W>
96  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
97  {
98  using Backend::Native;
99  using Backend::native;
100 
101  Dune::MatrixAdapter<Native<M>,
102  Native<V>,
103  Native<W>> opa(native(A));
104  Preconditioner<Native<M>,
105  Native<V>,
106  Native<W>,
107  1> prec(native(A), 3, 1.0);
108  Solver<Native<V>> solver(opa, prec, reduction, maxiter, verbose);
109  Dune::InverseOperatorResult stat;
110  solver.apply(native(z), native(r), stat);
111  res.converged = stat.converged;
112  res.iterations = stat.iterations;
113  res.elapsed = stat.elapsed;
114  res.reduction = stat.reduction;
115  res.conv_rate = stat.conv_rate;
116  }
117 
118  private:
119  unsigned maxiter;
120  int verbose;
121  };
122 
123  template<template<typename> class Solver>
125  : public SequentialNorm, public LinearResultStorage
126  {
127  public:
133  explicit ISTLBackend_SEQ_ILU0 (unsigned maxiter_=5000, int verbose_=1)
134  : maxiter(maxiter_), verbose(verbose_)
135  {}
143  template<class M, class V, class W>
144  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
145  {
146  using Backend::Native;
147  using Backend::native;
148  Dune::MatrixAdapter<Native<M>,
149  Native<V>,
150  Native<W>> opa(native(A));
151  Dune::SeqILU0<Native<M>,
152  Native<V>,
153  Native<W>
154  > ilu0(native(A), 1.0);
155  Solver<Native<V>> solver(opa, ilu0, reduction, maxiter, verbose);
156  Dune::InverseOperatorResult stat;
157  solver.apply(native(z), native(r), stat);
158  res.converged = stat.converged;
159  res.iterations = stat.iterations;
160  res.elapsed = stat.elapsed;
161  res.reduction = stat.reduction;
162  res.conv_rate = stat.conv_rate;
163  }
164  private:
165  unsigned maxiter;
166  int verbose;
167  };
168 
169  template<template<typename> class Solver>
171  : public SequentialNorm, public LinearResultStorage
172  {
173  public:
180  ISTLBackend_SEQ_ILUn (int n, double w, unsigned maxiter_=5000, int verbose_=1)
181  : n_(n), w_(w), maxiter(maxiter_), verbose(verbose_)
182  {}
190  template<class M, class V, class W>
191  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
192  {
193  using Backend::Native;
194  using Backend::native;
195  Dune::MatrixAdapter<Native<M>,
196  Native<V>,
197  Native<W>
198  > opa(native(A));
199  Dune::SeqILUn<Native<M>,
200  Native<V>,
201  Native<W>
202  > ilun(native(A), n_, w_);
203  Solver<Native<V>> solver(opa, ilun, reduction, maxiter, verbose);
204  Dune::InverseOperatorResult stat;
205  solver.apply(native(z), native(r), stat);
206  res.converged = stat.converged;
207  res.iterations = stat.iterations;
208  res.elapsed = stat.elapsed;
209  res.reduction = stat.reduction;
210  res.conv_rate = stat.conv_rate;
211  }
212  private:
213  int n_;
214  double w_;
215 
216  unsigned maxiter;
217  int verbose;
218  };
219 
222 
227  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>
228  {
229  public:
234  explicit ISTLBackend_SEQ_LOOP_Jac (unsigned maxiter_=5000, int verbose_=1)
235  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::LoopSolver>(maxiter_, verbose_)
236  {}
237  };
238 
243  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>
244  {
245  public:
250  explicit ISTLBackend_SEQ_BCGS_Jac (unsigned maxiter_=5000, int verbose_=1)
251  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::BiCGSTABSolver>(maxiter_, verbose_)
252  {}
253  };
254 
259  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>
260  {
261  public:
267  explicit ISTLBackend_SEQ_BCGS_SSOR (unsigned maxiter_=5000, int verbose_=1)
268  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::BiCGSTABSolver>(maxiter_, verbose_)
269  {}
270  };
271 
276  : public ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>
277  {
278  public:
284  explicit ISTLBackend_SEQ_BCGS_ILU0 (unsigned maxiter_=5000, int verbose_=1)
285  : ISTLBackend_SEQ_ILU0<Dune::BiCGSTABSolver>(maxiter_, verbose_)
286  {}
287  };
288 
293  : public ISTLBackend_SEQ_ILU0<Dune::CGSolver>
294  {
295  public:
301  explicit ISTLBackend_SEQ_CG_ILU0 (unsigned maxiter_=5000, int verbose_=1)
302  : ISTLBackend_SEQ_ILU0<Dune::CGSolver>(maxiter_, verbose_)
303  {}
304  };
305 
308  : public ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>
309  {
310  public:
319  explicit ISTLBackend_SEQ_BCGS_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
320  : ISTLBackend_SEQ_ILUn<Dune::BiCGSTABSolver>(n_, w_, maxiter_, verbose_)
321  {}
322  };
323 
326  : public ISTLBackend_SEQ_ILUn<Dune::CGSolver>
327  {
328  public:
337  explicit ISTLBackend_SEQ_CG_ILUn (int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
338  : ISTLBackend_SEQ_ILUn<Dune::CGSolver>(n_, w_, maxiter_, verbose_)
339  {}
340  };
341 
346  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>
347  {
348  public:
354  explicit ISTLBackend_SEQ_CG_SSOR (unsigned maxiter_=5000, int verbose_=1)
355  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::CGSolver>(maxiter_, verbose_)
356  {}
357  };
358 
363  : public ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>
364  {
365  public:
371  explicit ISTLBackend_SEQ_MINRES_SSOR (unsigned maxiter_=5000, int verbose_=1)
372  : ISTLBackend_SEQ_Base<Dune::SeqSSOR, Dune::MINRESSolver>(maxiter_, verbose_)
373  {}
374  };
375 
380  : public ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>
381  {
382  public:
387  explicit ISTLBackend_SEQ_CG_Jac (unsigned maxiter_=5000, int verbose_=1)
388  : ISTLBackend_SEQ_Base<Dune::SeqJac, Dune::CGSolver>(maxiter_, verbose_)
389  {}
390  };
391 
392 #if HAVE_SUPERLU || DOXYGEN
393 
397  : public SequentialNorm, public LinearResultStorage
398  {
399  public:
404  explicit ISTLBackend_SEQ_SuperLU (int verbose_=1)
405  : verbose(verbose_)
406  {}
407 
408 
414  ISTLBackend_SEQ_SuperLU (int maxiter, int verbose_)
415  : verbose(verbose_)
416  {}
417 
425  template<class M, class V, class W>
426  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
427  {
428  using Backend::Native;
429  using Backend::native;
430  using ISTLM = Native<M>;
431  Dune::SuperLU<ISTLM> solver(native(A), verbose);
432  Dune::InverseOperatorResult stat;
433  solver.apply(native(z), native(r), stat);
434  res.converged = stat.converged;
435  res.iterations = stat.iterations;
436  res.elapsed = stat.elapsed;
437  res.reduction = stat.reduction;
438  res.conv_rate = stat.conv_rate;
439  }
440 
441  private:
442  int verbose;
443  };
444 #endif // HAVE_SUPERLU || DOXYGEN
445 
446 #if HAVE_UMFPACK || DOXYGEN
447 
451  : public SequentialNorm, public LinearResultStorage
452  {
453  public:
458  explicit ISTLBackend_SEQ_UMFPack (int verbose_=1)
459  : verbose(verbose_)
460  {}
461 
462 
468  ISTLBackend_SEQ_UMFPack (int maxiter, int verbose_)
469  : verbose(verbose_)
470  {}
471 
479  template<class M, class V, class W>
480  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
481  {
482  using Backend::native;
483  using ISTLM = Backend::Native<M>;
484  Dune::UMFPack<ISTLM> solver(native(A), verbose);
485  Dune::InverseOperatorResult stat;
486  solver.apply(native(z), native(r), stat);
487  res.converged = stat.converged;
488  res.iterations = stat.iterations;
489  res.elapsed = stat.elapsed;
490  res.reduction = stat.reduction;
491  res.conv_rate = stat.conv_rate;
492  }
493 
494  private:
495  int verbose;
496  };
497 #endif // HAVE_UMFPACK || DOXYGEN
498 
501  : public SequentialNorm, public LinearResultStorage
502  {
503  public:
507  {}
508 
516  template<class M, class V, class W>
517  void apply(M& A, V& z, W& r, typename W::ElementType reduction)
518  {
519  using Backend::Native;
520  Dune::SeqJac<Native<M>,
521  Native<V>,
522  Native<W>
523  > jac(Backend::native(A),1,1.0);
524  jac.pre(z,r);
525  jac.apply(z,r);
526  jac.post(z);
527  res.converged = true;
528  res.iterations = 1;
529  res.elapsed = 0.0;
530  res.reduction = reduction;
531  res.conv_rate = reduction; // pow(reduction,1.0/1)
532  }
533  };
534 
536 
542  {
547  double tprepare;
549  int levels;
551  double tsolve;
553  double tsetup;
558  };
559 
560  template<class GO, template<class,class,class,int> class Preconditioner, template<class> class Solver,
561  bool skipBlocksizeCheck = false>
563  {
564  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
565  typedef typename GO::Traits::Jacobian M;
566  typedef Backend::Native<M> MatrixType;
567  typedef typename GO::Traits::Domain V;
568  typedef Backend::Native<V> VectorType;
569  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
570  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
571  typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
572  typedef Dune::Amg::AMG<Operator,VectorType,Smoother> AMG;
573  typedef Dune::Amg::Parameters Parameters;
574 
575  public:
576  ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1,
577  bool reuse_=false, bool usesuperlu_=true)
578  : maxiter(maxiter_), params(15,2000), verbose(verbose_),
579  reuse(reuse_), firstapply(true), usesuperlu(usesuperlu_)
580  {
581  params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
582  params.setDebugLevel(verbose_);
583 #if !HAVE_SUPERLU
584  if (usesuperlu == true)
585  {
586  std::cout << "WARNING: You are using AMG without SuperLU!"
587  << " Please consider installing SuperLU,"
588  << " or set the usesuperlu flag to false"
589  << " to suppress this warning." << std::endl;
590  }
591 #endif
592  }
593 
598  void setparams(Parameters params_)
599  {
600  params = params_;
601  }
602 
607  typename V::ElementType norm (const V& v) const
608  {
609  return Backend::native(v).two_norm();
610  }
611 
619  void apply(M& A, V& z, V& r, typename V::ElementType reduction)
620  {
621  Timer watch;
622  MatrixType& mat = Backend::native(A);
623  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
624  Dune::Amg::FirstDiagonal> > Criterion;
625  SmootherArgs smootherArgs;
626  smootherArgs.iterations = 1;
627  smootherArgs.relaxationFactor = 1;
628 
629  Criterion criterion(params);
630  Operator oop(mat);
631  //only construct a new AMG if the matrix changes
632  if (reuse==false || firstapply==true){
633  amg.reset(new AMG(oop, criterion, smootherArgs));
634  firstapply = false;
635  stats.tsetup = watch.elapsed();
636  stats.levels = amg->maxlevels();
637  stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
638  }
639  watch.reset();
640  Dune::InverseOperatorResult stat;
641 
642  Solver<VectorType> solver(oop,*amg,reduction,maxiter,verbose);
643  solver.apply(Backend::native(z),Backend::native(r),stat);
644  stats.tsolve= watch.elapsed();
645  res.converged = stat.converged;
646  res.iterations = stat.iterations;
647  res.elapsed = stat.elapsed;
648  res.reduction = stat.reduction;
649  res.conv_rate = stat.conv_rate;
650  }
651 
652 
658  {
659  return stats;
660  }
661 
662  private:
663  unsigned maxiter;
664  Parameters params;
665  int verbose;
666  bool reuse;
667  bool firstapply;
668  bool usesuperlu;
669  std::shared_ptr<AMG> amg;
670  ISTLAMGStatistics stats;
671  };
672 
675 
681  template<class GO>
683  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
684  {
685 
686  public:
695  ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
696  bool reuse_=false, bool usesuperlu_=true)
697  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::CGSolver>
698  (maxiter_, verbose_, reuse_, usesuperlu_)
699  {}
700  };
701 
707  template<class GO>
709  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
710  {
711 
712  public:
721  ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
722  bool reuse_=false, bool usesuperlu_=true)
723  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::BiCGSTABSolver>
724  (maxiter_, verbose_, reuse_, usesuperlu_)
725  {}
726  };
727 
733  template<class GO>
735  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
736  {
737 
738  public:
747  ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
748  bool reuse_=false, bool usesuperlu_=true)
749  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::BiCGSTABSolver>
750  (maxiter_, verbose_, reuse_, usesuperlu_)
751  {}
752  };
753 
759  template<class GO>
761  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
762  {
763 
764  public:
773  ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1,
774  bool reuse_=false, bool usesuperlu_=true)
775  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSSOR, Dune::LoopSolver>
776  (maxiter_, verbose_, reuse_, usesuperlu_)
777  {}
778  };
779 
785  template<class GO>
787  : public ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
788  {
789 
790  public:
799  ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1,
800  bool reuse_=false, bool usesuperlu_=true)
801  : ISTLBackend_SEQ_AMG<GO, Dune::SeqSOR, Dune::LoopSolver>
802  (maxiter_, verbose_, reuse_, usesuperlu_)
803  {}
804  };
805 
812  : public SequentialNorm, public LinearResultStorage
813  {
814  public :
815 
822  explicit ISTLBackend_SEQ_GMRES_ILU0(int restart_ = 200, int maxiter_ = 5000, int verbose_ = 1)
823  : restart(restart_), maxiter(maxiter_), verbose(verbose_)
824  {}
825 
833  template<class M, class V, class W>
834  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType>::real_type reduction)
835  {
836  using Backend::Native;
837  using Backend::native;
838  Dune::MatrixAdapter<
839  Native<M>,
840  Native<V>,
841  Native<W>
842  > opa(native(A));
843  Dune::SeqILU0<
844  Native<M>,
845  Native<V>,
846  Native<W>
847  > ilu0(native(A), 1.0);
848  Dune::RestartedGMResSolver<Native<V>> solver(opa,ilu0,reduction,restart,maxiter,verbose);
849  Dune::InverseOperatorResult stat;
850  solver.apply(native(z), native(r), stat);
851  res.converged = stat.converged;
852  res.iterations = stat.iterations;
853  res.elapsed = stat.elapsed;
854  res.reduction = stat.reduction;
855  res.conv_rate = stat.conv_rate;
856  }
857 
858  private :
859  int restart, maxiter, verbose;
860  };
861 
864 
865  } // namespace PDELab
866 } // namespace Dune
867 
868 #endif
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:96
ISTLBackend_SEQ_LOOP_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:234
ISTLBackend_SEQ_BCGS_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:284
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:426
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: istl/seqistlsolverbackend.hh:657
Backend for sequential BiCGSTAB solver with Jacobi preconditioner.
Definition: istl/seqistlsolverbackend.hh:242
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:144
virtual void apply(const X &x, Y &y) const
Definition: istl/seqistlsolverbackend.hh:48
Definition: solver.hh:42
double tsetup
The time needed for building the AMG hierarchy (coarsening).
Definition: istl/seqistlsolverbackend.hh:553
void setparams(Parameters params_)
set AMG parameters
Definition: istl/seqistlsolverbackend.hh:598
Linear solver backend for Restarted GMRes preconditioned with ILU(0)
Definition: istl/seqistlsolverbackend.hh:811
Sequential BiCGSTAB solver preconditioned with AMG smoothed by SOR.
Definition: istl/seqistlsolverbackend.hh:734
RFType reduction
Definition: solver.hh:35
Backend for sequential loop solver with Jacobi preconditioner.
Definition: istl/seqistlsolverbackend.hh:226
Sequential BiCGStab solver with ILU0 preconditioner.
Definition: istl/seqistlsolverbackend.hh:307
double tprepare
The needed for computing the parallel information and for adapting the linear system.
Definition: istl/seqistlsolverbackend.hh:547
Solver backend using UMFPack as a direct solver.
Definition: istl/seqistlsolverbackend.hh:450
Backend using a MINRes solver preconditioned by SSOR.
Definition: istl/seqistlsolverbackend.hh:362
ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: istl/seqistlsolverbackend.hh:773
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:480
ISTLBackend_SEQ_BCGS_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:267
Sequential conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: istl/seqistlsolverbackend.hh:682
unsigned int iterations
Definition: solver.hh:33
ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: istl/seqistlsolverbackend.hh:576
ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: istl/seqistlsolverbackend.hh:721
ISTLBackend_SEQ_UMFPack(int maxiter, int verbose_)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:468
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: istl/seqistlsolverbackend.hh:500
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:834
ISTLBackend_SEQ_CG_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:354
Definition: istl/seqistlsolverbackend.hh:73
ISTLBackend_SEQ_BCGS_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:319
ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: istl/seqistlsolverbackend.hh:799
Definition: istl/seqistlsolverbackend.hh:562
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const
Definition: istl/seqistlsolverbackend.hh:54
bool directCoarseLevelSolver
True if a direct solver was used on the coarset level.
Definition: istl/seqistlsolverbackend.hh:557
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:198
Definition: istl/seqistlsolverbackend.hh:42
bool converged
Definition: solver.hh:32
Solver backend using SuperLU as a direct solver.
Definition: istl/seqistlsolverbackend.hh:396
Sequential Loop solver preconditioned with AMG smoothed by SSOR.
Definition: istl/seqistlsolverbackend.hh:760
ISTLBackend_SEQ_CG_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:387
RFType conv_rate
Definition: solver.hh:36
double elapsed
Definition: solver.hh:34
ISTLBackend_SEQ_CG_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:301
ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: istl/seqistlsolverbackend.hh:695
ISTLBackend_SEQ_UMFPack(int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:458
ISTLBackend_SEQ_ExplicitDiagonal()
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:506
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:517
ISTLBackend_SEQ_SuperLU(int maxiter, int verbose_)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:414
Definition: istl/seqistlsolverbackend.hh:170
Backend for sequential BiCGSTAB solver with ILU0 preconditioner.
Definition: istl/seqistlsolverbackend.hh:275
Definition: adaptivity.hh:27
ISTLBackend_SEQ_GMRES_ILU0(int restart_=200, int maxiter_=5000, int verbose_=1)
make linear solver object
Definition: istl/seqistlsolverbackend.hh:822
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:182
Backend for sequential BiCGSTAB solver with SSOR preconditioner.
Definition: istl/seqistlsolverbackend.hh:258
X domain_type
Definition: istl/seqistlsolverbackend.hh:38
Definition: istl/seqistlsolverbackend.hh:35
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/seqistlsolverbackend.hh:607
ISTLBackend_SEQ_BCGS_Jac(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:250
ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:82
int levels
the number of levels in the AMG hierarchy.
Definition: istl/seqistlsolverbackend.hh:549
Backend for sequential conjugate gradient solver with SSOR preconditioner.
Definition: istl/seqistlsolverbackend.hh:345
int iterations
The number of iterations performed until convergence was reached.
Definition: istl/seqistlsolverbackend.hh:555
ISTLBackend_SEQ_MINRES_SSOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:371
ISTLBackend_SEQ_CG_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:337
Sequential congute gradient solver with ILU0 preconditioner.
Definition: istl/seqistlsolverbackend.hh:325
ISTLBackend_SEQ_SuperLU(int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:404
Class providing some statistics of the AMG solver.
Definition: istl/seqistlsolverbackend.hh:541
ISTLBackend_SEQ_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:133
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:52
Definition: istl/seqistlsolverbackend.hh:124
X::field_type field_type
Definition: istl/seqistlsolverbackend.hh:40
Sequential BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: istl/seqistlsolverbackend.hh:708
Definition: solver.hh:16
Sequential Loop solver preconditioned with AMG smoothed by SOR.
Definition: istl/seqistlsolverbackend.hh:786
ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: istl/seqistlsolverbackend.hh:747
void apply(M &A, V &z, W &r, typename W::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:191
Y range_type
Definition: istl/seqistlsolverbackend.hh:39
VTKWriter & w
Definition: function.hh:1101
double tsolve
The time spent in solving the system (without building the hierarchy.
Definition: istl/seqistlsolverbackend.hh:551
void apply(M &A, V &z, V &r, typename V::ElementType reduction)
solve the given linear system
Definition: istl/seqistlsolverbackend.hh:619
ISTLBackend_SEQ_ILUn(int n, double w, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/seqistlsolverbackend.hh:180
OnTheFlyOperator(GOS &gos_)
Definition: istl/seqistlsolverbackend.hh:44
Backend for conjugate gradient solver with Jacobi preconditioner.
Definition: istl/seqistlsolverbackend.hh:379
Backend for sequential conjugate gradient solver with ILU0 preconditioner.
Definition: istl/seqistlsolverbackend.hh:292