4#ifndef DUNE_PDELAB_BACKEND_ISTL_MATRIXFREE_BACKENDS_HH
5#define DUNE_PDELAB_BACKEND_ISTL_MATRIXFREE_BACKENDS_HH
7#include<dune/grid/common/rangegenerators.hh>
8#include<dune/istl/preconditioner.hh>
9#include<dune/istl/operators.hh>
10#include<dune/istl/solvers.hh>
28 template <
typename JacobianLOP,
typename BlockOffDiagonalLOP,
typename Gr
idFunctionSpace>
38 template<
typename GFSU,
typename GFSV,
typename LOP,
39 typename MB,
typename DF,
typename RF,
typename JF,
40 typename CU,
typename CV>
62 template<
class GO,
class PrecGO,
63 template<
class>
class Solver>
68 using V =
typename GO::Traits::Domain;
69 using W =
typename GO::Traits::Range;
75 unsigned maxiter=5000,
int verbose=1)
76 : opa_(go), seqprec_(precgo)
77 , maxiter_(maxiter), verbose_(verbose)
86 typename PrecGO::Traits::LocalAssembler::LocalOperator>
::value
88 "If you use the BlockSORPreconditioner you need to use FastDGGridOperator!");
94 "In order to use the matrix-free solvers your grid operator must follow the new\n"
95 "local operator interface for nonlinear jacobian applys, where the current\n"
96 "solution and the linearization point can have different types. This is best\n"
97 "explained by an example.\n\n"
98 "old: void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const X& z, const LFSV& lfsv, Y& y) const {}\n"
99 "new: void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const Z& z, const LFSV& lfsv, Y& y) const {}\n\n"
100 "Note that in the new interface the type of z is Z and doesn't have to be X.\n"
101 "You need to adjust all the nonlinear jacobian apply methods in your local operator\n"
102 "in a similar way.");
105 "In order to use the matrix-free solvers your grid operator must follow the new\n"
106 "local operator interface for nonlinear jacobian applys, where the current\n"
107 "solution and the linearization point can have different types. This is best\n"
108 "explained by an example.\n\n"
109 "old: void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const X& z, const LFSV& lfsv, Y& y) const {}\n"
110 "new: void jacobian_apply_volume (const EG& eg, const LFSU& lfsu, const X& x, const Z& z, const LFSV& lfsv, Y& y) const {}\n\n"
111 "Note that in the new interface the type of z is Z and doesn't have to be X.\n"
112 "You need to adjust all the nonlinear jacobian apply methods in your local operator\n"
113 "in a similar way.");
116 void apply(V& z, W& r,
typename V::ElementType reduction)
118 Solver<V> solver(opa_,seqprec_,reduction,maxiter_,verbose_);
119 Dune::InverseOperatorResult stat;
120 solver.apply(z, r, stat);
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
constexpr auto hasOldLOPInterface(T &t) -> typename std::enable_if<(decltype(hasOldOrNewJacobianApplyVolume(t))::value &&!decltype(hasNewJacobianApplyVolume(t))::value)||(decltype(hasOldOrNewJacobianApplyVolumePostSkeleton(t))::value &&!decltype(hasNewJacobianApplyVolumePostSkeleton(t))::value)||(decltype(hasOldOrNewJacobianApplySkeleton(t))::value &&!decltype(hasNewJacobianApplySkeleton(t))::value)||(decltype(hasOldOrNewJacobianApplyBoundary(t))::value &&!decltype(hasNewJacobianApplyBoundary(t))::value), std::true_type >::type
Definition: checklopinterface.hh:153
Definition: backends.hh:26
Definition: backends.hh:36
Sequential matrix-free solver backend.
Definition: backends.hh:67
ISTLBackend_SEQ_MatrixFree_Base(const GO &go, const PrecGO &precgo, unsigned maxiter=5000, int verbose=1)
Definition: backends.hh:74
void setLinearizationPoint(const V &u)
Definition: backends.hh:130
static constexpr bool isMatrixFree
Definition: backends.hh:136
void apply(V &z, W &r, typename V::ElementType reduction)
Definition: backends.hh:116
Definition: backends.hh:149
ISTLBackend_SEQ_BCGS_SOR(unsigned maxiter_=5000, int verbose_=1)
make a linear solver object
Definition: backends.hh:156
Turn a matrix-free Jacobi-type local preconditioner to a SOR one.
Definition: blocksorpreconditioner.hh:41
Turn a grid operator that represents a preconditioner into an ISTL preconditioner.
Definition: gridoperatorpreconditioner.hh:21
void setLinearizationPoint(const Domain &u)
Definition: gridoperatorpreconditioner.hh:39
void setLinearizationPoint(const X &u)
Definition: seqistlsolverbackend.hh:61
Definition: seqistlsolverbackend.hh:211
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
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:63
A grid function space.
Definition: gridfunctionspace.hh:191
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139