dune-pdelab  2.5-dev
ovlpistlsolverbackend.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_OVLPISTLSOLVERBACKEND_HH
4 #define DUNE_PDELAB_BACKEND_ISTL_OVLPISTLSOLVERBACKEND_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 
26 
27 namespace Dune {
28  namespace PDELab {
29 
33 
34  //========================================================
35  // Generic support for overlapping grids
36  // (need to be used with appropriate constraints)
37  //========================================================
38 
39  // operator that resets result to zero at constrained DOFS
40  template<class CC, class M, class X, class Y>
42  : public Dune::AssembledLinearOperator<M,X,Y>
43  {
44  public:
46  typedef M matrix_type;
47  typedef X domain_type;
48  typedef Y range_type;
49  typedef typename X::ElementType field_type;
50 
51  OverlappingOperator (const CC& cc_, const M& A)
52  : cc(cc_), _A_(A)
53  {}
54 
56  virtual void apply (const domain_type& x, range_type& y) const override
57  {
58  using Backend::native;
59  native(_A_).mv(native(x),native(y));
61  }
62 
64  virtual void applyscaleadd (field_type alpha, const domain_type& x, range_type& y) const override
65  {
66  using Backend::native;
67  native(_A_).usmv(alpha,native(x),native(y));
69  }
70 
71  SolverCategory::Category category() const override
72  {
73  return SolverCategory::overlapping;
74  }
75 
77  virtual const M& getmat () const override
78  {
79  return _A_;
80  }
81 
82  private:
83  const CC& cc;
84  const M& _A_;
85  };
86 
87  // new scalar product assuming at least overlap 1
88  // uses unique partitioning of nodes for parallelization
89  template<class GFS, class X>
91  : public Dune::ScalarProduct<X>
92  {
93  public:
95  typedef X domain_type;
96  typedef typename X::ElementType field_type;
97 
100  OverlappingScalarProduct (const GFS& gfs_, const ISTL::ParallelHelper<GFS>& helper_)
101  : gfs(gfs_), helper(helper_)
102  {}
103 
104 
109  virtual field_type dot (const X& x, const X& y) override
110  {
111  // do local scalar product on unique partition
112  field_type sum = helper.disjointDot(x,y);
113 
114  // do global communication
115  return gfs.gridView().comm().sum(sum);
116  }
117 
121  virtual double norm (const X& x) override
122  {
123  return sqrt(static_cast<double>(this->dot(x,x)));
124  }
125 
126  SolverCategory::Category category() const override
127  {
128  return SolverCategory::overlapping;
129  }
130 
131  private:
132  const GFS& gfs;
133  const ISTL::ParallelHelper<GFS>& helper;
134  };
135 
136  // wrapped sequential preconditioner
137  template<class CC, class GFS, class P>
139  : public Dune::Preconditioner<Dune::PDELab::Backend::Vector<GFS,typename P::domain_type::field_type>,
140  Dune::PDELab::Backend::Vector<GFS,typename P::range_type::field_type>>
141  {
142  public:
147 
149  OverlappingWrappedPreconditioner (const GFS& gfs_, P& prec_, const CC& cc_,
150  const ISTL::ParallelHelper<GFS>& helper_)
151  : gfs(gfs_), prec(prec_), cc(cc_), helper(helper_)
152  {}
153 
157  virtual void pre (domain_type& x, range_type& b) override
158  {
159  prec.pre(x,b);
160  }
161 
165  virtual void apply (domain_type& v, const range_type& d) override
166  {
167  range_type dd(d);
168  set_constrained_dofs(cc,0.0,dd);
169  prec.apply(Backend::native(v),Backend::native(dd));
171  if (gfs.gridView().comm().size()>1)
172  gfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
173  }
174 
175  SolverCategory::Category category() const override
176  {
177  return SolverCategory::overlapping;
178  }
179 
183  virtual void post (domain_type& x) override
184  {
185  prec.post(Backend::native(x));
186  }
187 
188  private:
189  const GFS& gfs;
190  P& prec;
191  const CC& cc;
192  const ISTL::ParallelHelper<GFS>& helper;
193  };
194 
195 
196 #if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
197  // exact subdomain solves with UMFPack as preconditioner
198  template<class GFS, class M, class X, class Y>
199  class UMFPackSubdomainSolver : public Dune::Preconditioner<X,Y>
200  {
201  typedef Backend::Native<M> ISTLM;
202 
203  public:
205  typedef X domain_type;
207  typedef Y range_type;
209  typedef typename X::ElementType field_type;
210 
217  UMFPackSubdomainSolver (const GFS& gfs_, const M& A_)
218  : gfs(gfs_), solver(Backend::native(A_),false) // this does the decomposition
219  {}
220 
224  virtual void pre (X& x, Y& b) override {}
225 
229  virtual void apply (X& v, const Y& d) override
230  {
231  Dune::InverseOperatorResult stat;
232  Y b(d); // need copy, since solver overwrites right hand side
233  solver.apply(Backend::native(v),Backend::native(b),stat);
234  if (gfs.gridView().comm().size()>1)
235  {
236  AddDataHandle<GFS,X> adddh(gfs,v);
237  gfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
238  }
239  }
240 
241  SolverCategory::Category category() const override
242  {
243  return SolverCategory::overlapping;
244  }
245 
249  virtual void post (X& x) override {}
250 
251  private:
252  const GFS& gfs;
253  Dune::UMFPack<ISTLM> solver;
254  };
255 #endif
256 
257 #if HAVE_SUPERLU
258  // exact subdomain solves with SuperLU as preconditioner
259  template<class GFS, class M, class X, class Y>
260  class SuperLUSubdomainSolver : public Dune::Preconditioner<X,Y>
261  {
262  typedef Backend::Native<M> ISTLM;
263 
264  public:
266  typedef X domain_type;
268  typedef Y range_type;
270  typedef typename X::ElementType field_type;
271 
278  SuperLUSubdomainSolver (const GFS& gfs_, const M& A_)
279  : gfs(gfs_), solver(Backend::native(A_),false) // this does the decomposition
280  {}
281 
285  virtual void pre (X& x, Y& b) override {}
286 
290  virtual void apply (X& v, const Y& d) override
291  {
292  Dune::InverseOperatorResult stat;
293  Y b(d); // need copy, since solver overwrites right hand side
294  solver.apply(Backend::native(v),Backend::native(b),stat);
295  if (gfs.gridView().comm().size()>1)
296  {
297  AddDataHandle<GFS,X> adddh(gfs,v);
298  gfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
299  }
300  }
301 
302  SolverCategory::Category category() const override
303  {
304  return SolverCategory::overlapping;
305  }
306 
310  virtual void post (X& x) override {}
311 
312  private:
313  const GFS& gfs;
314  Dune::SuperLU<ISTLM> solver;
315  };
316 
317  // exact subdomain solves with SuperLU as preconditioner
318  template<class GFS, class M, class X, class Y>
319  class RestrictedSuperLUSubdomainSolver : public Dune::Preconditioner<X,Y>
320  {
321  typedef typename M::Container ISTLM;
322 
323  public:
325  typedef X domain_type;
327  typedef Y range_type;
329  typedef typename X::ElementType field_type;
330 
338  RestrictedSuperLUSubdomainSolver (const GFS& gfs_, const M& A_,
339  const ISTL::ParallelHelper<GFS>& helper_)
340  : gfs(gfs_), solver(Backend::native(A_),false), helper(helper_) // this does the decomposition
341  {}
342 
346  virtual void pre (X& x, Y& b) override {}
347 
351  virtual void apply (X& v, const Y& d) override
352  {
353  using Backend::native;
354  Dune::InverseOperatorResult stat;
355  Y b(d); // need copy, since solver overwrites right hand side
356  solver.apply(native(v),native(b),stat);
357  if (gfs.gridView().comm().size()>1)
358  {
359  helper.maskForeignDOFs(native(v));
360  AddDataHandle<GFS,X> adddh(gfs,v);
361  gfs.gridView().communicate(adddh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
362  }
363  }
364 
365  SolverCategory::Category category() const override
366  {
367  return SolverCategory::overlapping;
368  }
369 
373  virtual void post (X& x) override {}
374 
375  private:
376  const GFS& gfs;
377  Dune::SuperLU<ISTLM> solver;
378  const ISTL::ParallelHelper<GFS>& helper;
379  };
380 #endif
381 
382  template<typename GFS>
384  {
385  public:
387  : gfs(gfs_), helper(gfs_)
388  {}
389 
394  template<typename X>
395  typename X::ElementType dot (const X& x, const X& y) const
396  {
397  // do local scalar product on unique partition
398  typename X::ElementType sum = helper.disjointDot(x,y);
399 
400  // do global communication
401  return gfs.gridView().comm().sum(sum);
402  }
403 
407  template<typename X>
408  typename Dune::template FieldTraits<typename X::ElementType >::real_type norm (const X& x) const
409  {
410  using namespace std;
411  return sqrt(static_cast<double>(this->dot(x,x)));
412  }
413 
415  {
416  return helper;
417  }
418 
419  // need also non-const version;
420  ISTL::ParallelHelper<GFS>& parallelHelper() // P.B.: needed for createIndexSetAndProjectForAMG
421  {
422  return helper;
423  }
424 
425  private:
426  const GFS& gfs;
428  };
429 
430 
431  template<typename GFS, typename X>
433  : public ScalarProduct<X>
434  {
435  public:
436  SolverCategory::Category category() const override
437  {
438  return SolverCategory::overlapping;
439  }
440 
442  : implementation(implementation_)
443  {}
444 
445  virtual typename X::Container::field_type dot(const X& x, const X& y) override
446  {
447  return implementation.dot(x,y);
448  }
449 
450  virtual typename X::Container::field_type norm (const X& x) override
451  {
452  using namespace std;
453  return sqrt(static_cast<double>(this->dot(x,x)));
454  }
455 
456  private:
457  const OVLPScalarProductImplementation<GFS>& implementation;
458  };
459 
460  template<class GFS, class C,
461  template<class,class,class,int> class Preconditioner,
462  template<class> class Solver>
465  {
466  public:
475  ISTLBackend_OVLP_Base (const GFS& gfs_, const C& c_, unsigned maxiter_=5000,
476  int steps_=5, int verbose_=1)
477  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), steps(steps_), verbose(verbose_)
478  {}
479 
487  template<class M, class V, class W>
488  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
489  {
490  using Backend::Native;
491  using Backend::native;
492  typedef OverlappingOperator<C,M,V,W> POP;
493  POP pop(c,A);
494  typedef OVLPScalarProduct<GFS,V> PSP;
495  PSP psp(*this);
496  typedef Preconditioner<
497  Native<M>,
498  Native<V>,
499  Native<W>,
500  1
501  > SeqPrec;
502  SeqPrec seqprec(native(A),steps,1.0);
504  WPREC wprec(gfs,seqprec,c,this->parallelHelper());
505  int verb=0;
506  if (gfs.gridView().comm().rank()==0) verb=verbose;
507  Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
508  Dune::InverseOperatorResult stat;
509  solver.apply(z,r,stat);
510  res.converged = stat.converged;
511  res.iterations = stat.iterations;
512  res.elapsed = stat.elapsed;
513  res.reduction = stat.reduction;
514  res.conv_rate = stat.conv_rate;
515  }
516  private:
517  const GFS& gfs;
518  const C& c;
519  unsigned maxiter;
520  int steps;
521  int verbose;
522  };
523 
524  // Base class for ILU0 as preconditioner
525  template<class GFS, class C,
526  template<class> class Solver>
529  {
530  public:
538  ISTLBackend_OVLP_ILU0_Base (const GFS& gfs_, const C& c_, unsigned maxiter_=5000, int verbose_=1)
539  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), verbose(verbose_)
540  {}
541 
549  template<class M, class V, class W>
550  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
551  {
552  using Backend::Native;
553  using Backend::native;
554  typedef OverlappingOperator<C,M,V,W> POP;
555  POP pop(c,A);
556  typedef OVLPScalarProduct<GFS,V> PSP;
557  PSP psp(*this);
558  typedef SeqILU0<
559  Native<M>,
560  Native<V>,
561  Native<W>,
562  1
563  > SeqPrec;
564  SeqPrec seqprec(native(A),1.0);
566  WPREC wprec(gfs,seqprec,c,this->parallelHelper());
567  int verb=0;
568  if (gfs.gridView().comm().rank()==0) verb=verbose;
569  Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
570  Dune::InverseOperatorResult stat;
571  solver.apply(z,r,stat);
572  res.converged = stat.converged;
573  res.iterations = stat.iterations;
574  res.elapsed = stat.elapsed;
575  res.reduction = stat.reduction;
576  res.conv_rate = stat.conv_rate;
577  }
578  private:
579  const GFS& gfs;
580  const C& c;
581  unsigned maxiter;
582  int steps;
583  int verbose;
584  };
585 
586  // Base class for ILUn as preconditioner
587  template<class GFS, class C,
588  template<class> class Solver>
591  {
592  public:
601  ISTLBackend_OVLP_ILUn_Base (const GFS& gfs_, const C& c_, int n_=1, unsigned maxiter_=5000, int verbose_=1)
602  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), n(n_), maxiter(maxiter_), verbose(verbose_)
603  {}
604 
612  template<class M, class V, class W>
613  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
614  {
615  using Backend::Native;
616  using Backend::native;
617  typedef OverlappingOperator<C,M,V,W> POP;
618  POP pop(c,A);
619  typedef OVLPScalarProduct<GFS,V> PSP;
620  PSP psp(*this);
621  typedef SeqILUn<
622  Native<M>,
623  Native<V>,
624  Native<W>,
625  1
626  > SeqPrec;
627  SeqPrec seqprec(native(A),n,1.0);
629  WPREC wprec(gfs,seqprec,c,this->parallelHelper());
630  int verb=0;
631  if (gfs.gridView().comm().rank()==0) verb=verbose;
632  Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
633  Dune::InverseOperatorResult stat;
634  solver.apply(z,r,stat);
635  res.converged = stat.converged;
636  res.iterations = stat.iterations;
637  res.elapsed = stat.elapsed;
638  res.reduction = stat.reduction;
639  res.conv_rate = stat.conv_rate;
640  }
641  private:
642  const GFS& gfs;
643  const C& c;
644  int n;
645  unsigned maxiter;
646  int steps;
647  int verbose;
648  };
649 
652 
658  template<class GFS, class CC>
660  : public ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::BiCGSTABSolver>
661  {
662  public:
671  ISTLBackend_OVLP_BCGS_SSORk (const GFS& gfs, const CC& cc, unsigned maxiter=5000,
672  int steps=5, int verbose=1)
673  : ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::BiCGSTABSolver>(gfs, cc, maxiter, steps, verbose)
674  {}
675  };
681  template<class GFS, class CC>
683  : public ISTLBackend_OVLP_ILU0_Base<GFS,CC,Dune::BiCGSTABSolver>
684  {
685  public:
693  ISTLBackend_OVLP_BCGS_ILU0 (const GFS& gfs, const CC& cc, unsigned maxiter=5000, int verbose=1)
694  : ISTLBackend_OVLP_ILU0_Base<GFS,CC,Dune::BiCGSTABSolver>(gfs, cc, maxiter, verbose)
695  {}
696  };
702  template<class GFS, class CC>
704  : public ISTLBackend_OVLP_ILUn_Base<GFS,CC,Dune::BiCGSTABSolver>
705  {
706  public:
715  ISTLBackend_OVLP_BCGS_ILUn (const GFS& gfs, const CC& cc, int n=1, unsigned maxiter=5000, int verbose=1)
716  : ISTLBackend_OVLP_ILUn_Base<GFS,CC,Dune::BiCGSTABSolver>(gfs, cc, n, maxiter, verbose)
717  {}
718  };
724  template<class GFS, class CC>
726  : public ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::CGSolver>
727  {
728  public:
737  ISTLBackend_OVLP_CG_SSORk (const GFS& gfs, const CC& cc, unsigned maxiter=5000,
738  int steps=5, int verbose=1)
739  : ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::CGSolver>(gfs, cc, maxiter, steps, verbose)
740  {}
741  };
742 
748  template<class GFS, class CC>
751  {
752  public:
760  ISTLBackend_OVLP_GMRES_ILU0 (const GFS& gfs_, const CC& cc_, unsigned maxiter_=5000, int verbose_=1,
761  int restart_ = 20)
762  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), cc(cc_), maxiter(maxiter_), verbose(verbose_),
763  restart(restart_)
764  {}
765 
772  template<class M, class V, class W>
773  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
774  {
775  using Backend::Native;
776  using Backend::native;
777  typedef OverlappingOperator<CC,M,V,W> POP;
778  POP pop(cc,A);
779  typedef OVLPScalarProduct<GFS,V> PSP;
780  PSP psp(*this);
781  typedef SeqILU0<
782  Native<M>,
783  Native<V>,
784  Native<W>,
785  1
786  > SeqPrec;
787  SeqPrec seqprec(native(A),1.0);
789  WPREC wprec(gfs,seqprec,cc,this->parallelHelper());
790  int verb=0;
791  if (gfs.gridView().comm().rank()==0) verb=verbose;
792  RestartedGMResSolver<V> solver(pop,psp,wprec,reduction,restart,maxiter,verb);
793  Dune::InverseOperatorResult stat;
794  solver.apply(z,r,stat);
795  res.converged = stat.converged;
796  res.iterations = stat.iterations;
797  res.elapsed = stat.elapsed;
798  res.reduction = stat.reduction;
799  res.conv_rate = stat.conv_rate;
800  }
801 
802  private:
803  const GFS& gfs;
804  const CC& cc;
805  unsigned maxiter;
806  int steps;
807  int verbose;
808  int restart;
809  };
810 
812 
813  template<class GFS, class C, template<typename> class Solver>
816  {
817  public:
825  ISTLBackend_OVLP_SuperLU_Base (const GFS& gfs_, const C& c_, unsigned maxiter_=5000,
826  int verbose_=1)
827  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), verbose(verbose_)
828  {}
829 
837  template<class M, class V, class W>
838  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
839  {
840  typedef OverlappingOperator<C,M,V,W> POP;
841  POP pop(c,A);
842  typedef OVLPScalarProduct<GFS,V> PSP;
843  PSP psp(*this);
844 #if HAVE_SUPERLU
845  typedef SuperLUSubdomainSolver<GFS,M,V,W> PREC;
846  PREC prec(gfs,A);
847  int verb=0;
848  if (gfs.gridView().comm().rank()==0) verb=verbose;
849  Solver<V> solver(pop,psp,prec,reduction,maxiter,verb);
850  Dune::InverseOperatorResult stat;
851  solver.apply(z,r,stat);
852  res.converged = stat.converged;
853  res.iterations = stat.iterations;
854  res.elapsed = stat.elapsed;
855  res.reduction = stat.reduction;
856  res.conv_rate = stat.conv_rate;
857 #else
858  std::cout << "No superLU support, please install and configure it." << std::endl;
859 #endif
860  }
861 
862  private:
863  const GFS& gfs;
864  const C& c;
865  unsigned maxiter;
866  int verbose;
867  };
868 
870 
871  template<class GFS, class C, template<typename> class Solver>
874  {
875  public:
883  ISTLBackend_OVLP_UMFPack_Base (const GFS& gfs_, const C& c_, unsigned maxiter_=5000,
884  int verbose_=1)
885  : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), verbose(verbose_)
886  {}
887 
895  template<class M, class V, class W>
896  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
897  {
898  typedef OverlappingOperator<C,M,V,W> POP;
899  POP pop(c,A);
900  typedef OVLPScalarProduct<GFS,V> PSP;
901  PSP psp(*this);
902 #if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
904  PREC prec(gfs,A);
905  int verb=0;
906  if (gfs.gridView().comm().rank()==0) verb=verbose;
907  Solver<V> solver(pop,psp,prec,reduction,maxiter,verb);
908  Dune::InverseOperatorResult stat;
909  solver.apply(z,r,stat);
910  res.converged = stat.converged;
911  res.iterations = stat.iterations;
912  res.elapsed = stat.elapsed;
913  res.reduction = stat.reduction;
914  res.conv_rate = stat.conv_rate;
915 #else
916  std::cout << "No UMFPack support, please install and configure it." << std::endl;
917 #endif
918  }
919 
920  private:
921  const GFS& gfs;
922  const C& c;
923  unsigned maxiter;
924  int verbose;
925  };
926 
929 
934  template<class GFS, class CC>
936  : public ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::BiCGSTABSolver>
937  {
938  public:
939 
947  ISTLBackend_OVLP_BCGS_SuperLU (const GFS& gfs_, const CC& cc_, unsigned maxiter_=5000,
948  int verbose_=1)
949  : ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::BiCGSTABSolver>(gfs_,cc_,maxiter_,verbose_)
950  {}
951  };
952 
958  template<class GFS, class CC>
960  : public ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::CGSolver>
961  {
962  public:
963 
971  ISTLBackend_OVLP_CG_SuperLU (const GFS& gfs_, const CC& cc_,
972  unsigned maxiter_=5000,
973  int verbose_=1)
974  : ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::CGSolver>(gfs_,cc_,maxiter_,verbose_)
975  {}
976  };
977 
983  template<class GFS, class CC>
985  : public ISTLBackend_OVLP_UMFPack_Base<GFS,CC,Dune::CGSolver>
986  {
987  public:
988 
996  ISTLBackend_OVLP_CG_UMFPack (const GFS& gfs_, const CC& cc_,
997  unsigned maxiter_=5000,
998  int verbose_=1)
999  : ISTLBackend_OVLP_UMFPack_Base<GFS,CC,Dune::CGSolver>(gfs_,cc_,maxiter_,verbose_)
1000  {}
1001  };
1002 
1003 
1007  template<class GFS>
1009  : public LinearResultStorage
1010  {
1011  public:
1016  explicit ISTLBackend_OVLP_ExplicitDiagonal (const GFS& gfs_)
1017  : gfs(gfs_)
1018  {}
1019 
1021  : gfs(other_.gfs)
1022  {}
1023 
1028  template<class V>
1029  typename V::ElementType norm(const V& v) const
1030  {
1031  static_assert
1033  "ISTLBackend_OVLP_ExplicitDiagonal::norm() should not be "
1034  "neccessary, so we skipped the implementation. If you have a "
1035  "scenario where you need it, please implement it or report back to "
1036  "us.");
1037  }
1038 
1046  template<class M, class V, class W>
1047  void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
1048  {
1049  using Backend::Native;
1050  using Backend::native;
1051  Dune::SeqJac<
1052  Native<M>,
1053  Native<V>,
1054  Native<W>
1055  > jac(native(A),1,1.0);
1056  jac.pre(native(z),native(r));
1057  jac.apply(native(z),native(r));
1058  jac.post(native(z));
1059  if (gfs.gridView().comm().size()>1)
1060  {
1061  CopyDataHandle<GFS,V> copydh(gfs,z);
1062  gfs.gridView().communicate(copydh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
1063  }
1064  res.converged = true;
1065  res.iterations = 1;
1066  res.elapsed = 0.0;
1067  res.reduction = static_cast<double>(reduction);
1068  res.conv_rate = static_cast<double>(reduction); // pow(reduction,1.0/1)
1069  }
1070 
1071  private:
1072  const GFS& gfs;
1073  };
1075 
1076  template<class GO, int s, template<class,class,class,int> class Preconditioner,
1077  template<class> class Solver>
1079  {
1080  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1082  typedef typename GO::Traits::Jacobian M;
1083  typedef Backend::Native<M> MatrixType;
1084  typedef typename GO::Traits::Domain V;
1085  typedef Backend::Native<V> VectorType;
1087 #if HAVE_MPI
1088  typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
1089  typedef Dune::BlockPreconditioner<VectorType,VectorType,Comm,Smoother> ParSmoother;
1090  typedef Dune::OverlappingSchwarzOperator<MatrixType,VectorType,VectorType,Comm> Operator;
1091 #else
1092  typedef Preconditioner<MatrixType,VectorType,VectorType,1> ParSmoother;
1093  typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
1094 #endif
1095  typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
1096  typedef Dune::Amg::AMG<Operator,VectorType,ParSmoother,Comm> AMG;
1097 
1098  typedef typename V::ElementType RF;
1099 
1100  public:
1101 
1105  typedef Dune::Amg::Parameters Parameters;
1106 
1107  public:
1108  ISTLBackend_AMG(const GFS& gfs_, unsigned maxiter_=5000,
1109  int verbose_=1, bool reuse_=false,
1110  bool usesuperlu_=true)
1111  : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), params(15,2000),
1112  verbose(verbose_), reuse(reuse_), firstapply(true),
1113  usesuperlu(usesuperlu_)
1114  {
1115  params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
1116  params.setDebugLevel(verbose_);
1117 #if !HAVE_SUPERLU
1118  if (gfs.gridView().comm().rank() == 0 && usesuperlu == true)
1119  {
1120  std::cout << "WARNING: You are using AMG without SuperLU!"
1121  << " Please consider installing SuperLU,"
1122  << " or set the usesuperlu flag to false"
1123  << " to suppress this warning." << std::endl;
1124  }
1125 #endif
1126  }
1127 
1132  void setParameters(const Parameters& params_)
1133  {
1134  params = params_;
1135  }
1136 
1144  const Parameters& parameters() const
1145  {
1146  return params;
1147  }
1148 
1150  void setReuse(bool reuse_)
1151  {
1152  reuse = reuse_;
1153  }
1154 
1156  bool getReuse() const
1157  {
1158  return reuse;
1159  }
1160 
1165  typename V::ElementType norm (const V& v) const
1166  {
1167  typedef OverlappingScalarProduct<GFS,V> PSP;
1168  PSP psp(gfs,phelper);
1169  return psp.norm(v);
1170  }
1171 
1179  void apply(M& A, V& z, V& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
1180  {
1181  Timer watch;
1182  Comm oocc(gfs.gridView().comm());
1183  MatrixType& mat=Backend::native(A);
1184  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
1185  Dune::Amg::FirstDiagonal> > Criterion;
1186 #if HAVE_MPI
1187  phelper.createIndexSetAndProjectForAMG(A, oocc);
1188  Operator oop(mat, oocc);
1189  Dune::OverlappingSchwarzScalarProduct<VectorType,Comm> sp(oocc);
1190 #else
1191  Operator oop(mat);
1192  Dune::SeqScalarProduct<VectorType> sp;
1193 #endif
1194  SmootherArgs smootherArgs;
1195  smootherArgs.iterations = 1;
1196  smootherArgs.relaxationFactor = 1;
1197  Criterion criterion(params);
1198  stats.tprepare=watch.elapsed();
1199  watch.reset();
1200 
1201  int verb=0;
1202  if (gfs.gridView().comm().rank()==0) verb=verbose;
1203  //only construct a new AMG if the matrix changes
1204  if (reuse==false || firstapply==true){
1205  amg.reset(new AMG(oop, criterion, smootherArgs, oocc));
1206  firstapply = false;
1207  stats.tsetup = watch.elapsed();
1208  stats.levels = amg->maxlevels();
1209  stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
1210  }
1211  watch.reset();
1212  Solver<VectorType> solver(oop,sp,*amg,RF(reduction),maxiter,verb);
1213  Dune::InverseOperatorResult stat;
1214 
1215  solver.apply(Backend::native(z),Backend::native(r),stat);
1216  stats.tsolve= watch.elapsed();
1217  res.converged = stat.converged;
1218  res.iterations = stat.iterations;
1219  res.elapsed = stat.elapsed;
1220  res.reduction = stat.reduction;
1221  res.conv_rate = stat.conv_rate;
1222  }
1223 
1229  {
1230  return stats;
1231  }
1232 
1233  private:
1234  const GFS& gfs;
1235  PHELPER phelper;
1236  unsigned maxiter;
1237  Parameters params;
1238  int verbose;
1239  bool reuse;
1240  bool firstapply;
1241  bool usesuperlu;
1242  shared_ptr<AMG> amg;
1243  ISTLAMGStatistics stats;
1244  };
1245 
1248 
1255  template<class GO, int s=96>
1257  : public ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::CGSolver>
1258  {
1259  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1260  public:
1270  ISTLBackend_CG_AMG_SSOR(const GFS& gfs_, unsigned maxiter_=5000,
1271  int verbose_=1, bool reuse_=false,
1272  bool usesuperlu_=true)
1273  : ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::CGSolver>
1274  (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
1275  {}
1276  };
1277 
1284  template<class GO, int s=96>
1286  : public ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>
1287  {
1288  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1289  public:
1299  ISTLBackend_BCGS_AMG_SSOR(const GFS& gfs_, unsigned maxiter_=5000,
1300  int verbose_=1, bool reuse_=false,
1301  bool usesuperlu_=true)
1302  : ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>
1303  (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
1304  {}
1305  };
1306 
1313  template<class GO, int s=96>
1315  : public ISTLBackend_AMG<GO, s, Dune::SeqILU0, Dune::BiCGSTABSolver>
1316  {
1317  typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1318  public:
1328  ISTLBackend_BCGS_AMG_ILU0(const GFS& gfs_, unsigned maxiter_=5000,
1329  int verbose_=1, bool reuse_=false,
1330  bool usesuperlu_=true)
1331  : ISTLBackend_AMG<GO, s, Dune::SeqILU0, Dune::BiCGSTABSolver>
1332  (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
1333  {}
1334  };
1335 
1337 
1339 
1340  } // namespace PDELab
1341 } // namespace Dune
1342 
1343 #endif // DUNE_PDELAB_BACKEND_ISTL_OVLPISTLSOLVERBACKEND_HH
Dune::PDELab::CopyDataHandle
Definition: genericdatahandle.hh:728
Dune::PDELab::OVLPScalarProductImplementation::parallelHelper
const ISTL::ParallelHelper< GFS > & parallelHelper() const
Definition: ovlpistlsolverbackend.hh:414
bcrsmatrix.hh
Dune::PDELab::ISTLBackend_AMG::getReuse
bool getReuse() const
Return whether the AMG is reused during call to apply()
Definition: ovlpistlsolverbackend.hh:1156
Dune::PDELab::ISTLBackend_BCGS_AMG_ILU0
Overlapping parallel BiCGStab solver preconditioned with AMG smoothed by ILU0.
Definition: ovlpistlsolverbackend.hh:1314
Dune::PDELab::ISTLBackend_AMG
Definition: ovlpistlsolverbackend.hh:1078
Dune::PDELab::OVLPScalarProduct::OVLPScalarProduct
OVLPScalarProduct(const OVLPScalarProductImplementation< GFS > &implementation_)
Definition: ovlpistlsolverbackend.hh:441
Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR::ISTLBackend_BCGS_AMG_SSOR
ISTLBackend_BCGS_AMG_SSOR(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1299
Dune::PDELab::OverlappingScalarProduct::OverlappingScalarProduct
OverlappingScalarProduct(const GFS &gfs_, const ISTL::ParallelHelper< GFS > &helper_)
Constructor needs to know the grid function space.
Definition: ovlpistlsolverbackend.hh:100
Dune::PDELab::UMFPackSubdomainSolver::UMFPackSubdomainSolver
UMFPackSubdomainSolver(const GFS &gfs_, const M &A_)
Constructor.
Definition: ovlpistlsolverbackend.hh:217
Dune::PDELab::OverlappingOperator::range_type
Y range_type
Definition: ovlpistlsolverbackend.hh:48
Dune::PDELab::ISTLBackend_OVLP_BCGS_SSORk
Overlapping parallel BiCGStab solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:659
Dune::PDELab::ISTLBackend_OVLP_CG_SSORk::ISTLBackend_OVLP_CG_SSORk
ISTLBackend_OVLP_CG_SSORk(const GFS &gfs, const CC &cc, unsigned maxiter=5000, int steps=5, int verbose=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:737
Dune::PDELab::ISTLBackend_OVLP_CG_UMFPack
Overlapping parallel CG solver with UMFPack preconditioner.
Definition: ovlpistlsolverbackend.hh:984
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal::norm
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: ovlpistlsolverbackend.hh:1029
Dune::PDELab::ISTLBackend_AMG::parameters
const Parameters & parameters() const
Get the parameters describing the behaviuour of AMG.
Definition: ovlpistlsolverbackend.hh:1144
Dune::PDELab::ISTLBackend_OVLP_UMFPack_Base
Definition: ovlpistlsolverbackend.hh:872
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal::ISTLBackend_OVLP_ExplicitDiagonal
ISTLBackend_OVLP_ExplicitDiagonal(const ISTLBackend_OVLP_ExplicitDiagonal &other_)
Definition: ovlpistlsolverbackend.hh:1020
Dune::PDELab::ISTLBackend_OVLP_Base::ISTLBackend_OVLP_Base
ISTLBackend_OVLP_Base(const GFS &gfs_, const C &c_, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:475
Dune::PDELab::ISTLBackend_OVLP_UMFPack_Base::ISTLBackend_OVLP_UMFPack_Base
ISTLBackend_OVLP_UMFPack_Base(const GFS &gfs_, const C &c_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:883
Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR
Overlapping parallel BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: ovlpistlsolverbackend.hh:1285
Dune::PDELab::ISTLBackend_OVLP_ILUn_Base::ISTLBackend_OVLP_ILUn_Base
ISTLBackend_OVLP_ILUn_Base(const GFS &gfs_, const C &c_, int n_=1, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:601
parallelhelper.hh
Dune::PDELab::ISTLBackend_OVLP_BCGS_ILUn::ISTLBackend_OVLP_BCGS_ILUn
ISTLBackend_OVLP_BCGS_ILUn(const GFS &gfs, const CC &cc, int n=1, unsigned maxiter=5000, int verbose=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:715
Dune::PDELab::Backend::native
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:192
Dune::PDELab::ISTLBackend_AMG::ISTLBackend_AMG
ISTLBackend_AMG(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: ovlpistlsolverbackend.hh:1108
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::OVLPScalarProductImplementation
Definition: ovlpistlsolverbackend.hh:383
Dune::PDELab::OverlappingOperator::apply
virtual void apply(const domain_type &x, range_type &y) const override
apply operator to x:
Definition: ovlpistlsolverbackend.hh:56
genericdatahandle.hh
Dune::PDELab::ISTLBackend_OVLP_GMRES_ILU0
Overlapping parallel restarted GMRes solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:749
Dune::PDELab::OverlappingOperator::OverlappingOperator
OverlappingOperator(const CC &cc_, const M &A)
Definition: ovlpistlsolverbackend.hh:51
Dune::PDELab::OverlappingOperator::matrix_type
M matrix_type
export types
Definition: ovlpistlsolverbackend.hh:46
Dune::PDELab::OverlappingOperator::domain_type
X domain_type
Definition: ovlpistlsolverbackend.hh:47
Dune::PDELab::ISTLBackend_OVLP_CG_SuperLU
Overlapping parallel CG solver with SuperLU preconditioner.
Definition: ovlpistlsolverbackend.hh:959
Dune::PDELab::OverlappingOperator::getmat
virtual const M & getmat() const override
get matrix via *
Definition: ovlpistlsolverbackend.hh:77
Dune::PDELab::OverlappingOperator::category
SolverCategory::Category category() const override
Definition: ovlpistlsolverbackend.hh:71
Dune::PDELab::ISTLBackend_CG_AMG_SSOR::ISTLBackend_CG_AMG_SSOR
ISTLBackend_CG_AMG_SSOR(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1270
Dune::PDELab::ISTLBackend_OVLP_ILUn_Base::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:613
Dune::PDELab::ISTLBackend_AMG::Parameters
Dune::Amg::Parameters Parameters
Parameters object to customize matrix hierachy building.
Definition: ovlpistlsolverbackend.hh:1105
Dune::PDELab::ISTL::ParallelHelper
Definition: parallelhelper.hh:50
Dune::PDELab::ISTLBackend_OVLP_Base
Definition: ovlpistlsolverbackend.hh:463
Dune::PDELab::ISTLBackend_OVLP_GMRES_ILU0::ISTLBackend_OVLP_GMRES_ILU0
ISTLBackend_OVLP_GMRES_ILU0(const GFS &gfs_, const CC &cc_, unsigned maxiter_=5000, int verbose_=1, int restart_=20)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:760
Dune::PDELab::OVLPScalarProduct::norm
virtual X::Container::field_type norm(const X &x) override
Definition: ovlpistlsolverbackend.hh:450
vector.hh
Dune::PDELab::ISTLBackend_OVLP_CG_SSORk
Overlapping parallel CGS solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:725
Dune::PDELab::OverlappingWrappedPreconditioner::post
virtual void post(domain_type &x) override
Clean up.
Definition: ovlpistlsolverbackend.hh:183
Dune::PDELab::OverlappingOperator::field_type
X::ElementType field_type
Definition: ovlpistlsolverbackend.hh:49
Dune::PDELab::ISTLBackend_OVLP_Base::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:488
Dune::PDELab::ISTLBackend_OVLP_GMRES_ILU0::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:773
Dune::PDELab::OverlappingWrappedPreconditioner::apply
virtual void apply(domain_type &v, const range_type &d) override
Apply the preconditioner.
Definition: ovlpistlsolverbackend.hh:165
Dune::PDELab::OverlappingScalarProduct
Definition: ovlpistlsolverbackend.hh:90
Dune::PDELab::OverlappingScalarProduct::norm
virtual double norm(const X &x) override
Norm of a right-hand side vector. The vector must be consistent on the interior+border partition.
Definition: ovlpistlsolverbackend.hh:121
Dune::PDELab::OVLPScalarProductImplementation::OVLPScalarProductImplementation
OVLPScalarProductImplementation(const GFS &gfs_)
Definition: ovlpistlsolverbackend.hh:386
Dune::PDELab::LinearResultStorage
Definition: solver.hh:42
Dune::PDELab::ISTLBackend_OVLP_SuperLU_Base::ISTLBackend_OVLP_SuperLU_Base
ISTLBackend_OVLP_SuperLU_Base(const GFS &gfs_, const C &c_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:825
Dune::PDELab::OverlappingScalarProduct::category
SolverCategory::Category category() const override
Definition: ovlpistlsolverbackend.hh:126
Dune::PDELab::ISTLBackend_AMG::setReuse
void setReuse(bool reuse_)
Set whether the AMG should be reused again during call to apply().
Definition: ovlpistlsolverbackend.hh:1150
Dune::PDELab::OVLPScalarProductImplementation::norm
Dune::template FieldTraits< typename X::ElementType >::real_type norm(const X &x) const
Norm of a right-hand side vector. The vector must be consistent on the interior+border partition.
Definition: ovlpistlsolverbackend.hh:408
Dune::PDELab::ISTLBackend_OVLP_BCGS_ILU0
Overlapping parallel BiCGStab solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:682
Dune::PDELab::ISTLBackend_OVLP_ILU0_Base::ISTLBackend_OVLP_ILU0_Base
ISTLBackend_OVLP_ILU0_Base(const GFS &gfs_, const C &c_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:538
Dune::PDELab::ISTLBackend_OVLP_CG_UMFPack::ISTLBackend_OVLP_CG_UMFPack
ISTLBackend_OVLP_CG_UMFPack(const GFS &gfs_, const CC &cc_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:996
Dune::PDELab::OVLPScalarProduct::dot
virtual X::Container::field_type dot(const X &x, const X &y) override
Definition: ovlpistlsolverbackend.hh:445
Dune::PDELab::OverlappingScalarProduct::field_type
X::ElementType field_type
Definition: ovlpistlsolverbackend.hh:96
Dune::PDELab::UMFPackSubdomainSolver::pre
virtual void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: ovlpistlsolverbackend.hh:224
Dune::PDELab::UMFPackSubdomainSolver::apply
virtual void apply(X &v, const Y &d) override
Apply the precondioner.
Definition: ovlpistlsolverbackend.hh:229
Dune::PDELab::OverlappingWrappedPreconditioner::OverlappingWrappedPreconditioner
OverlappingWrappedPreconditioner(const GFS &gfs_, P &prec_, const CC &cc_, const ISTL::ParallelHelper< GFS > &helper_)
Constructor.
Definition: ovlpistlsolverbackend.hh:149
Dune::PDELab::ISTLBackend_OVLP_ILUn_Base
Definition: ovlpistlsolverbackend.hh:589
Dune::PDELab::UMFPackSubdomainSolver::range_type
Y range_type
The range type of the preconditioner.
Definition: ovlpistlsolverbackend.hh:207
Dune::PDELab::OVLPScalarProduct
Definition: ovlpistlsolverbackend.hh:432
Dune::PDELab::ISTLBackend_OVLP_BCGS_ILU0::ISTLBackend_OVLP_BCGS_ILU0
ISTLBackend_OVLP_BCGS_ILU0(const GFS &gfs, const CC &cc, unsigned maxiter=5000, int verbose=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:693
Dune::PDELab::ISTLBackend_OVLP_BCGS_SuperLU::ISTLBackend_OVLP_BCGS_SuperLU
ISTLBackend_OVLP_BCGS_SuperLU(const GFS &gfs_, const CC &cc_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:947
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: ovlpistlsolverbackend.hh:1008
Dune::PDELab::ISTLBackend_AMG::setParameters
void setParameters(const Parameters &params_)
set AMG parameters
Definition: ovlpistlsolverbackend.hh:1132
constraints.hh
Dune::PDELab::OverlappingScalarProduct::dot
virtual field_type dot(const X &x, const X &y) override
Dot product of two vectors. It is assumed that the vectors are consistent on the interior+border part...
Definition: ovlpistlsolverbackend.hh:109
seqistlsolverbackend.hh
Dune::PDELab::UMFPackSubdomainSolver::field_type
X::ElementType field_type
The field type of the preconditioner.
Definition: ovlpistlsolverbackend.hh:209
Dune::PDELab::ISTLBackend_OVLP_BCGS_SSORk::ISTLBackend_OVLP_BCGS_SSORk
ISTLBackend_OVLP_BCGS_SSORk(const GFS &gfs, const CC &cc, unsigned maxiter=5000, int steps=5, int verbose=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:671
Dune::PDELab::UMFPackSubdomainSolver
Definition: ovlpistlsolverbackend.hh:199
Dune::PDELab::OverlappingWrappedPreconditioner::range_type
Dune::PDELab::Backend::Vector< GFS, typename P::range_type::field_type > range_type
The range type of the preconditioner.
Definition: ovlpistlsolverbackend.hh:146
Dune::PDELab::OverlappingOperator
Definition: ovlpistlsolverbackend.hh:41
Dune::PDELab::UMFPackSubdomainSolver::category
SolverCategory::Category category() const override
Definition: ovlpistlsolverbackend.hh:241
Dune::PDELab::ISTLBackend_OVLP_ILU0_Base::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:550
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:1047
Dune::PDELab::OVLPScalarProduct::category
SolverCategory::Category category() const override
Definition: ovlpistlsolverbackend.hh:436
Dune::PDELab::ISTLBackend_OVLP_ILU0_Base
Definition: ovlpistlsolverbackend.hh:527
Dune::PDELab::UMFPackSubdomainSolver::domain_type
X domain_type
The domain type of the preconditioner.
Definition: ovlpistlsolverbackend.hh:205
Dune::PDELab::OVLPScalarProductImplementation::dot
X::ElementType dot(const X &x, const X &y) const
Dot product of two vectors. It is assumed that the vectors are consistent on the interior+border part...
Definition: ovlpistlsolverbackend.hh:395
Dune::PDELab::ISTLBackend_OVLP_SuperLU_Base
Definition: ovlpistlsolverbackend.hh:814
Dune::PDELab::ISTLBackend_OVLP_CG_SuperLU::ISTLBackend_OVLP_CG_SuperLU
ISTLBackend_OVLP_CG_SuperLU(const GFS &gfs_, const CC &cc_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:971
Dune::PDELab::ISTLBackend_OVLP_UMFPack_Base::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:896
Dune::PDELab::ISTL::CommSelector::type
Dune::Amg::SequentialInformation type
Definition: parallelhelper.hh:429
Dune::PDELab::ISTLBackend_CG_AMG_SSOR
Overlapping parallel conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: ovlpistlsolverbackend.hh:1256
Dune::PDELab::ISTLAMGStatistics
Class providing some statistics of the AMG solver.
Definition: seqistlsolverbackend.hh:544
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal::ISTLBackend_OVLP_ExplicitDiagonal
ISTLBackend_OVLP_ExplicitDiagonal(const GFS &gfs_)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:1016
Dune::PDELab::ISTLBackend_AMG::apply
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:1179
Dune::PDELab::OverlappingScalarProduct::domain_type
X domain_type
export types
Definition: ovlpistlsolverbackend.hh:95
Dune::PDELab::set_constrained_dofs
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: constraints.hh:798
Dune::PDELab::AddDataHandle
Definition: genericdatahandle.hh:665
s
const std::string s
Definition: function.hh:830
Dune::PDELab::OverlappingOperator::applyscaleadd
virtual void applyscaleadd(field_type alpha, const domain_type &x, range_type &y) const override
apply operator to x, scale and add:
Definition: ovlpistlsolverbackend.hh:64
Dune::PDELab::OverlappingWrappedPreconditioner::domain_type
Dune::PDELab::Backend::Vector< GFS, typename P::domain_type::field_type > domain_type
The domain type of the preconditioner.
Definition: ovlpistlsolverbackend.hh:144
Dune::PDELab::OverlappingWrappedPreconditioner
Definition: ovlpistlsolverbackend.hh:138
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::Backend::Vector
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
Dune::PDELab::ISTLBackend_OVLP_BCGS_SuperLU
Overlapping parallel BiCGStab solver with SuperLU preconditioner.
Definition: ovlpistlsolverbackend.hh:935
Dune::PDELab::ISTLBackend_BCGS_AMG_ILU0::ISTLBackend_BCGS_AMG_ILU0
ISTLBackend_BCGS_AMG_ILU0(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1328
Dune::PDELab::OVLPScalarProductImplementation::parallelHelper
ISTL::ParallelHelper< GFS > & parallelHelper()
Definition: ovlpistlsolverbackend.hh:420
Dune::PDELab::OverlappingWrappedPreconditioner::category
SolverCategory::Category category() const override
Definition: ovlpistlsolverbackend.hh:175
Dune::PDELab::ISTLBackend_AMG::statistics
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: ovlpistlsolverbackend.hh:1228
Dune::PDELab::ISTLBackend_AMG::norm
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: ovlpistlsolverbackend.hh:1165
Dune::PDELab::ISTLBackend_OVLP_BCGS_ILUn
Overlapping parallel BiCGStab solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:703
Dune::PDELab::ISTLBackend_OVLP_SuperLU_Base::apply
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlpistlsolverbackend.hh:838
Dune::PDELab::UMFPackSubdomainSolver::post
virtual void post(X &x) override
Clean up.
Definition: ovlpistlsolverbackend.hh:249
Dune::PDELab::OverlappingWrappedPreconditioner::pre
virtual void pre(domain_type &x, range_type &b) override
Prepare the preconditioner.
Definition: ovlpistlsolverbackend.hh:157
Dune::PDELab::Backend::Native
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