dune-pdelab  2.4.1
istl/novlpistlsolverbackend.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_NOVLPISTLSOLVERBACKEND_HH
4 #define DUNE_NOVLPISTLSOLVERBACKEND_HH
5 
6 #include <cstddef>
7 
8 #include <dune/common/deprecated.hh>
9 #include <dune/common/parallel/mpihelper.hh>
10 
11 #include <dune/grid/common/gridenums.hh>
12 
13 #include <dune/istl/io.hh>
14 #include <dune/istl/operators.hh>
15 #include <dune/istl/owneroverlapcopy.hh>
16 #include <dune/istl/paamg/amg.hh>
17 #include <dune/istl/paamg/pinfo.hh>
18 #include <dune/istl/preconditioners.hh>
19 #include <dune/istl/scalarproducts.hh>
20 #include <dune/istl/solvercategory.hh>
21 #include <dune/istl/solvers.hh>
22 #include <dune/istl/superlu.hh>
23 
31 
32 namespace Dune {
33  namespace PDELab {
34 
38 
39  //========================================================
40  // Generic support for nonoverlapping grids
41  //========================================================
42 
44 
52  template<typename GFS, typename M, typename X, typename Y>
54  : public Dune::AssembledLinearOperator<M,X,Y>
55  {
56  public:
64  typedef typename X::field_type field_type;
65 
66  //redefine the category, that is the only difference
67  enum {category=Dune::SolverCategory::nonoverlapping};
68 
70 
80  NonoverlappingOperator (const GFS& gfs_, const M& A)
81  : gfs(gfs_), _A_(A)
82  { }
83 
85 
90  virtual void apply (const X& x, Y& y) const
91  {
92  using Backend::native;
93  // apply local operator; now we have sum y_p = sequential y
94  native(_A_).mv(native(x),native(y));
95 
96  // accumulate y on border
98  if (gfs.gridView().comm().size()>1)
99  gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
100  }
101 
103 
108  virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const
109  {
110  using Backend::native;
111  // apply local operator; now we have sum y_p = sequential y
112  native(_A_).usmv(alpha,native(x),native(y));
113 
114  // accumulate y on border
116  if (gfs.gridView().comm().size()>1)
117  gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
118  }
119 
121  virtual const M& getmat () const
122  {
123  return _A_;
124  }
125 
126  private:
127  const GFS& gfs;
128  const M& _A_;
129  };
130 
131  // parallel scalar product assuming no overlap
132  template<class GFS, class X>
133  class NonoverlappingScalarProduct : public Dune::ScalarProduct<X>
134  {
135  public:
137  typedef X domain_type;
138  typedef typename X::ElementType field_type;
139 
141  enum {category=Dune::SolverCategory::nonoverlapping};
142 
145  NonoverlappingScalarProduct (const GFS& gfs_, const istl::ParallelHelper<GFS>& helper_)
146  : gfs(gfs_), helper(helper_)
147  {}
148 
153  virtual field_type dot (const X& x, const X& y)
154  {
155  // do local scalar product on unique partition
156  field_type sum = helper.disjointDot(x,y);
157 
158  // do global communication
159  return gfs.gridView().comm().sum(sum);
160  }
161 
165  virtual double norm (const X& x)
166  {
167  return sqrt(static_cast<double>(this->dot(x,x)));
168  }
169 
172  void make_consistent (X& x) const
173  {
175  if (gfs.gridView().comm().size()>1)
176  gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
177  }
178 
179  private:
180  const GFS& gfs;
181  const istl::ParallelHelper<GFS>& helper;
182  };
183 
184  // parallel Richardson preconditioner
185  template<class GFS, class X, class Y>
186  class NonoverlappingRichardson : public Dune::Preconditioner<X,Y>
187  {
188  public:
190  typedef X domain_type;
192  typedef Y range_type;
194  typedef typename X::ElementType field_type;
195 
196  // define the category
197  enum {
199  category=Dune::SolverCategory::nonoverlapping
200  };
201 
203  NonoverlappingRichardson (const GFS& gfs_, const istl::ParallelHelper<GFS>& helper_)
204  : gfs(gfs_), helper(helper_)
205  {
206  }
207 
211  virtual void pre (X& x, Y& b) {}
212 
216  virtual void apply (X& v, const Y& d)
217  {
218  v = d;
219  }
220 
224  virtual void post (X& x) {}
225 
226  private:
227  const GFS& gfs;
228  const istl::ParallelHelper<GFS>& helper;
229  };
230 
232 
244  template<typename A, typename X, typename Y>
246  : public Dune::Preconditioner<X,Y>
247  {
248 
250 
251  Diagonal _inverse_diagonal;
252 
253  public:
255 
259  typedef X domain_type;
261 
265  typedef Y range_type;
267  typedef typename X::ElementType field_type;
268 
269  enum {
271  category=Dune::SolverCategory::nonoverlapping
272  };
273 
275 
286  template<typename GFS>
287  NonoverlappingJacobi(const GFS& gfs, const A &m)
288  : _inverse_diagonal(m)
289  {
290  // make the diagonal consistent...
291  typename istl::BlockMatrixDiagonal<A>::template AddMatrixElementVectorDataHandle<GFS> addDH(gfs, _inverse_diagonal);
292  gfs.gridView().communicate(addDH,
293  InteriorBorder_InteriorBorder_Interface,
294  ForwardCommunication);
295 
296  // ... and then invert it
297  _inverse_diagonal.invert();
298 
299  }
300 
302  virtual void pre (X& x, Y& b) {}
303 
305  /*
306  * For this preconditioner, this method works with both consistent and
307  * inconsistent vectors: if d is consistent, v will be consistent, if d
308  * is inconsistent, v will be inconsistent.
309  */
310  virtual void apply (X& v, const Y& d)
311  {
312  _inverse_diagonal.mv(d,v);
313  }
314 
316  virtual void post (X& x) {}
317  };
318 
321 
323  template<class GFS>
325  {
327 
328  public:
335  explicit ISTLBackend_NOVLP_CG_NOPREC (const GFS& gfs_,
336  unsigned maxiter_=5000,
337  int verbose_=1)
338  : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
339  {}
340 
345  template<class V>
346  typename V::ElementType norm (const V& v) const
347  {
348  V x(v); // make a copy because it has to be made consistent
350  PSP psp(gfs,phelper);
351  psp.make_consistent(x);
352  return psp.norm(x);
353  }
354 
362  template<class M, class V, class W>
363  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
364  {
366  POP pop(gfs,A);
368  PSP psp(gfs,phelper);
370  PRICH prich(gfs,phelper);
371  int verb=0;
372  if (gfs.gridView().comm().rank()==0) verb=verbose;
373  Dune::CGSolver<V> solver(pop,psp,prich,reduction,maxiter,verb);
374  Dune::InverseOperatorResult stat;
375  solver.apply(z,r,stat);
376  res.converged = stat.converged;
377  res.iterations = stat.iterations;
378  res.elapsed = stat.elapsed;
379  res.reduction = stat.reduction;
380  res.conv_rate = stat.conv_rate;
381  }
382 
385  {
386  return res;
387  }
388 
389  private:
390  const GFS& gfs;
391  PHELPER phelper;
393  unsigned maxiter;
394  int verbose;
395  };
396 
398  template<class GFS>
400  {
402 
403  const GFS& gfs;
404  PHELPER phelper;
406  unsigned maxiter;
407  int verbose;
408 
409  public:
411 
416  explicit ISTLBackend_NOVLP_CG_Jacobi(const GFS& gfs_,
417  unsigned maxiter_ = 5000,
418  int verbose_ = 1) :
419  gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
420  {}
421 
423 
428  template<class V>
429  typename V::ElementType norm (const V& v) const
430  {
431  V x(v); // make a copy because it has to be made consistent
433  PSP psp(gfs,phelper);
434  psp.make_consistent(x);
435  return psp.norm(x);
436  }
437 
439 
451  template<class M, class V, class W>
452  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
453  {
455  POP pop(gfs,A);
457  PSP psp(gfs,phelper);
458 
459  typedef NonoverlappingJacobi<M,V,W> PPre;
460  PPre ppre(gfs,Backend::native(A));
461 
462  int verb=0;
463  if (gfs.gridView().comm().rank()==0) verb=verbose;
464  CGSolver<V> solver(pop,psp,ppre,reduction,maxiter,verb);
465  InverseOperatorResult stat;
466  solver.apply(z,r,stat);
467  res.converged = stat.converged;
468  res.iterations = stat.iterations;
469  res.elapsed = stat.elapsed;
470  res.reduction = stat.reduction;
471  res.conv_rate = stat.conv_rate;
472  }
473 
476  { return res; }
477  };
478 
480  template<class GFS>
482  {
484 
485  public:
492  explicit ISTLBackend_NOVLP_BCGS_NOPREC (const GFS& gfs_, unsigned maxiter_=5000, int verbose_=1)
493  : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
494  {}
495 
500  template<class V>
501  typename V::ElementType norm (const V& v) const
502  {
503  V x(v); // make a copy because it has to be made consistent
505  PSP psp(gfs,phelper);
506  psp.make_consistent(x);
507  return psp.norm(x);
508  }
509 
517  template<class M, class V, class W>
518  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
519  {
521  POP pop(gfs,A);
523  PSP psp(gfs,phelper);
525  PRICH prich(gfs,phelper);
526  int verb=0;
527  if (gfs.gridView().comm().rank()==0) verb=verbose;
528  Dune::BiCGSTABSolver<V> solver(pop,psp,prich,reduction,maxiter,verb);
529  Dune::InverseOperatorResult stat;
530  solver.apply(z,r,stat);
531  res.converged = stat.converged;
532  res.iterations = stat.iterations;
533  res.elapsed = stat.elapsed;
534  res.reduction = stat.reduction;
535  res.conv_rate = stat.conv_rate;
536  }
537 
540  {
541  return res;
542  }
543 
544  private:
545  const GFS& gfs;
546  PHELPER phelper;
548  unsigned maxiter;
549  int verbose;
550  };
551 
552 
554  template<class GFS>
556  {
558 
559  public:
566  explicit ISTLBackend_NOVLP_BCGS_Jacobi (const GFS& gfs_, unsigned maxiter_=5000, int verbose_=1)
567  : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
568  {}
569 
574  template<class V>
575  typename V::ElementType norm (const V& v) const
576  {
577  V x(v); // make a copy because it has to be made consistent
579  PSP psp(gfs,phelper);
580  psp.make_consistent(x);
581  return psp.norm(x);
582  }
583 
591  template<class M, class V, class W>
592  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
593  {
595  POP pop(gfs,A);
597  PSP psp(gfs,phelper);
598 
599  typedef NonoverlappingJacobi<M,V,W> PPre;
600  PPre ppre(gfs,A);
601 
602  int verb=0;
603  if (gfs.gridView().comm().rank()==0) verb=verbose;
604  Dune::BiCGSTABSolver<V> solver(pop,psp,ppre,reduction,maxiter,verb);
605  Dune::InverseOperatorResult stat;
606  solver.apply(z,r,stat);
607  res.converged = stat.converged;
608  res.iterations = stat.iterations;
609  res.elapsed = stat.elapsed;
610  res.reduction = stat.reduction;
611  res.conv_rate = stat.conv_rate;
612  }
613 
616  {
617  return res;
618  }
619 
620  private:
621  const GFS& gfs;
622  PHELPER phelper;
624  unsigned maxiter;
625  int verbose;
626  };
627 
629  template<typename GFS>
631  {
633 
634  const GFS& gfs;
635  PHELPER phelper;
637 
638  public:
644  explicit ISTLBackend_NOVLP_ExplicitDiagonal(const GFS& gfs_)
645  : gfs(gfs_), phelper(gfs)
646  {}
647 
652  template<class V>
653  typename V::ElementType norm (const V& v) const
654  {
656  V x(v); // make a copy because it has to be made consistent
657  PSP psp(gfs,phelper);
658  psp.make_consistent(x);
659  return psp.norm(x);
660  }
661 
669  template<class M, class V, class W>
670  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
671  {
672  Dune::SeqJac<M,V,W> jac(A,1,1.0);
673  jac.pre(z,r);
674  jac.apply(z,r);
675  jac.post(z);
676  if (gfs.gridView().comm().size()>1)
677  {
679  gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
680  }
681  res.converged = true;
682  res.iterations = 1;
683  res.elapsed = 0.0;
684  res.reduction = reduction;
685  res.conv_rate = reduction; // pow(reduction,1.0/1)
686  }
687 
690  {
691  return res;
692  }
693  };
695 
696 
697  template<class GO,
698  template<class,class,class,int> class Preconditioner,
699  template<class> class Solver>
701  {
702  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
704 
705  public:
713  explicit ISTLBackend_NOVLP_BASE_PREC (const GO& grid_operator, unsigned maxiter_ = 5000, unsigned steps_ = 5, int verbose_ = 1)
714  : _grid_operator(grid_operator)
715  , gfs(grid_operator.trialGridFunctionSpace())
716  , phelper(gfs,verbose_)
717  , maxiter(maxiter_)
718  , steps(steps_)
719  , verbose(verbose_)
720  {}
721 
726  template<class Vector>
727  typename Vector::ElementType norm (const Vector& v) const
728  {
729  Vector x(v); // make a copy because it has to be made consistent
731  PSP psp(gfs,phelper);
732  psp.make_consistent(x);
733  return psp.norm(x);
734  }
735 
743  template<class M, class V, class W>
744  void apply(M& A, V& z, W& r, typename V::ElementType reduction)
745  {
746  using MatrixType = Backend::Native<M>;
747  MatrixType& mat = Backend::native(A);
748  using VectorType = Backend::Native<W>;
749 #if HAVE_MPI
751  _grid_operator.make_consistent(A);
752  istl::assertParallelUG(gfs.gridView().comm());
753  Comm oocc(gfs.gridView().comm(),Dune::SolverCategory::nonoverlapping);
754  phelper.createIndexSetAndProjectForAMG(mat, oocc);
755  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
756  Smoother smoother(mat, steps, 1.0);
757  typedef Dune::NonoverlappingSchwarzScalarProduct<VectorType,Comm> PSP;
758  PSP psp(oocc);
759  typedef Dune::NonoverlappingSchwarzOperator<MatrixType,VectorType,VectorType,Comm> Operator;
760  Operator oop(mat,oocc);
761  typedef Dune::NonoverlappingBlockPreconditioner<Comm, Smoother> ParSmoother;
762  ParSmoother parsmoother(smoother, oocc);
763 #else
764  typedef Preconditioner<MatrixType,VectorType,VectorType,1> ParSmoother;
765  ParSmoother parsmoother(mat, steps, 1.0);
766  typedef Dune::SeqScalarProduct<VectorType> PSP;
767  PSP psp;
768  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
769  Operator oop(mat);
770 #endif
771  int verb=0;
772  if (gfs.gridView().comm().rank()==0) verb=verbose;
773  Solver<VectorType> solver(oop,psp,parsmoother,reduction,maxiter,verb);
774  Dune::InverseOperatorResult stat;
775  //make r consistent
776  if (gfs.gridView().comm().size()>1){
778  gfs.gridView().communicate(adddh,
779  Dune::InteriorBorder_InteriorBorder_Interface,
780  Dune::ForwardCommunication);
781  }
782 
783  solver.apply(z,r,stat);
784  res.converged = stat.converged;
785  res.iterations = stat.iterations;
786  res.elapsed = stat.elapsed;
787  res.reduction = stat.reduction;
788  res.conv_rate = stat.conv_rate;
789  }
790 
793  {
794  return res;
795  }
796 
797  private:
798  const GO& _grid_operator;
799  const GFS& gfs;
800  PHELPER phelper;
802  unsigned maxiter;
803  unsigned steps;
804  int verbose;
805  };
806 
809 
821  template<class GO>
823  : public ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::BiCGSTABSolver>
824  {
825 
826  public:
834  explicit ISTLBackend_NOVLP_BCGS_SSORk (const GO& grid_operator, unsigned maxiter_=5000,
835  int steps_=5, int verbose_=1)
836  : ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::BiCGSTABSolver>(grid_operator, maxiter_, steps_, verbose_)
837  {}
838  };
839 
843  template<class GO>
845  : public ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::CGSolver>
846  {
847 
848  public:
856  explicit ISTLBackend_NOVLP_CG_SSORk (const GO& grid_operator, unsigned maxiter_=5000,
857  int steps_=5, int verbose_=1)
858  : ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::CGSolver>(grid_operator, maxiter_, steps_, verbose_)
859  {}
860  };
863 
864  template<class GO,int s, template<class,class,class,int> class Preconditioner,
865  template<class> class Solver>
867  {
868  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
869  typedef typename istl::ParallelHelper<GFS> PHELPER;
870  typedef typename GO::Traits::Jacobian M;
871  typedef Backend::Native<M> MatrixType;
872  typedef typename GO::Traits::Domain V;
873  typedef Backend::Native<V> VectorType;
875 #if HAVE_MPI
876  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
877  typedef Dune::NonoverlappingBlockPreconditioner<Comm,Smoother> ParSmoother;
878  typedef Dune::NonoverlappingSchwarzOperator<MatrixType,VectorType,VectorType,Comm> Operator;
879 #else
880  typedef Preconditioner<MatrixType,VectorType,VectorType,1> ParSmoother;
881  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
882 #endif
883  typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
884  typedef Dune::Amg::AMG<Operator,VectorType,ParSmoother,Comm> AMG;
885  typedef Dune::Amg::Parameters Parameters;
886 
887  public:
888  ISTLBackend_AMG_NOVLP(const GO& grid_operator, unsigned maxiter_=5000,
889  int verbose_=1, bool reuse_=false,
890  bool usesuperlu_=true)
891  : _grid_operator(grid_operator)
892  , gfs(grid_operator.trialGridFunctionSpace())
893  , phelper(gfs,verbose_)
894  , maxiter(maxiter_)
895  , params(15,2000,1.2,1.6,Dune::Amg::atOnceAccu)
896  , verbose(verbose_)
897  , reuse(reuse_)
898  , firstapply(true)
899  , usesuperlu(usesuperlu_)
900  {
901  params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
902  params.setDebugLevel(verbose_);
903 #if !HAVE_SUPERLU
904  if (phelper.rank() == 0 && usesuperlu == true)
905  {
906  std::cout << "WARNING: You are using AMG without SuperLU!"
907  << " Please consider installing SuperLU,"
908  << " or set the usesuperlu flag to false"
909  << " to suppress this warning." << std::endl;
910  }
911 #endif
912  }
913 
918  void setParameters(const Parameters& params_)
919  {
920  params = params_;
921  }
922 
923  void setparams(Parameters params_) DUNE_DEPRECATED_MSG("setparams() is deprecated, use setParameters() instead")
924  {
925  params = params_;
926  }
927 
935  const Parameters& parameters() const
936  {
937  return params;
938  }
939 
944  typename V::ElementType norm (const V& v) const
945  {
946  V x(v); // make a copy because it has to be made consistent
948  PSP psp(gfs,phelper);
949  psp.make_consistent(x);
950  return psp.norm(x);
951  }
952 
953  void apply(M& A, V& z, V& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
954  {
955  Timer watch;
956  MatrixType& mat = Backend::native(A);
957  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
958  Dune::Amg::FirstDiagonal> > Criterion;
959 #if HAVE_MPI
960  Comm oocc(gfs.gridView().comm(),Dune::SolverCategory::nonoverlapping);
961  _grid_operator.make_consistent(A);
962  phelper.createIndexSetAndProjectForAMG(A, oocc);
963  Dune::NonoverlappingSchwarzScalarProduct<VectorType,Comm> sp(oocc);
964  Operator oop(mat, oocc);
965 #else
966  Comm oocc(gfs.gridView().comm());
967  Operator oop(mat);
968  Dune::SeqScalarProduct<VectorType> sp;
969 #endif
970  SmootherArgs smootherArgs;
971  smootherArgs.iterations = 1;
972  smootherArgs.relaxationFactor = 1;
973  //use noAccu or atOnceAccu
974  Criterion criterion(params);
975  stats.tprepare=watch.elapsed();
976  watch.reset();
977 
978  int verb=0;
979  if (gfs.gridView().comm().rank()==0) verb=verbose;
980  //only construct a new AMG if the matrix changes
981  if (reuse==false || firstapply==true){
982  amg.reset(new AMG(oop, criterion, smootherArgs, oocc));
983  firstapply = false;
984  stats.tsetup = watch.elapsed();
985  stats.levels = amg->maxlevels();
986  stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
987  }
988 
989  Dune::InverseOperatorResult stat;
990  // make r consistent
991  if (gfs.gridView().comm().size()>1) {
993  gfs.gridView().communicate(adddh,
994  Dune::InteriorBorder_InteriorBorder_Interface,
995  Dune::ForwardCommunication);
996  }
997  watch.reset();
998  Solver<VectorType> solver(oop,sp,*amg,reduction,maxiter,verb);
999  solver.apply(Backend::native(z),Backend::native(r),stat);
1000  stats.tsolve= watch.elapsed();
1001  res.converged = stat.converged;
1002  res.iterations = stat.iterations;
1003  res.elapsed = stat.elapsed;
1004  res.reduction = stat.reduction;
1005  res.conv_rate = stat.conv_rate;
1006  }
1007 
1013  {
1014  return stats;
1015  }
1016 
1017  private:
1018  const GO& _grid_operator;
1019  const GFS& gfs;
1020  PHELPER phelper;
1021  unsigned maxiter;
1022  Parameters params;
1023  int verbose;
1024  bool reuse;
1025  bool firstapply;
1026  bool usesuperlu;
1027  std::shared_ptr<AMG> amg;
1028  ISTLAMGStatistics stats;
1029  };
1030 
1033 
1046  template<class GO, int s=96>
1048  : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::CGSolver>
1049  {
1050 
1051  public:
1052  ISTLBackend_NOVLP_CG_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1053  int verbose_=1, bool reuse_=false,
1054  bool usesuperlu_=true)
1055  : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::CGSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1056  {}
1057  };
1058 
1071  template<class GO, int s=96>
1073  : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>
1074  {
1075 
1076  public:
1077  ISTLBackend_NOVLP_BCGS_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1078  int verbose_=1, bool reuse_=false,
1079  bool usesuperlu_=true)
1080  : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1081  {}
1082  };
1083 
1096  template<class GO, int s=96>
1098  : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::LoopSolver>
1099  {
1100 
1101  public:
1102  ISTLBackend_NOVLP_LS_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1103  int verbose_=1, bool reuse_=false,
1104  bool usesuperlu_=true)
1105  : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::LoopSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1106  {}
1107  };
1110 
1111  } // namespace PDELab
1112 } // namespace Dune
1113 
1114 #endif
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:199
void make_consistent(X &x) const
make additive vector consistent
Definition: istl/novlpistlsolverbackend.hh:172
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const
apply operator to x, scale and add:
Definition: istl/novlpistlsolverbackend.hh:108
NonoverlappingScalarProduct(const GFS &gfs_, const istl::ParallelHelper< GFS > &helper_)
Constructor needs to know the grid function space.
Definition: istl/novlpistlsolverbackend.hh:145
Backend::Native< X > domain_type
export type of vectors the matrix is applied to
Definition: istl/novlpistlsolverbackend.hh:60
NonoverlappingJacobi(const GFS &gfs, const A &m)
Constructor.
Definition: istl/novlpistlsolverbackend.hh:287
ISTLBackend_NOVLP_CG_SSORk(const GO &grid_operator, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:856
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:944
virtual const M & getmat() const
extract the matrix
Definition: istl/novlpistlsolverbackend.hh:121
Definition: solver.hh:42
ISTLBackend_NOVLP_CG_NOPREC(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:335
parallel non-overlapping Jacobi preconditioner
Definition: istl/novlpistlsolverbackend.hh:245
Nonoverlapping parallel BiCGSTAB solver preconditioned by block SSOR.
Definition: istl/novlpistlsolverbackend.hh:822
virtual void apply(const X &x, Y &y) const
apply operator
Definition: istl/novlpistlsolverbackend.hh:90
const LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:475
Definition: parallelhelper.hh:51
virtual field_type dot(const X &x, const X &y)
Dot product of two vectors. It is assumed that the vectors are consistent on the interior+border part...
Definition: istl/novlpistlsolverbackend.hh:153
ISTLBackend_NOVLP_CG_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: istl/novlpistlsolverbackend.hh:1052
virtual double norm(const X &x)
Norm of a right-hand side vector. The vector must be consistent on the interior+border partition...
Definition: istl/novlpistlsolverbackend.hh:165
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:183
Nonoverlapping parallel CG solver without preconditioner.
Definition: istl/novlpistlsolverbackend.hh:324
X domain_type
The domain type of the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:190
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/novlpistlsolverbackend.hh:518
void mv(const X &x, Y &y) const
Definition: blockmatrixdiagonal.hh:239
unsigned int iterations
Definition: solver.hh:33
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:211
Nonoverlapping parallel LoopSolver preconditioned with AMG smoothed by SSOR.
Definition: istl/novlpistlsolverbackend.hh:1097
Definition: genericdatahandle.hh:622
Backend::Native< M > matrix_type
export type of matrix
Definition: istl/novlpistlsolverbackend.hh:58
double elapsed
Definition: solver.hh:34
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:429
ISTLBackend_NOVLP_BCGS_SSORk(const GO &grid_operator, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:834
ISTLBackend_AMG_NOVLP(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: istl/novlpistlsolverbackend.hh:888
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:384
Definition: adaptivity.hh:27
Nonoverlapping parallel BiCGStab solver without preconditioner.
Definition: istl/novlpistlsolverbackend.hh:481
virtual void post(X &x)
Clean up.
Definition: istl/novlpistlsolverbackend.hh:316
X::ElementType field_type
The field type of the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:267
Definition: istl/novlpistlsolverbackend.hh:186
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/novlpistlsolverbackend.hh:363
Vector::ElementType norm(const Vector &v) const
Compute global norm of a vector.
Definition: istl/novlpistlsolverbackend.hh:727
bool converged
Definition: solver.hh:32
Dune::Amg::SequentialInformation type
Definition: parallelhelper.hh:415
NonoverlappingRichardson(const GFS &gfs_, const istl::ParallelHelper< GFS > &helper_)
Constructor.
Definition: istl/novlpistlsolverbackend.hh:203
void setparams(Parameters params_)
Definition: istl/novlpistlsolverbackend.hh:923
Definition: istl/novlpistlsolverbackend.hh:67
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: istl/novlpistlsolverbackend.hh:630
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: istl/novlpistlsolverbackend.hh:310
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:615
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:792
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:346
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:689
Backend::Native< Y > range_type
export type of result vectors
Definition: istl/novlpistlsolverbackend.hh:62
ISTLBackend_NOVLP_ExplicitDiagonal(const GFS &gfs_)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:644
Nonoverlapping parallel CG solver with Jacobi preconditioner.
Definition: istl/novlpistlsolverbackend.hh:399
Class providing some statistics of the AMG solver.
Definition: istl/seqistlsolverbackend.hh:541
ISTLBackend_NOVLP_CG_Jacobi(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:416
X domain_type
The domain type of the operator.
Definition: istl/novlpistlsolverbackend.hh:259
virtual void post(X &x)
Clean up.
Definition: istl/novlpistlsolverbackend.hh:224
void invert()
Definition: blockmatrixdiagonal.hh:233
X::ElementType field_type
Definition: istl/novlpistlsolverbackend.hh:138
RFType reduction
Definition: solver.hh:35
void setParameters(const Parameters &params_)
set AMG parameters
Definition: istl/novlpistlsolverbackend.hh:918
Definition: blockmatrixdiagonal.hh:214
Y range_type
The range type of the operator.
Definition: istl/novlpistlsolverbackend.hh:265
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:575
Nonoverlapping parallel BiCGStab solver with Jacobi preconditioner.
Definition: istl/novlpistlsolverbackend.hh:555
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/novlpistlsolverbackend.hh:670
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: istl/novlpistlsolverbackend.hh:1012
Definition: istl/novlpistlsolverbackend.hh:700
X::ElementType field_type
The field type of the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:194
X domain_type
export types
Definition: istl/novlpistlsolverbackend.hh:137
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:302
Y range_type
The range type of the preconditioner.
Definition: istl/novlpistlsolverbackend.hh:192
Operator for the non-overlapping parallel case.
Definition: istl/novlpistlsolverbackend.hh:53
Nonoverlapping parallel CG solver preconditioned with AMG smoothed by SSOR.
Definition: istl/novlpistlsolverbackend.hh:1047
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:113
X::field_type field_type
export type of the entries for x
Definition: istl/novlpistlsolverbackend.hh:64
const Parameters & parameters() const
Get the parameters describing the behaviuour of AMG.
Definition: istl/novlpistlsolverbackend.hh:935
const std::string s
Definition: function.hh:1102
ISTLBackend_NOVLP_BASE_PREC(const GO &grid_operator, unsigned maxiter_=5000, unsigned steps_=5, int verbose_=1)
Constructor.
Definition: istl/novlpistlsolverbackend.hh:713
Nonoverlapping parallel CG solver preconditioned by block SSOR.
Definition: istl/novlpistlsolverbackend.hh:844
void apply(M &A, V &z, W &r, typename V::ElementType reduction)
Solve the given linear system.
Definition: istl/novlpistlsolverbackend.hh:744
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
Definition: istl/novlpistlsolverbackend.hh:953
ISTLBackend_NOVLP_BCGS_Jacobi(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:566
void assertParallelUG(T comm)
Definition: parallelhelper.hh:431
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:501
ISTLBackend_NOVLP_BCGS_NOPREC(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: istl/novlpistlsolverbackend.hh:492
Nonoverlapping parallel BiCGStab solver preconditioned with AMG smoothed by SSOR. ...
Definition: istl/novlpistlsolverbackend.hh:1072
NonoverlappingOperator(const GFS &gfs_, const M &A)
Construct a non-overlapping operator.
Definition: istl/novlpistlsolverbackend.hh:80
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: istl/novlpistlsolverbackend.hh:216
RFType conv_rate
Definition: solver.hh:36
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: istl/novlpistlsolverbackend.hh:539
ISTLBackend_NOVLP_BCGS_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: istl/novlpistlsolverbackend.hh:1077
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/novlpistlsolverbackend.hh:592
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: istl/novlpistlsolverbackend.hh:653
ISTLBackend_NOVLP_LS_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: istl/novlpistlsolverbackend.hh:1102
Definition: istl/novlpistlsolverbackend.hh:866
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: istl/novlpistlsolverbackend.hh:452
Definition: istl/novlpistlsolverbackend.hh:133