2 #ifndef DUNE_PDELAB_L2_HH
3 #define DUNE_PDELAB_L2_HH
7 #include<dune/common/exceptions.hh>
8 #include<dune/common/fvector.hh>
9 #include<dune/common/static_assert.hh>
11 #include<dune/geometry/type.hh>
12 #include<dune/geometry/referenceelements.hh>
13 #include<dune/geometry/quadraturerules.hh>
15 #include <dune/localfunctions/common/interfaceswitch.hh>
46 L2 (
int intorder_=2,
double scaling=1.0)
52 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
53 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
56 typedef FiniteElementInterfaceSwitch<
57 typename LFSU::Traits::FiniteElementType
59 typedef BasisInterfaceSwitch<
60 typename FESwitch::Basis
64 typedef typename BasisSwitch::DomainField DF;
65 typedef typename BasisSwitch::RangeField RF;
66 typedef typename BasisSwitch::Range RangeType;
68 typedef typename LFSU::Traits::SizeType size_type;
71 const int dim = EG::Geometry::dimension;
74 Dune::GeometryType gt = eg.geometry().type();
75 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
78 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
81 std::vector<RangeType> phi(lfsu.size());
82 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(it->position(),phi);
86 for (size_type i=0; i<lfsu.size(); i++)
87 u += RF(x(lfsu,i)*phi[i]);
90 RF factor = _scaling * it->weight() * eg.geometry().integrationElement(it->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::dimension;
119 Dune::GeometryType gt = eg.geometry().type();
120 const Dune::QuadratureRule<DF,dim>& rule = Dune::QuadratureRules<DF,dim>::rule(gt,intorder);
123 for (
typename Dune::QuadratureRule<DF,dim>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
126 std::vector<RangeType> phi(lfsu.size());
127 FESwitch::basis(lfsu.finiteElement()).evaluateFunction(it->position(),phi);
130 RF factor = _scaling * it->weight() * eg.geometry().integrationElement(it->position());
131 for (size_type j=0; j<lfsu.size(); j++)
132 for (size_type i=0; i<lfsu.size(); i++)
133 mat.accumulate(lfsv,i,lfsu,j, phi[j]*phi[i]*factor);
139 const double _scaling;
161 : scalar_operator(intorder_)
165 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename R>
166 void alpha_volume (
const EG&
eg,
const LFSU& lfsu,
const X& x,
const LFSV& lfsv, R& r)
const
168 for(
unsigned int i=0; i<LFSU::CHILDREN; ++i)
169 scalar_operator.
alpha_volume(eg,lfsu.child(i),x,lfsv.child(i),r);
173 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
177 for(
unsigned int i=0; i<LFSU::CHILDREN; ++i)
L2(int intorder_=2, double scaling=1.0)
Definition: l2.hh:46
Default flags for all local operators.
Definition: flags.hh:18
static const int dim
Definition: adaptivity.hh:82
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:174
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:53
PowerL2(int intorder_=2)
Definition: l2.hh:160
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
const EG & eg
Definition: common/constraints.hh:277
void alpha_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
Definition: l2.hh:166
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: l2.hh:98