2 #ifndef DUNE_PDELAB_L2_HH
3 #define DUNE_PDELAB_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>
45 L2 (
int intorder_=2,
double scaling=1.0)
51 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
52 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
55 typedef FiniteElementInterfaceSwitch<
56 typename LFSU::Traits::FiniteElementType
58 typedef BasisInterfaceSwitch<
59 typename FESwitch::Basis
63 typedef typename BasisSwitch::DomainField DF;
64 typedef typename BasisSwitch::RangeField RF;
65 typedef typename BasisSwitch::Range RangeType;
67 typedef typename LFSU::Traits::SizeType size_type;
70 const int dim = EG::Geometry::mydimension;
73 const auto& geometry = eg.geometry();
74 Dune::GeometryType gt = geometry.type();
75 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
78 for (
const auto& qp : rule)
81 std::vector<RangeType> phi(lfsu.size());
82 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
86 for (size_type i=0; i<lfsu.size(); i++)
87 u += RF(x(lfsu,i)*phi[i]);
90 RF factor = _scaling * qp.weight() * geometry.integrationElement(qp.position());
91 for (size_type i=0; i<lfsu.size(); i++)
92 r.accumulate(lfsv,i, u*phi[i]*factor);
97 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
102 typedef FiniteElementInterfaceSwitch<
103 typename LFSU::Traits::FiniteElementType
105 typedef BasisInterfaceSwitch<
106 typename FESwitch::Basis
110 typedef typename BasisSwitch::DomainField DF;
111 typedef typename BasisSwitch::RangeField RF;
112 typedef typename BasisSwitch::Range RangeType;
113 typedef typename LFSU::Traits::SizeType size_type;
116 const int dim = EG::Geometry::mydimension;
119 const auto& geometry = eg.geometry();
120 Dune::GeometryType gt = geometry.type();
121 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
124 for (
const auto& qp : rule)
127 std::vector<RangeType> phi(lfsu.size());
128 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(qp.position(),phi);
131 RF factor = _scaling * qp.weight() * geometry.integrationElement(qp.position());
132 for (size_type j=0; j<lfsu.size(); j++)
133 for (size_type i=0; i<lfsu.size(); i++)
134 mat.accumulate(lfsv,i,lfsu,j, phi[j]*phi[i]*factor);
140 const double _scaling;
162 : scalar_operator(intorder_)
166 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
167 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
169 for(
unsigned int i=0; i<LFSU::CHILDREN; ++i)
170 scalar_operator.
alpha_volume(eg,lfsu.child(i),x,lfsv.child(i),r);
174 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
178 for(
unsigned int i=0; i<LFSU::CHILDREN; ++i)
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
L2(int intorder_=2, double scaling=1.0)
Definition: l2.hh:45
PowerL2(int intorder_=2)
Definition: l2.hh:161
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:167
Definition: adaptivity.hh:27
static const int dim
Definition: adaptivity.hh:83
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:52
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:175
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:98
const EG & eg
Definition: constraints.hh:280
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89