1#ifndef DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
2#define DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
9#include <dune/common/fvector.hh>
29 template<
class LocalFunction,
class BasisFunctionSet >
30 struct LocalFunctionWrapper
40 LocalFunctionWrapper (
const LocalFunction &lf,
const BasisFunctionSet &bset ) : lf_( lf ) {}
43 void evaluate (
const Arg &x,
typename Traits::RangeType &y )
const
48 typename Traits::RangeType operator()(
const Arg &x)
const
50 typename Traits::RangeType y;
56 const LocalFunction &lf_;
60 template<
class LocalFunction,
class Entity,
class ShapeFunctionSet,
class Transformation >
61 struct LocalFunctionWrapper< LocalFunction, TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > >
63 typedef TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > BasisFunctionSetType;
73 LocalFunctionWrapper (
const LocalFunction &lf,
const BasisFunctionSetType &bset )
74 : lf_( lf ), geometry_( lf_.entity().geometry() ), bset_( bset ) {}
77 void evaluate (
const Arg &x,
typename Traits::RangeType &y )
const
79 typename Traits::RangeType help;
80 lf_.evaluate( x, help );
81 typename Transformation::InverseTransformationType transf( geometry_, x );
82 y = transf.apply( help );
85 typename Traits::RangeType operator() (
const Arg &x )
const
87 typename Traits::RangeType y;
93 const LocalFunction &lf_;
94 const typename Entity::Geometry geometry_;
95 const BasisFunctionSetType &bset_;
104 template<
class Space,
class LocalInterpolation,
bool scalarSFS >
110 template<
class Space,
class LocalInterpolation >
120 typedef typename BasisFunctionSetType::FunctionSpaceType FunctionSpaceType;
122 typedef typename FunctionSpaceType::RangeType RangeType;
123 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
124 static const int dimRange = FunctionSpaceType::dimRange;
126 typedef std::size_t size_type;
128 template<
class LocalFunction >
129 using LocalFunctionWrapper = Impl::LocalFunctionWrapper< LocalFunction, BasisFunctionSetType >;
133 : basisFunctionSet_()
134 , localInterpolation_()
140 : basisFunctionSet_( basisFunctionSet ),
141 localInterpolation_( localInterpolation )
146 basisFunctionSet_ = other.basisFunctionSet_;
156 template<
class LocalFunction,
class Dof >
159 LocalFunctionWrapper< LocalFunction > wrapper(
localFunction, basisFunctionSet() );
160 localInterpolation().interpolate( wrapper, localDofVector );
163 template<
class LocalFunction,
class Dof,
class Allocator >
169 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
178 assert( localInterpolation_.has_value() );
179 return *localInterpolation_;
183 BasisFunctionSetType basisFunctionSet_;
184 std::optional< LocalInterpolationType > localInterpolation_;
189 template <
int dimRange>
190 struct RangeConverter
192 RangeConverter ( std::size_t range ) : range_( range ) {}
195 FieldVector< T, 1 > operator() (
const FieldVector< T, dimRange > &
in )
const
200 template<
class T,
int j >
201 FieldMatrix< T, 1, j > operator() (
const FieldMatrix< T, dimRange, j > &
in )
const
204 FieldMatrix<T,1,j> ret;
212 template <
class DofVector,
class DofAlignment>
213 struct SubDofVectorWrapper
214 :
public SubDofVector< DofVector, DofAlignment >
216 typedef SubDofVector< DofVector, DofAlignment > BaseType;
218 SubDofVectorWrapper( DofVector&
dofs,
int coordinate,
const DofAlignment &dofAlignment )
224 void resize(
const unsigned int) {}
232 template<
class Space,
class LocalInterpolation >
242 typedef typename BasisFunctionSetType::FunctionSpaceType FunctionSpaceType;
245 typedef typename FunctionSpaceType::RangeType RangeType;
246 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
247 static const int dimRange = FunctionSpaceType::dimRange;
248 static const int dimR = Space::FunctionSpaceType::dimRange;
250 typedef std::size_t size_type;
256 : basisFunctionSet_(),
257 localInterpolation_(),
263 : basisFunctionSet_( basisFunctionSet ),
264 localInterpolation_( localInterpolation ),
265 dofAlignment_( basisFunctionSet_ )
269 : basisFunctionSet_( other.basisFunctionSet_ ),
270 localInterpolation_( other.localInterpolation_ ),
271 dofAlignment_( basisFunctionSet_ )
276 basisFunctionSet_ = other.basisFunctionSet_;
286 template<
class LocalFunction,
class Vector>
291 std::fill(localDofVector.begin(),localDofVector.end(),0);
292 for( std::size_t i = 0; i < dimR; ++i )
294 Impl::SubDofVectorWrapper< Vector, DofAlignmentType > subLdv( localDofVector, i, dofAlignment_ );
297 subLdv, PriorityTag<42>()
302 template<
class LocalFunction,
class Dof,
class Allocator >
308 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
317 assert( localInterpolation_.has_value() );
318 return *localInterpolation_;
322 template<
class LocalFunction,
class Vector>
325 std::declval<LocalInterpolationType>().interpolate(
328 localInterpolation().interpolate(
localFunction, localDofVector );
330 template<
class LocalFunction,
class Vector>
331 void operator() (
const LocalFunction &
localFunction, Vector &localDofVector, PriorityTag<0> )
const
333 std::vector<typename Vector::value_type> tmp(basisFunctionSet_.size());
335 for (
unsigned int i=0;i<tmp.size();++i)
336 localDofVector[i] = tmp[i];
338 BasisFunctionSetType basisFunctionSet_;
339 std::optional< LocalInterpolationType > localInterpolation_;
340 DofAlignmentType dofAlignment_;
343 template <
class DiscreteFunctionSpace >
350 typedef std::vector< typename DiscreteFunctionSpace::RangeFieldType >
363 template<
class LocalFunction,
class Dof >
370 template<
class LocalFunction,
class Dof,
class Allocator >
378 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
386 template<
class LocalFunction,
class LocalDofVector >
389 const int size =
dofs.size();
391 tmpDofs.resize( size );
395 for(
int i=0; i<size; ++i )
397 dofs[ i ] = tmpDofs[ i ];
Definition: bindguard.hh:11
static GridFunctionView< GF > localFunction(const GF &gf)
Definition: gridfunctionview.hh:118
LocalFunctionConverter< HostLocalFunction, Converter, __InstationaryFunction::HoldCopy > localFunctionConverter(HostLocalFunction hostLocalFunction, const Converter &converter=Converter())
Definition: converter.hh:199
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
Definition: common/localcontribution.hh:14
interface for local functions
Definition: localfunction.hh:77
FunctionSpaceType::DomainType DomainType
type of domain vectors, i.e., type of coordinates
Definition: localfunction.hh:105
FunctionSpaceType::RangeType RangeType
type of range vectors, i.e., type of function values
Definition: localfunction.hh:107
Implementation of DofAlignment.
Definition: basisfunctionset/vectorial.hh:141
Definition: common/localinterpolation.hh:22
Definition: common/localinterpolation.hh:74
void unbind()
Definition: common/localinterpolation.hh:91
const InterpolationImplType & interpolation() const
Definition: common/localinterpolation.hh:104
void bind(const EntityType &entity)
Definition: common/localinterpolation.hh:86
std::optional< InterpolationImplType > interpolation_
Definition: common/localinterpolation.hh:116
Definition: localfiniteelement/interpolation.hh:105
Definition: localfiniteelement/interpolation.hh:112
LocalFiniteElementInterpolation()
Definition: localfiniteelement/interpolation.hh:132
const LocalInterpolationType & localInterpolation() const
Definition: localfiniteelement/interpolation.hh:176
LocalFiniteElementInterpolation(const BasisFunctionSetType &basisFunctionSet, const LocalInterpolationType &localInterpolation=LocalInterpolationType())
Definition: localfiniteelement/interpolation.hh:138
LocalInterpolation LocalInterpolationType
Definition: localfiniteelement/interpolation.hh:117
Space::BasisFunctionSetType BasisFunctionSetType
Definition: localfiniteelement/interpolation.hh:116
BasisFunctionSetType basisFunctionSet() const
Definition: localfiniteelement/interpolation.hh:175
void unbind()
Definition: localfiniteelement/interpolation.hh:151
ThisType & operator=(const ThisType &other)
Definition: localfiniteelement/interpolation.hh:144
Definition: localfiniteelement/interpolation.hh:234
void unbind()
Definition: localfiniteelement/interpolation.hh:281
ThisType & operator=(const ThisType &other)
Definition: localfiniteelement/interpolation.hh:274
Space::BasisFunctionSetType BasisFunctionSetType
Definition: localfiniteelement/interpolation.hh:238
LocalFiniteElementInterpolation(const ThisType &other)
Definition: localfiniteelement/interpolation.hh:268
LocalFiniteElementInterpolation()
Definition: localfiniteelement/interpolation.hh:255
LocalInterpolation LocalInterpolationType
Definition: localfiniteelement/interpolation.hh:239
const LocalInterpolationType & localInterpolation() const
Definition: localfiniteelement/interpolation.hh:315
BasisFunctionSetType basisFunctionSet() const
Definition: localfiniteelement/interpolation.hh:314
LocalFiniteElementInterpolation(const BasisFunctionSetType &basisFunctionSet, const LocalInterpolationType &localInterpolation=LocalInterpolationType())
Definition: localfiniteelement/interpolation.hh:261
Definition: localfiniteelement/interpolation.hh:346
ThreadSafeValue< TemporarayDofVectorType > tmpDofs_
Definition: localfiniteelement/interpolation.hh:352
const InterpolationImplType & interpolation() const
Definition: common/localinterpolation.hh:104
void operator()(const LocalFunction &localFunction, std::vector< Dof > &localDofVector) const
Definition: localfiniteelement/interpolation.hh:364
std::vector< typename DiscreteFunctionSpace::RangeFieldType > TemporarayDofVectorType
Definition: localfiniteelement/interpolation.hh:351
LocalFEInterpolationWrapper(const DiscreteFunctionSpace &space)
Definition: localfiniteelement/interpolation.hh:356