dune-pdelab  2.0.0
l2volumefunctional.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_LOCALOPERATOR_L2VOLUMEFUNCTIONAL_HH
4 #define DUNE_PDELAB_LOCALOPERATOR_L2VOLUMEFUNCTIONAL_HH
5 
6 #include <vector>
7 
8 #include <dune/common/exceptions.hh>
9 #include <dune/common/fvector.hh>
10 
11 #include <dune/geometry/type.hh>
12 #include <dune/geometry/quadraturerules.hh>
13 
14 #include <dune/localfunctions/common/interfaceswitch.hh>
15 
19 
20 namespace Dune {
21  namespace PDELab {
25 
35  template<typename F>
38  {
39  public:
40  // residual assembly flags
41  enum { doLambdaVolume = true };
42 
47  L2VolumeFunctional (const F& f, unsigned int quadOrder)
48  : f_(f), quadOrder_(quadOrder)
49  {}
50 
51  // volume integral depending only on test functions
52  template<typename EG, typename LFSV, typename R>
53  void lambda_volume (const EG& eg, const LFSV& lfsv, R& r) const
54  {
55  // domain and range field type
56  typedef FiniteElementInterfaceSwitch<typename LFSV::Traits::FiniteElementType> FESwitch;
57  typedef BasisInterfaceSwitch<typename FESwitch::Basis> BasisSwitch;
58 
59  typedef typename BasisSwitch::DomainField DF;
60  typedef typename BasisSwitch::RangeField RF;
61  typedef typename BasisSwitch::Range Range;
62 
63  // dimensions
64  static const int dimLocal = EG::Geometry::mydimension;
65 
66  // select quadrature rule
67  Dune::GeometryType gt = eg.geometry().type();
68  const Dune::QuadratureRule<DF,dimLocal>& rule =
69  Dune::QuadratureRules<DF,dimLocal>::rule(gt,quadOrder_);
70 
71  // loop over quadrature points
72  for (typename Dune::QuadratureRule<DF,dimLocal>::const_iterator it =
73  rule.begin(); it!=rule.end(); ++it)
74  {
75  // evaluate shape functions
76  std::vector<Range> phi(lfsv.size());
77  FESwitch::basis(lfsv.finiteElement()).
78  evaluateFunction(it->position(),phi);
79 
80  // evaluate right hand side parameter function
81  typename F::Traits::RangeType y(0.0);
82  f_.evaluate(eg.entity(),it->position(),y);
83 
84  // integrate f
85  RF factor = r.weight() * it->weight() * eg.geometry().integrationElement(it->position());
86  for (size_t i=0; i<lfsv.size(); i++)
87  r.rawAccumulate(lfsv,i,y*phi[i]*factor);
88  }
89  }
90 
91  protected:
92  const F& f_;
93 
94  // Quadrature rule order
95  unsigned int quadOrder_;
96  };
97 
99  } // namespace PDELab
100 } // namespace Dune
101 
102 #endif
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: l2volumefunctional.hh:53
Definition: l2volumefunctional.hh:41
L2VolumeFunctional(const F &f, unsigned int quadOrder)
Constructor.
Definition: l2volumefunctional.hh:47
Default flags for all local operators.
Definition: flags.hh:18
unsigned int quadOrder_
Definition: l2volumefunctional.hh:95
A local operator that tests a function against a test function space defined on the entire grid...
Definition: l2volumefunctional.hh:36
const F & f_
Definition: l2volumefunctional.hh:92
const F & f
Definition: common/constraints.hh:145
const EG & eg
Definition: common/constraints.hh:277