2 #ifndef DUNE_PDELAB_LOCALOPERATOR_L2_HH
3 #define DUNE_PDELAB_LOCALOPERATOR_L2_HH
7 #include <dune/common/fvector.hh>
9 #include <dune/localfunctions/common/interfaceswitch.hh>
39 : _intorderadd(intorderadd)
44 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
45 void alpha_volume(
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
48 using FESwitch = FiniteElementInterfaceSwitch<
49 typename LFSU::Traits::FiniteElementType>;
50 using BasisSwitch = BasisInterfaceSwitch<
51 typename FESwitch::Basis>;
54 using RF =
typename BasisSwitch::RangeField;
55 using RangeType =
typename BasisSwitch::Range;
56 using size_type =
typename LFSU::Traits::SizeType;
59 auto geo = eg.geometry();
62 std::vector<RangeType> phi(lfsu.size());
65 auto intorder = 2*FESwitch::basis(lfsu.finiteElement()).order() + _intorderadd;
71 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
75 for (size_type i=0; i<lfsu.size(); i++)
76 u += RF(x(lfsu,i)*phi[i]);
79 auto factor = _scaling * qp.weight() * geo.integrationElement(qp.position());
80 for (size_type i=0; i<lfsu.size(); i++)
81 r.accumulate(lfsv,i, u*phi[i]*factor);
86 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename Y>
93 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
94 void jacobian_volume(
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, M & mat)
const
97 using FESwitch = FiniteElementInterfaceSwitch<
98 typename LFSU::Traits::FiniteElementType>;
99 using BasisSwitch = BasisInterfaceSwitch<
100 typename FESwitch::Basis>;
103 using RangeType =
typename BasisSwitch::Range;
104 using size_type =
typename LFSU::Traits::SizeType;
107 auto geo = eg.geometry();
110 std::vector<RangeType> phi(lfsu.size());
113 auto intorder = 2*FESwitch::basis(lfsu.finiteElement()).order() + _intorderadd;
119 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
122 auto factor = _scaling * qp.weight() * geo.integrationElement(qp.position());
123 for (size_type j=0; j<lfsu.size(); j++)
124 for (size_type i=0; i<lfsu.size(); i++)
125 mat.accumulate(lfsv,i,lfsu,j, phi[j]*phi[i]*factor);
164 L2 (
int intorderadd = 0,
double scaling = 1.0)
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
QuadratureRuleWrapper< QuadratureRule< typename Geometry::ctype, Geometry::mydimension > > quadratureRule(const Geometry &geo, std::size_t order, QuadratureType::Enum quadrature_type=QuadratureType::GaussLegendre)
Returns a quadrature rule for the given geometry.
Definition: quadraturerules.hh:117
Block diagonal extension of scalar local operator.
Definition: blockdiagonal.hh:78
Default flags for all local operators.
Definition: flags.hh:19
Default class for additional methods in instationary local operators.
Definition: idefault.hh:90
@ doPatternVolume
Definition: l2.hh:33
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:45
@ doAlphaVolume
Definition: l2.hh:36
ScalarL2(int intorderadd, double scaling)
Definition: l2.hh:38
void jacobian_apply_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
Definition: l2.hh:87
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:94
L2(int intorderadd=0, double scaling=1.0)
Constructs a new L2 operator.
Definition: l2.hh:164
sparsity pattern generator
Definition: pattern.hh:14