dune-pdelab  2.4-dev
darcy_FEM.hh
Go to the documentation of this file.
1 #ifndef DarcyVelocityFromHeadFEM_HH
2 #define DarcyVelocityFromHeadFEM_HH
3 
5 
13 template<typename P, typename T, typename X>
16  Dune::PDELab::GridFunctionTraits<
17  typename T::Traits::GridViewType,
18  typename T::Traits::FiniteElementType::Traits::LocalBasisType
19  ::Traits::RangeFieldType,
20  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
21  ::dimDomain,
22  Dune::FieldVector<
23  typename T::Traits::FiniteElementType::Traits
24  ::LocalBasisType::Traits::RangeFieldType,
25  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
26  ::dimDomain> >,
27  DarcyVelocityFromHeadFEM<P,T,X> >
28 {
29  typedef T GFS;
30  typedef typename GFS::Traits::FiniteElementType::Traits::
31  LocalBasisType::Traits LBTraits;
32 
35  typedef typename X::template LocalView<LFSCache> LView;
36 
37 public:
39  typename GFS::Traits::GridViewType,
40  typename LBTraits::RangeFieldType,
41  LBTraits::dimDomain,
42  Dune::FieldVector<
43  typename LBTraits::RangeFieldType,
44  LBTraits::dimDomain> > Traits;
45 
46 private:
48  Traits,
50 
51 public:
57  DarcyVelocityFromHeadFEM (const P& p, const GFS& gfs, X& x_)
58  : pp(stackobject_to_shared_ptr(p)),
59  pgfs(stackobject_to_shared_ptr(gfs)),
60  pxg(stackobject_to_shared_ptr(x_)),
61  lfs(pgfs),
62  lfs_cache(lfs),
63  lview(*pxg)
64  {}
65 
66  // Evaluate
67  inline void evaluate (const typename Traits::ElementType& e,
68  const typename Traits::DomainType& x,
69  typename Traits::RangeType& y) const
70  {
71  // get and bind local functions space
72  lfs.bind(e);
73  lfs_cache.update();
74 
75  // get local coefficients
76  std::vector<typename Traits::RangeFieldType> xl(lfs.size());
77  lview.bind(lfs_cache);
78  lview.read(xl);
79  lview.unbind();
80 
81  // get Jacobian of geometry
82  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
83  JgeoIT(e.geometry().jacobianInverseTransposed(x));
84 
85  // get local Jacobians/gradients of the shape functions
86  std::vector<typename LBTraits::JacobianType> J(lfs.size());
87  lfs.finiteElement().localBasis().evaluateJacobian(x,J);
88 
89  typename Traits::RangeType gradphi;
90  typename Traits::RangeType minusgrad(0);
91  for(unsigned int i = 0; i < lfs.size(); ++i) {
92  // compute global gradient of shape function i
93  gradphi = 0;
94  JgeoIT.umv(J[i][0], gradphi);
95 
96  // sum up global gradients, weighting them with the appropriate coeff
97  minusgrad.axpy(-xl[i], gradphi);
98  }
99 
100  // multiply with permeability tensor
101  typedef typename Traits::DomainFieldType DF;
102  const int dim = LBTraits::dimDomain;
103  const Dune::FieldVector<DF,dim>
104  inside_cell_center_local = Dune::ReferenceElements<DF,dim>::
105  general(e.type()).position(0,0);
106  typename P::Traits::PermTensorType A(pp->A(e,inside_cell_center_local));
107  A.mv(minusgrad,y);
108  }
109 
111  inline const typename Traits::GridViewType& getGridView () const
112  {
113  return pgfs->gridView();
114  }
115 
116 private:
117  Dune::shared_ptr<const GFS> pgfs;
118  Dune::shared_ptr<X> pxg;
119  Dune::shared_ptr<const P> pp;
120  mutable LFS lfs;
121  mutable LFSCache lfs_cache;
122  mutable LView lview;
123 };
124 
125 #endif
Dune::PDELab::GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, Dune::FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: darcy_FEM.hh:44
Provide Darcy velocity as a vector-valued grid function.
Definition: darcy_FEM.hh:14
const E & e
Definition: interpolate.hh:172
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: darcy_FEM.hh:67
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
Traits::IndexContainer::size_type size() const
get current size
Definition: localfunctionspace.hh:206
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: darcy_FEM.hh:111
DarcyVelocityFromHeadFEM(const P &p, const GFS &gfs, X &x_)
Construct a DarcyVelocityFromHeadFEM.
Definition: darcy_FEM.hh:57
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
void update()
Definition: lfsindexcache.hh:300
traits class holding the function signature, same as in local function
Definition: function.hh:175
static const int dim
Definition: adaptivity.hh:83
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
const P & p
Definition: constraints.hh:146
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186