dune-pdelab  2.5-dev
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 
20 
21 namespace Dune {
22  namespace PDELab {
26 
36  template<typename F>
39  {
40  public:
41  // residual assembly flags
42  enum { doLambdaVolume = true };
43 
49  L2VolumeFunctional (const F& f, unsigned int quadOrder)
50  : f_(f), quadOrder_(quadOrder)
51  {}
52 
53  // volume integral depending only on test functions
54  template<typename EG, typename LFSV, typename R>
55  void lambda_volume (const EG& eg, const LFSV& lfsv, R& r) const
56  {
57  // Define types
58  using FESwitch = FiniteElementInterfaceSwitch<typename LFSV::Traits::FiniteElementType>;
59  using BasisSwitch = BasisInterfaceSwitch<typename FESwitch::Basis>;
60  using Range = typename BasisSwitch::Range;
61 
62  // Reference to cell
63  const auto& cell = eg.entity();
64 
65  // get geometry
66  auto geo = eg.geometry();
67 
68  // Initialize vectors outside for loop
69  std::vector<Range> phi(lfsv.size());
70  typename F::Traits::RangeType y(0.0);
71 
72  // loop over quadrature points
73  for (const auto& ip : quadratureRule(geo,quadOrder_))
74  {
75  // evaluate shape functions
76  FESwitch::basis(lfsv.finiteElement()).
77  evaluateFunction(ip.position(),phi);
78 
79  // evaluate right hand side parameter function
80  f_.evaluate(cell,ip.position(),y);
81 
82  // integrate f
83  auto factor = r.weight() * ip.weight() * geo.integrationElement(ip.position());
84  for (size_t i=0; i<lfsv.size(); i++)
85  r.rawAccumulate(lfsv,i,y*phi[i]*factor);
86  }
87  }
88 
89  protected:
90  const F& f_;
91 
92  // Quadrature rule order
93  unsigned int quadOrder_;
94  };
95 
97  } // namespace PDELab
98 } // namespace Dune
99 
100 #endif // DUNE_PDELAB_LOCALOPERATOR_L2VOLUMEFUNCTIONAL_HH
Dune::PDELab::L2VolumeFunctional
A local operator that tests a function against a test function space defined on the entire grid.
Definition: l2volumefunctional.hh:37
Dune::PDELab::L2VolumeFunctional::lambda_volume
void lambda_volume(const EG &eg, const LFSV &lfsv, R &r) const
Definition: l2volumefunctional.hh:55
flags.hh
Dune::PDELab::LocalOperatorDefaultFlags
Default flags for all local operators.
Definition: flags.hh:18
Dune::PDELab::L2VolumeFunctional::quadOrder_
unsigned int quadOrder_
Definition: l2volumefunctional.hh:93
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
idefault.hh
Dune::PDELab::L2VolumeFunctional::L2VolumeFunctional
L2VolumeFunctional(const F &f, unsigned int quadOrder)
Constructor.
Definition: l2volumefunctional.hh:49
Dune::PDELab::L2VolumeFunctional::doLambdaVolume
@ doLambdaVolume
Definition: l2volumefunctional.hh:42
Dune::PDELab::quadratureRule
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
pattern.hh
quadraturerules.hh
Dune::PDELab::L2VolumeFunctional::f_
const F & f_
Definition: l2volumefunctional.hh:90