1#ifndef DUNE_FEM_SPACE_SHAPEFUNCTIONSET_CACHING_HH
2#define DUNE_FEM_SPACE_SHAPEFUNCTIONSET_CACHING_HH
9#include <dune/geometry/quadraturerules.hh>
30 template<
class ShapeFunctionSet >
54 static const int pointSetId = detail::SelectPointSetId< ShapeFunctionSetType >::value;
59 , shapeFunctionSet_( shapeFunctionSet )
68 return shapeFunctionSet_.
order();
73 return shapeFunctionSet_.
size();
76 template<
class Po
int,
class Functor >
82 template<
class Quadrature,
class Functor >
85 const bool cacheable = std::is_convertible< Quadrature, CachingInterface >::value;
89 template<
class Po
int,
class Functor >
95 template<
class Quadrature,
class Functor >
98 const bool cacheable = std::is_convertible< Quadrature, CachingInterface >::value;
102 template<
class Po
int,
class Functor >
105 return shapeFunctionSet_.
hessianEach( x, functor );
108 GeometryType
type ()
const {
return type_; }
110 template <
class QuadratureType >
113 return ReturnCache< QuadratureType, std::is_base_of< CachingInterface, QuadratureType >::value > ::
114 ranges( *
this, quadrature, rangeCaches_, localRangeCache_ );
117 template <
class QuadratureType >
120 return ReturnCache< QuadratureType, std::is_base_of< CachingInterface, QuadratureType >::value > ::
121 jacobians( *
this, quadrature, jacobianCaches_, localJacobianCache_ );
128 template<
class Quad,
bool cacheable >
132 ranges(
const ThisType& shapeFunctionSet,
143 const unsigned int nop = quad.nop();
144 const unsigned int size = shapeFunctionSet.size();
147 storage.resize(
size * nop );
150 for(
unsigned int qp = 0 ; qp < nop; ++ qp )
152 const int cacheQp = quad.cachingPoint( qp );
154 shapeFunctionSet.evaluateEach( quad[ qp ], funztor );
160 jacobians(
const ThisType& shapeFunctionSet,
171 const unsigned int nop = quad.nop();
172 const unsigned int size = shapeFunctionSet.size();
175 storage.resize(
size * nop );
178 for(
unsigned int qp = 0 ; qp < nop; ++ qp )
180 const int cacheQp = quad.cachingPoint( qp );
181 AssignFunctor< JacobianRangeType* > funztor( data + ( cacheQp *
size ) );
182 shapeFunctionSet.jacobianEach( quad[ qp ], funztor );
188 template<
class Quad >
189 struct ReturnCache< Quad, true >
192 ranges(
const ThisType& shapeFunctionSet,
197 return cache[ quad.id() ];
201 jacobians(
const ThisType& shapeFunctionSet,
206 return cache[ quad.id() ];
211 template<
class Quadrature,
class Functor >
212 void evaluateEach (
const Quadrature &quadrature, std::size_t pt, Functor functor,
213 std::integral_constant< bool, false > )
const
218 template<
class Quadrature,
class Functor >
219 void evaluateEach (
const Quadrature &quadrature, std::size_t pt, Functor functor,
220 std::integral_constant< bool, true > )
const;
222 template<
class Quadrature,
class Functor >
223 void jacobianEach (
const Quadrature &quadrature, std::size_t pt, Functor functor,
224 std::integral_constant< bool, false > )
const
229 template<
class Quadrature,
class Functor >
230 void jacobianEach (
const Quadrature &quadrature, std::size_t pt, Functor functor,
231 std::integral_constant< bool, true > )
const;
234 void cacheQuadrature( std::size_t
id, std::size_t codim, std::size_t
size );
236 template<
class Po
intVector >
237 void cachePoints ( std::size_t
id,
const PointVector &points );
240 ShapeFunctionSet shapeFunctionSet_;
253 template<
class ShapeFunctionSet >
260 template<
class ShapeFunctionSet >
261 template<
class Quadrature,
class Functor >
264 std::integral_constant< bool, true > )
const
266 assert( (quadrature.
id() < rangeCaches_.size()) && !rangeCaches_[ quadrature.
id() ].empty() );
267 const RangeType *cache = rangeCaches_[ quadrature.
id() ].data();
269 const unsigned int numShapeFunctions = size();
270 const unsigned int cpt = quadrature.cachingPoint( pt );
278 if constexpr ( quadPointSetId == pointSetId )
281 if( quadrature.isInterpolationQuadrature(numShapeFunctions) )
285 assert( quadPointSetId >= 0 );
286 assert( pointSetId >= 0 );
291 const unsigned int point = quadrature.interpolationPoint( pt );
293 for(
unsigned int i = 0; i < numShapeFunctions; ++i )
294 assert( (cache[ cpt*numShapeFunctions + i ] - RangeType(i==point)).two_norm() < 1e-8 ) ;
297 functor( point, RangeType(1) );
302 for(
unsigned int i = 0; i < numShapeFunctions; ++i )
303 functor( i, cache[ cpt*numShapeFunctions + i ] );
307 template<
class ShapeFunctionSet >
308 template<
class Quadrature,
class Functor >
309 inline void CachingShapeFunctionSet< ShapeFunctionSet >
310 ::jacobianEach (
const Quadrature &quadrature, std::size_t pt, Functor functor,
311 std::integral_constant< bool, true > )
const
313 assert( (quadrature.id() < jacobianCaches_.size()) && !jacobianCaches_[ quadrature.id() ].empty() );
314 const JacobianRangeType *cache = jacobianCaches_[ quadrature.id() ].data();
316 const unsigned int numShapeFunctions = size();
317 const unsigned int cpt = quadrature.cachingPoint( pt );
318 for(
unsigned int i = 0; i < numShapeFunctions; ++i )
319 functor( i, cache[ cpt*numShapeFunctions + i ] );
323 template<
class ShapeFunctionSet >
324 inline void CachingShapeFunctionSet< ShapeFunctionSet >
325 ::cacheQuadrature( std::size_t
id, std::size_t codim, std::size_t size )
329 DUNE_THROW(SingleThreadModeError,
"CachingShapeFunctionSet::cacheQuadrature: only call in single thread mode!");
332 if(
id >= rangeCaches_.size() )
334 rangeCaches_.resize(
id+1, RangeVectorType() );
335 jacobianCaches_.resize(
id+1, JacobianRangeVectorType() );
338 assert( rangeCaches_[
id ].empty() == jacobianCaches_[
id ].empty() );
340 if( rangeCaches_[
id ].empty() )
342 typedef typename FunctionSpaceType::DomainFieldType ctype;
343 const int dim = FunctionSpaceType::dimDomain;
347 cachePoints(
id, PointProvider< ctype, dim, 0 >::getPoints(
id, type_ ) );
351 cachePoints(
id, PointProvider< ctype, dim, 1 >::getPoints(
id, type_ ) );
355 DUNE_THROW( NotImplemented,
"Caching for codim > 1 not implemented." );
361 template<
class ShapeFunctionSet >
362 template<
class Po
intVector >
363 inline void CachingShapeFunctionSet< ShapeFunctionSet >
364 ::cachePoints ( std::size_t
id,
const PointVector &points )
366 const unsigned int numShapeFunctions = size();
367 const unsigned int numPoints = points.size();
369 RangeVectorType& ranges = rangeCaches_[ id ];
370 ranges.resize( numShapeFunctions * numPoints );
372 JacobianRangeVectorType& jacobians = jacobianCaches_[ id ];
373 jacobians.resize( numShapeFunctions * numPoints );
375 if( ranges.empty() || jacobians.empty() )
376 DUNE_THROW( OutOfMemoryError,
"Unable to allocate shape function set caches." );
378 for(
unsigned int pt = 0; pt < numPoints; ++pt )
380 evaluateEach( points[ pt ], AssignFunctor< RangeType * >( ranges.data() + pt*numShapeFunctions ) );
381 jacobianEach( points[ pt ], AssignFunctor< JacobianRangeType * >( jacobians.data() + pt*numShapeFunctions ) );
Definition: bindguard.hh:11
detail::SelectPointSetId< Quadrature, -Dune::QuadratureType::size > SelectQuadraturePointSetId
Definition: quadrature.hh:541
Definition: explicitfieldvector.hh:75
Definition: misc/functor.hh:31
static bool singleThreadMode()
returns true if program is operating on one thread currently
Definition: mpimanager.hh:436
static void unregisterStorage(StorageInterface &storage)
Definition: registry.hh:85
static void registerStorage(StorageInterface &storage)
Definition: registry.hh:70
Definition: registry.hh:30
wrapper for a (Quadrature,int) pair
Definition: quadrature.hh:43
unsigned int index() const
Definition: quadrature.hh:70
const QuadratureType & quadrature() const
Definition: quadrature.hh:65
size_t id() const
obtain the identifier of the integration point list
Definition: quadrature.hh:327
actual interface class for quadratures
Definition: quadrature.hh:405
A vector valued function space.
Definition: functionspace.hh:60
Definition: caching.hh:33
void evaluateEach(const QuadraturePointWrapper< Quadrature > &x, Functor functor) const
Definition: caching.hh:83
ShapeFunctionSet::HessianRangeType HessianRangeType
Definition: caching.hh:44
const ThisType & impl() const
Definition: caching.hh:125
void jacobianEach(const Point &x, Functor functor) const
Definition: caching.hh:90
const RangeVectorType & rangeCache(const QuadratureType &quadrature) const
Definition: caching.hh:111
std::vector< JacobianRangeType > JacobianRangeVectorType
Definition: caching.hh:47
void jacobianEach(const QuadraturePointWrapper< Quadrature > &x, Functor functor) const
Definition: caching.hh:96
ShapeFunctionSet ShapeFunctionSetType
Definition: caching.hh:37
std::vector< RangeType > RangeVectorType
Definition: caching.hh:46
std::vector< JacobianRangeVectorType > JacobianCacheVectorType
Definition: caching.hh:50
void hessianEach(const Point &x, Functor functor) const
Definition: caching.hh:103
ShapeFunctionSet::RangeType RangeType
Definition: caching.hh:42
CachingShapeFunctionSet(const GeometryType &type, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
Definition: caching.hh:56
ShapeFunctionSet::DomainType DomainType
Definition: caching.hh:41
GeometryType type() const
Definition: caching.hh:108
ShapeFunctionSet::JacobianRangeType JacobianRangeType
Definition: caching.hh:43
std::vector< RangeVectorType > RangeCacheVectorType
Definition: caching.hh:49
ShapeFunctionSet::FunctionSpaceType FunctionSpaceType
Definition: caching.hh:39
~CachingShapeFunctionSet()
Definition: caching.hh:254
void evaluateEach(const Point &x, Functor functor) const
Definition: caching.hh:77
std::size_t size() const
Definition: caching.hh:71
int order() const
Definition: caching.hh:66
const ThisType & scalarShapeFunctionSet() const
Definition: caching.hh:124
const JacobianRangeVectorType & jacobianCache(const QuadratureType &quadrature) const
Definition: caching.hh:118
static const int pointSetId
Definition: caching.hh:54
Interface class for shape function sets.
Definition: shapefunctionset/shapefunctionset.hh:33
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: shapefunctionset/shapefunctionset.hh:43
FunctionSpaceType::RangeType RangeType
range type
Definition: shapefunctionset/shapefunctionset.hh:41
void hessianEach(const Point &x, Functor functor) const
evalute hessian of each shape function
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
std::size_t size() const
return number of shape functions
int order() const
return order of shape functions
FunctionSpaceType::DomainType DomainType
domain type
Definition: shapefunctionset/shapefunctionset.hh:39
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function