3#ifndef DUNE_PDELAB_FUNCTION_DISCRETEGRIDVIEWFUNCTION_HH
4#define DUNE_PDELAB_FUNCTION_DISCRETEGRIDVIEWFUNCTION_HH
12#include <dune/common/exceptions.hh>
13#include <dune/common/shared_ptr.hh>
14#include <dune/common/fvector.hh>
16#include <dune/localfunctions/common/interfaceswitch.hh>
23#include <dune/functions/gridfunctions/gridviewfunction.hh>
27 template<
typename T,
int N,
typename R2>
28 struct common_type<
Dune::FieldVector<T,N>, R2>
30 using type = Dune::FieldVector<typename std::common_type<T,R2>::type,N>;
33 template<
typename T,
int N,
typename R2>
34 struct common_type<
Dune::FieldVector<T,N>, Dune::FieldVector<R2,N>>
36 using type = Dune::FieldVector<typename std::common_type<T,R2>::type,N>;
43template<
typename Signature,
typename E,
template<
class>
class D,
int B,
46 Functions::Imp::GridFunctionTraits<
47 typename DiscreteGridViewFunctionTraits<Signature,E,D,B,diffOrder-1>::DerivativeSignature
51template<
typename Signature,
typename E,
template<
class>
class D,
int B>
53 Functions::Imp::GridFunctionTraits<Signature,E,D,B>
70template<
typename GFS,
typename V,
int diffOrder = 0>
74 using GridView =
typename GFS::Traits::GridView;
75 using EntitySet = Functions::GridViewEntitySet<GridView, 0>;
77 using Domain =
typename EntitySet::GlobalCoordinate;
78 using LocalBasisTraits =
typename GFS::Traits::FiniteElementMap::Traits::FiniteElement::Traits::LocalBasisType::Traits;
81 using ElementaryRange =
typename std::common_type<LocalBasisRange, VectorRange>::type;
84 using Element =
typename EntitySet::Element;
88 Functions::DefaultDerivativeTraits, 16, diffOrder>;
90 using Range =
typename Traits::Range;
100 using XView =
typename Vector::template ConstLocalView<LFSCache>;
110 LocalFunction(
const std::shared_ptr<const GridFunctionSpace> gfs,
const std::shared_ptr<const Vector> v)
146 DUNE_THROW(InvalidStateException,
"can't call localContext on unbound DiscreteGridViewFunction::LocalFunction");
187 return evaluate<diffOrder>(coord);
192 typename std::enable_if<(dOrder > 2),
194 evaluate(
const Domain& coord)
const
196 if (diffOrder > 2) DUNE_THROW(NotImplemented,
197 "Derivatives are only implemented up to degree 2");
201 typename std::enable_if<dOrder == 0,
203 evaluate(
const Domain& coord)
const
206 auto&
basis =
lfs_.finiteElement().localBasis();
216 typename std::enable_if<dOrder == 1,
218 evaluate(
const Domain& coord)
const
222 const typename Element::Geometry::JacobianInverseTransposed
223 JgeoIT =
element_->geometry().jacobianInverseTransposed(coord);
226 lfs_.finiteElement().localBasis().evaluateJacobian(coord,
yb_);
230 for(std::size_t i = 0; i <
yb_.size(); ++i) {
231 assert(gradphi.size() ==
yb_[i].size());
232 for(std::size_t j = 0; j < gradphi.size(); ++j) {
235 JgeoIT.mv(
yb_[i][j], gradphi[j]);
240 r[j].axpy(
xl_[i], gradphi[j]);
247 typename std::enable_if<dOrder == 2,
249 evaluate(
const Domain& coord)
const
253 if (!
element_->geometry().affine())
254 DUNE_THROW(NotImplemented,
"Due to missing features in the Geometry interface, "
255 "the computation of higher derivatives (>=2) works only for affine transformations.");
257 const typename Element::Geometry::JacobianInverseTransposed
258 JgeoIT =
element_->geometry().jacobianInverseTransposed(coord);
262 static const unsigned int dim = GridView::dimensionworld;
268 std::array<std::size_t, dim> directions;
269 for(std::size_t i = 0; i <
dim; ++i) {
274 for(std::size_t j = i+1; j <
dim; ++j) {
276 lfs_.finiteElement().localBasis().partial(directions,coord,
yb_);
277 assert(
yb_.size() == 1);
278 for(std::size_t n = 0; n <
yb_.size(); ++n) {
280 r[i][j] +=
xl_[i] *
yb_[j];
289 for(std::size_t i = 0; i <
dim; ++i)
290 for(std::size_t j = i; j <
dim; ++j)
291 r[i][j] *= JgeoIT[i][j] * JgeoIT[i][j];
297 const std::shared_ptr<const GridFunctionSpace>
pgfs_;
298 const std::shared_ptr<const Vector>
v_;
302 mutable std::vector<ElementaryRange>
xl_;
303 mutable std::vector<Range>
yb_;
308 : pgfs_(stackobject_to_shared_ptr(gfs)),v_(stackobject_to_shared_ptr(v))
345 DUNE_THROW(NotImplemented,
"not implemented");
372 return pgfs_->gridView();
377 const std::shared_ptr<const GridFunctionSpace> pgfs_;
378 const std::shared_ptr<const Vector> v_;
static const int dim
Definition: adaptivity.hh:84
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::FieldVector< typename std::common_type< T, R2 >::type, N > type
Definition: discretegridviewfunction.hh:30
Dune::FieldVector< typename std::common_type< T, R2 >::type, N > type
Definition: discretegridviewfunction.hh:36
Definition: discretegridviewfunction.hh:49
A discrete function defined over a GridFunctionSpace.
Definition: discretegridviewfunction.hh:72
friend DiscreteGridViewFunction< GFS, V, diffOrder+1 > derivative(const DiscreteGridViewFunction &t)
Definition: discretegridviewfunction.hh:348
const GridFunctionSpace & gridFunctionSpace() const
Definition: discretegridviewfunction.hh:320
DiscreteGridViewFunction(const GridFunctionSpace &gfs, const Vector &v)
Definition: discretegridviewfunction.hh:307
typename Traits::Range Range
Definition: discretegridviewfunction.hh:90
Functions::GridViewEntitySet< GridView, 0 > EntitySet
Definition: discretegridviewfunction.hh:75
typename EntitySet::GlobalCoordinate Domain
Definition: discretegridviewfunction.hh:77
typename GFS::Traits::GridView GridView
Definition: discretegridviewfunction.hh:74
typename EntitySet::LocalCoordinate LocalDomain
Definition: discretegridviewfunction.hh:83
typename std::common_type< LocalBasisRange, VectorRange >::type ElementaryRange
Definition: discretegridviewfunction.hh:81
typename V::ElementType VectorRange
Definition: discretegridviewfunction.hh:80
const V & dofs() const
Definition: discretegridviewfunction.hh:325
GFS Basis
Definition: discretegridviewfunction.hh:92
typename EntitySet::Element Element
Definition: discretegridviewfunction.hh:84
GFS GridFunctionSpace
Definition: discretegridviewfunction.hh:93
const Basis & basis() const
Definition: discretegridviewfunction.hh:316
EntitySet entitySet() const
Get associated EntitySet.
Definition: discretegridviewfunction.hh:370
DiscreteGridViewFunction(std::shared_ptr< const GridFunctionSpace > pgfs, std::shared_ptr< const Vector > v)
Definition: discretegridviewfunction.hh:311
auto dofsStorage() const
Returns storage object of the dof storage vector.
Definition: discretegridviewfunction.hh:337
Range operator()(const Domain &x) const
Definition: discretegridviewfunction.hh:343
friend LocalFunction localFunction(const DiscreteGridViewFunction &t)
Get local function of wrapped function.
Definition: discretegridviewfunction.hh:362
typename GFS::Traits::FiniteElementMap::Traits::FiniteElement::Traits::LocalBasisType::Traits LocalBasisTraits
Definition: discretegridviewfunction.hh:78
auto gridFunctionSpaceStorage() const
Returns storage object of the grid function space.
Definition: discretegridviewfunction.hh:331
V Vector
Definition: discretegridviewfunction.hh:94
typename LocalBasisTraits::RangeType LocalBasisRange
Definition: discretegridviewfunction.hh:79
Definition: discretegridviewfunction.hh:97
LocalFunction(const std::shared_ptr< const GridFunctionSpace > gfs, const std::shared_ptr< const Vector > v)
Definition: discretegridviewfunction.hh:110
XView x_view_
Definition: discretegridviewfunction.hh:301
const std::shared_ptr< const GridFunctionSpace > pgfs_
Definition: discretegridviewfunction.hh:297
LFS lfs_
Definition: discretegridviewfunction.hh:299
std::size_t size_type
Definition: discretegridviewfunction.hh:108
void bind(const Element &element)
Bind LocalFunction to grid element.
Definition: discretegridviewfunction.hh:127
GlobalFunction::Element Element
Definition: discretegridviewfunction.hh:107
const Element & localContext() const
Definition: discretegridviewfunction.hh:142
const std::shared_ptr< const Vector > v_
Definition: discretegridviewfunction.hh:298
void unbind()
Definition: discretegridviewfunction.hh:137
std::vector< ElementaryRange > xl_
Definition: discretegridviewfunction.hh:302
LocalDomain Domain
Definition: discretegridviewfunction.hh:105
GlobalFunction::Range Range
Definition: discretegridviewfunction.hh:106
friend DiscreteGridViewFunction< GFS, V, diffOrder+1 >::LocalFunction derivative(const LocalFunction &t)
free function to obtain the derivative of a LocalFunction
Definition: discretegridviewfunction.hh:166
const Element * element_
Definition: discretegridviewfunction.hh:304
LFSCache lfs_cache_
Definition: discretegridviewfunction.hh:300
Range operator()(const Domain &coord)
Evaluate LocalFunction at bound element.
Definition: discretegridviewfunction.hh:185
std::vector< Range > yb_
Definition: discretegridviewfunction.hh:303
void update()
Definition: lfsindexcache.hh:304
Definition: lfsindexcache.hh:979