3#ifndef DUNE_PDELAB_LOCALOPERATOR_POINTDIAGONALWRAPPER_HH
4#define DUNE_PDELAB_LOCALOPERATOR_POINTDIAGONALWRAPPER_HH
13 template <
typename AccumulationView>
18 : _view(view), _diagonal(diagonal)
21 template <
typename LFSU,
typename I,
typename LFSV,
typename J,
typename Value>
24 if (_diagonal && i == j){
25 _view.accumulate(lfsu, i,
value);
30 AccumulationView& _view;
51 template <
class LocalOperator>
67 static constexpr bool isLinear = LocalOperator::isLinear;
81 : _localOperator(localOperator)
87 : _localOperator(other._localOperator)
91 template<
typename LFSU,
typename LFSV,
typename LocalPattern>
93 LocalPattern& pattern)
const
95 _localOperator.pattern_volume(lfsu, lfsv, pattern);
98 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
99 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
102 _localOperator.jacobian_volume(eg, lfsu, x, lfsv, view);
105 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
107 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
108 const LFSU& lfsu_n,
const X& x_n,
const LFSV& lfsv_n,
109 R& r_s, R& r_n)
const
113 _localOperator.jacobian_skeleton(
ig,
116 view_ss, view_other, view_other, view_other);
119 template<
typename IG,
typename LFSU,
typename X,
typename LFSV,
typename R>
121 const LFSU& lfsu_s,
const X& x_s,
const LFSV& lfsv_s,
125 _localOperator.jacobian_boundary(
ig, lfsu_s, x_s, lfsv_s, view);
130 const LocalOperator& _localOperator;
138 template <
typename LocalOperator,
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
147 localOperator.alpha_volume(eg, lfsu, x, lfsv, y);
150 auto entitySet = lfsu.gridFunctionSpace().entitySet();
151 std::size_t intersectionIndex = 0;
152 for (
const auto& is : intersections(lfsu.gridFunctionSpace().gridView(), eg.entity()))
155 typename std::remove_reference<
156 typename std::remove_const<
160 >
ig(is, intersectionIndex++);
162 auto intersectionType = std::get<0>(intersectionData);
165 switch (intersectionType){
167 localOperator.alpha_skeleton(
ig, lfsu, x, lfsv, lfsu, x, lfsv, y, y);
170 localOperator.alpha_boundary(
ig, lfsu, x, lfsv, y);
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)
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 assembleLocalPointDiagonal(const LocalOperator &localOperator, const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y)
A function for assembling the point diagonal of a single block.
Definition: pointdiagonalwrapper.hh:139
Wrap intersection.
Definition: geometrywrapper.hh:57
Default flags for all local operators.
Definition: flags.hh:19
Definition: pointdiagonalwrapper.hh:15
void accumulate(const LFSU &lfsu, I i, const LFSV &lfsv, J j, Value value)
Definition: pointdiagonalwrapper.hh:22
PointDiagonalAccumulationViewWrapper(AccumulationView &view, bool diagonal)
Definition: pointdiagonalwrapper.hh:17
A local operator that accumulates the point diagonal.
Definition: pointdiagonalwrapper.hh:54
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, R &r_s, R &r_n) const
Definition: pointdiagonalwrapper.hh:106
PointDiagonalLocalOperatorWrapper(const LocalOperator &localOperator)
Construct new instance of class.
Definition: pointdiagonalwrapper.hh:80
static constexpr bool isLinear
Definition: pointdiagonalwrapper.hh:67
static constexpr bool doAlphaVolume
Definition: pointdiagonalwrapper.hh:62
static constexpr bool doPatternVolume
Definition: pointdiagonalwrapper.hh:58
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: pointdiagonalwrapper.hh:99
PointDiagonalLocalOperatorWrapper(const PointDiagonalLocalOperatorWrapper &other)
Copy constructor.
Definition: pointdiagonalwrapper.hh:86
static constexpr bool doAlphaBoundary
Definition: pointdiagonalwrapper.hh:64
void pattern_volume(const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
Definition: pointdiagonalwrapper.hh:92
static constexpr bool doAlphaSkeleton
Definition: pointdiagonalwrapper.hh:63
void alpha_boundary(const IG &ig, const LFSU &lfsu_s, const X &x_s, const LFSV &lfsv_s, R &r_s) const
Definition: pointdiagonalwrapper.hh:120
static constexpr bool doSkeletonTwoSided
Definition: pointdiagonalwrapper.hh:74
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139