4 #ifndef DUNE_PDELAB_LOCALOPERATOR_BLOCKDIAGONALWRAPPER_HH
5 #define DUNE_PDELAB_LOCALOPERATOR_BLOCKDIAGONALWRAPPER_HH
17 template <
typename AccumulationView>
23 : _view(view), _diagonal(diagonal)
28 return _view.weight();
41 template <
typename LFSU,
typename I,
typename Value>
45 _view.accumulate(lfsu, i,
value);
49 template <
typename LFSU,
typename I,
typename LFSV,
typename J,
typename Value>
53 _view.accumulate(lfsu, i, lfsv, j,
value);
58 AccumulationView& _view;
79 template <
class LocalOperator>
95 static constexpr
bool isLinear = LocalOperator::isLinear;
109 : _localOperator(localOperator)
114 : _localOperator(other._localOperator)
118 template<
typename LFSU,
typename LFSV,
typename LocalPattern>
120 LocalPattern& pattern)
const
122 _localOperator.pattern_volume(lfsu, lfsv, pattern);
125 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename MAT>
132 _localOperator.jacobian_volume(eg,lfsu,x,lfsv,mat);
135 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename MAT>
137 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
138 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
139 MAT& mat_ss, MAT& mat_sn,
140 MAT& mat_ns, MAT& mat_nn)
const
148 _localOperator.jacobian_skeleton(
ig, lfsu_s, x_s, lfsv_s, lfsu_n, x_n, lfsv_n, view_ss, view_other, view_other, view_other);
151 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename MAT>
153 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
156 _localOperator.jacobian_boundary(
ig, lfsu_s, x_s, lfsv_s, mat_ss);
160 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
165 template<
typename EG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
172 template<
typename IG,
typename LFSU,
typename Z,
typename LFSV,
typename Y>
174 const LFSU& lfsu_s,
const Z& z_s,
const LFSV& lfsv_s,
175 const LFSU& lfsu_n,
const Z& z_n,
const LFSV& lfsv_n,
176 Y& y_s, Y& y_n)
const
183 Z z_zero(z_n.size(), 0.0);
186 template<
typename IG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
188 const LFSU& lfsu_s,
const X& x_s,
const Z& z_s,
const LFSV& lfsv_s,
189 const LFSU& lfsu_n,
const X& x_n,
const Z& z_n,
const LFSV& lfsv_n,
190 Y& y_s, Y& y_n)
const
197 Z z_zero(z_n.size(), 0.0);
198 Dune::PDELab::impl::jacobianApplySkeleton(_localOperator,
ig, lfsu_s, x_s, z_s, lfsv_s, lfsu_n, x_n, z_zero, lfsv_n, view_s, view_n);
201 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
203 const LFSU& lfsu_s,
const X& z_s,
const LFSV& lfsv_s,
209 template<
typename IG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
211 const LFSU& lfsu_s,
const X& x_s,
const Z& z_s,
const LFSV& lfsv_s,
220 const LocalOperator& _localOperator;
225 template <
typename LocalOperator,
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename MAT>
234 localOperator.jacobian_volume(eg, lfsu, x, lfsv, mat);
237 auto entitySet = lfsu.gridFunctionSpace().entitySet();
238 std::size_t intersectionIndex = 0;
239 for (
const auto& is : intersections(lfsu.gridFunctionSpace().gridView(), eg.entity()))
243 auto intersectionType = std::get<0>(intersectionData);
246 switch (intersectionType){
248 localOperator.jacobian_skeleton(
ig, lfsu, x, lfsv, lfsu, x, lfsv, mat, mat, mat, mat);
251 localOperator.jacobian_boundary(
ig, lfsu, x, lfsv, mat);
260 template<
typename LocalOperator,
typename EG,
typename LFSU,
typename X,
typename Z,
typename LFSV,
typename Y>
270 if (LocalOperator::isLinear)
276 auto entitySet = lfsu.gridFunctionSpace().entitySet();
277 std::size_t intersectionIndex = 0;
278 for (
const auto& is : intersections(lfsu.gridFunctionSpace().gridView(), eg.entity()))
282 auto intersectionType = std::get<0>(intersectionData);
285 switch (intersectionType){
287 if (LocalOperator::isLinear)
290 Dune::PDELab::impl::jacobianApplySkeleton(localOperator,
ig, lfsu, x, z, lfsv, lfsu, x, z, lfsv, y, y);
293 if (LocalOperator::isLinear)
const IG & ig
Definition: constraints.hh:149
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
@ boundary
domain boundary intersection (neighbor() == false && boundary() == true)
@ skeleton
skeleton intersection (neighbor() == true && boundary() == false)
void assembleLocalDiagonalBlock(const LocalOperator &localOperator, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, MAT &mat)
A function for assembling a single diagonal block.
Definition: blockdiagonalwrapper.hh:226
std::tuple< IntersectionType, typename EntitySet::Element > classifyIntersection(const EntitySet &entity_set, const Intersection &is)
Classifies the type of an intersection wrt to the passed EntitySet.
Definition: intersectiontype.hh:37
void applyLocalDiagonalBlock(const LocalOperator &localOperator, const EG &eg, const LFSU &lfsu, const X &x, const Z &z, const LFSV &lfsv, Y &y)
A function for applying a single diagonal block.
Definition: blockdiagonalwrapper.hh:261
std::enable_if_t< LOP::isLinear > jacobianApplySkeleton(const LOP &lop, const IG &ig, const LFSU &lfsu_s, const X &z_s, const LFSV &lfsv_s, const LFSU &lfsu_n, const X &z_n, const LFSV &lfsv_n, Y &y_s, Y &y_n)
Definition: jacobianapplyhelper.hh:64
std::enable_if_t< LOP::isLinear > jacobianApplyVolume(const LOP &lop, const EG &eg, const LFSU &lfsu, const X &z, const LFSV &lfsv, Y &y)
Definition: jacobianapplyhelper.hh:23
std::enable_if_t< LOP::isLinear > jacobianApplyBoundary(const LOP &lop, const IG &ig, const LFSU &lfsu_s, const X &z_s, const LFSV &lfsv_s, Y &y_s)
Definition: jacobianapplyhelper.hh:109
Wrap intersection.
Definition: geometrywrapper.hh:57
Definition: blockdiagonalwrapper.hh:19
auto data()
Definition: blockdiagonalwrapper.hh:31
const auto data() const
Definition: blockdiagonalwrapper.hh:36
const auto weight()
Definition: blockdiagonalwrapper.hh:26
void accumulate(const LFSU &lfsu, I i, const LFSV &lfsv, J j, Value value)
Definition: blockdiagonalwrapper.hh:50
BlockDiagonalAccumulationViewWrapper(AccumulationView &view, bool diagonal)
Definition: blockdiagonalwrapper.hh:22
void accumulate(const LFSU &lfsu, I i, Value value)
Definition: blockdiagonalwrapper.hh:42
A local operator that accumulates the block diagonal.
Definition: blockdiagonalwrapper.hh:82
void jacobian_apply_boundary(const IG &ig, const LFSU &lfsu_s, const X &z_s, const LFSV &lfsv_s, Y &y_s) const
Definition: blockdiagonalwrapper.hh:202
static constexpr bool doAlphaBoundary
Definition: blockdiagonalwrapper.hh:92
BlockDiagonalLocalOperatorWrapper(const LocalOperator &localOperator)
Construct new instance of class.
Definition: blockdiagonalwrapper.hh:108
BlockDiagonalLocalOperatorWrapper(const BlockDiagonalLocalOperatorWrapper &other)
Copy constructor.
Definition: blockdiagonalwrapper.hh:113
static constexpr bool doSkeletonTwoSided
Definition: blockdiagonalwrapper.hh:102
void jacobian_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, MAT &mat_ss) const
Definition: blockdiagonalwrapper.hh:152
void jacobian_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, MAT &mat_ss, MAT &mat_sn, MAT &mat_ns, MAT &mat_nn) const
Definition: blockdiagonalwrapper.hh:136
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const Z &z, const LFSV &lfsv, Y &y) const
Definition: blockdiagonalwrapper.hh:166
static constexpr bool doAlphaVolume
Definition: blockdiagonalwrapper.hh:90
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, MAT &mat) const
Definition: blockdiagonalwrapper.hh:126
static constexpr bool isLinear
Definition: blockdiagonalwrapper.hh:95
static constexpr bool doPatternVolume
Definition: blockdiagonalwrapper.hh:86
void jacobian_apply_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const Z &z_s, const LFSV &lfsv_s, Y &y_s) const
Definition: blockdiagonalwrapper.hh:210
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
Definition: blockdiagonalwrapper.hh:173
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
Definition: blockdiagonalwrapper.hh:187
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &z, const LFSV &lfsv, Y &y) const
Definition: blockdiagonalwrapper.hh:161
static constexpr bool doAlphaSkeleton
Definition: blockdiagonalwrapper.hh:91
void pattern_volume(const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
Definition: blockdiagonalwrapper.hh:119
Default flags for all local operators.
Definition: flags.hh:19
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139