DUNE-FEM (unstable)

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