dune-fem 2.8.0
Loading...
Searching...
No Matches
localfiniteelement/interpolation.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
2#define DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
3
4#include <cstddef>
5
6#include <utility>
7#include <vector>
8
9#include <dune/common/fvector.hh>
10
12
18
19namespace Dune
20{
21
22 namespace Fem
23 {
24
25 namespace Impl
26 {
27 // LocalFunctionWrapper
28 // --------------------
29 template< class LocalFunction, class BasisFunctionSet >
30 struct LocalFunctionWrapper
31 {
32 struct Traits
33 {
34 typedef typename LocalFunction::DomainType DomainType;
35 typedef typename LocalFunction::RangeType RangeType;
36 };
37 typedef typename LocalFunction::DomainType DomainType;
38 typedef typename LocalFunction::RangeType RangeType;
39
40 LocalFunctionWrapper ( const LocalFunction &lf, const BasisFunctionSet &bset ) : lf_( lf ) {}
41
42 template< class Arg >
43 void evaluate ( const Arg &x, typename Traits::RangeType &y ) const
44 {
45 lf_.evaluate( x, y );
46 }
47 template< class Arg >
48 typename Traits::RangeType operator()(const Arg &x) const
49 {
50 typename Traits::RangeType y;
51 evaluate(x,y);
52 return y;
53 }
54
55 private:
56 const LocalFunction &lf_;
57 };
58 // LocalFunctionWrapper
59 // --------------------
60 template< class LocalFunction, class Entity, class ShapeFunctionSet, class Transformation >
61 struct LocalFunctionWrapper< LocalFunction, TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > >
62 {
63 typedef TransformedBasisFunctionSet< Entity, ShapeFunctionSet, Transformation > BasisFunctionSetType;
64
65 struct Traits
66 {
67 typedef typename LocalFunction::DomainType DomainType;
68 typedef typename LocalFunction::RangeType RangeType;
69 };
70 typedef typename LocalFunction::DomainType DomainType;
71 typedef typename LocalFunction::RangeType RangeType;
72
73 LocalFunctionWrapper ( const LocalFunction &lf, const BasisFunctionSetType &bset )
74 : lf_( lf ), geometry_( lf_.entity().geometry() ), bset_( bset ) {}
75
76 template< class Arg >
77 void evaluate ( const Arg &x, typename Traits::RangeType &y ) const
78 {
79 typename Traits::RangeType help;
80 lf_.evaluate( x, help );
81 typename Transformation::InverseTransformationType transf( geometry_, x );
82 y = transf.apply( help );
83 }
84 template< class Arg >
85 typename Traits::RangeType operator() ( const Arg &x ) const
86 {
87 typename Traits::RangeType y;
88 evaluate(x,y);
89 return y;
90 }
91
92 private:
93 const LocalFunction &lf_;
94 const typename Entity::Geometry geometry_;
95 const BasisFunctionSetType &bset_;
96 };
97
98 } // namespace Impl
99
100
101 // LocalFiniteElementInterpolation
102 // -------------------------------
103
104 template< class Space, class LocalInterpolation, bool scalarSFS >
106
107 // a vector valued shape basis function set taken from
108 // dune-localfunction - the vector value 'blow up' is not yet supported
109 // for this case
110 template< class Space, class LocalInterpolation >
112 {
114
115 public:
116 typedef typename Space::BasisFunctionSetType BasisFunctionSetType;
118
119 private:
120 typedef typename BasisFunctionSetType::FunctionSpaceType FunctionSpaceType;
121
122 typedef typename FunctionSpaceType::RangeType RangeType;
123 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
124 static const int dimRange = FunctionSpaceType::dimRange;
125
126 typedef std::size_t size_type;
127
128 template< class LocalFunction >
129 using LocalFunctionWrapper = Impl::LocalFunctionWrapper< LocalFunction, BasisFunctionSetType >;
130
131 public:
133 : basisFunctionSet_()
134 , localInterpolation_()
135 {
136 }
137
138 explicit LocalFiniteElementInterpolation ( const BasisFunctionSetType &basisFunctionSet,
139 const LocalInterpolationType &localInterpolation = LocalInterpolationType() )
140 : basisFunctionSet_( basisFunctionSet ),
141 localInterpolation_( localInterpolation )
142 {}
143
144 ThisType& operator=( const ThisType& other )
145 {
146 basisFunctionSet_ = other.basisFunctionSet_;
147 localInterpolation_.emplace( other.localInterpolation() );
148 return *this;
149 }
150
151 void unbind()
152 {
153 // basisFunctionSet_ = BasisFunctionSetType();
154 }
155
156 template< class LocalFunction, class Dof >
157 void operator() ( const LocalFunction &localFunction, std::vector< Dof > &localDofVector ) const
158 {
159 LocalFunctionWrapper< LocalFunction > wrapper( localFunction, basisFunctionSet() );
160 localInterpolation().interpolate( wrapper, localDofVector );
161 }
162
163 template< class LocalFunction, class Dof, class Allocator >
164 void operator() ( const LocalFunction &localFunction, Dune::DynamicVector< Dof, Allocator > &localDofVector ) const
165 {
166 (*this)(localFunction,localDofVector.container() );
167 }
168
169 template< class LocalFunction, class DiscreteFunction, template< class > class Assembly >
170 void operator() ( const LocalFunction &localFunction, LocalContribution< DiscreteFunction, Assembly > &localContribution ) const
171 {
172 (*this)(localFunction,localContribution.localDofVector());
173 }
174
175 BasisFunctionSetType basisFunctionSet () const { return basisFunctionSet_; }
177 {
178 assert( localInterpolation_.has_value() );
179 return *localInterpolation_;
180 }
181
182 private:
183 BasisFunctionSetType basisFunctionSet_;
184 std::optional< LocalInterpolationType > localInterpolation_;
185 };
186
187 namespace Impl
188 {
189 template <int dimRange>
190 struct RangeConverter
191 {
192 RangeConverter ( std::size_t range ) : range_( range ) {}
193
194 template< class T >
195 FieldVector< T, 1 > operator() ( const FieldVector< T, dimRange > &in ) const
196 {
197 return in[ range_ ];
198 }
199
200 template< class T, int j >
201 FieldMatrix< T, 1, j > operator() ( const FieldMatrix< T, dimRange, j > &in ) const
202 {
203 // return in[ range_ ]; // implicit conversion fails
204 FieldMatrix<T,1,j> ret;
205 ret[0] = in[range_];
206 return ret;
207 }
208
209 protected:
210 std::size_t range_;
211 };
212 template <class DofVector, class DofAlignment>
213 struct SubDofVectorWrapper
214 : public SubDofVector< DofVector, DofAlignment >
215 {
216 typedef SubDofVector< DofVector, DofAlignment > BaseType;
217
218 SubDofVectorWrapper( DofVector& dofs, int coordinate, const DofAlignment &dofAlignment )
219 : BaseType( dofs, coordinate, dofAlignment )
220 {}
221
223 void clear() {}
224 void resize( const unsigned int) {}
225 };
226
227 }
228
229 // a scalar valued shape basis function set taken from
230 // dune-localfunction - for this the vector value 'blow up' is supported
231 // as with other spaces
232 template< class Space, class LocalInterpolation >
234 {
236
237 public:
238 typedef typename Space::BasisFunctionSetType BasisFunctionSetType;
240
241 private:
242 typedef typename BasisFunctionSetType::FunctionSpaceType FunctionSpaceType;
243 // typedef typename Space::FunctionSpaceType FunctionSpaceType;
244
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;
249
250 typedef std::size_t size_type;
251
253
254 public:
256 : basisFunctionSet_(),
257 localInterpolation_(),
258 dofAlignment_()
259 {}
260
261 explicit LocalFiniteElementInterpolation ( const BasisFunctionSetType &basisFunctionSet,
262 const LocalInterpolationType &localInterpolation = LocalInterpolationType() )
263 : basisFunctionSet_( basisFunctionSet ),
264 localInterpolation_( localInterpolation ),
265 dofAlignment_( basisFunctionSet_ )
266 {}
267
269 : basisFunctionSet_( other.basisFunctionSet_ ),
270 localInterpolation_( other.localInterpolation_ ),
271 dofAlignment_( basisFunctionSet_ )
272 {}
273
274 ThisType& operator=( const ThisType& other )
275 {
276 basisFunctionSet_ = other.basisFunctionSet_;
277 localInterpolation_.emplace( other.localInterpolation() );
278 return *this;
279 }
280
281 void unbind()
282 {
283 // basisFunctionSet_ = BasisFunctionSetType();
284 }
285
286 template< class LocalFunction, class Vector>
287 void operator() ( const LocalFunction &localFunction, Vector &localDofVector ) const
288 {
289 // clear dofs before something is adedd
290 // localDofVector.clear(); // does not exist on DynVector so use 'fill' instead
291 std::fill(localDofVector.begin(),localDofVector.end(),0);
292 for( std::size_t i = 0; i < dimR; ++i )
293 {
294 Impl::SubDofVectorWrapper< Vector, DofAlignmentType > subLdv( localDofVector, i, dofAlignment_ );
295 (*this)(
296 localFunctionConverter( localFunction, Impl::RangeConverter<dimR>(i) ),
297 subLdv, PriorityTag<42>()
298 );
299 }
300 }
301
302 template< class LocalFunction, class Dof, class Allocator >
303 void operator() ( const LocalFunction &localFunction, Dune::DynamicVector< Dof, Allocator > &localDofVector ) const
304 {
305 (*this)(localFunction,localDofVector.container() );
306 }
307
308 template< class LocalFunction, class DiscreteFunction, template< class > class Assembly >
309 void operator() ( const LocalFunction &localFunction, LocalContribution< DiscreteFunction, Assembly > &localContribution ) const
310 {
311 (*this)(localFunction,localContribution.localDofVector());
312 }
313
314 BasisFunctionSetType basisFunctionSet () const { return basisFunctionSet_; }
316 {
317 assert( localInterpolation_.has_value() );
318 return *localInterpolation_;
319 }
320
321 private:
322 template< class LocalFunction, class Vector>
323 auto operator() ( const LocalFunction &localFunction, Vector &localDofVector, PriorityTag<1> ) const
324 -> void_t< decltype(
325 std::declval<LocalInterpolationType>().interpolate(
326 localFunction, localDofVector)) >
327 {
328 localInterpolation().interpolate( localFunction, localDofVector );
329 }
330 template< class LocalFunction, class Vector>
331 void operator() ( const LocalFunction &localFunction, Vector &localDofVector, PriorityTag<0> ) const
332 {
333 std::vector<typename Vector::value_type> tmp(basisFunctionSet_.size());
334 localInterpolation().interpolate( localFunction, tmp);
335 for (unsigned int i=0;i<tmp.size();++i)
336 localDofVector[i] = tmp[i];
337 }
338 BasisFunctionSetType basisFunctionSet_;
339 std::optional< LocalInterpolationType > localInterpolation_;
340 DofAlignmentType dofAlignment_;
341 };
342
343 template <class DiscreteFunctionSpace >
345 : public LocalInterpolationWrapper< DiscreteFunctionSpace >
346 {
348 protected:
350 typedef std::vector< typename DiscreteFunctionSpace::RangeFieldType >
353
355 public:
357 : BaseType( space )
358 {}
359
360 using BaseType :: bind;
361 using BaseType :: unbind;
362
363 template< class LocalFunction, class Dof >
364 void operator() ( const LocalFunction &localFunction, std::vector< Dof > &localDofVector ) const
365 {
366 // interpolation only works for std::vector storage
367 interpolation()( localFunction, localDofVector );
368 }
369
370 template< class LocalFunction, class Dof, class Allocator >
371 void operator() ( const LocalFunction &localFunction, Dune::DynamicVector< Dof, Allocator > &localDofVector ) const
372 {
373 // DynamicVector internally stores a std::vector
374 (*this)(localFunction, localDofVector.container() );
375 }
376
377
378 template< class LocalFunction, class DiscreteFunction, template< class > class Assembly >
380 {
381 // LocalContribution internally stores a std::vector
382 (*this)(localFunction, localContribution.localDofVector());
383 }
384
386 template< class LocalFunction, class LocalDofVector >
387 void operator () ( const LocalFunction &localFunction, LocalDofVector &dofs ) const
388 {
389 const int size = dofs.size();
391 tmpDofs.resize( size );
392 (*this)(localFunction, tmpDofs );
393
394 // copy to dof vector
395 for( int i=0; i<size; ++i )
396 {
397 dofs[ i ] = tmpDofs[ i ];
398 }
399 }
400 };
401
402 } // namespace Fem
403
404} // namespace Dune
405
406#endif // #ifndef DUNE_FEM_SPACE_LOCALFINITEELEMENT_INTERPOLATION_HH
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
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
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
discrete function space