1#ifndef DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
2#define DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
12#include <dune/fem/misc/threads/threadsafevalue.hh>
14#include <dune/fem/function/common/localcontribution.hh>
15#include <dune/fem/space/basisfunctionset/transformed.hh>
16#include <dune/fem/space/basisfunctionset/vectorial.hh>
17#include <dune/fem/space/combinedspace/interpolation.hh>
18#include <dune/fem/space/common/localinterpolation.hh>
30 template<
class LocalFunction,
class BasisFunctionSet >
31 struct LocalFunctionWrapper
41 LocalFunctionWrapper (
const LocalFunction &lf,
const BasisFunctionSet &bset ) : lf_( lf ) {}
44 void evaluate (
const Arg &x,
typename Traits::RangeType &y )
const
49 typename Traits::RangeType operator()(
const Arg &x)
const
51 typename Traits::RangeType y;
57 const LocalFunction &lf_;
61 template<
class LocalFunction,
class Entity,
class ShapeFunctionSet,
class Transformation >
62 struct LocalFunctionWrapper< LocalFunction, TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > >
64 typedef TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > BasisFunctionSetType;
74 LocalFunctionWrapper (
const LocalFunction &lf,
const BasisFunctionSetType &bset )
75 : lf_( lf ), bset_( bset )
81 void evaluate (
const Arg &x,
typename Traits::RangeType &y )
const
83 typename Traits::RangeType help;
84 lf_.evaluate( x, help );
87 typename Transformation::InverseTransformationType transf( lf_.geometry(), x );
88 y = transf.apply( help );
91 typename Traits::RangeType operator() (
const Arg &x )
const
93 typename Traits::RangeType y;
99 const LocalFunction &lf_;
100 const BasisFunctionSetType &bset_;
109 template<
class Space,
class LocalInterpolation,
bool scalarSFS >
110 class LocalFiniteElementInterpolation;
115 template<
class Space,
class LocalInterpolation >
116 class LocalFiniteElementInterpolation<Space,LocalInterpolation,false>
118 typedef LocalFiniteElementInterpolation< Space, LocalInterpolation,false > ThisType;
121 typedef typename Space::BasisFunctionSetType BasisFunctionSetType;
122 typedef LocalInterpolation LocalInterpolationType;
131 typedef std::size_t size_type;
133 template<
class LocalFunction >
134 using LocalFunctionWrapper = Impl::LocalFunctionWrapper< LocalFunction, BasisFunctionSetType >;
137 LocalFiniteElementInterpolation()
138 : basisFunctionSet_()
139 , localInterpolation_()
143 explicit LocalFiniteElementInterpolation (
const BasisFunctionSetType &basisFunctionSet,
144 const LocalInterpolationType &localInterpolation = LocalInterpolationType() )
145 : basisFunctionSet_( basisFunctionSet ),
146 localInterpolation_( localInterpolation )
149 ThisType& operator=(
const ThisType& other )
151 basisFunctionSet_ = other.basisFunctionSet_;
152 localInterpolation_.emplace( other.localInterpolation() );
161 template<
class LocalFunction,
class Dof >
162 void operator() (
const LocalFunction &localFunction, std::vector< Dof > &localDofVector )
const
164 LocalFunctionWrapper< LocalFunction > wrapper( localFunction, basisFunctionSet() );
165 localInterpolation().interpolate( wrapper, localDofVector );
168 template<
class LocalFunction,
class Dof,
class Allocator >
171 (*this)(localFunction,localDofVector.container() );
174 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
175 void operator() (
const LocalFunction &localFunction, LocalContribution< DiscreteFunction, Assembly > &localContribution )
const
177 (*this)(localFunction,localContribution.localDofVector());
180 BasisFunctionSetType basisFunctionSet ()
const {
return basisFunctionSet_; }
181 const LocalInterpolationType &localInterpolation ()
const
183 assert( localInterpolation_.has_value() );
184 return *localInterpolation_;
188 BasisFunctionSetType basisFunctionSet_;
189 std::optional< LocalInterpolationType > localInterpolation_;
194 template <
int dimRange>
195 struct RangeConverter
197 RangeConverter ( std::size_t range ) : range_( range ) {}
200 FieldVector< T, 1 > operator() (
const FieldVector< T, dimRange > &in )
const
205 template<
class T,
int j >
206 FieldMatrix< T, 1, j > operator() (
const FieldMatrix< T, dimRange, j > &in )
const
209 FieldMatrix<T,1,j> ret;
217 template <
class DofVector,
class DofAlignment>
218 struct SubDofVectorWrapper
219 :
public SubDofVector< DofVector, DofAlignment >
221 typedef SubDofVector< DofVector, DofAlignment > BaseType;
223 SubDofVectorWrapper( DofVector& dofs,
int coordinate,
const DofAlignment &dofAlignment )
224 : BaseType( dofs, coordinate, dofAlignment )
229 void resize(
const unsigned int) {}
237 template<
class Space,
class LocalInterpolation >
238 class LocalFiniteElementInterpolation<Space,LocalInterpolation,true>
240 typedef LocalFiniteElementInterpolation< Space, LocalInterpolation,true > ThisType;
243 typedef typename Space::BasisFunctionSetType BasisFunctionSetType;
244 typedef LocalInterpolation LocalInterpolationType;
253 static const int dimR = Space::FunctionSpaceType::dimRange;
255 typedef std::size_t size_type;
257 typedef VerticalDofAlignment< BasisFunctionSetType, RangeType> DofAlignmentType;
260 LocalFiniteElementInterpolation ()
261 : basisFunctionSet_(),
262 localInterpolation_(),
266 explicit LocalFiniteElementInterpolation (
const BasisFunctionSetType &basisFunctionSet,
267 const LocalInterpolationType &localInterpolation = LocalInterpolationType() )
268 : basisFunctionSet_( basisFunctionSet ),
269 localInterpolation_( localInterpolation ),
270 dofAlignment_( basisFunctionSet_ )
273 LocalFiniteElementInterpolation (
const ThisType& other )
274 : basisFunctionSet_( other.basisFunctionSet_ ),
275 localInterpolation_( other.localInterpolation_ ),
276 dofAlignment_( basisFunctionSet_ )
279 ThisType& operator=(
const ThisType& other )
281 basisFunctionSet_ = other.basisFunctionSet_;
282 localInterpolation_.emplace( other.localInterpolation() );
291 template<
class LocalFunction,
class Vector>
292 void operator() (
const LocalFunction &localFunction, Vector &localDofVector )
const
296 std::fill(localDofVector.begin(),localDofVector.end(),0);
297 for( std::size_t i = 0; i < dimR; ++i )
299 Impl::SubDofVectorWrapper< Vector, DofAlignmentType > subLdv( localDofVector, i, dofAlignment_ );
301 localFunctionConverter( localFunction, Impl::RangeConverter<dimR>(i) ),
302 subLdv, PriorityTag<42>()
307 template<
class LocalFunction,
class Dof,
class Allocator >
310 (*this)(localFunction,localDofVector.container() );
313 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
314 void operator() (
const LocalFunction &localFunction, LocalContribution< DiscreteFunction, Assembly > &localContribution )
const
316 (*this)(localFunction,localContribution.localDofVector());
319 BasisFunctionSetType basisFunctionSet ()
const {
return basisFunctionSet_; }
320 const LocalInterpolationType &localInterpolation ()
const
322 assert( localInterpolation_.has_value() );
323 return *localInterpolation_;
327 template<
class LocalFunction,
class Vector>
328 auto operator() (
const LocalFunction &localFunction, Vector &localDofVector, PriorityTag<1> )
const
330 std::declval<LocalInterpolationType>().interpolate(
331 localFunction, localDofVector)) >
333 localInterpolation().interpolate( localFunction, localDofVector );
335 template<
class LocalFunction,
class Vector>
336 void operator() (
const LocalFunction &localFunction, Vector &localDofVector, PriorityTag<0> )
const
338 std::vector<typename Vector::value_type> tmp(basisFunctionSet_.size());
339 localInterpolation().interpolate( localFunction, tmp);
340 for (
unsigned int i=0;i<tmp.size();++i)
341 localDofVector[i] = tmp[i];
343 BasisFunctionSetType basisFunctionSet_;
344 std::optional< LocalInterpolationType > localInterpolation_;
345 DofAlignmentType dofAlignment_;
348 template <
class DiscreteFunctionSpace >
349 class LocalFEInterpolationWrapper
350 :
public LocalInterpolationWrapper< DiscreteFunctionSpace >
352 typedef LocalInterpolationWrapper< DiscreteFunctionSpace > BaseType;
354 using BaseType :: interpolation_;
355 typedef std::vector< typename DiscreteFunctionSpace::RangeFieldType >
356 TemporarayDofVectorType;
359 using BaseType::interpolation;
365 using BaseType :: bind;
366 using BaseType :: unbind;
368 template<
class LocalFunction,
class Dof >
369 void operator() (
const LocalFunction &localFunction, std::vector< Dof > &localDofVector )
const
372 interpolation()( localFunction, localDofVector );
375 template<
class LocalFunction,
class Dof,
class Allocator >
379 (*this)(localFunction, localDofVector.container() );
383 template<
class LocalFunction,
class DiscreteFunction,
template<
class >
class Assembly >
384 void operator() (
const LocalFunction &localFunction, LocalContribution< DiscreteFunction, Assembly > &localContribution )
const
387 (*this)(localFunction, localContribution.localDofVector());
391 template<
class LocalFunction,
class LocalDofVector >
392 void operator () (
const LocalFunction &localFunction, LocalDofVector &dofs )
const
394 const int size = dofs.size();
395 TemporarayDofVectorType& tmpDofs = *tmpDofs_;
396 tmpDofs.resize(
size );
397 (*this)(localFunction, tmpDofs );
400 for(
int i=0; i<
size; ++i )
402 dofs[ i ] = tmpDofs[ i ];
Construct a vector with a dynamic size.
Definition: dynvector.hh:59
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceType::DomainType DomainType
type of domain vectors, i.e., type of coordinates
Definition: localfunction.hh:108
FunctionSpaceType::RangeType RangeType
type of range vectors, i.e., type of function values
Definition: localfunction.hh:110
Implements a vector constructed from a given type representing a field and a compile-time given size.
typename Impl::voider< Types... >::type void_t
Is void for all valid input types. The workhorse for C++11 SFINAE-techniques.
Definition: typetraits.hh:40
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75