3 #ifndef DUNE_PDELAB_BACKEND_ISTL_MATRIXFREE_BLOCKSORPRECONDITIONER_HH
4 #define DUNE_PDELAB_BACKEND_ISTL_MATRIXFREE_BLOCKSORPRECONDITIONER_HH
37 template<
typename JacobianLOP,
typename BlockOffDiagonalLOP,
typename Gr
idFunctionSpace>
42 using value_type =
typename GridFunctionSpace::Traits::GridView::Grid::ctype;
56 static constexpr
bool isLinear = JacobianLOP::isLinear;
59 BlockOffDiagonalLOP& blockOffDiagonalLOP,
61 const double omega=1.0)
62 : _jacobianlop(jacobianlop)
63 , _blockOffDiagonalLOP(blockOffDiagonalLOP)
65 , _a_i(gridFunctionSpace.ordering().maxLocalSize())
66 , _b_i(gridFunctionSpace.ordering().maxLocalSize())
75 return _jacobianlop.requireSetup();
79 _jacobianlop.setRequireSetup(v);
83 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
84 void alpha_volume(
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, Y& y)
const
86 _jacobianlop.alpha_volume(eg,lfsu,x,lfsv,y);
90 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
92 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
93 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
98 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
103 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
105 const LFSU& lfsu,
const X& x,
106 const LFSV& lfsv, Y& y)
const
128 template<
typename EG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
135 template<
typename IG,
typename LFSU,
typename Z,
typename LFSV,
typename Y>
137 const LFSU& lfsu_s,
const Z& z_s,
const LFSV& lfsv_s,
138 const LFSU& lfsu_n,
const Z& z_n,
const LFSV& lfsv_n,
139 Y& y_s, Y& y_n)
const
195 W _a_i_view(_a_i, y_s.weight());
196 if (
ig.inside().partitionType() == Dune::InteriorEntity){
202 _blockOffDiagonalLOP.jacobian_apply_skeleton(
ig, lfsu_s, y_s, lfsv_s, lfsu_n, y_n, lfsv_s, _a_i_view, _a_i_view);
207 template<
typename IG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
209 const LFSU& lfsu_s,
const X& x_s,
const Z& z_s,
const LFSV& lfsv_s,
210 const LFSU& lfsu_n,
const X& x_n,
const Z& z_n,
const LFSV& lfsv_n,
211 Y& y_s, Y& y_n)
const
214 W _a_i_view(_a_i, y_s.weight());
215 if (
ig.inside().partitionType() == Dune::InteriorEntity)
216 _blockOffDiagonalLOP.jacobian_apply_skeleton(
ig, lfsu_s, x_s, y_s, lfsv_s, lfsu_n, x_n, y_n, lfsv_s, _a_i_view, _a_i_view);
220 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
222 const LFSU& lfsu,
const X& x,
223 const LFSV& lfsv, Y& y)
const
228 std::transform(_a_i.
data(),
232 std::minus<value_type>{});
237 W _b_i_view(_b_i, y.weight());
238 _jacobianlop.jacobian_apply_volume(eg, lfsu, _a_i, lfsv, _b_i_view);
242 std::transform(y.data(),
249 return (1-_omega)*a - _omega*b;
254 template<
typename EG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
271 std::minus<value_type>{});
274 W _b_i_view(_b_i,y.weight());
275 _jacobianlop.jacobian_apply_volume(eg, lfsu, x, _a_i, lfsv, _b_i_view);
278 std::transform(y.data(),
285 return (1-_omega)*a - _omega*b;
291 JacobianLOP& _jacobianlop;
292 BlockOffDiagonalLOP& _blockOffDiagonalLOP;
const IG & ig
Definition: constraints.hh:149
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Turn a matrix-free Jacobi-type local preconditioner to a SOR one.
Definition: blocksorpreconditioner.hh:41
void jacobian_apply_skeleton(const IG &ig, const LFSU &lfsu_s, const X &x_s, const Z &z_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &x_n, const Z &z_n, const LFSV &lfsv_n, Y &y_s, Y &y_n) const
Gather off-block-diagonals in Gauss-Seidel process of linearized operator.
Definition: blocksorpreconditioner.hh:208
static constexpr bool doAlphaVolume
Definition: blocksorpreconditioner.hh:52
void setRequireSetup(bool v)
Definition: blocksorpreconditioner.hh:77
static constexpr bool doAlphaVolumePostSkeleton
Definition: blocksorpreconditioner.hh:54
BlockSORPreconditionerLocalOperator(JacobianLOP &jacobianlop, BlockOffDiagonalLOP &blockOffDiagonalLOP, const GridFunctionSpace &gridFunctionSpace, const double omega=1.0)
Definition: blocksorpreconditioner.hh:58
void jacobian_apply_volume_post_skeleton(const EG &eg, const LFSU &lfsu, const X &x, const Z &z, const LFSV &lfsv, Y &y) const
apply preconditioner after skeleton terms, linearized version
Definition: blocksorpreconditioner.hh:255
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Prepare underlying diagonal block preconditioner.
Definition: blocksorpreconditioner.hh:84
void alpha_skeleton(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &x_n, const LFSV &lfsv_n, Y &y_s, Y &y_n) const
Provide this method, but it actually does not nothing.
Definition: blocksorpreconditioner.hh:91
bool requireSetup()
Definition: blocksorpreconditioner.hh:73
static constexpr bool doSkeletonTwoSided
Definition: blocksorpreconditioner.hh:55
void jacobian_apply_volume_post_skeleton(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Apply preconditioner after skeleton terms, linear version.
Definition: blocksorpreconditioner.hh:221
static constexpr bool doPatternVolumePostSkeleton
Definition: blocksorpreconditioner.hh:51
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Linear operator application, volume terms.
Definition: blocksorpreconditioner.hh:104
static constexpr bool doAlphaSkeleton
Definition: blocksorpreconditioner.hh:53
static constexpr bool isLinear
Definition: blocksorpreconditioner.hh:56
static constexpr bool doPatternVolume
Definition: blocksorpreconditioner.hh:49
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const Z &z, const LFSV &lfsv, Y &y) const
linearized operator application, volume terms
Definition: blocksorpreconditioner.hh:129
void alpha_volume_post_skeleton(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Provide this method, but it actually does nothing.
Definition: blocksorpreconditioner.hh:99
static constexpr bool doPatternSkeleton
Definition: blocksorpreconditioner.hh:50
void jacobian_apply_skeleton(const IG &ig, const LFSU &lfsu_s, const Z &z_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const Z &z_n, const LFSV &lfsv_n, Y &y_s, Y &y_n) const
Gather off-block-diagonals in Gauss-Seidel process of linear operator.
Definition: blocksorpreconditioner.hh:136
A grid function space.
Definition: gridfunctionspace.hh:191
size_type size() const
The size of the container.
Definition: localvector.hh:317
WeightedVectorAccumulationView< LocalVector > WeightedAccumulationView
An accumulate-only view of this container that automatically applies a weight to all contributions.
Definition: localvector.hh:211
auto data()
Access underlying container.
Definition: localvector.hh:244
Default flags for all local operators.
Definition: flags.hh:19
sparsity pattern generator
Definition: pattern.hh:14