3 #ifndef DUNE_PDELAB_FUNCTION_DISCRETEGRIDVIEWFUNCTION_HH 4 #define DUNE_PDELAB_FUNCTION_DISCRETEGRIDVIEWFUNCTION_HH 10 #include <type_traits> 12 #include <dune/common/exceptions.hh> 13 #include <dune/common/fvector.hh> 15 #include <dune/localfunctions/common/interfaceswitch.hh> 22 #include <dune/functions/gridfunctions/gridviewfunction.hh> 26 template<
typename T,
int N,
typename R2>
27 struct common_type<
Dune::FieldVector<T,N>, R2>
29 using type = Dune::FieldVector<typename std::common_type<T,R2>::type,N>;
32 template<
typename T,
int N,
typename R2>
33 struct common_type<
Dune::FieldVector<T,N>, Dune::FieldVector<R2,N>>
35 using type = Dune::FieldVector<typename std::common_type<T,R2>::type,N>;
42 template<
typename Signature,
typename E,
template<
class>
class D,
int B,
45 Functions::Imp::GridFunctionTraits<
46 typename DiscreteGridViewFunctionTraits<Signature,E,D,B,diffOrder-1>::DerivativeSignature
50 template<
typename Signature,
typename E,
template<
class>
class D,
int B>
52 Functions::Imp::GridFunctionTraits<Signature,E,D,B>
69 template<
typename GFS,
typename V,
int diffOrder = 0>
73 using GridView =
typename GFS::Traits::GridView;
74 using EntitySet = Functions::GridViewEntitySet<GridView, 0>;
76 using Domain =
typename EntitySet::GlobalCoordinate;
77 using LocalBasisTraits =
typename GFS::Traits::FiniteElementMap::Traits::FiniteElement::Traits::LocalBasisType::Traits;
80 using ElementaryRange =
typename std::common_type<LocalBasisRange, VectorRange>::type;
83 using Element =
typename EntitySet::Element;
87 Functions::DefaultDerivativeTraits, 16, diffOrder>;
89 using Range =
typename Traits::Range;
99 using XView =
typename Vector::template ConstLocalView<LFSCache>;
109 LocalFunction(
const shared_ptr<const GridFunctionSpace> gfs,
const shared_ptr<const Vector> v)
115 , xl_(pgfs_->maxLocalSize())
116 , yb_(pgfs_->maxLocalSize())
131 x_view_.bind(lfs_cache_);
145 DUNE_THROW(InvalidStateException,
"can't call localContext on unbound DiscreteGridViewFunction::LocalFunction");
186 return evaluate<diffOrder>(coord);
191 typename std::enable_if<(dOrder > 2),
193 evaluate(
const Domain& coord)
const 195 if (diffOrder > 2) DUNE_THROW(NotImplemented,
196 "Derivatives are only implemented up to degree 2");
200 typename std::enable_if<dOrder == 0,
202 evaluate(
const Domain& coord)
const 205 auto& basis = lfs_.finiteElement().localBasis();
206 basis.evaluateFunction(coord,yb_);
207 for (
size_type i = 0; i < yb_.size(); ++i)
209 r.axpy(xl_[i],yb_[i]);
215 typename std::enable_if<dOrder == 1,
217 evaluate(
const Domain& coord)
const 221 const typename Element::Geometry::JacobianInverseTransposed
222 JgeoIT = element_->geometry().jacobianInverseTransposed(coord);
225 lfs_.finiteElement().localBasis().evaluateJacobian(coord,yb_);
229 for(std::size_t i = 0; i < yb_.size(); ++i) {
230 assert(gradphi.size() == yb_[i].size());
231 for(std::size_t j = 0; j < gradphi.size(); ++j) {
234 JgeoIT.mv(yb_[i][j], gradphi[j]);
239 r[j].axpy(xl_[i], gradphi[j]);
246 typename std::enable_if<dOrder == 2,
248 evaluate(
const Domain& coord)
const 252 if (! element_->geometry().affine())
253 DUNE_THROW(NotImplemented,
"Due to missing features in the Geometry interface, " 254 "the computation of higher derivatives (>=2) works only for affine transformations.");
256 const typename Element::Geometry::JacobianInverseTransposed
257 JgeoIT = element_->geometry().jacobianInverseTransposed(coord);
261 static const unsigned int dim = GridView::dimensionworld;
267 std::array<std::size_t, dim> directions;
268 for(std::size_t i = 0; i <
dim; ++i) {
273 for(std::size_t j = i+1; j <
dim; ++j) {
275 lfs_.finiteElement().localBasis().partial(directions,coord,yb_);
276 assert( yb_.size() == 1);
277 for(std::size_t n = 0; n < yb_.size(); ++n) {
279 r[i][j] += xl_[i] * yb_[j];
288 for(std::size_t i = 0; i <
dim; ++i)
289 for(std::size_t j = i; j <
dim; ++j)
290 r[i][j] *= JgeoIT[i][j] * JgeoIT[i][j];
296 const shared_ptr<const GridFunctionSpace>
pgfs_;
297 const shared_ptr<const Vector>
v_;
301 mutable std::vector<ElementaryRange>
xl_;
302 mutable std::vector<Range>
yb_;
307 : pgfs_(stackobject_to_shared_ptr(gfs)),v_(stackobject_to_shared_ptr(v))
332 DUNE_THROW(NotImplemented,
"not implemented");
359 return pgfs_->gridView();
364 const shared_ptr<const GridFunctionSpace> pgfs_;
365 const shared_ptr<const Vector> v_;
372 #endif // DUNE_PDELAB_FUNCTION_DISCRETEGRIDVIEWFUNCTION_HH const GridFunctionSpace & gridFunctionSpace() const
Definition: discretegridviewfunction.hh:319
static const int dim
Definition: adaptivity.hh:84
const Element * element_
Definition: discretegridviewfunction.hh:303
LocalFunction(const shared_ptr< const GridFunctionSpace > gfs, const shared_ptr< const Vector > v)
Definition: discretegridviewfunction.hh:109
GlobalFunction::Range Range
Definition: discretegridviewfunction.hh:105
LocalDomain Domain
Definition: discretegridviewfunction.hh:104
void unbind()
Definition: discretegridviewfunction.hh:136
Dune::FieldVector< typename std::common_type< T, R2 >::type, N > type
Definition: discretegridviewfunction.hh:35
const shared_ptr< const GridFunctionSpace > pgfs_
Definition: discretegridviewfunction.hh:296
const V & dofs() const
Definition: discretegridviewfunction.hh:324
friend DiscreteGridViewFunction< GFS, V, diffOrder+1 >::LocalFunction derivative(const LocalFunction &t)
free function to obtain the derivative of a LocalFunction
Definition: discretegridviewfunction.hh:165
DiscreteGridViewFunction(std::shared_ptr< const GridFunctionSpace > pgfs, std::shared_ptr< const Vector > v)
Definition: discretegridviewfunction.hh:310
friend LocalFunction localFunction(const DiscreteGridViewFunction &t)
Get local function of wrapped function.
Definition: discretegridviewfunction.hh:349
const Basis & basis() const
Definition: discretegridviewfunction.hh:315
typename GFS::Traits::GridView GridView
Definition: discretegridviewfunction.hh:73
friend DiscreteGridViewFunction< GFS, V, diffOrder+1 > derivative(const DiscreteGridViewFunction &t)
Definition: discretegridviewfunction.hh:335
LFS lfs_
Definition: discretegridviewfunction.hh:298
typename EntitySet::Element Element
Definition: discretegridviewfunction.hh:83
EntitySet entitySet() const
Get associated EntitySet.
Definition: discretegridviewfunction.hh:357
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
void bind(const Element &element)
Bind LocalFunction to grid element.
Definition: discretegridviewfunction.hh:126
GlobalFunction::Element Element
Definition: discretegridviewfunction.hh:106
A discrete function defined over a GridFunctionSpace.
Definition: discretegridviewfunction.hh:70
const Element & localContext() const
Definition: discretegridviewfunction.hh:141
typename Traits::Range Range
Definition: discretegridviewfunction.hh:89
XView x_view_
Definition: discretegridviewfunction.hh:300
std::size_t size_type
Definition: discretegridviewfunction.hh:107
std::vector< ElementaryRange > xl_
Definition: discretegridviewfunction.hh:301
typename EntitySet::GlobalCoordinate Domain
Definition: discretegridviewfunction.hh:76
typename LocalBasisTraits::RangeType LocalBasisRange
Definition: discretegridviewfunction.hh:78
std::vector< Range > yb_
Definition: discretegridviewfunction.hh:302
DiscreteGridViewFunction(const GridFunctionSpace &gfs, const Vector &v)
Definition: discretegridviewfunction.hh:306
GFS Basis
Definition: discretegridviewfunction.hh:91
GFS GridFunctionSpace
Definition: discretegridviewfunction.hh:92
Functions::GridViewEntitySet< GridView, 0 > EntitySet
Definition: discretegridviewfunction.hh:74
typename EntitySet::LocalCoordinate LocalDomain
Definition: discretegridviewfunction.hh:82
LFSCache lfs_cache_
Definition: discretegridviewfunction.hh:299
const shared_ptr< const Vector > v_
Definition: discretegridviewfunction.hh:297
typename V::ElementType VectorRange
Definition: discretegridviewfunction.hh:79
Range operator()(const Domain &coord)
Evaluate LocalFunction at bound element.
Definition: discretegridviewfunction.hh:184
typename GFS::Traits::FiniteElementMap::Traits::FiniteElement::Traits::LocalBasisType::Traits LocalBasisTraits
Definition: discretegridviewfunction.hh:77
Dune::FieldVector< typename std::common_type< T, R2 >::type, N > type
Definition: discretegridviewfunction.hh:29
V Vector
Definition: discretegridviewfunction.hh:93
typename std::common_type< LocalBasisRange, VectorRange >::type ElementaryRange
Definition: discretegridviewfunction.hh:80
Definition: discretegridviewfunction.hh:95
Definition: discretegridviewfunction.hh:44