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