Go to the documentation of this file.
3 #ifndef DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
4 #define DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
6 #include <dune/common/deprecated.hh>
7 #include <dune/common/parallel/mpihelper.hh>
9 #include <dune/istl/owneroverlapcopy.hh>
10 #include <dune/istl/solvercategory.hh>
11 #include <dune/istl/operators.hh>
12 #include <dune/istl/solvers.hh>
13 #include <dune/istl/preconditioners.hh>
14 #include <dune/istl/scalarproducts.hh>
15 #include <dune/istl/paamg/amg.hh>
16 #include <dune/istl/paamg/pinfo.hh>
17 #include <dune/istl/io.hh>
18 #include <dune/istl/superlu.hh>
19 #include <dune/istl/umfpack.hh>
44 template<
typename X,
typename Y,
typename GO>
51 static constexpr
bool isLinear = GO::LocalAssembler::isLinear();
66 virtual void apply (
const X& x, Y& y)
const override
70 go.jacobian_apply(x,y);
73 DUNE_THROW(Dune::InvalidStateException,
"You seem to apply a nonlinear operator without setting the linearization point first!");
74 go.jacobian_apply(*u_, x, y);
83 go.jacobian_apply(x,temp);
86 DUNE_THROW(Dune::InvalidStateException,
"You seem to apply a nonlinear operator without setting the linearization point first!");
87 go.jacobian_apply(*u_, x, temp);
92 SolverCategory::Category
category()
const override
94 return SolverCategory::sequential;
111 template<
template<
class>
class Solver>
122 : maxiter(maxiter_), verbose(verbose_)
132 template<
class M,
class V,
class W>
133 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
138 Dune::MatrixAdapter<Native<M>,
140 Native<W>> opa(
native(A));
141 Dune::Richardson<Native<V>,Native<W> > prec(0.7);
142 Solver<Native<V> > solver(opa, prec, reduction, maxiter, verbose);
143 Dune::InverseOperatorResult stat;
157 template<
class GO,
template<
class>
class Solver>
161 using V =
typename GO::Traits::Domain;
162 using W =
typename GO::Traits::Range;
182 void apply(V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
184 Dune::Richardson<V,W> prec(1.0);
185 Solver<V> solver(opa_, prec, reduction, maxiter_, verbose_);
186 Dune::InverseOperatorResult stat;
187 solver.apply(z, r, stat);
207 template<
template<
class,
class,
class,
int>
class Preconditioner,
208 template<
class>
class Solver>
219 : maxiter(maxiter_), verbose(verbose_), preconditioner_steps(preconditioner_steps_)
231 template<
class M,
class V,
class W>
232 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
237 Dune::MatrixAdapter<Native<M>,
239 Native<W>> opa(
native(A));
240 Preconditioner<Native<M>,
243 1> prec(
native(A), preconditioner_steps, 1.0);
244 Solver<Native<V>> solver(opa, prec, reduction, maxiter, verbose);
245 Dune::InverseOperatorResult stat;
257 unsigned preconditioner_steps;
260 template<
template<
typename>
class Solver>
271 : maxiter(maxiter_), verbose(verbose_)
280 template<
class M,
class V,
class W>
281 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
285 Dune::MatrixAdapter<Native<M>,
287 Native<W>> opa(
native(A));
288 Dune::SeqILU<Native<M>,
292 Solver<Native<V>> solver(opa, ilu0, reduction, maxiter, verbose);
293 Dune::InverseOperatorResult stat;
306 template<
template<
typename>
class Solver>
318 : n_(n), w_(
w), maxiter(maxiter_), verbose(verbose_)
327 template<
class M,
class V,
class W>
328 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
332 Dune::MatrixAdapter<Native<M>,
336 Dune::SeqILU<Native<M>,
339 > ilun(
native(A), n_, w_,
false);
340 Solver<Native<V>> solver(opa, ilun, reduction, maxiter, verbose);
341 Dune::InverseOperatorResult stat;
551 #if HAVE_SUPERLU || DOXYGEN
584 template<
class M,
class V,
class W>
585 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
589 using ISTLM = Native<M>;
590 Dune::SuperLU<ISTLM> solver(
native(A), verbose);
591 Dune::InverseOperatorResult stat;
603 #endif // HAVE_SUPERLU || DOXYGEN
605 #if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
638 template<
class M,
class V,
class W>
639 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
643 Dune::UMFPack<ISTLM> solver(
native(A), verbose);
644 Dune::InverseOperatorResult stat;
656 #endif // HAVE_SUITESPARSE_UMFPACK || DOXYGEN
675 template<
class M,
class V,
class W>
676 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType >::real_type reduction)
679 Dune::SeqJac<Native<M>,
719 template<
class GO,
template<
class,
class,
class,
int>
class Preconditioner,
template<
class>
class Solver,
720 bool skipBlocksizeCheck =
false>
723 typedef typename GO::Traits::TrialGridFunctionSpace GFS;
724 typedef typename GO::Traits::Jacobian M;
726 typedef typename GO::Traits::Domain V;
728 typedef Preconditioner<MatrixType,VectorType,VectorType,1> Smoother;
729 typedef Dune::MatrixAdapter<MatrixType,VectorType,VectorType> Operator;
730 typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
731 typedef Dune::Amg::AMG<Operator,VectorType,Smoother> AMG;
736 bool reuse_=
false,
bool usesuperlu_=
true)
737 : maxiter(maxiter_), params(15,2000), verbose(verbose_),
738 reuse(reuse_), firstapply(true), usesuperlu(usesuperlu_)
740 params.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
741 params.setDebugLevel(verbose_);
743 if (usesuperlu ==
true)
745 std::cout <<
"WARNING: You are using AMG without SuperLU!"
746 <<
" Please consider installing SuperLU,"
747 <<
" or set the usesuperlu flag to false"
748 <<
" to suppress this warning." << std::endl;
778 typename V::ElementType
norm (
const V& v)
const
790 void apply(M& A, V& z, V& r,
typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
794 typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<MatrixType,
795 Dune::Amg::FirstDiagonal> > Criterion;
796 SmootherArgs smootherArgs;
797 smootherArgs.iterations = 1;
798 smootherArgs.relaxationFactor = 1;
800 Criterion criterion(params);
802 if (reuse==
false || firstapply==
true){
803 oop.reset(
new Operator(mat));
804 amg.reset(
new AMG(*oop, criterion, smootherArgs));
806 stats.
tsetup = watch.elapsed();
807 stats.
levels = amg->maxlevels();
811 Dune::InverseOperatorResult stat;
813 Solver<VectorType> solver(*oop,*amg,reduction,maxiter,verbose);
815 stats.
tsolve= watch.elapsed();
840 std::shared_ptr<Operator> oop;
841 std::shared_ptr<AMG> amg;
868 bool reuse_=
false,
bool usesuperlu_=
true)
870 (maxiter_, verbose_, reuse_, usesuperlu_)
894 bool reuse_=
false,
bool usesuperlu_=
true)
896 (maxiter_, verbose_, reuse_, usesuperlu_)
920 bool reuse_=
false,
bool usesuperlu_=
true)
922 (maxiter_, verbose_, reuse_, usesuperlu_)
946 bool reuse_=
false,
bool usesuperlu_=
true)
948 (maxiter_, verbose_, reuse_, usesuperlu_)
972 bool reuse_=
false,
bool usesuperlu_=
true)
974 (maxiter_, verbose_, reuse_, usesuperlu_)
995 : restart(restart_), maxiter(maxiter_), verbose(verbose_)
1005 template<
class M,
class V,
class W>
1006 void apply(M& A, V& z, W& r,
typename Dune::template FieldTraits<typename W::ElementType>::real_type reduction)
1010 Dune::MatrixAdapter<
1020 Dune::RestartedGMResSolver<Native<V>> solver(opa,ilu0,reduction,restart,maxiter,verbose);
1021 Dune::InverseOperatorResult stat;
1031 int restart, maxiter, verbose;
1060 #endif // DUNE_PDELAB_BACKEND_ISTL_SEQISTLSOLVERBACKEND_HH
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:585
void apply(V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:182
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:790
ISTLBackend_SEQ_AMG(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: seqistlsolverbackend.hh:735
RFType reduction
Definition: solver.hh:35
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: seqistlsolverbackend.hh:659
ISTLBackend_SEQ_MatrixFree_BCGS_Richardson(const GO &go_, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:1048
Backend for conjugate gradient solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:538
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:328
ISTLBackend_SEQ_BCGS_Richardson(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:393
Definition: seqistlsolverbackend.hh:307
bool directCoarseLevelSolver
True if a direct solver was used on the coarset level.
Definition: seqistlsolverbackend.hh:716
OnTheFlyOperator(const GO &go_)
Definition: seqistlsolverbackend.hh:54
Backend for sequential conjugate gradient solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:504
RFType conv_rate
Definition: solver.hh:36
bool getReuse() const
Return whether the AMG is reused during call to apply()
Definition: seqistlsolverbackend.hh:769
Backend for sequential conjugate gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:451
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
Solver backend using UMFPack as a direct solver.
Definition: seqistlsolverbackend.hh:609
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:63
virtual void apply(const X &x, Y &y) const override
Definition: seqistlsolverbackend.hh:66
Definition: seqistlsolverbackend.hh:45
ISTLBackend_SEQ_CG_Jac(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:546
ISTLBackend_SEQ_CG_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:460
Sequential Loop solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:932
void setparams(Parameters params_)
set AMG parameters
Definition: seqistlsolverbackend.hh:757
Backend for sequential BiCGSTAB solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:417
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
bool converged
Definition: solver.hh:32
ISTLBackend_SEQ_GMRES_ILU0(int restart_=200, int maxiter_=5000, int verbose_=1)
make linear solver object
Definition: seqistlsolverbackend.hh:994
Sequential congute gradient solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:484
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:232
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:192
Class providing some statistics of the AMG solver.
Definition: seqistlsolverbackend.hh:700
Definition: seqistlsolverbackend.hh:158
unsigned int iterations
Definition: solver.hh:33
Linear solver backend for Restarted GMRes preconditioned with ILU(0)
Definition: seqistlsolverbackend.hh:983
ISTLBackend_SEQ_UMFPack(int maxiter, int verbose_)
make a linear solver object
Definition: seqistlsolverbackend.hh:627
ISTLBackend_SEQ_MatrixFree_Richardson(const GO &go, unsigned maxiter=5000, int verbose=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:169
Backend for sequential BiCGSTAB solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:401
Definition: recipe-operator-splitting.cc:107
ISTLBackend_SEQ_ILUn(int n, double w, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:317
ISTLBackend_SEQ_BCGS_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:478
X::field_type field_type
Definition: seqistlsolverbackend.hh:50
ISTLBackend_SEQ_BCGS_Jac(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:409
ISTLBackend_SEQ_CG_SSOR(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:513
ISTLBackend_SEQ_BCGS_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:443
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:676
int levels
the number of levels in the AMG hierarchy.
Definition: seqistlsolverbackend.hh:708
Sequential conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:854
X domain_type
Definition: seqistlsolverbackend.hh:48
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:281
Y range_type
Definition: seqistlsolverbackend.hh:49
static constexpr bool isLinear
Definition: seqistlsolverbackend.hh:51
ISTLBackend_SEQ_LS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:971
Backend for sequential loop solver with Jacobi preconditioner.
Definition: seqistlsolverbackend.hh:367
double tsolve
The time spent in solving the system (without building the hierarchy.
Definition: seqistlsolverbackend.hh:710
double tsetup
The time needed for building the AMG hierarchy (coarsening).
Definition: seqistlsolverbackend.hh:712
void setLinearizationPoint(const V &u)
Set position of jacobian, ust be called before apply() for nonlinear problems.
Definition: seqistlsolverbackend.hh:196
ISTLBackend_SEQ_ExplicitDiagonal()
make a linear solver object
Definition: seqistlsolverbackend.hh:665
ISTLBackend_SEQ_Richardson(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:121
ISTLBackend_SEQ_MINRES_SSOR(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:530
Definition: seqistlsolverbackend.hh:261
ISTLBackend_SEQ_BCGS_AMG_SOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:919
ISTLBackend_SEQ_LS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:945
ISTLBackend_SEQ_ILU0(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:270
ISTLBackend_SEQ_CG_ILUn(int n_, double w_=1.0, unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:496
Definition: seqistlsolverbackend.hh:112
Sequential Loop solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:958
Backend using a MINRes solver preconditioned by SSOR.
Definition: seqistlsolverbackend.hh:521
SolverCategory::Category category() const override
Definition: seqistlsolverbackend.hh:92
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:639
ISTLBackend_SEQ_BCGS_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:893
int iterations
The number of iterations performed until convergence was reached.
Definition: seqistlsolverbackend.hh:714
void setLinearizationPoint(const X &u)
Definition: seqistlsolverbackend.hh:61
Sequential BiCGSTAB solver preconditioned with AMG smoothed by SOR.
Definition: seqistlsolverbackend.hh:906
ISTLBackend_SEQ_SuperLU(int maxiter, int verbose_)
make a linear solver object
Definition: seqistlsolverbackend.hh:573
ISTLBackend_SEQ_SuperLU(int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:563
Definition: seqistlsolverbackend.hh:721
Definition: seqistlsolverbackend.hh:209
Backend for sequential BiCGSTAB solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:434
double elapsed
Definition: solver.hh:34
Definition: seqistlsolverbackend.hh:1039
double tprepare
The needed for computing the parallel information and for adapting the linear system.
Definition: seqistlsolverbackend.hh:706
Solver backend using SuperLU as a direct solver.
Definition: seqistlsolverbackend.hh:555
ISTLBackend_SEQ_UMFPack(int verbose_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:617
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:133
ISTLBackend_SEQ_CG_AMG_SSOR(unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Constructor.
Definition: seqistlsolverbackend.hh:867
void apply(M &A, V &z, W &r, typename Dune::template FieldTraits< typename W::ElementType >::real_type reduction)
solve the given linear system
Definition: seqistlsolverbackend.hh:1006
Sequential BiCGStab solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:880
Sequential BiCGStab solver with ILU0 preconditioner.
Definition: seqistlsolverbackend.hh:466
ISTLBackend_SEQ_BCGS_SSOR(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:426
VTKWriter & w
Definition: function.hh:842
ISTLBackend_SEQ_Base(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:218
ISTLBackend_SEQ_LOOP_Jac(unsigned maxiter_=5000, int verbose_=1, unsigned preconditioner_steps_=1)
make a linear solver object
Definition: seqistlsolverbackend.hh:375
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const override
Definition: seqistlsolverbackend.hh:78
const ISTLAMGStatistics & statistics() const
Get statistics of the AMG solver (no of levels, timings).
Definition: seqistlsolverbackend.hh:828
void setReuse(bool reuse_)
Set whether the AMG should be reused again during call to apply().
Definition: seqistlsolverbackend.hh:763
Backend for sequential BiCGSTAB solver with Richardson precondition (equivalent to no preconditioner ...
Definition: seqistlsolverbackend.hh:385
V::ElementType norm(const V &v) const
compute global norm of a vector
Definition: seqistlsolverbackend.hh:778