DUNE PDELab (git)

monomialbasis.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
4// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5#ifndef DUNE_MONOMIALBASIS_HH
6#define DUNE_MONOMIALBASIS_HH
7
8#include <vector>
9
12
13#include <dune/geometry/type.hh>
14#include <dune/geometry/topologyfactory.hh>
15
16#include <dune/localfunctions/utility/field.hh>
17#include <dune/localfunctions/utility/multiindex.hh>
18#include <dune/localfunctions/utility/tensor.hh>
19
20namespace Dune
21{
22 /************************************************
23 * Classes for evaluating ''Monomials'' on any order
24 * for all reference element type.
25 * For a simplex topology these are the normal
26 * monomials for cube topologies the bimonomials.
27 * The construction follows the construction of the
28 * generic geometries using tensor products for
29 * prism generation and duffy transform for pyramid
30 * construction.
31 * A derivative argument can be applied, in which case
32 * all derivatives up to the desired order are
33 * evaluated. Note that for higher order derivatives
34 * only the ''lower'' part of the symmetric tensor
35 * is evaluated, e.g., passing derivative equal to 2
36 * to the class will provide the vector
37 * (d/dxdx p, d/dxydx p, d/dydy p,
38 * d/dx p, d/dy p, p)
39 * Important:
40 * So far the computation of the derivatives has not
41 * been fully implemented for general pyramid
42 * construction, i.e., in the case where a pyramid is
43 * build over a non simplex base geometry.
44 *
45 * Central classes:
46 * 1) template< GeometryType::Id geometryId, class F >
47 * class MonomialBasisImpl;
48 * Implementation of the monomial evaluation for
49 * a given topology and field type.
50 * The method evaluate fills a F* vector
51 * 2) template< GeometryType::Id geometryId, class F >
52 * class MonomialBasis
53 * The base class for the static monomial evaluation
54 * providing additional evaluate methods including
55 * one taking std::vector<F>.
56 * 3) template< int dim, class F >
57 * class VirtualMonomialBasis
58 * Virtualization of the MonomialBasis.
59 * 4) template< int dim, class F >
60 * struct MonomialBasisFactory;
61 * A factory class for the VirtualMonomialBasis
62 * 5) template< int dim, class F >
63 * struct MonomialBasisProvider
64 * A singleton container for the virtual monomial
65 * basis
66 ************************************************/
67
68 // Internal Forward Declarations
69 // -----------------------------
70
71 template< GeometryType::Id geometryId >
72 class MonomialBasisSize;
73
74 template< GeometryType::Id geometryId, class F >
75 class MonomialBasis;
76
77
78
79 // MonomialBasisSize
80 // -----------------
81
82 template< GeometryType::Id geometryId >
83 class MonomialBasisSize
84 {
85 typedef MonomialBasisSize< geometryId > This;
86
87 public:
88 static This &instance ()
89 {
90 static This _instance;
91 return _instance;
92 }
93
94 unsigned int maxOrder_;
95
96 // sizes_[ k ]: number of basis functions of exactly order k
97 mutable unsigned int *sizes_;
98
99 // numBaseFunctions_[ k ] = sizes_[ 0 ] + ... + sizes_[ k ]
100 mutable unsigned int *numBaseFunctions_;
101
102 MonomialBasisSize ()
103 : maxOrder_( 0 ),
104 sizes_( 0 ),
105 numBaseFunctions_( 0 )
106 {
107 computeSizes( 2 );
108 }
109
110 ~MonomialBasisSize ()
111 {
112 delete[] sizes_;
113 delete[] numBaseFunctions_;
114 }
115
116 unsigned int operator() ( const unsigned int order ) const
117 {
118 return numBaseFunctions_[ order ];
119 }
120
121 unsigned int maxOrder() const
122 {
123 return maxOrder_;
124 }
125
126 void computeSizes ( unsigned int order )
127 {
128 if (order <= maxOrder_)
129 return;
130
131 maxOrder_ = order;
132
133 delete[] sizes_;
134 delete[] numBaseFunctions_;
135 sizes_ = new unsigned int[ order+1 ];
136 numBaseFunctions_ = new unsigned int[ order+1 ];
137
138 constexpr GeometryType geometry = geometryId;
139 constexpr auto dim = geometry.dim();
140
141 sizes_[ 0 ] = 1;
142 for( unsigned int k = 1; k <= order; ++k )
143 sizes_[ k ] = 0;
144
145 std::fill(numBaseFunctions_, numBaseFunctions_+order+1, 1);
146
147 for( int codim=dim-1; codim>=0; codim--)
148 {
149 if (Impl::isPrism(geometry.id(),dim,codim))
150 {
151 for( unsigned int k = 1; k <= order; ++k )
152 {
153 sizes_[ k ] = numBaseFunctions_[ k ] + k*sizes_[ k ];
154 numBaseFunctions_[ k ] = numBaseFunctions_[ k-1 ] + sizes_[ k ];
155 }
156 }
157 else
158 {
159 for( unsigned int k = 1; k <= order; ++k )
160 {
161 sizes_[ k ] = numBaseFunctions_[ k ];
162 numBaseFunctions_[ k ] = numBaseFunctions_[ k-1 ] + sizes_[ k ];
163 }
164 }
165 }
166 }
167 };
168
169
170
171 // MonomialBasisHelper
172 // -------------------
173
174
175 template< int mydim, int dim, class F >
176 struct MonomialBasisHelper
177 {
178 typedef MonomialBasisSize< GeometryTypes::simplex(mydim).toId() > MySize;
179 typedef MonomialBasisSize< GeometryTypes::simplex(dim).toId() > Size;
180
181 static void copy ( const unsigned int deriv, F *&wit, F *&rit,
182 const unsigned int numBaseFunctions, const F &z )
183 {
184 // n(d,k) = size<k>[d];
185 MySize &mySize = MySize::instance();
186 Size &size = Size::instance();
187
188 const F *const rend = rit + size( deriv )*numBaseFunctions;
189 for( ; rit != rend; )
190 {
191 F *prit = rit;
192
193 *wit = z * *rit;
194 ++rit, ++wit;
195
196 for( unsigned d = 1; d <= deriv; ++d )
197 {
198 #ifndef NDEBUG
199 const F *const derivEnd = rit + mySize.sizes_[ d ];
200 #endif
201
202 {
203 const F *const drend = rit + mySize.sizes_[ d ] - mySize.sizes_[ d-1 ];
204 for( ; rit != drend ; ++rit, ++wit )
205 *wit = z * *rit;
206 }
207
208 for (unsigned int j=1; j<d; ++j)
209 {
210 const F *const drend = rit + mySize.sizes_[ d-j ] - mySize.sizes_[ d-j-1 ];
211 for( ; rit != drend ; ++prit, ++rit, ++wit )
212 *wit = F(j) * *prit + z * *rit;
213 }
214 *wit = F(d) * *prit + z * *rit;
215 ++prit, ++rit, ++wit;
216 assert(derivEnd == rit);
217 rit += size.sizes_[d] - mySize.sizes_[d];
218 prit += size.sizes_[d-1] - mySize.sizes_[d-1];
219 const F *const emptyWitEnd = wit + size.sizes_[d] - mySize.sizes_[d];
220 for ( ; wit != emptyWitEnd; ++wit )
221 *wit = Zero<F>();
222 }
223 }
224 }
225 };
226
227
228
229 // MonomialBasisImpl
230 // -----------------
231
232 template< GeometryType::Id geometryId, class F>
233 class MonomialBasisImpl
234 {
235 public:
236 typedef F Field;
237
238 static constexpr GeometryType geometry = geometryId;
239
240 static const unsigned int dimDomain = geometry.dim();
241
242 typedef FieldVector< Field, dimDomain > DomainVector;
243
244 private:
245 friend class MonomialBasis< geometryId, Field >;
246
247 MonomialBasisImpl ()
248 {}
249
250 template< int dimD >
251 void evaluate ( const unsigned int deriv, const unsigned int order,
252 const FieldVector< Field, dimD > &x,
253 const unsigned int block, const unsigned int *const offsets,
254 Field *const values ) const
255 {
256 //start with vertex
257 *values = Unity< F >();
258 F *const end = values + block;
259 for( Field *it = values+1 ; it != end; ++it )
260 *it = Zero< F >();
261
263
264 if constexpr ( geometry == gt)
265 return;
266 else
267 evaluate<gt,dimD>(deriv, order, x, block, offsets, values );
268 }
269
270 template<GeometryType::Id baseGeometryId, int dimD >
271 void evaluate ( const unsigned int deriv, const unsigned int order,
272 const FieldVector< Field, dimD > &x,
273 const unsigned int block, const unsigned int *const offsets,
274 Field *const values ) const
275 {
276
277 static constexpr GeometryType baseGeometry = baseGeometryId;
278
279 auto constexpr isPrismatic = geometry.isPrismatic(baseGeometry.dim());
280
281 // compute
282 typedef MonomialBasisHelper< baseGeometry.dim() + 1, dimD, Field > Helper;
283 typedef MonomialBasisSize<baseGeometryId> BaseSize;
284
285 const BaseSize &size = BaseSize::instance();
286 const_cast<BaseSize&>(size).computeSizes(order);
287
288 const Field &z = x[ baseGeometry.dim() ];
289
290 Field *row0 = values;
291 for( unsigned int k = 1; k <= order; ++k )
292 {
293 Field *row1 = values + block*offsets[ k-1 ];
294 Field *wit = row1 + block*size.sizes_[ k ];
295 if constexpr ( isPrismatic )
296 Helper::copy( deriv, wit, row1, k*size.sizes_[ k ], z );
297 Helper::copy( deriv, wit, row0, size( k-1 ), z );
298 row0 = row1;
299 }
300
301 // stop if desired dimension is reached
302 if constexpr( baseGeometry.dim() == dimDomain-1)
303 return;
304 else
305 {
306 constexpr GeometryType nextGeometry = isPrismatic ? GeometryTypes::prismaticExtension(baseGeometry)
307 : GeometryTypes::conicalExtension(baseGeometry);
308
309 evaluate<nextGeometry.toId(),dimD>(deriv, order, x, block, offsets, values );
310 }
311 }
312
313 void integrate ( const unsigned int order,
314 const unsigned int *const offsets,
315 Field *const values ) const
316 {
317 //start with vertex
318 values[ 0 ] = Unity< Field >();
319 static constexpr GeometryType gt = GeometryTypes::vertex;
320
321 if constexpr ( geometry == gt)
322 return;
323 else
324 integrate<gt>(order, offsets, values);
325 }
326
327 template<GeometryType::Id baseGeometryId>
328 void integrate ( const unsigned int order,
329 const unsigned int *const offsets,
330 Field *const values) const
331 {
332 static constexpr GeometryType baseGeometry = baseGeometryId;
333
334 auto constexpr isPrismatic = geometry.isPrismatic(baseGeometry.dim());
335
336 // decide which kind of integration should be performed
337 if constexpr ( isPrismatic )
338 integratePrismatic<baseGeometry>(order, offsets, values);
339 else
340 integrateConical<baseGeometry>(order, offsets, values);
341
342 // stop if the desired dimension is reached
343 if constexpr( baseGeometry.dim() == dimDomain-1)
344 return;
345 else
346 {
347 static constexpr GeometryType nextGeometry = (isPrismatic ? GeometryTypes::prismaticExtension(baseGeometry)
348 : GeometryTypes::conicalExtension(baseGeometry));
349
350 integrate<nextGeometry.toId()>(order, offsets, values);
351 }
352
353 }
354
355 template<GeometryType::Id baseGeometryId>
356 void integratePrismatic ( const unsigned int order,
357 const unsigned int *const offsets,
358 Field *const values ) const
359 {
360 typedef MonomialBasisSize<baseGeometryId> BaseSize;
361 static const BaseSize &size = BaseSize::instance();
362 const unsigned int *const baseSizes = size.sizes_;
363
364 static constexpr GeometryType baseGeometry = baseGeometryId;
365 static constexpr GeometryType nextGeometry = GeometryTypes::prismaticExtension(baseGeometry);
366
367 typedef MonomialBasisSize<nextGeometry.toId()> Size;
368 static const Size &mySize = Size::instance();
369
370 Field *row0 = values;
371 for( unsigned int k = 1; k <= order; ++k )
372 {
373 Field *const row1begin = values + offsets[ k-1 ];
374 Field *const row1End = row1begin + mySize.sizes_[ k ];
375 assert( (unsigned int)(row1End - values) <= offsets[ k ] );
376
377 Field *row1 = row1begin;
378 Field *it = row1begin + baseSizes[ k ];
379 for( unsigned int j = 1; j <= k; ++j )
380 {
381 Field *const end = it + baseSizes[ k ];
382 assert( (unsigned int)(end - values) <= offsets[ k ] );
383 for( ; it != end; ++row1, ++it )
384 *it = (Field( j ) / Field( j+1 )) * (*row1);
385 }
386 for( ; it != row1End; ++row0, ++it )
387 *it = (Field( k ) / Field( k+1 )) * (*row0);
388 row0 = row1;
389 }
390 }
391
392
393 template<GeometryType::Id baseGeometryId>
394 void integrateConical ( const unsigned int order,
395 const unsigned int *const offsets,
396 Field *const values) const
397 {
398 typedef MonomialBasisSize<baseGeometryId> BaseSize;
399 static const BaseSize &size = BaseSize::instance();
400 const unsigned int *const baseSizes = size.sizes_;
401
402 static constexpr GeometryType baseGeometry = baseGeometryId;
403
404 {
405 Field *const col0End = values + baseSizes[ 0 ];
406 for( Field *it = values; it != col0End; ++it )
407 *it *= Field( 1 ) / Field( int(baseGeometry.dim()+1) );
408 }
409
410 Field *row0 = values;
411 for( unsigned int k = 1; k <= order; ++k )
412 {
413 const Field factor = (Field( 1 ) / Field( k + baseGeometry.dim()+1));
414
415 Field *const row1 = values+offsets[ k-1 ];
416 Field *const col0End = row1 + baseSizes[ k ];
417 Field *it = row1;
418 for( ; it != col0End; ++it )
419 *it *= factor;
420 for( unsigned int i = 1; i <= k; ++i )
421 {
422 Field *const end = it + baseSizes[ k-i ];
423 assert( (unsigned int)(end - values) <= offsets[ k ] );
424 for( ; it != end; ++row0, ++it )
425 *it = (*row0) * (Field( i ) * factor);
426 }
427 row0 = row1;
428 }
429 }
430
431 };
432
433
434 // MonomialBasis
435 // -------------
436
437 template< GeometryType::Id geometryId, class F >
438 class MonomialBasis
439 : public MonomialBasisImpl< geometryId, F >
440 {
441 static constexpr GeometryType geometry = geometryId;
442 typedef MonomialBasis< geometryId, F > This;
443 typedef MonomialBasisImpl< geometryId, F > Base;
444
445 public:
446 static const unsigned int dimension = Base::dimDomain;
447 static const unsigned int dimRange = 1;
448
449 typedef typename Base::Field Field;
450
451 typedef typename Base::DomainVector DomainVector;
452
453 typedef Dune::FieldVector<Field,dimRange> RangeVector;
454
455 typedef MonomialBasisSize<geometryId> Size;
456
457 MonomialBasis (unsigned int order)
458 : Base(),
459 order_(order),
460 size_(Size::instance())
461 {
462 assert(order<=1024); // avoid wrapping of unsigned int (0-1) order=1024 is quite high...)
463 }
464
465 const unsigned int *sizes ( unsigned int order ) const
466 {
467 size_.computeSizes( order );
468 return size_.numBaseFunctions_;
469 }
470
471 const unsigned int *sizes () const
472 {
473 return sizes( order_ );
474 }
475
476 unsigned int size () const
477 {
478 size_.computeSizes( order_ );
479 return size_( order_ );
480 }
481
482 unsigned int derivSize ( const unsigned int deriv ) const
483 {
484 MonomialBasisSize< GeometryTypes::simplex(dimension).toId() >::instance().computeSizes( deriv );
485 return MonomialBasisSize< GeometryTypes::simplex(dimension).toId() >::instance() ( deriv );
486 }
487
488 unsigned int order () const
489 {
490 return order_ ;
491 }
492
493 unsigned int topologyId ( ) const
494 {
495 return geometry.id();
496 }
497
498 void evaluate ( const unsigned int deriv, const DomainVector &x,
499 Field *const values ) const
500 {
501 Base::evaluate( deriv, order_, x, derivSize( deriv ), sizes( order_ ), values );
502 }
503
504 template <unsigned int deriv>
505 void evaluate ( const DomainVector &x,
506 Field *const values ) const
507 {
508 evaluate( deriv, x, values );
509 }
510
511 template<unsigned int deriv, class Vector >
512 void evaluate ( const DomainVector &x,
513 Vector &values ) const
514 {
515 evaluate<deriv>(x,&(values[0]));
516 }
517 template<unsigned int deriv, DerivativeLayoutNS::DerivativeLayout layout >
518 void evaluate ( const DomainVector &x,
519 Derivatives<Field,dimension,1,deriv,layout> *values ) const
520 {
521 evaluate<deriv>(x,&(values->block()));
522 }
523 template< unsigned int deriv >
524 void evaluate ( const DomainVector &x,
525 FieldVector<Field,Derivatives<Field,dimension,1,deriv,DerivativeLayoutNS::value>::size> *values ) const
526 {
527 evaluate(0,x,&(values[0][0]));
528 }
529
530 template<class Vector >
531 void evaluate ( const DomainVector &x,
532 Vector &values ) const
533 {
534 evaluate<0>(x,&(values[0]));
535 }
536
537 template< class DVector, class RVector >
538 void evaluate ( const DVector &x, RVector &values ) const
539 {
540 assert( DVector::dimension == dimension);
541 DomainVector bx;
542 for( int d = 0; d < dimension; ++d )
543 field_cast( x[ d ], bx[ d ] );
544 evaluate<0>( bx, values );
545 }
546
547 void integrate ( Field *const values ) const
548 {
549 Base::integrate( order_, sizes( order_ ), values );
550 }
551 template <class Vector>
552 void integrate ( Vector &values ) const
553 {
554 integrate( &(values[ 0 ]) );
555 }
556 private:
557 MonomialBasis(const This&);
558 This& operator=(const This&);
559 unsigned int order_;
560 Size &size_;
561 };
562
563
564
565 // StdMonomialBasis
566 // ----------------
567
568 template< int dim,class F >
569 class StandardMonomialBasis
570 : public MonomialBasis< GeometryTypes::simplex(dim).toId() , F >
571 {
572 typedef StandardMonomialBasis< dim, F > This;
573 typedef MonomialBasis< GeometryTypes::simplex(dim).toId(), F > Base;
574
575 public:
576 static constexpr GeometryType geometry = GeometryTypes::simplex(dim);
577 static const int dimension = dim;
578
579 StandardMonomialBasis ( unsigned int order )
580 : Base( order )
581 {}
582 };
583
584
585
586 // StandardBiMonomialBasis
587 // -----------------------
588
589 template< int dim, class F >
590 class StandardBiMonomialBasis
591 : public MonomialBasis< GeometryTypes::cube(dim).toId() , F >
592 {
593 typedef StandardBiMonomialBasis< dim, F > This;
594 typedef MonomialBasis< GeometryTypes::cube(dim).toId() , F > Base;
595
596 public:
597 static constexpr GeometryType geometry = GeometryTypes::cube(dim);
598 static const int dimension = dim;
599
600 StandardBiMonomialBasis ( unsigned int order )
601 : Base( order )
602 {}
603 };
604
605 // -----------------------------------------------------------
606 // -----------------------------------------------------------
607 // VirtualMonomialBasis
608 // -------------------
609
610 template< int dim, class F >
611 class VirtualMonomialBasis
612 {
613 typedef VirtualMonomialBasis< dim, F > This;
614
615 public:
616 typedef F Field;
617 typedef F StorageField;
618 static const int dimension = dim;
619 static const unsigned int dimRange = 1;
620
621 typedef FieldVector<Field,dimension> DomainVector;
622 typedef FieldVector<Field,dimRange> RangeVector;
623
624 explicit VirtualMonomialBasis(const GeometryType& gt,
625 unsigned int order)
626 : order_(order), geometry_(gt) {}
627
628 virtual ~VirtualMonomialBasis() {}
629
630 virtual const unsigned int *sizes ( ) const = 0;
631
632 unsigned int size ( ) const
633 {
634 return sizes( )[ order_ ];
635 }
636
637 unsigned int order () const
638 {
639 return order_;
640 }
641
642 GeometryType type() const
643 {
644 return geometry_;
645 }
646
647 virtual void evaluate ( const unsigned int deriv, const DomainVector &x,
648 Field *const values ) const = 0;
649 template < unsigned int deriv >
650 void evaluate ( const DomainVector &x,
651 Field *const values ) const
652 {
653 evaluate( deriv, x, values );
654 }
655 template < unsigned int deriv, int size >
656 void evaluate ( const DomainVector &x,
657 Dune::FieldVector<Field,size> *const values ) const
658 {
659 evaluate( deriv, x, &(values[0][0]) );
660 }
661 template<unsigned int deriv, DerivativeLayoutNS::DerivativeLayout layout >
662 void evaluate ( const DomainVector &x,
663 Derivatives<Field,dimension,1,deriv,layout> *values ) const
664 {
665 evaluate<deriv>(x,&(values->block()));
666 }
667 template <unsigned int deriv, class Vector>
668 void evaluate ( const DomainVector &x,
669 Vector &values ) const
670 {
671 evaluate<deriv>( x, &(values[ 0 ]) );
672 }
673 template< class Vector >
674 void evaluate ( const DomainVector &x,
675 Vector &values ) const
676 {
677 evaluate<0>(x,values);
678 }
679 template< class DVector, class RVector >
680 void evaluate ( const DVector &x, RVector &values ) const
681 {
682 assert( DVector::dimension == dimension);
683 DomainVector bx;
684 for( int d = 0; d < dimension; ++d )
685 field_cast( x[ d ], bx[ d ] );
686 evaluate<0>( bx, values );
687 }
688 template< unsigned int deriv, class DVector, class RVector >
689 void evaluate ( const DVector &x, RVector &values ) const
690 {
691 assert( DVector::dimension == dimension);
692 DomainVector bx;
693 for( int d = 0; d < dimension; ++d )
694 field_cast( x[ d ], bx[ d ] );
695 evaluate<deriv>( bx, values );
696 }
697
698 virtual void integrate ( Field *const values ) const = 0;
699 template <class Vector>
700 void integrate ( Vector &values ) const
701 {
702 integrate( &(values[ 0 ]) );
703 }
704 protected:
705 unsigned int order_;
706 GeometryType geometry_;
707 };
708
709 template< GeometryType::Id geometryId, class F >
710 class VirtualMonomialBasisImpl
711 : public VirtualMonomialBasis< static_cast<GeometryType>(geometryId).dim(), F >
712 {
713 static constexpr GeometryType geometry = geometryId;
714 typedef VirtualMonomialBasis< geometry.dim(), F > Base;
715 typedef VirtualMonomialBasisImpl< geometryId, F > This;
716
717 public:
718 typedef typename Base::Field Field;
719 typedef typename Base::DomainVector DomainVector;
720
721 VirtualMonomialBasisImpl(unsigned int order)
722 : Base(geometry,order), basis_(order)
723 {}
724
725 const unsigned int *sizes ( ) const
726 {
727 return basis_.sizes(order_);
728 }
729
730 void evaluate ( const unsigned int deriv, const DomainVector &x,
731 Field *const values ) const
732 {
733 basis_.evaluate(deriv,x,values);
734 }
735
736 void integrate ( Field *const values ) const
737 {
738 basis_.integrate(values);
739 }
740
741 private:
742 MonomialBasis<geometryId,Field> basis_;
743 using Base::order_;
744 };
745
746 // MonomialBasisFactory
747 // --------------------
748
749 template< int dim, class F >
750 struct MonomialBasisFactory
751 {
752 static const unsigned int dimension = dim;
753 typedef F StorageField;
754
755 typedef unsigned int Key;
756 typedef const VirtualMonomialBasis< dimension, F > Object;
757
758 template < int dd, class FF >
759 struct EvaluationBasisFactory
760 {
761 typedef MonomialBasisFactory<dd,FF> Type;
762 };
763
764 template< GeometryType::Id geometryId >
765 static Object* create ( const Key &order )
766 {
767 return new VirtualMonomialBasisImpl< geometryId, StorageField >( order );
768 }
769 static void release( Object *object ) { delete object; }
770 };
771
772
773
774 // MonomialBasisProvider
775 // ---------------------
776
777 template< int dim, class SF >
778 struct MonomialBasisProvider
779 : public TopologySingletonFactory< MonomialBasisFactory< dim, SF > >
780 {
781 static const unsigned int dimension = dim;
782 typedef SF StorageField;
783 template < int dd, class FF >
784 struct EvaluationBasisFactory
785 {
786 typedef MonomialBasisProvider<dd,FF> Type;
787 };
788 };
789
790}
791
792#endif
vector space out of a tensor product of fields.
Definition: fvector.hh:91
constexpr Id toId() const
Create an Id representation of this GeometryType.
Definition: type.hh:210
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Implements a vector constructed from a given type representing a field and a compile-time given size.
bool gt(const T &first, const T &second, typename EpsilonType< T >::Type epsilon)
test if first greater than second
Definition: float_cmp.cc:158
constexpr GeometryType cube(unsigned int dim)
Returns a GeometryType representing a hypercube of dimension dim.
Definition: type.hh:462
constexpr GeometryType simplex(unsigned int dim)
Returns a GeometryType representing a simplex of dimension dim.
Definition: type.hh:453
constexpr GeometryType vertex
GeometryType representing a vertex.
Definition: type.hh:492
constexpr GeometryType prismaticExtension(const GeometryType &gt)
Return GeometryType of a prismatic construction with gt as base
Definition: type.hh:483
constexpr GeometryType conicalExtension(const GeometryType &gt)
Return GeometryType of a conical construction with gt as base
Definition: type.hh:477
Dune namespace.
Definition: alignedallocator.hh:13
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:159
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
A unique label for each type of element that can occur in a grid.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 24, 23:30, 2024)