2 #ifndef DUNE_PDELAB_LOCALOPERATOR_L2_HH 3 #define DUNE_PDELAB_LOCALOPERATOR_L2_HH 7 #include<dune/common/exceptions.hh> 8 #include<dune/common/fvector.hh> 10 #include<dune/geometry/type.hh> 11 #include<dune/geometry/referenceelements.hh> 12 #include<dune/geometry/quadraturerules.hh> 14 #include<dune/localfunctions/common/interfaceswitch.hh> 48 L2 (
int intorder_=2,
double scaling=1.0)
54 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
55 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 58 using FESwitch = FiniteElementInterfaceSwitch<
59 typename LFSU::Traits::FiniteElementType>;
60 using BasisSwitch = BasisInterfaceSwitch<
61 typename FESwitch::Basis>;
64 using RF =
typename BasisSwitch::RangeField;
65 using RangeType =
typename BasisSwitch::Range;
66 using size_type =
typename LFSU::Traits::SizeType;
69 auto geo = eg.geometry();
72 std::vector<RangeType> phi(lfsu.size());
78 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
82 for (size_type i=0; i<lfsu.size(); i++)
83 u += RF(x(lfsu,i)*phi[i]);
86 auto factor = _scaling * qp.weight() * geo.integrationElement(qp.position());
87 for (size_type i=0; i<lfsu.size(); i++)
88 r.accumulate(lfsv,i, u*phi[i]*factor);
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,
98 using FESwitch = FiniteElementInterfaceSwitch<
99 typename LFSU::Traits::FiniteElementType>;
100 using BasisSwitch = BasisInterfaceSwitch<
101 typename FESwitch::Basis>;
104 using RangeType =
typename BasisSwitch::Range;
105 using size_type =
typename LFSU::Traits::SizeType;
108 auto geo = eg.geometry();
111 std::vector<RangeType> phi(lfsu.size());
117 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
120 auto factor = _scaling * qp.weight() * geo.integrationElement(qp.position());
121 for (size_type j=0; j<lfsu.size(); j++)
122 for (size_type i=0; i<lfsu.size(); i++)
123 mat.accumulate(lfsv,i,lfsu,j, phi[j]*phi[i]*factor);
129 const double _scaling;
151 : scalar_operator(intorder_)
155 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
156 void alpha_volume (
const EG& eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const 158 for(
unsigned int i=0; i<TypeTree::degree(lfsu); ++i)
159 scalar_operator.alpha_volume(eg,lfsu.child(i),x,lfsv.child(i),r);
163 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
167 for(
unsigned int i=0; i<TypeTree::degree(lfsu); ++i)
168 scalar_operator.jacobian_volume(eg,lfsu.child(i),x,lfsv.child(i),mat);
179 #endif // DUNE_PDELAB_LOCALOPERATOR_L2_HH void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:55
PowerL2(int intorder_=2)
Definition: l2.hh:150
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:94
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:156
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
Implements linear and nonlinear versions of jacobian_apply_volume() based on alpha_volume() ...
Definition: numericaljacobianapply.hh:32
sparsity pattern generator
Definition: pattern.hh:13
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
L2(int intorder_=2, double scaling=1.0)
Definition: l2.hh:48
Default flags for all local operators.
Definition: flags.hh:18
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:164