3#ifndef DUNE_PDELAB_BACKEND_ISTL_OVLPISTLSOLVERBACKEND_HH
4#define DUNE_PDELAB_BACKEND_ISTL_OVLPISTLSOLVERBACKEND_HH
10#include <dune/istl/solvercategory.hh>
16#include <dune/istl/paamg/pinfo.hh>
20#include <dune/pdelab/constraints/common/constraints.hh>
21#include <dune/pdelab/gridfunctionspace/genericdatahandle.hh>
22#include <dune/pdelab/backend/istl/vector.hh>
23#include <dune/pdelab/backend/istl/bcrsmatrix.hh>
24#include <dune/pdelab/backend/istl/parallelhelper.hh>
25#include <dune/pdelab/backend/istl/seqistlsolverbackend.hh>
40 template<
class CC,
class M,
class X,
class Y>
41 class OverlappingOperator
46 typedef M matrix_type;
47 typedef X domain_type;
49 typedef typename X::ElementType field_type;
51 OverlappingOperator (
const CC& cc_,
const M& A)
56 virtual void apply (
const domain_type& x, range_type& y)
const override
58 using Backend::native;
59 native(_A_).mv(native(x),native(y));
64 virtual void applyscaleadd (field_type alpha,
const domain_type& x, range_type& y)
const override
66 using Backend::native;
67 native(_A_).usmv(alpha,native(x),native(y));
77 virtual const M& getmat ()
const override
89 template<
class GFS,
class X>
90 class OverlappingScalarProduct
95 typedef X domain_type;
96 typedef typename X::ElementType field_type;
100 OverlappingScalarProduct (
const GFS& gfs_,
const ISTL::ParallelHelper<GFS>& helper_)
101 : gfs(gfs_), helper(helper_)
109 virtual field_type
dot (
const X& x,
const X& y)
const override
112 field_type sum = helper.disjointDot(x,y);
115 return gfs.gridView().comm().sum(sum);
121 virtual double norm (
const X& x)
const override
123 return sqrt(
static_cast<double>(this->
dot(x,x)));
133 const ISTL::ParallelHelper<GFS>& helper;
137 template<
class CC,
class GFS,
class P>
138 class OverlappingWrappedPreconditioner
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>>
144 using domain_type = Dune::PDELab::Backend::Vector<GFS,typename P::domain_type::field_type>;
146 using range_type = Dune::PDELab::Backend::Vector<GFS,typename P::range_type::field_type>;
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_)
157 virtual void pre (domain_type& x, range_type& b)
override
165 virtual void apply (domain_type& v,
const range_type& d)
override
169 prec.apply(Backend::native(v),Backend::native(dd));
170 Dune::PDELab::AddDataHandle<GFS,domain_type> adddh(gfs,v);
171 if (gfs.gridView().comm().size()>1)
183 virtual void post (domain_type& x)
override
185 prec.post(Backend::native(x));
192 const ISTL::ParallelHelper<GFS>& helper;
196#if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
198 template<
class GFS,
class M,
class X,
class Y>
201 typedef Backend::Native<M> ISTLM;
205 typedef X domain_type;
207 typedef Y range_type;
209 typedef typename X::ElementType field_type;
217 UMFPackSubdomainSolver (
const GFS& gfs_,
const M& A_)
218 : gfs(gfs_), solver(Backend::native(A_),false)
224 virtual void pre (X& x, Y& b)
override {}
229 virtual void apply (X& v,
const Y& d)
override
233 solver.apply(Backend::native(v),Backend::native(b),stat);
234 if (gfs.gridView().comm().size()>1)
236 AddDataHandle<GFS,X> adddh(gfs,v);
249 virtual void post (X& x)
override {}
259 template<
class GFS,
class M,
class X,
class Y>
262 typedef Backend::Native<M> ISTLM;
266 typedef X domain_type;
268 typedef Y range_type;
270 typedef typename X::ElementType field_type;
278 SuperLUSubdomainSolver (
const GFS& gfs_,
const M& A_)
279 : gfs(gfs_), solver(Backend::native(A_),false)
285 virtual void pre (X& x, Y& b)
override {}
290 virtual void apply (X& v,
const Y& d)
override
294 solver.apply(Backend::native(v),Backend::native(b),stat);
295 if (gfs.gridView().comm().size()>1)
297 AddDataHandle<GFS,X> adddh(gfs,v);
310 virtual void post (X& x)
override {}
318 template<
class GFS,
class M,
class X,
class Y>
321 typedef typename M::Container ISTLM;
325 typedef X domain_type;
327 typedef Y range_type;
329 typedef typename X::ElementType field_type;
338 RestrictedSuperLUSubdomainSolver (
const GFS& gfs_,
const M& A_,
339 const ISTL::ParallelHelper<GFS>& helper_)
340 : gfs(gfs_), solver(Backend::native(A_),false), helper(helper_)
346 virtual void pre (X& x, Y& b)
override {}
351 virtual void apply (X& v,
const Y& d)
override
353 using Backend::native;
356 solver.apply(native(v),native(b),stat);
357 if (gfs.gridView().comm().size()>1)
359 helper.maskForeignDOFs(native(v));
360 AddDataHandle<GFS,X> adddh(gfs,v);
373 virtual void post (X& x)
override {}
378 const ISTL::ParallelHelper<GFS>& helper;
382 template<
typename GFS>
383 class OVLPScalarProductImplementation
386 OVLPScalarProductImplementation(
const GFS& gfs_)
387 : gfs(gfs_), helper(gfs_)
395 typename X::ElementType
dot (
const X& x,
const X& y)
const
398 typename X::ElementType sum = helper.disjointDot(x,y);
401 return gfs.gridView().comm().sum(sum);
408 typename Dune::template FieldTraits<typename X::ElementType >::real_type norm (
const X& x)
const
411 return sqrt(
static_cast<double>(this->
dot(x,x)));
414 const ISTL::ParallelHelper<GFS>& parallelHelper()
const
420 ISTL::ParallelHelper<GFS>& parallelHelper()
427 ISTL::ParallelHelper<GFS> helper;
431 template<
typename GFS,
typename X>
432 class OVLPScalarProduct
433 :
public ScalarProduct<X>
441 OVLPScalarProduct(
const OVLPScalarProductImplementation<GFS>& implementation_)
442 : implementation(implementation_)
445 virtual typename X::Container::field_type
dot(
const X& x,
const X& y)
const override
447 return implementation.dot(x,y);
450 virtual typename X::Container::field_type norm (
const X& x)
const override
453 return sqrt(
static_cast<double>(this->
dot(x,x)));
457 const OVLPScalarProductImplementation<GFS>& implementation;
460 template<
class GFS,
class C,
461 template<
class,
class,
class,
int>
class Preconditioner,
462 template<
class>
class Solver>
463 class ISTLBackend_OVLP_Base
464 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
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_)
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)
490 using Backend::Native;
491 using Backend::native;
492 typedef OverlappingOperator<C,M,V,W> POP;
494 typedef OVLPScalarProduct<GFS,V> PSP;
496 typedef Preconditioner<
502 SeqPrec seqprec(native(A),steps,1.0);
503 typedef OverlappingWrappedPreconditioner<C,GFS,SeqPrec> WPREC;
504 WPREC wprec(gfs,seqprec,c,this->parallelHelper());
506 if (gfs.gridView().comm().rank()==0) verb=verbose;
507 Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
509 solver.apply(z,r,stat);
525 template<
class GFS,
class C,
526 template<
class>
class Solver>
527 class ISTLBackend_OVLP_ILU0_Base
528 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
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_)
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)
552 using Backend::Native;
553 using Backend::native;
554 typedef OverlappingOperator<C,M,V,W> POP;
556 typedef OVLPScalarProduct<GFS,V> PSP;
564 SeqPrec seqprec(native(A),1.0);
565 typedef OverlappingWrappedPreconditioner<C,GFS,SeqPrec> WPREC;
566 WPREC wprec(gfs,seqprec,c,this->parallelHelper());
568 if (gfs.gridView().comm().rank()==0) verb=verbose;
569 Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
571 solver.apply(z,r,stat);
587 template<
class GFS,
class C,
588 template<
class>
class Solver>
589 class ISTLBackend_OVLP_ILUn_Base
590 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
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_)
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)
615 using Backend::Native;
616 using Backend::native;
617 typedef OverlappingOperator<C,M,V,W> POP;
619 typedef OVLPScalarProduct<GFS,V> PSP;
627 SeqPrec seqprec(native(A),n,1.0);
628 typedef OverlappingWrappedPreconditioner<C,GFS,SeqPrec> WPREC;
629 WPREC wprec(gfs,seqprec,c,this->parallelHelper());
631 if (gfs.gridView().comm().rank()==0) verb=verbose;
632 Solver<V> solver(pop,psp,wprec,reduction,maxiter,verb);
634 solver.apply(z,r,stat);
658 template<
class GFS,
class CC>
660 :
public ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::BiCGSTABSolver>
672 int steps=5,
int verbose=1)
681 template<
class GFS,
class CC>
683 :
public ISTLBackend_OVLP_ILU0_Base<GFS,CC,Dune::BiCGSTABSolver>
702 template<
class GFS,
class CC>
704 :
public ISTLBackend_OVLP_ILUn_Base<GFS,CC,Dune::BiCGSTABSolver>
716 : ISTLBackend_OVLP_ILUn_Base<GFS,CC,
Dune::
BiCGSTABSolver>(gfs, cc, n, maxiter, verbose)
724 template<
class GFS,
class CC>
726 :
public ISTLBackend_OVLP_Base<GFS,CC,Dune::SeqSSOR, Dune::CGSolver>
738 int steps=5,
int verbose=1)
748 template<
class GFS,
class CC>
750 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
762 : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), cc(cc_), maxiter(maxiter_), verbose(verbose_),
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)
775 using Backend::Native;
776 using Backend::native;
777 typedef OverlappingOperator<CC,M,V,W> POP;
779 typedef OVLPScalarProduct<GFS,V> PSP;
787 SeqPrec seqprec(native(A),1.0);
788 typedef OverlappingWrappedPreconditioner<CC,GFS,SeqPrec> WPREC;
789 WPREC wprec(gfs,seqprec,cc,this->parallelHelper());
791 if (gfs.gridView().comm().rank()==0) verb=verbose;
794 solver.
apply(z,r,stat);
813 template<
class GFS,
class C,
template<
typename>
class Solver>
814 class ISTLBackend_OVLP_SuperLU_Base
815 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
825 ISTLBackend_OVLP_SuperLU_Base (
const GFS& gfs_,
const C& c_,
unsigned maxiter_=5000,
827 : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), verbose(verbose_)
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)
840 typedef OverlappingOperator<C,M,V,W> POP;
842 typedef OVLPScalarProduct<GFS,V> PSP;
845 typedef SuperLUSubdomainSolver<GFS,M,V,W> PREC;
848 if (gfs.gridView().comm().rank()==0) verb=verbose;
849 Solver<V> solver(pop,psp,prec,reduction,maxiter,verb);
851 solver.apply(z,r,stat);
858 std::cout <<
"No superLU support, please install and configure it." << std::endl;
871 template<
class GFS,
class C,
template<
typename>
class Solver>
872 class ISTLBackend_OVLP_UMFPack_Base
873 :
public OVLPScalarProductImplementation<GFS>,
public LinearResultStorage
883 ISTLBackend_OVLP_UMFPack_Base (
const GFS& gfs_,
const C& c_,
unsigned maxiter_=5000,
885 : OVLPScalarProductImplementation<GFS>(gfs_), gfs(gfs_), c(c_), maxiter(maxiter_), verbose(verbose_)
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)
898 typedef OverlappingOperator<C,M,V,W> POP;
900 typedef OVLPScalarProduct<GFS,V> PSP;
902#if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
903 typedef UMFPackSubdomainSolver<GFS,M,V,W> PREC;
906 if (gfs.gridView().comm().rank()==0) verb=verbose;
907 Solver<V> solver(pop,psp,prec,reduction,maxiter,verb);
909 solver.apply(z,r,stat);
916 std::cout <<
"No UMFPack support, please install and configure it." << std::endl;
934 template<
class GFS,
class CC>
936 :
public ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::BiCGSTABSolver>
949 : ISTLBackend_OVLP_SuperLU_Base<GFS,CC,
Dune::
BiCGSTABSolver>(gfs_,cc_,maxiter_,verbose_)
958 template<
class GFS,
class CC>
960 :
public ISTLBackend_OVLP_SuperLU_Base<GFS,CC,Dune::CGSolver>
972 unsigned maxiter_=5000,
974 : ISTLBackend_OVLP_SuperLU_Base<GFS,CC,
Dune::
CGSolver>(gfs_,cc_,maxiter_,verbose_)
983 template<
class GFS,
class CC>
985 :
public ISTLBackend_OVLP_UMFPack_Base<GFS,CC,Dune::CGSolver>
997 unsigned maxiter_=5000,
999 : ISTLBackend_OVLP_UMFPack_Base<GFS,CC,
Dune::
CGSolver>(gfs_,cc_,maxiter_,verbose_)
1009 :
public LinearResultStorage
1029 typename V::ElementType
norm(
const V& v)
const
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 "
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)
1049 using Backend::Native;
1050 using Backend::native;
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)
1061 CopyDataHandle<GFS,V> copydh(gfs,z);
1064 res.converged =
true;
1067 res.reduction =
static_cast<double>(reduction);
1068 res.conv_rate =
static_cast<double>(reduction);
1076 template<
class GO,
int s,
template<
class,
class,
class,
int>
class Preconditioner,
1077 template<
class>
class Solver>
1078 class ISTLBackend_AMG :
public LinearResultStorage
1080 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1081 typedef ISTL::ParallelHelper<GFS> PHELPER;
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;
1086 typedef typename ISTL::CommSelector<s,Dune::MPIHelper::isFake>::type Comm;
1098 typedef typename V::ElementType RF;
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_)
1115 params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
1116 params.setDebugLevel(verbose_);
1118 if (gfs.gridView().comm().rank() == 0 && usesuperlu ==
true)
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;
1132 void setParameters(
const Parameters& params_)
1150 void setReuse(
bool reuse_)
1156 bool getReuse()
const
1165 typename V::ElementType norm (
const V& v)
const
1167 typedef OverlappingScalarProduct<GFS,V> PSP;
1168 PSP psp(gfs,phelper);
1179 void apply(M& A, V& z, V& r,
typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
1182 Comm oocc(gfs.gridView().comm());
1183 MatrixType& mat=Backend::native(A);
1187 phelper.createIndexSetAndProjectForAMG(A, oocc);
1188 Operator oop(mat, oocc);
1194 SmootherArgs smootherArgs;
1195 smootherArgs.iterations = 1;
1196 smootherArgs.relaxationFactor = 1;
1197 Criterion criterion(params);
1198 stats.tprepare=watch.elapsed();
1202 if (gfs.gridView().comm().rank()==0) verb=verbose;
1204 if (reuse==
false || firstapply==
true){
1205 amg.reset(
new AMG(oop, criterion, smootherArgs, oocc));
1207 stats.tsetup = watch.elapsed();
1208 stats.levels = amg->maxlevels();
1209 stats.directCoarseLevelSolver=amg->usesDirectCoarseLevelSolver();
1212 Solver<VectorType> solver(oop,sp,*amg,RF(reduction),maxiter,verb);
1215 solver.apply(Backend::native(z),Backend::native(r),stat);
1216 stats.tsolve= watch.
elapsed();
1228 const ISTLAMGStatistics& statistics()
const
1242 std::shared_ptr<AMG> amg;
1243 ISTLAMGStatistics stats;
1255 template<
class GO,
int s=96>
1257 :
public ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::CGSolver>
1259 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1271 int verbose_=1,
bool reuse_=
false,
1272 bool usesuperlu_=
true)
1274 (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
1284 template<
class GO,
int s=96>
1286 :
public ISTLBackend_AMG<GO, s, Dune::SeqSSOR, Dune::BiCGSTABSolver>
1288 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1300 int verbose_=1,
bool reuse_=
false,
1301 bool usesuperlu_=
true)
1303 (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
1313 template<
class GO,
int s=96>
1315 :
public ISTLBackend_AMG<GO, s, Dune::SeqILU, Dune::BiCGSTABSolver>
1317 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
1329 int verbose_=1,
bool reuse_=
false,
1330 bool usesuperlu_=
true)
1332 (gfs_, maxiter_, verbose_, reuse_, usesuperlu_)
Parallel algebraic multigrid based on agglomeration.
Definition: amg.hh:62
The criterion describing the stop criteria for the coarsening process.
Definition: matrixhierarchy.hh:283
Norm that uses only the [0][0] entry of the block to determine couplings.
Definition: aggregates.hh:451
All parameters for AMG.
Definition: parameters.hh:391
Criterion taking advantage of symmetric matrices.
Definition: aggregates.hh:515
A linear operator exporting itself in matrix form.
Definition: operators.hh:107
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition: solvers.hh:417
Block parallel preconditioner.
Definition: schwarz.hh:273
conjugate gradient method
Definition: solvers.hh:189
Adapter to turn a matrix into a linear operator.
Definition: operators.hh:135
An overlapping Schwarz operator.
Definition: schwarz.hh:76
Scalar product for overlapping Schwarz methods.
Definition: scalarproducts.hh:183
Overlapping parallel BiCGStab solver preconditioned with AMG smoothed by ILU0.
Definition: ovlpistlsolverbackend.hh:1316
ISTLBackend_BCGS_AMG_ILU0(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1328
Overlapping parallel BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: ovlpistlsolverbackend.hh:1287
ISTLBackend_BCGS_AMG_SSOR(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1299
Overlapping parallel conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: ovlpistlsolverbackend.hh:1258
ISTLBackend_CG_AMG_SSOR(const GFS &gfs_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: ovlpistlsolverbackend.hh:1270
Overlapping parallel BiCGStab solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:684
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
Overlapping parallel BiCGStab solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:705
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
Overlapping parallel BiCGStab solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:661
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
Overlapping parallel BiCGStab solver with SuperLU preconditioner.
Definition: ovlpistlsolverbackend.hh:937
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
Overlapping parallel CGS solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:727
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
Overlapping parallel CG solver with SuperLU preconditioner.
Definition: ovlpistlsolverbackend.hh:961
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
Overlapping parallel CG solver with UMFPack preconditioner.
Definition: ovlpistlsolverbackend.hh:986
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
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: ovlpistlsolverbackend.hh:1010
ISTLBackend_OVLP_ExplicitDiagonal(const GFS &gfs_)
make a linear solver object
Definition: ovlpistlsolverbackend.hh:1016
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: ovlpistlsolverbackend.hh:1029
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
Overlapping parallel restarted GMRes solver with ILU0 preconditioner.
Definition: ovlpistlsolverbackend.hh:751
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
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
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
virtual void post(Dune::PDELab::Backend::Vector< GFS, P::domain_type::field_type > &x)=0
Clean up.
virtual void pre(Dune::PDELab::Backend::Vector< GFS, P::domain_type::field_type > &x, Dune::PDELab::Backend::Vector< GFS, P::range_type::field_type > &b)=0
Prepare the preconditioner.
implements the Generalized Minimal Residual (GMRes) method
Definition: solvers.hh:812
virtual void apply(X &x, Y &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:895
Base class for scalar product and norm computation.
Definition: scalarproducts.hh:48
Sequential ILU preconditioner.
Definition: preconditioners.hh:500
The sequential jacobian preconditioner.
Definition: preconditioners.hh:390
Sequential SSOR preconditioner.
Definition: preconditioners.hh:138
Default implementation for the scalar case.
Definition: scalarproducts.hh:153
Definition: recipe-operator-splitting.cc:108
Definition of the DUNE_DEPRECATED macro for the case that config.h is not available.
Some generic functions for pretty printing vectors and matrices.
Implementations of the inverse operator interface.
auto dot(const A &a, const B &b) -> typename std::enable_if<!IsVector< A >::value &&!std::is_same< typename FieldTraits< A >::field_type, typename FieldTraits< A >::real_type > ::value, decltype(conj(a) *b)>::type
computes the dot product for fundamental data types according to Petsc's VectDot function: dot(a,...
Definition: dotproduct.hh:40
decltype(auto) apply(F &&f, ArgTuple &&args)
Apply function with arguments given as tuple.
Definition: apply.hh:46
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: constraints.hh:796
Helpers for dealing with MPI.
Dune namespace.
Definition: alignedallocator.hh:14
Define general, extensible interface for operators. The available implementation wraps a matrix.
Classes providing communication interfaces for overlapping Schwarz methods.
Define general preconditioner interface.
Define base class for scalar product and norm.
template which always yields a false value
Definition: typetraits.hh:126
The default class for the smoother arguments.
Definition: smoother.hh:37
Statistics about the application of an inverse operator.
Definition: solver.hh:46
double elapsed
Elapsed time in seconds.
Definition: solver.hh:80
int iterations
Number of iterations.
Definition: solver.hh:65
double reduction
Reduction achieved: .
Definition: solver.hh:68
double conv_rate
Convergence rate (average reduction per step)
Definition: solver.hh:74
bool converged
True if convergence criterion has been met.
Definition: solver.hh:71
Category
Definition: solvercategory.hh:21
@ overlapping
Category for overlapping solvers.
Definition: solvercategory.hh:27
Classes for using SuperLU with ISTL matrices.