DUNE-FEM (unstable)

genericbasefunctions.hh
1#ifndef DUNE_FEM_SPACE_LAGRANGE_GENERICBASEFUNCTIONS_HH
2#define DUNE_FEM_SPACE_LAGRANGE_GENERICBASEFUNCTIONS_HH
3
5
6#include "genericlagrangepoints.hh"
7
8namespace Dune
9{
10
11 namespace Fem
12 {
13
14 template< class FunctionSpace, class GeometryType, unsigned int order >
15 class GenericLagrangeBaseFunction;
16
17
18 template< class FunctionSpace, unsigned int order >
19 class GenericLagrangeBaseFunction< FunctionSpace, PointGeometry, order >
20 {
21 typedef GenericLagrangeBaseFunction< FunctionSpace, PointGeometry, order > ThisType;
22
23 static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
24
25 public:
26 typedef FunctionSpace FunctionSpaceType;
27
28 typedef PointGeometry GeometryType;
29
30 static constexpr unsigned int polynomialOrder = order;
31
32 typedef GenericLagrangePoint< GeometryType, polynomialOrder >
33 LagrangePointType;
34 static const unsigned int numBaseFunctions = LagrangePointType :: numLagrangePoints;
35
36 typedef typename FunctionSpaceType :: DomainType DomainType;
37 typedef typename FunctionSpaceType :: RangeType RangeType;
38
39 typedef typename FunctionSpaceType :: DomainFieldType DomainFieldType;
40 typedef typename FunctionSpaceType :: RangeFieldType RangeFieldType;
41 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
42
43 private:
44 const LagrangePointType lagrangePoint_;
45
46 public:
47 explicit GenericLagrangeBaseFunction ( unsigned int baseNum )
48 : lagrangePoint_( baseNum )
49 {}
50
51 template< class LocalDofCoordinateType, class LocalCoordinateType >
52 inline static void evaluate ( LocalDofCoordinateType &dofCoordinate,
53 const FieldVector< int, 0 > &diffVariable,
54 DomainFieldType factor,
55 const LocalCoordinateType &x,
56 RangeType &phi )
57 {
58 phi[ 0 ] = 1;
59 }
60
61 template< class LocalDofCoordinateType, class LocalCoordinateType >
62 inline static void evaluate ( LocalDofCoordinateType &dofCoordinate,
63 const FieldVector< int, 1 > &diffVariable,
64 DomainFieldType factor,
65 const LocalCoordinateType &x,
66 RangeType &phi )
67 {
68 phi[ 0 ] = 0;
69 }
70
71 template< class LocalDofCoordinateType, class LocalCoordinateType >
72 inline static void evaluate ( LocalDofCoordinateType &dofCoordinate,
73 const FieldVector< int, 2 > &diffVariable,
74 DomainFieldType factor,
75 const LocalCoordinateType &x,
76 RangeType &phi )
77 {
78 phi[ 0 ] = 0;
79 }
80
81 template< int diffOrder >
82 inline void
83 evaluate ( const FieldVector< int, diffOrder > &diffVariable,
84 const DomainType &x,
85 RangeType &phi ) const
86 {
87 const LocalCoordinate< GeometryType, DomainFieldType > xlocal( x );
88 LagrangePointType point( lagrangePoint_ );
89 evaluate( point.dofCoordinate_, diffVariable, 1, xlocal, phi );
90 }
91 };
92
93
94 template< class FunctionSpace, class BaseGeometryType >
95 class GenericLagrangeBaseFunction< FunctionSpace, PyramidGeometry< BaseGeometryType >, 0 >
96 {
97 typedef GenericLagrangeBaseFunction< FunctionSpace, PyramidGeometry< BaseGeometryType >, 0 > ThisType;
98
99 static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
100
101 public:
102 typedef FunctionSpace FunctionSpaceType;
103
104 typedef PyramidGeometry< BaseGeometryType > GeometryType;
105
106 enum { polynomialOrder = 0 };
107
108 typedef GenericLagrangePoint< GeometryType, polynomialOrder > LagrangePointType;
109 static const unsigned int numBaseFunctions = LagrangePointType :: numLagrangePoints;
110
111
112 typedef typename FunctionSpaceType :: DomainType DomainType;
113 typedef typename FunctionSpaceType :: RangeType RangeType;
114
115 typedef typename FunctionSpaceType :: DomainFieldType DomainFieldType;
116 typedef typename FunctionSpaceType :: RangeFieldType RangeFieldType;
117 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
118
119 private:
120 const LagrangePointType lagrangePoint_;
121
122 public:
123 explicit GenericLagrangeBaseFunction ( unsigned int baseNum )
124 : lagrangePoint_( baseNum )
125 {}
126
127 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
128 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
129 const FieldVector< int, 0 > &diffVariable,
130 DomainFieldType factor,
131 const LocalCoordinateType &x,
132 RangeType &phi )
133 {
134 phi[ 0 ] = 1;
135 }
136
137 template< class LocalDofCoordinateType, class LocalCoordinateType >
138 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
139 const FieldVector< int, 0 > &diffVariable,
140 DomainFieldType factor,
141 const LocalCoordinateType &x,
142 RangeType &phi )
143 {
144 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
145 }
146
147 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
148 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
149 const FieldVector< int, 1 > &diffVariable,
150 DomainFieldType factor,
151 const LocalCoordinateType &x,
152 RangeType &phi )
153 {
154 phi[ 0 ] = 0;
155 }
156
157 template< class LocalDofCoordinateType, class LocalCoordinateType >
158 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
159 const FieldVector< int, 1 > &diffVariable,
160 DomainFieldType factor,
161 const LocalCoordinateType &x,
162 RangeType &phi )
163 {
164 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
165 }
166
167 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
168 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
169 const FieldVector< int, 2 > &diffVariable,
170 DomainFieldType factor,
171 const LocalCoordinateType &x,
172 RangeType &phi )
173 {
174 phi[ 0 ] = 0;
175 }
176
177 template< class LocalDofCoordinateType, class LocalCoordinateType >
178 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
179 const FieldVector< int, 2 > &diffVariable,
180 DomainFieldType factor,
181 const LocalCoordinateType &x,
182 RangeType &phi )
183 {
184 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
185 }
186
187 template< int diffOrder >
188 void evaluate ( const FieldVector< int, diffOrder > &diffVariable,
189 const DomainType &x,
190 RangeType &phi ) const
191 {
192 const LocalCoordinate< GeometryType, DomainFieldType > xlocal( x );
193 LagrangePointType point( lagrangePoint_ );
194 evaluate( point.dofCoordinate_, diffVariable, 1, xlocal, phi );
195 }
196 };
197
198
199 template< class FunctionSpace, class BaseGeometryType, unsigned int order >
200 class GenericLagrangeBaseFunction< FunctionSpace, PyramidGeometry< BaseGeometryType >, order >
201 {
202 typedef GenericLagrangeBaseFunction< FunctionSpace, PyramidGeometry< BaseGeometryType >, order > ThisType;
203
204 static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
205
206 public:
207 typedef FunctionSpace FunctionSpaceType;
208
209 typedef PyramidGeometry< BaseGeometryType > GeometryType;
210
211 static constexpr unsigned int polynomialOrder = order;
212
213 typedef GenericLagrangePoint< GeometryType, polynomialOrder >
214 LagrangePointType;
215 static const unsigned int numBaseFunctions = LagrangePointType :: numLagrangePoints;
216
217
218 typedef typename FunctionSpaceType :: DomainType DomainType;
219 typedef typename FunctionSpaceType :: RangeType RangeType;
220
221 typedef typename FunctionSpaceType :: DomainFieldType DomainFieldType;
222 typedef typename FunctionSpaceType :: RangeFieldType RangeFieldType;
223 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
224
225 private:
226 typedef GenericLagrangeBaseFunction
227 < FunctionSpaceType, BaseGeometryType, polynomialOrder >
228 DimensionReductionType;
229 typedef GenericLagrangeBaseFunction
230 < FunctionSpaceType, GeometryType, polynomialOrder - 1 >
231 OrderReductionType;
232
233 private:
234 const LagrangePointType lagrangePoint_;
235
236 public:
237 explicit GenericLagrangeBaseFunction ( unsigned int baseNum )
238 : lagrangePoint_( baseNum )
239 {}
240
241 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
242 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
243 const FieldVector< int, 0 > &diffVariable,
244 DomainFieldType factor,
245 const LocalCoordinateType &x,
246 RangeType &phi )
247 {
248 const DomainFieldType divisor = DomainFieldType( 1 ) / ((DomainFieldType)polynomialOrder);
249 const DomainFieldType myfactor = porder * divisor;
250 const RealType myshift = (porder - polynomialOrder) * divisor;
251
252 if( LagrangePointType :: useDimReduction( dofCoordinate ) )
253 {
254 DimensionReductionType::evaluate( dofCoordinate.base(), diffVariable, myfactor * factor, x.base(), phi );
255
256 const unsigned int height
257 = LagrangePointType :: height( dofCoordinate );
258 for( unsigned int i = 0; i < height; ++i )
259 {
260 ++(*dofCoordinate);
261 RangeType psi;
262 evaluate< porder >( dofCoordinate, diffVariable, factor, x, psi );
263 phi -= psi;
264 }
265 (*dofCoordinate) -= height;
266 }
267 else
268 {
269 --(*dofCoordinate);
270 OrderReductionType::template evaluate< porder >( dofCoordinate, diffVariable, factor, x, phi );
271 ++(*dofCoordinate);
272 phi[ 0 ] *= (polynomialOrder / ((RealType)(*dofCoordinate)))
273 * (myfactor * factor * (*x) - myshift);
274 }
275 }
276
277 template< class LocalDofCoordinateType, class LocalCoordinateType >
278 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
279 const FieldVector< int, 0 > &diffVariable,
280 DomainFieldType factor,
281 const LocalCoordinateType &x,
282 RangeType &phi )
283 {
284 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
285 }
286
287 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
288 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
289 const FieldVector< int, 1 > &diffVariable,
290 DomainFieldType factor,
291 const LocalCoordinateType &x,
292 RangeType &phi )
293 {
294 const DomainFieldType divisor = polynomialOrder;
295 const DomainFieldType myfactor = porder / divisor;
296 const RealType myshift = (porder - polynomialOrder) / divisor;
297
298 FieldVector< int, 0 > dv;
299
300 if( LagrangePointType :: useDimReduction( dofCoordinate ) )
301 {
302 if( (unsigned int)diffVariable[ 0 ] != LocalDofCoordinateType::index )
303 DimensionReductionType::evaluate( dofCoordinate.base(), diffVariable, myfactor * factor, x.base(), phi );
304 else
305 phi = 0;
306
307 const unsigned int height
308 = LagrangePointType :: height( dofCoordinate );
309 for( unsigned int i = 0; i < height; ++i )
310 {
311 ++(*dofCoordinate);
312 RangeType psi;
313 evaluate< porder >( dofCoordinate, diffVariable, factor, x, psi );
314 phi -= psi;
315 }
316 (*dofCoordinate) -= height;
317 }
318 else
319 {
320 --(*dofCoordinate);
321 OrderReductionType::template evaluate< porder >( dofCoordinate, diffVariable, factor, x, phi );
322 phi *= (myfactor * factor * (*x) - myshift);
323
324 if( (unsigned int)diffVariable[ 0 ] == LocalDofCoordinateType::index )
325 {
326 RangeType psi;
327 OrderReductionType::template evaluate< porder >( dofCoordinate, dv, factor, x, psi );
328 phi.axpy( myfactor * factor, psi );
329 }
330 ++(*dofCoordinate);
331 phi *= (polynomialOrder) / ((RealType)(*dofCoordinate));
332 }
333 }
334
335 template< class LocalDofCoordinateType, class LocalCoordinateType >
336 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
337 const FieldVector< int, 1 > &diffVariable,
338 DomainFieldType factor,
339 const LocalCoordinateType &x,
340 RangeType &phi )
341 {
342 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
343 }
344
345 template< unsigned int porder, class LocalDofCoordinateType, class LocalCoordinateType >
346 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
347 const FieldVector< int, 2 > &diffVariable,
348 DomainFieldType factor,
349 const LocalCoordinateType &x,
350 RangeType &phi )
351 {
352 const DomainFieldType divisor = polynomialOrder;
353 const DomainFieldType myfactor = porder / divisor;
354 const RealType myshift = (porder - polynomialOrder) / divisor;
355
356 FieldVector< int, 1 > dv0( diffVariable[ 0 ] );
357 FieldVector< int, 1 > dv1( diffVariable[ 1 ] );
358
359 if( LagrangePointType :: useDimReduction( dofCoordinate ) )
360 {
361 if( ((unsigned int)diffVariable[ 0 ] != LocalDofCoordinateType::index)
362 && ((unsigned int)diffVariable[ 1 ] != LocalDofCoordinateType::index) )
363 DimensionReductionType::evaluate( dofCoordinate.base(), diffVariable, myfactor * factor, x.base(), phi );
364 else
365 phi = 0;
366
367 const unsigned int height = LagrangePointType::height( dofCoordinate );
368 for( unsigned int i = 0; i < height; ++i )
369 {
370 ++(*dofCoordinate);
371 RangeType psi;
372 evaluate< porder >( dofCoordinate, diffVariable, factor, x, psi );
373 phi -= psi;
374 }
375 (*dofCoordinate) -= height;
376 }
377 else
378 {
379 RangeType psi;
380 --(*dofCoordinate);
381 OrderReductionType::template evaluate< porder >( dofCoordinate, diffVariable, factor, x, phi );
382 phi *= (myfactor * factor * (*x) - myshift);
383
384 if( (unsigned int)diffVariable[ 0 ] == LocalDofCoordinateType::index )
385 {
386 OrderReductionType::template evaluate< porder >( dofCoordinate, dv1, factor, x, psi );
387 phi.axpy( myfactor * factor, psi );
388 }
389
390 if( (unsigned int)diffVariable[ 1 ] == LocalDofCoordinateType::index )
391 {
392 OrderReductionType::template evaluate< porder >( dofCoordinate, dv0, factor, x, psi );
393 phi.axpy( myfactor * factor, psi );
394 }
395
396 ++(*dofCoordinate);
397 phi *= (polynomialOrder) / ((RealType)(*dofCoordinate));
398 }
399 }
400
401 template< class LocalDofCoordinateType, class LocalCoordinateType >
402 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
403 const FieldVector< int, 2 > &diffVariable,
404 DomainFieldType factor,
405 const LocalCoordinateType &x,
406 RangeType &phi )
407 {
408 return evaluate< polynomialOrder >( dofCoordinate, diffVariable, factor, x, phi );
409 }
410
411 template< int diffOrder >
412 void evaluate ( const FieldVector< int, diffOrder > &diffVariable,
413 const DomainType &x,
414 RangeType &phi ) const
415 {
416 const LocalCoordinate< GeometryType, DomainFieldType > xlocal( x );
417 LagrangePointType point( lagrangePoint_ );
418 evaluate( point.dofCoordinate_, diffVariable, 1, xlocal, phi );
419 }
420 };
421
422
423 template< class FunctionSpace, class FirstGeometryType, class SecondGeometryType, unsigned int order >
424 class GenericLagrangeBaseFunction< FunctionSpace, ProductGeometry< FirstGeometryType, SecondGeometryType >, order >
425 {
426 typedef GenericLagrangeBaseFunction< FunctionSpace, ProductGeometry< FirstGeometryType, SecondGeometryType >, order > ThisType;
427
428 static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
429
430 public:
431 typedef FunctionSpace FunctionSpaceType;
432
433 typedef ProductGeometry< FirstGeometryType, SecondGeometryType > GeometryType;
434
435 static constexpr unsigned int polynomialOrder = order;
436
437 typedef GenericLagrangePoint< GeometryType, polynomialOrder > LagrangePointType;
438 static const unsigned int numBaseFunctions = LagrangePointType :: numLagrangePoints;
439
440 typedef typename FunctionSpaceType :: DomainType DomainType;
441 typedef typename FunctionSpaceType :: RangeType RangeType;
442
443 typedef typename FunctionSpaceType :: DomainFieldType DomainFieldType;
444 typedef typename FunctionSpaceType :: RangeFieldType RangeFieldType;
445 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
446
447 private:
448 typedef GenericLagrangeBaseFunction
449 < FunctionSpaceType, FirstGeometryType, polynomialOrder >
450 FirstReductionType;
451 typedef GenericLagrangeBaseFunction
452 < FunctionSpaceType, SecondGeometryType, polynomialOrder >
453 SecondReductionType;
454
455 private:
456 const LagrangePointType lagrangePoint_;
457
458 public:
459 explicit GenericLagrangeBaseFunction ( unsigned int baseNum )
460 : lagrangePoint_( baseNum )
461 {}
462
463 template< class LocalDofCoordinateType, class LocalCoordinateType >
464 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
465 const FieldVector< int, 0 > &diffVariable,
466 DomainFieldType factor,
467 const LocalCoordinateType &x,
468 RangeType &phi )
469 {
470 RangeType psi;
471 FirstReductionType::evaluate( dofCoordinate.first(), diffVariable, factor, x.first(), phi );
472 SecondReductionType::evaluate( dofCoordinate.second(), diffVariable, factor, x.second(), psi );
473 phi[ 0 ] *= psi[ 0 ];
474 }
475
476 template< class LocalDofCoordinateType, class LocalCoordinateType >
477 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
478 const FieldVector< int, 1 > &diffVariable,
479 DomainFieldType factor,
480 const LocalCoordinateType &x,
481 RangeType &phi )
482 {
483 FieldVector< int, 0 > dv;
484 RangeType psi1, psi2;
485
486 FirstReductionType::evaluate( dofCoordinate.first(), diffVariable, factor, x.first(), psi1 );
487 SecondReductionType::evaluate( dofCoordinate.second(), dv, factor, x.second(), psi2 );
488 phi[ 0 ] = psi1[ 0 ] * psi2[ 0 ];
489
490 FirstReductionType::evaluate( dofCoordinate.first(), dv, factor, x.first(), psi1 );
491 SecondReductionType::evaluate( dofCoordinate.second(), diffVariable, factor, x.second(), psi2 );
492 phi[ 0 ] += psi1[ 0 ] * psi2[ 0 ];
493 }
494
495 template< class LocalDofCoordinateType, class LocalCoordinateType >
496 static void evaluate ( LocalDofCoordinateType &dofCoordinate,
497 const FieldVector< int, 2 > &diffVariable,
498 DomainFieldType factor,
499 const LocalCoordinateType &x,
500 RangeType &phi )
501 {
502 FieldVector< int, 0 > dv;
503 FieldVector< int, 1 > dv0( diffVariable[ 0 ] );
504 FieldVector< int, 1 > dv1( diffVariable[ 1 ] );
505 RangeType psi1, psi2;
506
507 FirstReductionType::evaluate( dofCoordinate.first(), diffVariable, factor, x.first(), psi1 );
508 SecondReductionType::evaluate( dofCoordinate.second(), dv, factor, x.second(), psi2 );
509 phi[ 0 ] = psi1[ 0 ] * psi2[ 0 ];
510
511 FirstReductionType::evaluate( dofCoordinate.first(), dv0, factor, x.first(), psi1 );
512 SecondReductionType::evaluate( dofCoordinate.second(), dv1, factor, x.second(), psi2 );
513 phi[ 0 ] += psi1[ 0 ] * psi2[ 0 ];
514
515 FirstReductionType::evaluate( dofCoordinate.first(), dv1, factor, x.first(), psi1 );
516 SecondReductionType::evaluate( dofCoordinate.second(), dv0, factor, x.second(), psi2 );
517 phi[ 0 ] += psi1[ 0 ] * psi2[ 0 ];
518
519 FirstReductionType::evaluate( dofCoordinate.first(), dv, factor, x.first(), psi1 );
520 SecondReductionType::evaluate( dofCoordinate.second(), diffVariable, factor, x.second(), psi2 );
521 phi[ 0 ] += psi1[ 0 ] * psi2[ 0 ];
522 }
523
524 template< int diffOrder >
525 void evaluate ( const FieldVector< int, diffOrder > &diffVariable,
526 const DomainType &x,
527 RangeType &phi ) const
528 {
529 const LocalCoordinate< GeometryType, DomainFieldType > xlocal( x );
530 LagrangePointType point( lagrangePoint_ );
531 evaluate( point.dofCoordinate_, diffVariable, 1, xlocal, phi );
532 }
533 };
534
535 } // namespace Fem
536
537} // namespace Dune
538
539#endif // #ifndef DUNE_FEM_SPACE_LAGRANGE_GENERICBASEFUNCTIONS_HH
FunctionSpaceTraits::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:60
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::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
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
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Implements a vector constructed from a given type representing a field and a compile-time given size.
Dune namespace.
Definition: alignedallocator.hh:13
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 21, 23:30, 2024)