dune-pdelab  2.4.1
darcyfem.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH
3 #define DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH
4 
5 #include <dune/common/fvector.hh>
6 #include <dune/geometry/referenceelements.hh>
7 
12 
14 
22 template<typename P, typename T, typename X>
25  Dune::PDELab::GridFunctionTraits<
26  typename T::Traits::GridViewType,
27  typename T::Traits::FiniteElementType::Traits::LocalBasisType
28  ::Traits::RangeFieldType,
29  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
30  ::dimDomain,
31  Dune::FieldVector<
32  typename T::Traits::FiniteElementType::Traits
33  ::LocalBasisType::Traits::RangeFieldType,
34  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
35  ::dimDomain> >,
36  DarcyVelocityFromHeadFEM<P,T,X> >
37 {
38  using GFS = T;
39  using LBTraits = typename GFS::Traits::FiniteElementType::
40  Traits::LocalBasisType::Traits;
41 
44  using LView = typename X::template LocalView<LFSCache>;
45 
46 public:
48  typename GFS::Traits::GridViewType,
49  typename LBTraits::RangeFieldType,
50  LBTraits::dimDomain,
51  Dune::FieldVector<
52  typename LBTraits::RangeFieldType,
53  LBTraits::dimDomain> >;
54 
55 private:
57  Traits,
59 
60 public:
66  DarcyVelocityFromHeadFEM (const P& p, const GFS& gfs, X& x_)
67  : pgfs(stackobject_to_shared_ptr(gfs))
68  , pxg(stackobject_to_shared_ptr(x_))
69  , pp(stackobject_to_shared_ptr(p))
70  , lfs(pgfs)
71  , lfs_cache(lfs)
72  , lview(*pxg)
73  {}
74 
75  // Evaluate
76  inline void evaluate (const typename Traits::ElementType& e,
77  const typename Traits::DomainType& x,
78  typename Traits::RangeType& y) const
79  {
80  // get and bind local functions space
81  lfs.bind(e);
82  lfs_cache.update();
83 
84  // get local coefficients
85  std::vector<typename Traits::RangeFieldType> xl(lfs.size());
86  lview.bind(lfs_cache);
87  lview.read(xl);
88  lview.unbind();
89 
90  // get geometry
91  auto geo = e.geometry();
92 
93  // get Jacobian of geometry
94  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
95  JgeoIT(geo.jacobianInverseTransposed(x));
96 
97  // get local Jacobians/gradients of the shape functions
98  std::vector<typename LBTraits::JacobianType> J(lfs.size());
99  lfs.finiteElement().localBasis().evaluateJacobian(x,J);
100 
101  typename Traits::RangeType gradphi;
102  typename Traits::RangeType minusgrad(0);
103  for(unsigned int i = 0; i < lfs.size(); ++i) {
104  // compute global gradient of shape function i
105  gradphi = 0;
106  JgeoIT.umv(J[i][0], gradphi);
107 
108  // sum up global gradients, weighting them with the appropriate coeff
109  minusgrad.axpy(-xl[i], gradphi);
110  }
111 
112  // multiply with permeability tensor
113  auto inside_cell_center_local = referenceElement(geo).position(0,0);
114  typename P::Traits::PermTensorType A(pp->A(e,inside_cell_center_local));
115  A.mv(minusgrad,y);
116  }
117 
119  inline const typename Traits::GridViewType& getGridView () const
120  {
121  return pgfs->gridView();
122  }
123 
124 private:
125  std::shared_ptr<const GFS> pgfs;
126  std::shared_ptr<X> pxg;
127  std::shared_ptr<const P> pp;
128  mutable LFS lfs;
129  mutable LFSCache lfs_cache;
130  mutable LView lview;
131 };
132 
133 #endif
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
ReferenceElementWrapper< ReferenceElement< typename Geometry::ctype, Geometry::mydimension > > referenceElement(const Geometry &geo)
Returns the reference element for the given geometry.
Definition: referenceelements.hh:144
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
DarcyVelocityFromHeadFEM(const P &p, const GFS &gfs, X &x_)
Construct a DarcyVelocityFromHeadFEM.
Definition: darcyfem.hh:66
const Entity & e
Definition: localfunctionspace.hh:111
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
Traits::IndexContainer::size_type size() const
get current size
Definition: localfunctionspace.hh:206
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: darcyfem.hh:76
const P & p
Definition: constraints.hh:147
Provide Darcy velocity as a vector-valued grid function.
Definition: darcyfem.hh:23
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: darcyfem.hh:119
void update()
Definition: lfsindexcache.hh:300
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
Dune::PDELab::GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, Dune::FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: darcyfem.hh:53
traits class holding the function signature, same as in local function
Definition: function.hh:175