2 #ifndef DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
3 #define DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
7 #include <dune/common/fmatrix.hh>
16 template<
class DomainType,
26 typedef typename LFS::Traits::FiniteElementType::
27 Traits::LocalBasisType::Traits::RangeType RangeType;
29 std::vector<RangeType> phi(lfs.size());
30 lfs.finiteElement().localBasis().evaluateFunction( location, phi );
33 valU = RangeFieldType(0.0);
34 for( std::size_t i=0; i<lfs.size(); i++ )
35 valU += xlocal( lfs, i ) * phi[i];
41 template<
class DomainType,
54 typedef typename LFS::Traits::FiniteElementType::
55 Traits::LocalBasisType::Traits::JacobianType JacobianType;
57 typedef typename LFS::Traits::SizeType size_type;
59 std::vector<JacobianType> js(lfs.size());
60 lfs.finiteElement().localBasis().evaluateJacobian( location, js );
62 enum {
dim = RangeType::dimension };
64 typename EG::Geometry::JacobianInverseTransposed jac;
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] );
74 for (size_type i=0; i<lfs.size(); i++)
75 gradu.axpy( xlocal(lfs,i), gradphi[i] );
82 #endif // DUNE_PDELAB_LOCALOPERATOR_EVAL_HH
void evalGradient(const DomainType &location, const EG &eg, const LFS &lfs, const LV &xlocal, RangeType &gradu)
Definition: eval.hh:47
static const int dim
Definition: adaptivity.hh:82
void evalFunction(const DomainType &location, const LFS &lfs, const LV &xlocal, RangeFieldType &valU)
Definition: eval.hh:20
const EG & eg
Definition: common/constraints.hh:277