dune-pdelab 2.7-git
Loading...
Searching...
No Matches
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_PDELAB_BACKEND_ISTL_NOVLPISTLSOLVERBACKEND_HH
4#define DUNE_PDELAB_BACKEND_ISTL_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
32namespace 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
67
77 NonoverlappingOperator (const GFS& gfs_, const M& A)
78 : gfs(gfs_), _A_(A)
79 { }
80
82
87 virtual void apply (const X& x, Y& y) const override
88 {
89 using Backend::native;
90 // apply local operator; now we have sum y_p = sequential y
91 native(_A_).mv(native(x),native(y));
92
93 // accumulate y on border
95 if (gfs.gridView().comm().size()>1)
96 gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
97 }
98
100
105 virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const override
106 {
107 using Backend::native;
108 // apply local operator; now we have sum y_p = sequential y
109 native(_A_).usmv(alpha,native(x),native(y));
110
111 // accumulate y on border
113 if (gfs.gridView().comm().size()>1)
114 gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
115 }
116
117 SolverCategory::Category category() const override
118 {
119 return SolverCategory::nonoverlapping;
120 }
121
123 virtual const M& getmat () const override
124 {
125 return _A_;
126 }
127
128 private:
129 const GFS& gfs;
130 const M& _A_;
131 };
132
133 // parallel scalar product assuming no overlap
134 template<class GFS, class X>
135 class NonoverlappingScalarProduct : public Dune::ScalarProduct<X>
136 {
137 public:
139 typedef X domain_type;
140 typedef typename X::ElementType field_type;
141
142 SolverCategory::Category category() const override
143 {
144 return SolverCategory::nonoverlapping;
145 }
146
150 : gfs(gfs_), helper(helper_)
151 {}
152
157 virtual field_type dot (const X& x, const X& y) const override
158 {
159 // do local scalar product on unique partition
160 field_type sum = helper.disjointDot(x,y);
161
162 // do global communication
163 return gfs.gridView().comm().sum(sum);
164 }
165
169 virtual double norm (const X& x) const override
170 {
171 return sqrt(static_cast<double>(this->dot(x,x)));
172 }
173
176 void make_consistent (X& x) const
177 {
179 if (gfs.gridView().comm().size()>1)
180 gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
181 }
182
183 private:
184 const GFS& gfs;
185 const ISTL::ParallelHelper<GFS>& helper;
186 };
187
188 // parallel Richardson preconditioner
189 template<class GFS, class X, class Y>
190 class NonoverlappingRichardson : public Dune::Preconditioner<X,Y>
191 {
192 public:
194 typedef X domain_type;
196 typedef Y range_type;
198 typedef typename X::ElementType field_type;
199
200 // define the category
201 SolverCategory::Category category() const override
202 {
203 return SolverCategory::nonoverlapping;
204 }
205
207 NonoverlappingRichardson (const GFS& gfs_, const ISTL::ParallelHelper<GFS>& helper_)
208 : gfs(gfs_), helper(helper_)
209 {
210 }
211
215 virtual void pre (X& x, Y& b) const override {}
216
220 virtual void apply (X& v, const Y& d) const override
221 {
222 v = d;
223 }
224
228 virtual void post (X& x) override {}
229
230 private:
231 const GFS& gfs;
232 const ISTL::ParallelHelper<GFS>& helper;
233 };
234
236
249 template<typename A, typename X, typename Y>
251 : public Dune::Preconditioner<X,Y>
252 {
253
255
256 Diagonal _inverse_diagonal;
257
258 public:
260
264 typedef X domain_type;
266
270 typedef Y range_type;
272 typedef typename X::ElementType field_type;
273
274 SolverCategory::Category category() const override
275 {
276 return SolverCategory::nonoverlapping;
277 }
278
280
291 template<typename GFS>
292 NonoverlappingJacobi(const GFS& gfs, const A &m)
293 : _inverse_diagonal(m)
294 {
295 // make the diagonal consistent...
296 typename ISTL::BlockMatrixDiagonal<A>::template AddMatrixElementVectorDataHandle<GFS> addDH(gfs, _inverse_diagonal);
297 gfs.gridView().communicate(addDH,
298 InteriorBorder_InteriorBorder_Interface,
299 ForwardCommunication);
300
301 // ... and then invert it
302 _inverse_diagonal.invert();
303
304 }
305
307 virtual void pre (X& x, Y& b) override {}
308
310 /*
311 * For this preconditioner, this method works with both consistent and
312 * inconsistent vectors: if d is consistent, v will be consistent, if d
313 * is inconsistent, v will be inconsistent.
314 */
315 virtual void apply (X& v, const Y& d) override
316 {
317 _inverse_diagonal.mv(d,v);
318 }
319
321 virtual void post (X& x) override {}
322 };
323
326
328 template<class GFS>
330 {
332
333 public:
340 explicit ISTLBackend_NOVLP_CG_NOPREC (const GFS& gfs_,
341 unsigned maxiter_=5000,
342 int verbose_=1)
343 : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
344 {}
345
350 template<class V>
351 typename V::ElementType norm (const V& v) const
352 {
353 V x(v); // make a copy because it has to be made consistent
355 PSP psp(gfs,phelper);
356 psp.make_consistent(x);
357 return psp.norm(x);
358 }
359
367 template<class M, class V, class W>
368 void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
369 {
371 POP pop(gfs,A);
373 PSP psp(gfs,phelper);
375 PRICH prich(gfs,phelper);
376 int verb=0;
377 if (gfs.gridView().comm().rank()==0) verb=verbose;
378 Dune::CGSolver<V> solver(pop,psp,prich,reduction,maxiter,verb);
379 Dune::InverseOperatorResult stat;
380 solver.apply(z,r,stat);
381 res.converged = stat.converged;
382 res.iterations = stat.iterations;
383 res.elapsed = stat.elapsed;
384 res.reduction = stat.reduction;
385 res.conv_rate = stat.conv_rate;
386 }
387
390 {
391 return res;
392 }
393
394 private:
395 const GFS& gfs;
396 PHELPER phelper;
398 unsigned maxiter;
399 int verbose;
400 };
401
403 template<class GFS>
405 {
407
408 const GFS& gfs;
409 PHELPER phelper;
411 unsigned maxiter;
412 int verbose;
413
414 public:
416
421 explicit ISTLBackend_NOVLP_CG_Jacobi(const GFS& gfs_,
422 unsigned maxiter_ = 5000,
423 int verbose_ = 1) :
424 gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
425 {}
426
428
433 template<class V>
434 typename V::ElementType norm (const V& v) const
435 {
436 V x(v); // make a copy because it has to be made consistent
438 PSP psp(gfs,phelper);
439 psp.make_consistent(x);
440 return psp.norm(x);
441 }
442
444
456 template<class M, class V, class W>
457 void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
458 {
460 POP pop(gfs,A);
462 PSP psp(gfs,phelper);
463
464 typedef NonoverlappingJacobi<M,V,W> PPre;
465 PPre ppre(gfs,Backend::native(A));
466
467 int verb=0;
468 if (gfs.gridView().comm().rank()==0) verb=verbose;
469 CGSolver<V> solver(pop,psp,ppre,reduction,maxiter,verb);
470 InverseOperatorResult stat;
471 solver.apply(z,r,stat);
472 res.converged = stat.converged;
473 res.iterations = stat.iterations;
474 res.elapsed = stat.elapsed;
475 res.reduction = stat.reduction;
476 res.conv_rate = stat.conv_rate;
477 }
478
481 { return res; }
482 };
483
485 template<class GFS>
487 {
489
490 public:
497 explicit ISTLBackend_NOVLP_BCGS_NOPREC (const GFS& gfs_, unsigned maxiter_=5000, int verbose_=1)
498 : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
499 {}
500
505 template<class V>
506 typename V::ElementType norm (const V& v) const
507 {
508 V x(v); // make a copy because it has to be made consistent
510 PSP psp(gfs,phelper);
511 psp.make_consistent(x);
512 return psp.norm(x);
513 }
514
522 template<class M, class V, class W>
523 void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
524 {
526 POP pop(gfs,A);
528 PSP psp(gfs,phelper);
530 PRICH prich(gfs,phelper);
531 int verb=0;
532 if (gfs.gridView().comm().rank()==0) verb=verbose;
533 Dune::BiCGSTABSolver<V> solver(pop,psp,prich,reduction,maxiter,verb);
534 Dune::InverseOperatorResult stat;
535 solver.apply(z,r,stat);
536 res.converged = stat.converged;
537 res.iterations = stat.iterations;
538 res.elapsed = stat.elapsed;
539 res.reduction = stat.reduction;
540 res.conv_rate = stat.conv_rate;
541 }
542
545 {
546 return res;
547 }
548
549 private:
550 const GFS& gfs;
551 PHELPER phelper;
553 unsigned maxiter;
554 int verbose;
555 };
556
557
559 template<class GFS>
561 {
563
564 public:
571 explicit ISTLBackend_NOVLP_BCGS_Jacobi (const GFS& gfs_, unsigned maxiter_=5000, int verbose_=1)
572 : gfs(gfs_), phelper(gfs,verbose_), maxiter(maxiter_), verbose(verbose_)
573 {}
574
579 template<class V>
580 typename V::ElementType norm (const V& v) const
581 {
582 V x(v); // make a copy because it has to be made consistent
584 PSP psp(gfs,phelper);
585 psp.make_consistent(x);
586 return psp.norm(x);
587 }
588
596 template<class M, class V, class W>
597 void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
598 {
600 POP pop(gfs,A);
602 PSP psp(gfs,phelper);
603
604 typedef NonoverlappingJacobi<M,V,W> PPre;
605 PPre ppre(gfs,A);
606
607 int verb=0;
608 if (gfs.gridView().comm().rank()==0) verb=verbose;
609 Dune::BiCGSTABSolver<V> solver(pop,psp,ppre,reduction,maxiter,verb);
610 Dune::InverseOperatorResult stat;
611 solver.apply(z,r,stat);
612 res.converged = stat.converged;
613 res.iterations = stat.iterations;
614 res.elapsed = stat.elapsed;
615 res.reduction = stat.reduction;
616 res.conv_rate = stat.conv_rate;
617 }
618
621 {
622 return res;
623 }
624
625 private:
626 const GFS& gfs;
627 PHELPER phelper;
629 unsigned maxiter;
630 int verbose;
631 };
632
634 template<typename GFS>
636 {
638
639 const GFS& gfs;
640 PHELPER phelper;
642
643 public:
649 explicit ISTLBackend_NOVLP_ExplicitDiagonal(const GFS& gfs_)
650 : gfs(gfs_), phelper(gfs)
651 {}
652
657 template<class V>
658 typename V::ElementType norm (const V& v) const
659 {
661 V x(v); // make a copy because it has to be made consistent
662 PSP psp(gfs,phelper);
663 psp.make_consistent(x);
664 return psp.norm(x);
665 }
666
674 template<class M, class V, class W>
675 void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
676 {
677 Dune::SeqJac<M,V,W> jac(A,1,1.0);
678 jac.pre(z,r);
679 jac.apply(z,r);
680 jac.post(z);
681 if (gfs.gridView().comm().size()>1)
682 {
684 gfs.gridView().communicate(adddh,Dune::InteriorBorder_InteriorBorder_Interface,Dune::ForwardCommunication);
685 }
686 res.converged = true;
687 res.iterations = 1;
688 res.elapsed = 0.0;
689 res.reduction = reduction;
690 res.conv_rate = reduction; // pow(reduction,1.0/1)
691 }
692
695 {
696 return res;
697 }
698 };
700
701
709 template<class GO,
710 template<class,class,class,int> class Preconditioner,
711 template<class> class Solver>
713 {
714 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
716
717 public:
725 explicit ISTLBackend_NOVLP_BASE_PREC (const GO& grid_operator, unsigned maxiter_ = 5000, unsigned steps_ = 5, int verbose_ = 1)
726 : _grid_operator(grid_operator)
727 , gfs(grid_operator.trialGridFunctionSpace())
728 , phelper(gfs,verbose_)
729 , maxiter(maxiter_)
730 , steps(steps_)
731 , verbose(verbose_)
732 {}
733
738 template<class Vector>
739 typename Vector::ElementType norm (const Vector& v) const
740 {
741 Vector x(v); // make a copy because it has to be made consistent
743 PSP psp(gfs,phelper);
744 psp.make_consistent(x);
745 return psp.norm(x);
746 }
747
755 template<class M, class V, class W>
756 void apply(M& A, V& z, W& r, typename V::ElementType reduction)
757 {
758 using MatrixType = Backend::Native<M>;
759 MatrixType& mat = Backend::native(A);
760 using VectorType = Backend::Native<W>;
761#if HAVE_MPI
763 _grid_operator.make_consistent(A);
764 ISTL::assertParallelUG(gfs.gridView().comm());
765 Comm oocc(gfs.gridView().comm(),Dune::SolverCategory::nonoverlapping);
766 phelper.createIndexSetAndProjectForAMG(mat, oocc);
767 typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
768 Smoother smoother(mat, steps, 1.0);
769 typedef Dune::NonoverlappingSchwarzScalarProduct<VectorType,Comm> PSP;
770 PSP psp(oocc);
771 typedef Dune::NonoverlappingSchwarzOperator<MatrixType,VectorType,VectorType,Comm> Operator;
772 Operator oop(mat,oocc);
773 typedef Dune::NonoverlappingBlockPreconditioner<Comm, Smoother> ParSmoother;
774 ParSmoother parsmoother(smoother, oocc);
775#else
776 typedef Preconditioner<MatrixType,VectorType,VectorType,1> ParSmoother;
777 ParSmoother parsmoother(mat, steps, 1.0);
778 typedef Dune::SeqScalarProduct<VectorType> PSP;
779 PSP psp;
780 typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
781 Operator oop(mat);
782#endif
783 int verb=0;
784 if (gfs.gridView().comm().rank()==0) verb=verbose;
785 Solver<VectorType> solver(oop,psp,parsmoother,reduction,maxiter,verb);
786 Dune::InverseOperatorResult stat;
787 //make r consistent
788 if (gfs.gridView().comm().size()>1){
790 gfs.gridView().communicate(adddh,
791 Dune::InteriorBorder_InteriorBorder_Interface,
792 Dune::ForwardCommunication);
793 }
794
795 solver.apply(z,r,stat);
796 res.converged = stat.converged;
797 res.iterations = stat.iterations;
798 res.elapsed = stat.elapsed;
799 res.reduction = stat.reduction;
800 res.conv_rate = stat.conv_rate;
801 }
802
805 {
806 return res;
807 }
808
809 private:
810 const GO& _grid_operator;
811 const GFS& gfs;
812 PHELPER phelper;
814 unsigned maxiter;
815 unsigned steps;
816 int verbose;
817 };
818
821
834 template<class GO>
836 : public ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::BiCGSTABSolver>
837 {
838
839 public:
847 explicit ISTLBackend_NOVLP_BCGS_SSORk (const GO& grid_operator, unsigned maxiter_=5000,
848 int steps_=5, int verbose_=1)
849 : ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::BiCGSTABSolver>(grid_operator, maxiter_, steps_, verbose_)
850 {}
851 };
852
859 template<class GO>
861 : public ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::CGSolver>
862 {
863
864 public:
872 explicit ISTLBackend_NOVLP_CG_SSORk (const GO& grid_operator, unsigned maxiter_=5000,
873 int steps_=5, int verbose_=1)
874 : ISTLBackend_NOVLP_BASE_PREC<GO,Dune::SeqSSOR, Dune::CGSolver>(grid_operator, maxiter_, steps_, verbose_)
875 {}
876 };
879
880 template<class GO,int s, template<class,class,class,int> class Preconditioner,
881 template<class> class Solver>
883 {
884 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
885 typedef typename ISTL::ParallelHelper<GFS> PHELPER;
886 typedef typename GO::Traits::Jacobian M;
887 typedef Backend::Native<M> MatrixType;
888 typedef typename GO::Traits::Domain V;
889 typedef Backend::Native<V> VectorType;
891#if HAVE_MPI
892 typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
893 typedef Dune::NonoverlappingBlockPreconditioner<Comm,Smoother> ParSmoother;
894 typedef Dune::NonoverlappingSchwarzOperator<MatrixType,VectorType,VectorType,Comm> Operator;
895#else
896 typedef Preconditioner<MatrixType,VectorType,VectorType,1> ParSmoother;
897 typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
898#endif
899 typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
900 typedef Dune::Amg::AMG<Operator,VectorType,ParSmoother,Comm> AMG;
901 typedef Dune::Amg::Parameters Parameters;
902
903 public:
904 ISTLBackend_AMG_NOVLP(const GO& grid_operator, unsigned maxiter_=5000,
905 int verbose_=1, bool reuse_=false,
906 bool usesuperlu_=true)
907 : _grid_operator(grid_operator)
908 , gfs(grid_operator.trialGridFunctionSpace())
909 , phelper(gfs,verbose_)
910 , maxiter(maxiter_)
911 , params(15,2000,1.2,1.6,Dune::Amg::atOnceAccu)
912 , verbose(verbose_)
913 , reuse(reuse_)
914 , firstapply(true)
915 , usesuperlu(usesuperlu_)
916 {
917 params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
918 params.setDebugLevel(verbose_);
919#if !HAVE_SUPERLU
920 if (phelper.rank() == 0 && usesuperlu == true)
921 {
922 std::cout << "WARNING: You are using AMG without SuperLU!"
923 << " Please consider installing SuperLU,"
924 << " or set the usesuperlu flag to false"
925 << " to suppress this warning." << std::endl;
926 }
927#endif
928 }
929
934 void setParameters(const Parameters& params_)
935 {
936 params = params_;
937 }
938
946 const Parameters& parameters() const
947 {
948 return params;
949 }
950
952 void setReuse(bool reuse_)
953 {
954 reuse = reuse_;
955 }
956
958 bool getReuse() const
959 {
960 return reuse;
961 }
962
963
968 typename V::ElementType norm (const V& v) const
969 {
970 V x(v); // make a copy because it has to be made consistent
972 PSP psp(gfs,phelper);
973 psp.make_consistent(x);
974 return psp.norm(x);
975 }
976
977 void apply(M& A, V& z, V& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
978 {
979 Timer watch;
980 MatrixType& mat = Backend::native(A);
981 typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
982 Dune::Amg::FirstDiagonal> > Criterion;
983#if HAVE_MPI
984 Comm oocc(gfs.gridView().comm(),Dune::SolverCategory::nonoverlapping);
985 _grid_operator.make_consistent(A);
986 phelper.createIndexSetAndProjectForAMG(A, oocc);
987 Dune::NonoverlappingSchwarzScalarProduct<VectorType,Comm> sp(oocc);
988 Operator oop(mat, oocc);
989#else
990 Comm oocc(gfs.gridView().comm());
991 Operator oop(mat);
992 Dune::SeqScalarProduct<VectorType> sp;
993#endif
994 SmootherArgs smootherArgs;
995 smootherArgs.iterations = 1;
996 smootherArgs.relaxationFactor = 1;
997 //use noAccu or atOnceAccu
998 Criterion criterion(params);
999 stats.tprepare=watch.elapsed();
1000 watch.reset();
1001
1002 int verb=0;
1003 if (gfs.gridView().comm().rank()==0) verb=verbose;
1004 //only construct a new AMG if the matrix changes
1005 if (reuse==false || firstapply==true){
1006 amg.reset(new AMG(oop, criterion, smootherArgs, oocc));
1007 firstapply = false;
1008 stats.tsetup = watch.elapsed();
1009 stats.levels = amg->maxlevels();
1010 stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
1011 }
1012
1013 Dune::InverseOperatorResult stat;
1014 // make r consistent
1015 if (gfs.gridView().comm().size()>1) {
1017 gfs.gridView().communicate(adddh,
1018 Dune::InteriorBorder_InteriorBorder_Interface,
1019 Dune::ForwardCommunication);
1020 }
1021 watch.reset();
1022 Solver<VectorType> solver(oop,sp,*amg,reduction,maxiter,verb);
1023 solver.apply(Backend::native(z),Backend::native(r),stat);
1024 stats.tsolve= watch.elapsed();
1025 res.converged = stat.converged;
1026 res.iterations = stat.iterations;
1027 res.elapsed = stat.elapsed;
1028 res.reduction = stat.reduction;
1029 res.conv_rate = stat.conv_rate;
1030 }
1031
1037 {
1038 return stats;
1039 }
1040
1041 private:
1042 const GO& _grid_operator;
1043 const GFS& gfs;
1044 PHELPER phelper;
1045 unsigned maxiter;
1046 Parameters params;
1047 int verbose;
1048 bool reuse;
1049 bool firstapply;
1050 bool usesuperlu;
1051 std::shared_ptr<AMG> amg;
1052 ISTLAMGStatistics stats;
1053 };
1054
1057
1072 template<class GO, int s=96>
1074 : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::CGSolver>
1075 {
1076
1077 public:
1078 ISTLBackend_NOVLP_CG_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1079 int verbose_=1, bool reuse_=false,
1080 bool usesuperlu_=true)
1081 : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::CGSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1082 {}
1083 };
1084
1099 template<class GO, int s=96>
1101 : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>
1102 {
1103
1104 public:
1105 ISTLBackend_NOVLP_BCGS_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1106 int verbose_=1, bool reuse_=false,
1107 bool usesuperlu_=true)
1108 : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1109 {}
1110 };
1111
1126 template<class GO, int s=96>
1128 : public ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::LoopSolver>
1129 {
1130
1131 public:
1132 ISTLBackend_NOVLP_LS_AMG_SSOR(const GO& grid_operator, unsigned maxiter_=5000,
1133 int verbose_=1, bool reuse_=false,
1134 bool usesuperlu_=true)
1135 : ISTLBackend_AMG_NOVLP<GO, s, Dune::SeqSSOR, Dune::LoopSolver>(grid_operator, maxiter_,verbose_,reuse_,usesuperlu_)
1136 {}
1137 };
1140
1141 } // namespace PDELab
1142} // namespace Dune
1143
1144#endif // DUNE_PDELAB_BACKEND_ISTL_NOVLPISTLSOLVERBACKEND_HH
const std::string s
Definition: function.hh:843
void assertParallelUG(T comm)
Definition: parallelhelper.hh:445
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:192
typename native_type< T >::type Native
Alias of the native container type associated with T or T itself if it is not a backend wrapper.
Definition: backend/interface.hh:176
Definition: blockmatrixdiagonal.hh:215
void invert()
Definition: blockmatrixdiagonal.hh:233
void mv(const X &x, Y &y) const
Definition: blockmatrixdiagonal.hh:239
Operator for the non-overlapping parallel case.
Definition: novlpistlsolverbackend.hh:55
NonoverlappingOperator(const GFS &gfs_, const M &A)
Construct a non-overlapping operator.
Definition: novlpistlsolverbackend.hh:77
virtual const M & getmat() const override
extract the matrix
Definition: novlpistlsolverbackend.hh:123
virtual void apply(const X &x, Y &y) const override
apply operator
Definition: novlpistlsolverbackend.hh:87
Backend::Native< X > domain_type
export type of vectors the matrix is applied to
Definition: novlpistlsolverbackend.hh:60
X::field_type field_type
export type of the entries for x
Definition: novlpistlsolverbackend.hh:64
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const override
apply operator to x, scale and add:
Definition: novlpistlsolverbackend.hh:105
Backend::Native< Y > range_type
export type of result vectors
Definition: novlpistlsolverbackend.hh:62
SolverCategory::Category category() const override
Definition: novlpistlsolverbackend.hh:117
Backend::Native< M > matrix_type
export type of matrix
Definition: novlpistlsolverbackend.hh:58
Definition: novlpistlsolverbackend.hh:136
X::ElementType field_type
Definition: novlpistlsolverbackend.hh:140
SolverCategory::Category category() const override
Definition: novlpistlsolverbackend.hh:142
X domain_type
export types
Definition: novlpistlsolverbackend.hh:139
NonoverlappingScalarProduct(const GFS &gfs_, const ISTL::ParallelHelper< GFS > &helper_)
Constructor needs to know the grid function space.
Definition: novlpistlsolverbackend.hh:149
virtual double norm(const X &x) const override
Norm of a right-hand side vector. The vector must be consistent on the interior+border partition.
Definition: novlpistlsolverbackend.hh:169
void make_consistent(X &x) const
make additive vector consistent
Definition: novlpistlsolverbackend.hh:176
virtual field_type dot(const X &x, const X &y) const override
Dot product of two vectors. It is assumed that the vectors are consistent on the interior+border part...
Definition: novlpistlsolverbackend.hh:157
Definition: novlpistlsolverbackend.hh:191
X::ElementType field_type
The field type of the preconditioner.
Definition: novlpistlsolverbackend.hh:198
NonoverlappingRichardson(const GFS &gfs_, const ISTL::ParallelHelper< GFS > &helper_)
Constructor.
Definition: novlpistlsolverbackend.hh:207
virtual void apply(X &v, const Y &d) const override
Apply the precondioner.
Definition: novlpistlsolverbackend.hh:220
X domain_type
The domain type of the preconditioner.
Definition: novlpistlsolverbackend.hh:194
Y range_type
The range type of the preconditioner.
Definition: novlpistlsolverbackend.hh:196
virtual void post(X &x) override
Clean up.
Definition: novlpistlsolverbackend.hh:228
virtual void pre(X &x, Y &b) const override
Prepare the preconditioner.
Definition: novlpistlsolverbackend.hh:215
SolverCategory::Category category() const override
Definition: novlpistlsolverbackend.hh:201
parallel non-overlapping Jacobi preconditioner
Definition: novlpistlsolverbackend.hh:252
virtual void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: novlpistlsolverbackend.hh:307
virtual void post(X &x) override
Clean up.
Definition: novlpistlsolverbackend.hh:321
virtual void apply(X &v, const Y &d) override
Apply the precondioner.
Definition: novlpistlsolverbackend.hh:315
SolverCategory::Category category() const override
Definition: novlpistlsolverbackend.hh:274
X::ElementType field_type
The field type of the preconditioner.
Definition: novlpistlsolverbackend.hh:272
X domain_type
The domain type of the operator.
Definition: novlpistlsolverbackend.hh:264
NonoverlappingJacobi(const GFS &gfs, const A &m)
Constructor.
Definition: novlpistlsolverbackend.hh:292
Y range_type
The range type of the operator.
Definition: novlpistlsolverbackend.hh:270
Nonoverlapping parallel CG solver without preconditioner.
Definition: novlpistlsolverbackend.hh:330
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: novlpistlsolverbackend.hh:368
ISTLBackend_NOVLP_CG_NOPREC(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:340
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:389
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:351
Nonoverlapping parallel CG solver with Jacobi preconditioner.
Definition: novlpistlsolverbackend.hh:405
const LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:480
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:434
ISTLBackend_NOVLP_CG_Jacobi(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:421
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: novlpistlsolverbackend.hh:457
Nonoverlapping parallel BiCGStab solver without preconditioner.
Definition: novlpistlsolverbackend.hh:487
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:544
ISTLBackend_NOVLP_BCGS_NOPREC(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:497
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: novlpistlsolverbackend.hh:523
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:506
Nonoverlapping parallel BiCGStab solver with Jacobi preconditioner.
Definition: novlpistlsolverbackend.hh:561
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:620
ISTLBackend_NOVLP_BCGS_Jacobi(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:571
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:580
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: novlpistlsolverbackend.hh:597
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: novlpistlsolverbackend.hh:636
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:658
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:694
ISTLBackend_NOVLP_ExplicitDiagonal(const GFS &gfs_)
make a linear solver object
Definition: novlpistlsolverbackend.hh:649
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: novlpistlsolverbackend.hh:675
Utility base class for preconditioned novlp backends.
Definition: novlpistlsolverbackend.hh:713
ISTLBackend_NOVLP_BASE_PREC(const GO &grid_operator, unsigned maxiter_=5000, unsigned steps_=5, int verbose_=1)
Constructor.
Definition: novlpistlsolverbackend.hh:725
const Dune::PDELab::LinearSolverResult< double > & result() const
Return access to result data.
Definition: novlpistlsolverbackend.hh:804
Vector::ElementType norm(const Vector &v) const
Compute global norm of a vector.
Definition: novlpistlsolverbackend.hh:739
void apply(M &A, V &z, W &r, typename V::ElementType reduction)
Solve the given linear system.
Definition: novlpistlsolverbackend.hh:756
Nonoverlapping parallel BiCGSTAB solver preconditioned by block SSOR.
Definition: novlpistlsolverbackend.hh:837
ISTLBackend_NOVLP_BCGS_SSORk(const GO &grid_operator, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:847
Nonoverlapping parallel CG solver preconditioned by block SSOR.
Definition: novlpistlsolverbackend.hh:862
ISTLBackend_NOVLP_CG_SSORk(const GO &grid_operator, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
make a linear solver object
Definition: novlpistlsolverbackend.hh:872
Definition: novlpistlsolverbackend.hh:883
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
Definition: novlpistlsolverbackend.hh:977
const Parameters & parameters() const
Get the parameters describing the behaviuour of AMG.
Definition: novlpistlsolverbackend.hh:946
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: novlpistlsolverbackend.hh:1036
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: novlpistlsolverbackend.hh:968
void setParameters(const Parameters &params_)
set AMG parameters
Definition: novlpistlsolverbackend.hh:934
ISTLBackend_AMG_NOVLP(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: novlpistlsolverbackend.hh:904
bool getReuse() const
Return whether the AMG is reused during call to apply()
Definition: novlpistlsolverbackend.hh:958
void setReuse(bool reuse_)
Set whether the AMG should be reused again during call to apply().
Definition: novlpistlsolverbackend.hh:952
Nonoverlapping parallel CG solver preconditioned with AMG smoothed by SSOR.
Definition: novlpistlsolverbackend.hh:1075
ISTLBackend_NOVLP_CG_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: novlpistlsolverbackend.hh:1078
Nonoverlapping parallel BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: novlpistlsolverbackend.hh:1102
ISTLBackend_NOVLP_BCGS_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: novlpistlsolverbackend.hh:1105
Nonoverlapping parallel LoopSolver preconditioned with AMG smoothed by SSOR.
Definition: novlpistlsolverbackend.hh:1129
ISTLBackend_NOVLP_LS_AMG_SSOR(const GO &grid_operator, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: novlpistlsolverbackend.hh:1132
Definition: parallelhelper.hh:51
RankIndex rank() const
Returns the MPI rank of this process.
Definition: parallelhelper.hh:248
void createIndexSetAndProjectForAMG(MatrixType &m, Comm &c)
Makes the matrix consistent and creates the parallel information for AMG.
Dune::Amg::SequentialInformation type
Definition: parallelhelper.hh:429
Class providing some statistics of the AMG solver.
Definition: seqistlsolverbackend.hh:701
double tprepare
The needed for computing the parallel information and for adapting the linear system.
Definition: seqistlsolverbackend.hh:706
double tsetup
The time needed for building the AMG hierarchy (coarsening).
Definition: seqistlsolverbackend.hh:712
double tsolve
The time spent in solving the system (without building the hierarchy.
Definition: seqistlsolverbackend.hh:710
bool directCoarseLevelSolver
True if a direct solver was used on the coarset level.
Definition: seqistlsolverbackend.hh:716
int levels
the number of levels in the AMG hierarchy.
Definition: seqistlsolverbackend.hh:708
Definition: solver.hh:31
RFType conv_rate
Definition: solver.hh:36
bool converged
Definition: solver.hh:32
RFType reduction
Definition: solver.hh:35
unsigned int iterations
Definition: solver.hh:33
double elapsed
Definition: solver.hh:34
Definition: solver.hh:54
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:63
Definition: genericdatahandle.hh:667
Definition: recipe-operator-splitting.cc:108