dune-pdelab  2.5-dev
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 
11 
13 
21 template<typename P, typename T, typename X>
24  Dune::PDELab::GridFunctionTraits<
25  typename T::Traits::GridViewType,
26  typename T::Traits::FiniteElementType::Traits::LocalBasisType
27  ::Traits::RangeFieldType,
28  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
29  ::dimDomain,
30  Dune::FieldVector<
31  typename T::Traits::FiniteElementType::Traits
32  ::LocalBasisType::Traits::RangeFieldType,
33  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
34  ::dimDomain> >,
35  DarcyVelocityFromHeadFEM<P,T,X> >
36 {
37  using GFS = T;
38  using LBTraits = typename GFS::Traits::FiniteElementType::
39  Traits::LocalBasisType::Traits;
40 
43  using LView = typename X::template LocalView<LFSCache>;
44 
45 public:
47  typename GFS::Traits::GridViewType,
48  typename LBTraits::RangeFieldType,
49  LBTraits::dimDomain,
50  Dune::FieldVector<
51  typename LBTraits::RangeFieldType,
52  LBTraits::dimDomain> >;
53 
54 private:
56  Traits,
58 
59 public:
65  DarcyVelocityFromHeadFEM (const P& p, const GFS& gfs, X& x_)
66  : pgfs(stackobject_to_shared_ptr(gfs))
67  , pxg(stackobject_to_shared_ptr(x_))
68  , pp(stackobject_to_shared_ptr(p))
69  , lfs(pgfs)
70  , lfs_cache(lfs)
71  , lview(*pxg)
72  {}
73 
74  // Evaluate
75  inline void evaluate (const typename Traits::ElementType& e,
76  const typename Traits::DomainType& x,
77  typename Traits::RangeType& y) const
78  {
79  // get and bind local functions space
80  lfs.bind(e);
81  lfs_cache.update();
82 
83  // get local coefficients
84  std::vector<typename Traits::RangeFieldType> xl(lfs.size());
85  lview.bind(lfs_cache);
86  lview.read(xl);
87  lview.unbind();
88 
89  // get geometry
90  auto geo = e.geometry();
91 
92  // get Jacobian of geometry
93  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
94  JgeoIT(geo.jacobianInverseTransposed(x));
95 
96  // get local Jacobians/gradients of the shape functions
97  std::vector<typename LBTraits::JacobianType> J(lfs.size());
98  lfs.finiteElement().localBasis().evaluateJacobian(x,J);
99 
100  typename Traits::RangeType gradphi;
101  typename Traits::RangeType minusgrad(0);
102  for(unsigned int i = 0; i < lfs.size(); ++i) {
103  // compute global gradient of shape function i
104  gradphi = 0;
105  JgeoIT.umv(J[i][0], gradphi);
106 
107  // sum up global gradients, weighting them with the appropriate coeff
108  minusgrad.axpy(-xl[i], gradphi);
109  }
110 
111  // multiply with permeability tensor
112  auto inside_cell_center_local = referenceElement(geo).position(0,0);
113  typename P::Traits::PermTensorType A(pp->A(e,inside_cell_center_local));
114  A.mv(minusgrad,y);
115  }
116 
118  inline const typename Traits::GridViewType& getGridView () const
119  {
120  return pgfs->gridView();
121  }
122 
123 private:
124  std::shared_ptr<const GFS> pgfs;
125  std::shared_ptr<X> pxg;
126  std::shared_ptr<const P> pp;
127  mutable LFS lfs;
128  mutable LFSCache lfs_cache;
129  mutable LView lview;
130 };
131 
132 #endif // DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: darcyfem.hh:75
Dune::PDELab::GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, Dune::FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: darcyfem.hh:52
const Entity & e
Definition: localfunctionspace.hh:120
leaf of a function tree
Definition: function.hh:298
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:49
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:118
traits class holding the function signature, same as in local function
Definition: function.hh:176
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: darcyfem.hh:118
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:115
void update()
Definition: lfsindexcache.hh:304
DarcyVelocityFromHeadFEM(const P &p, const GFS &gfs, X &x_)
Construct a DarcyVelocityFromHeadFEM.
Definition: darcyfem.hh:65
const P & p
Definition: constraints.hh:147
Provide Darcy velocity as a vector-valued grid function.
Definition: darcyfem.hh:22
Traits::IndexContainer::size_type size() const
get current size
Definition: localfunctionspace.hh:220