dune-pdelab  2.5-dev
eval.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
3 #define DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
4 
5 #include <vector>
6 
7 #include <dune/common/fmatrix.hh>
8 
9 namespace Dune {
10  namespace PDELab {
11 
16  template< class DomainType,
17  class LFS,
18  class LV,
19  class RangeFieldType>
20  inline void evalFunction( const DomainType& location, // expects element local coordinates!
21  const LFS& lfs,
22  const LV& xlocal,
23  RangeFieldType &valU
24  ){
25 
26  typedef typename LFS::Traits::FiniteElementType::
27  Traits::LocalBasisType::Traits::RangeType RangeType; // normally = Dune::FieldVector<double,1>
28 
29  std::vector<RangeType> phi(lfs.size());
30  lfs.finiteElement().localBasis().evaluateFunction( location, phi );
31 
32  // evaluate u
33  valU = RangeFieldType(0.0);
34  for( std::size_t i=0; i<lfs.size(); i++ )
35  valU += xlocal( lfs, i ) * phi[i];
36  }
37 
38 
39 
40 
41  template< class DomainType,
42  class EG,
43  class LFS,
44  class LV,
45  class RangeType
46  >
47  inline void evalGradient( const DomainType& location, // expects element local coordinates!
48  const EG& eg,
49  const LFS& lfs,
50  const LV& xlocal,
51  RangeType &gradu
52  ){
53 
54  typedef typename LFS::Traits::FiniteElementType::
55  Traits::LocalBasisType::Traits::JacobianType JacobianType;
56 
57  typedef typename LFS::Traits::SizeType size_type;
58 
59  std::vector<JacobianType> js(lfs.size());
60  lfs.finiteElement().localBasis().evaluateJacobian( location, js );
61 
62  enum { dim = RangeType::dimension };
63  // transformation
64  typename EG::Geometry::JacobianInverseTransposed jac;
65 
66  // transform gradients of shape functions to real element
67  jac = eg.geometry().jacobianInverseTransposed( location );
68  std::vector<RangeType> gradphi(lfs.size());
69  for (size_type i=0; i<lfs.size(); i++)
70  jac.mv( js[i][0], gradphi[i] );
71 
72  // compute gradient of u
73  gradu = RangeType(0);
74  for (size_type i=0; i<lfs.size(); i++)
75  gradu.axpy( xlocal(lfs,i), gradphi[i] );
76 
77  }
78 
79 
80  }
81 }
82 #endif // DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
static const int dim
Definition: adaptivity.hh:83
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void evalFunction(const DomainType &location, const LFS &lfs, const LV &xlocal, RangeFieldType &valU)
Definition: eval.hh:20
void evalGradient(const DomainType &location, const EG &eg, const LFS &lfs, const LV &xlocal, RangeType &gradu)
Definition: eval.hh:47