1#ifndef DUNE_FEM_SPACE_LAGRANGE_GENERICLAGRANGEPOINTS_HH
2#define DUNE_FEM_SPACE_LAGRANGE_GENERICLAGRANGEPOINTS_HH
11#include "genericgeometry.hh"
20 template<
class ImplType,
unsigned int order,
bool bottom = true >
21 class GenericLagrangePoint;
24 template<
unsigned int order,
bool bottom >
25 class GenericLagrangePoint< PointGeometry, order, bottom >
29 static const unsigned int dimension = GeometryType::dimension;
31 typedef LocalCoordinate< GeometryType, unsigned int > DofCoordinateType;
33 static const unsigned int polynomialOrder = order;
35 template<
class,
unsigned int,
bool >
36 friend class GenericLagrangePoint;
38 template<
class,
class,
unsigned int >
39 friend class GenericLagrangeBaseFunction;
42 typedef GenericLagrangePoint< GeometryType, polynomialOrder > ThisType;
45 static const unsigned int numLagrangePoints = 1;
48 DofCoordinateType dofCoordinate_;
51 template<
unsigned int codim >
55 GenericLagrangePoint (
unsigned int index )
57 dofCoordinate( index, dofCoordinate_ );
60 GenericLagrangePoint (
const ThisType &point )
61 : dofCoordinate_( point.dofCoordinate_ )
64 template<
class LocalCoordinateType >
65 static void dofSubEntity ( LocalCoordinateType &coordinate,
67 unsigned int &subEntity )
73 template<
class LocalCoordinateType >
74 static void dofSubEntity ( LocalCoordinateType &coordinate,
76 unsigned int &subEntity,
77 unsigned int &dofNumber )
84 void dofSubEntity (
unsigned int &codim,
85 unsigned int &subEntity )
87 dofSubEntity( dofCoordinate_, codim, subEntity );
90 void dofSubEntity (
unsigned int &codim,
91 unsigned int &subEntity,
92 unsigned int &dofNumber )
94 dofSubEntity( dofCoordinate_, codim, subEntity, dofNumber );
97 static unsigned int entityDofNumber (
unsigned int codim,
98 unsigned int subEntity,
105 template<
class LocalCoordinateType >
106 static unsigned int height ( LocalCoordinateType &coordinate )
108 return polynomialOrder;
111 unsigned int height ()
113 return height( dofCoordinate_ );
116 template<
class FieldType >
117 void local ( FieldVector< FieldType, dimension > &coordinate )
const
119 const FieldType factor = FieldType( 1 ) / FieldType( polynomialOrder );
120 for(
unsigned int i = 0; i < dimension; ++i )
121 coordinate[ i ] = factor * dofCoordinate_[ i ];
130 static unsigned int maxDofs (
unsigned int codim )
132 return ((codim == 0) ? 1 : 0);
142 static unsigned int numDofs (
unsigned int codim,
unsigned int subEntity )
144 return ((codim == 0) ? 1 : 0);
153 static unsigned int numDofs (
unsigned int codim )
155 return ((codim == 0) ? 1 : 0);
159 template<
class LocalCoordinateType >
160 static inline void dofCoordinate (
unsigned int index,
161 LocalCoordinateType &coordinate )
163 assert( index <= numLagrangePoints );
170 template<
unsigned int order,
bool bottom >
171 template<
unsigned int codim >
172 struct GenericLagrangePoint< PointGeometry, order, bottom >::Codim
174 static unsigned int maxDofs ()
176 return ((codim == 0) ? 1 : 0);
182 template<
class BaseGeometry,
bool bottom >
183 class GenericLagrangePoint< PyramidGeometry< BaseGeometry >, 0, bottom >
186 typedef BaseGeometry BaseGeometryType;
187 typedef PyramidGeometry< BaseGeometryType >
GeometryType;
189 static const unsigned int dimension = GeometryType::dimension;
191 typedef LocalCoordinate< GeometryType, unsigned int > DofCoordinateType;
193 static const unsigned int polynomialOrder = 0;
195 template<
class,
unsigned int,
bool >
196 friend class GenericLagrangePoint;
198 template<
class,
class,
unsigned int >
199 friend class GenericLagrangeBaseFunction;
202 typedef GenericLagrangePoint< GeometryType, polynomialOrder > ThisType;
205 static const unsigned int numLagrangePoints = 1;
208 DofCoordinateType dofCoordinate_;
211 template<
unsigned int codim >
215 GenericLagrangePoint (
unsigned int index )
217 dofCoordinate( index, dofCoordinate_ );
220 GenericLagrangePoint (
const ThisType &point )
221 : dofCoordinate_( point.dofCoordinate_ )
224 template<
class LocalCoordinateType >
225 static void dofSubEntity ( LocalCoordinateType &coordinate,
227 unsigned int &subEntity )
229 codim = (bottom ? 0 : dimension);
233 template<
class LocalCoordinateType >
234 static void dofSubEntity ( LocalCoordinateType &coordinate,
236 unsigned int &subEntity,
237 unsigned int &dofNumber )
239 codim = (bottom ? 0 : dimension);
244 void dofSubEntity (
unsigned int &codim,
245 unsigned int &subEntity )
247 dofSubEntity( dofCoordinate_, codim, subEntity );
250 void dofSubEntity (
unsigned int &codim,
251 unsigned int &subEntity,
252 unsigned int &dofNumber )
254 dofSubEntity( dofCoordinate_, codim, subEntity, dofNumber );
257 static unsigned int entityDofNumber (
unsigned int codim,
258 unsigned int subEntity,
264 template<
class LocalCoordinateType >
265 static unsigned int height ( LocalCoordinateType &coordinate )
267 return polynomialOrder;
270 unsigned int height ()
272 return height( dofCoordinate_ );
275 template<
class FieldType >
276 void local ( FieldVector< FieldType, dimension > &coordinate )
const
278 const FieldType factor = FieldType( 1 ) / FieldType( polynomialOrder );
279 for(
unsigned int i = 0; i < dimension; ++i )
280 coordinate[ i ] = factor * dofCoordinate_[ i ];
289 static unsigned int maxDofs (
unsigned int codim )
292 return ((codim == 0) ? 1 : 0);
294 return ((codim == dimension) ? 1 : 0);
304 static unsigned int numDofs (
unsigned int codim,
unsigned int subEntity )
307 return ((codim == 0) ? 1 : 0);
309 return ((codim == dimension) ? 1 : 0);
318 static unsigned int numDofs (
unsigned int codim )
321 return ((codim == 0) ? 1 : 0);
323 return ((codim == dimension) ? 1 : 0);
327 template<
class LocalCoordinateType >
328 static void dofCoordinate (
unsigned int index,
329 LocalCoordinateType &coordinate )
331 assert( index <= numLagrangePoints );
338 template<
class BaseGeometry,
bool bottom >
339 template<
unsigned int codim >
340 struct GenericLagrangePoint< PyramidGeometry< BaseGeometry >, 0, bottom >::Codim
342 static inline unsigned int maxDofs ()
345 return ((codim == 0) ? 1 : 0);
347 return ((codim == dimension) ? 1 : 0);
353 template<
class BaseGeometry,
unsigned int order,
bool bottom >
354 class GenericLagrangePoint< PyramidGeometry< BaseGeometry >, order, bottom >
357 typedef BaseGeometry BaseGeometryType;
358 typedef PyramidGeometry< BaseGeometryType >
GeometryType;
359 static const unsigned int dimension = GeometryType :: dimension;
361 typedef LocalCoordinate< GeometryType, unsigned int > DofCoordinateType;
363 static const unsigned int polynomialOrder = order;
365 template<
class,
unsigned int,
bool >
366 friend class GenericLagrangePoint;
368 template<
class,
class,
unsigned int >
369 friend class GenericLagrangeBaseFunction;
372 typedef GenericLagrangePoint< GeometryType, polynomialOrder > ThisType;
374 typedef GenericLagrangePoint<
GeometryType, polynomialOrder - 1,
false >
376 typedef GenericLagrangePoint< BaseGeometryType, polynomialOrder >
377 DimensionReductionType;
380 static const unsigned int numLagrangePoints
381 = DimensionReductionType::numLagrangePoints + OrderReductionType::numLagrangePoints;
384 DofCoordinateType dofCoordinate_;
387 template<
unsigned int codim >
391 GenericLagrangePoint (
unsigned int index )
393 dofCoordinate( index, dofCoordinate_ );
396 GenericLagrangePoint (
const ThisType &point )
397 : dofCoordinate_( point.dofCoordinate_ )
400 template<
class LocalCoordinateType >
401 static void dofSubEntity ( LocalCoordinateType &coordinate,
403 unsigned int &subEntity )
405 if( !useDimReduction( coordinate ) )
408 OrderReductionType::dofSubEntity( coordinate, codim, subEntity );
411 if( bottom && (codim > 0) )
412 subEntity += BaseGeometryType::numSubEntities( codim - 1 );
416 DimensionReductionType::dofSubEntity( coordinate.base(), codim, subEntity );
422 template<
class LocalCoordinateType >
423 static void dofSubEntity ( LocalCoordinateType &coordinate,
425 unsigned int &subEntity,
426 unsigned int &dofNumber )
428 if( !useDimReduction( coordinate ) )
431 OrderReductionType::template dofSubEntity( coordinate, codim, subEntity, dofNumber );
435 subEntity += (codim > 0 ? BaseGeometryType :: numSubEntities( codim - 1 ) : 0);
437 dofNumber += DimensionReductionType::numDofs( codim, subEntity );
441 DimensionReductionType::dofSubEntity( coordinate.base(), codim, subEntity, dofNumber );
448 void dofSubEntity (
unsigned int &codim,
449 unsigned int &subEntity )
451 dofSubEntity( dofCoordinate_, codim, subEntity );
454 void dofSubEntity (
unsigned int &codim,
455 unsigned int &subEntity,
456 unsigned int &dofNumber )
458 dofSubEntity( dofCoordinate_, codim, subEntity, dofNumber );
461 static unsigned int entityDofNumber (
unsigned int codim,
462 unsigned int subEntity,
468 return OrderReductionType::entityDofNumber( codim, subEntity, dof )
469 + DimensionReductionType::numLagrangePoints;
471 const unsigned int numBaseSubEntities
472 = BaseGeometryType::numSubEntities( codim - 1 );
473 if( subEntity >= numBaseSubEntities )
474 return OrderReductionType::entityDofNumber( codim, subEntity - numBaseSubEntities, dof )
475 + DimensionReductionType::numLagrangePoints;
477 return DimensionReductionType::entityDofNumber( codim - 1, subEntity, dof );
481 const unsigned int numBaseEntityDofs
482 = DimensionReductionType::numDofs( codim, subEntity );
483 if( dof >= numBaseEntityDofs )
484 return OrderReductionType::entityDofNumber( codim, subEntity, dof - numBaseEntityDofs )
485 + DimensionReductionType::numLagrangePoints;
487 return DimensionReductionType::entityDofNumber( codim, subEntity, dof );
491 template<
class LocalCoordinateType >
492 static unsigned int height ( LocalCoordinateType &coordinate )
494 if( !useDimReduction( coordinate ) ) {
496 unsigned int h = OrderReductionType :: height( coordinate );
500 return DimensionReductionType :: height( coordinate.base() );
503 unsigned int height ()
505 return height( dofCoordinate_ );
508 template<
class FieldType >
509 void local ( FieldVector< FieldType, dimension > &coordinate )
const
511 const FieldType factor = FieldType( 1 ) / FieldType( polynomialOrder );
512 for(
unsigned int i = 0; i < dimension; ++i )
513 coordinate[ i ] = factor * dofCoordinate_[ i ];
522 static unsigned int maxDofs (
unsigned int codim )
524 const unsigned int maxOrderDofs
525 = OrderReductionType::maxDofs( codim );
527 if( bottom && (codim == 0) )
530 const unsigned int maxDimDofs
531 = DimensionReductionType::maxDofs( bottom ? codim - 1 : codim );
533 return (bottom ? std::max( maxDimDofs, maxOrderDofs )
534 : maxDimDofs + maxOrderDofs);
544 static unsigned int numDofs (
unsigned int codim,
unsigned int subEntity )
548 if( bottom && (codim == 0) )
549 return OrderReductionType::numDofs( codim, subEntity );
551 const unsigned int numBaseSubEntities
552 = BaseGeometryType::numSubEntities( codim - 1 );
553 if( subEntity < numBaseSubEntities )
554 return DimensionReductionType::numDofs( codim - 1, subEntity );
556 return OrderReductionType::numDofs( codim, subEntity - numBaseSubEntities );
560 return DimensionReductionType::numDofs( codim, subEntity )
561 + OrderReductionType::numDofs( codim, subEntity );
571 static unsigned int numDofs (
unsigned int codim )
575 const unsigned int orderDofs = OrderReductionType::numDofs( codim );
577 return orderDofs + DimensionReductionType::numDofs( codim - 1 );
583 return DimensionReductionType::numDofs( codim )
584 + OrderReductionType::numDofs( codim );
588 template<
class LocalCoordinateType >
589 static bool useDimReduction (
const LocalCoordinateType &coordinate )
591 return (*coordinate == 0);
595 template<
class LocalCoordinateType >
596 static void dofCoordinate (
unsigned int index,
597 LocalCoordinateType &coordinate )
599 assert( index <= numLagrangePoints );
601 if( index < DimensionReductionType::numLagrangePoints )
604 DimensionReductionType :: dofCoordinate( index, coordinate.base() );
608 const int orderIndex = index - DimensionReductionType::numLagrangePoints;
609 OrderReductionType::dofCoordinate( orderIndex, coordinate );
617 template<
class BaseGeometry,
unsigned int order,
bool bottom >
618 template<
unsigned int codim >
619 struct GenericLagrangePoint< PyramidGeometry< BaseGeometry >, order, bottom >::Codim
621 static unsigned int maxDofs ()
623 const unsigned int maxOrderDofs
624 = OrderReductionType::template Codim< codim >::maxDofs();
626 const unsigned int maxDimDofs
627 = DimensionReductionType::template Codim< (bottom ? codim - 1 : codim) >::maxDofs();
630 return std::max( maxDimDofs, maxOrderDofs );
632 return maxDimDofs + maxOrderDofs;
638 template<
class FirstGeometry,
class SecondGeometry,
unsigned int order,
bool bottom >
639 class GenericLagrangePoint< ProductGeometry< FirstGeometry, SecondGeometry >, order, bottom >
642 typedef ProductGeometry< FirstGeometry, SecondGeometry >
GeometryType;
643 static const unsigned int dimension = GeometryType::dimension;
645 typedef LocalCoordinate< GeometryType, unsigned int > DofCoordinateType;
647 static const unsigned int polynomialOrder = order;
649 template<
class,
unsigned int,
bool >
650 friend class GenericLagrangePoint;
652 template<
class,
class,
unsigned int >
653 friend class GenericLagrangeBaseFunction;
656 typedef GenericLagrangePoint< GeometryType, polynomialOrder > ThisType;
658 typedef GenericLagrangePoint< FirstGeometry, polynomialOrder > FirstReductionType;
659 typedef GenericLagrangePoint< SecondGeometry, polynomialOrder > SecondReductionType;
662 static const unsigned int numLagrangePoints
663 = FirstReductionType::numLagrangePoints * SecondReductionType::numLagrangePoints;
666 DofCoordinateType dofCoordinate_;
669 template<
unsigned int codim,
unsigned int i >
670 struct CodimIterator;
673 template<
unsigned int codim >
677 GenericLagrangePoint (
unsigned int index )
679 dofCoordinate( index, dofCoordinate_ );
682 GenericLagrangePoint (
const ThisType &point )
683 : dofCoordinate_( point.dofCoordinate_ )
686 template<
class LocalCoordinateType >
687 static void dofSubEntity( LocalCoordinateType &coordinate,
689 unsigned int &subEntity )
691 unsigned int firstCodim, secondCodim;
692 unsigned int firstSubEntity, secondSubEntity;
694 FirstReductionType::dofSubEntity( coordinate.first(), firstCodim, firstSubEntity );
695 SecondReductionType::dofSubEntity( coordinate.second(), secondCodim, secondSubEntity );
697 codim = firstCodim + secondCodim;
700 for(
unsigned int i = 0; i < secondCodim; ++i )
701 subEntity += FirstGeometry::numSubEntities( codim-i ) * SecondGeometry::numSubEntities( i );
702 subEntity += firstSubEntity + secondSubEntity * FirstGeometry::numSubEntities( firstCodim );
705 template<
class LocalCoordinateType >
706 static void dofSubEntity( LocalCoordinateType &coordinate,
708 unsigned int &subEntity,
709 unsigned int &dofNumber )
711 unsigned int firstCodim, secondCodim;
712 unsigned int firstSubEntity, secondSubEntity;
713 unsigned int firstDofNumber, secondDofNumber;
715 FirstReductionType::dofSubEntity( coordinate.first(), firstCodim, firstSubEntity, firstDofNumber );
716 SecondReductionType::dofSubEntity( coordinate.second(), secondCodim, secondSubEntity, secondDofNumber );
718 codim = firstCodim + secondCodim;
721 for(
unsigned int i = 0; i < secondCodim; ++i )
722 subEntity += FirstGeometry::numSubEntities( codim-i ) * SecondGeometry::numSubEntities( i );
723 subEntity += firstSubEntity + secondSubEntity * FirstGeometry::numSubEntities( firstCodim );
725 dofNumber = firstDofNumber + secondDofNumber * FirstReductionType::numDofs( firstCodim, firstSubEntity );
728 void dofSubEntity (
unsigned int &codim,
unsigned int &subEntity )
730 dofSubEntity( dofCoordinate_, codim, subEntity );
733 void dofSubEntity (
unsigned int &codim,
734 unsigned int &subEntity,
735 unsigned int &dofNumber )
737 dofSubEntity( dofCoordinate_, codim, subEntity, dofNumber );
740 static unsigned int entityDofNumber (
unsigned int codim,
741 unsigned int subEntity,
742 unsigned int dofNumber )
744 unsigned int firstCodim = codim;
745 unsigned int secondCodim = 0;
746 for( ; secondCodim < codim; --firstCodim, ++secondCodim )
748 const unsigned int num
749 = FirstGeometry::numSubEntities( firstCodim ) * SecondGeometry::numSubEntities( secondCodim );
751 if( subEntity < num )
756 const unsigned int n = FirstGeometry::numSubEntities( firstCodim );
757 const unsigned int firstSubEntity = subEntity % n;
758 const unsigned int secondSubEntity = subEntity / n;
760 const unsigned int m = FirstReductionType::numDofs( firstCodim, firstSubEntity );
761 const unsigned int firstDofNumber = dofNumber % m;
762 const unsigned int secondDofNumber = dofNumber / m;
764 const unsigned int firstEntityDofNumber
765 = FirstReductionType::entityDofNumber( firstCodim, firstSubEntity, firstDofNumber );
766 const unsigned int secondEntityDofNumber
767 = SecondReductionType::entityDofNumber( secondCodim, secondSubEntity, secondDofNumber );
769 return firstEntityDofNumber + secondEntityDofNumber * FirstReductionType::numLagrangePoints;
772 template<
class LocalCoordinateType >
773 static unsigned int height ( LocalCoordinateType &coordinate )
775 const unsigned int firstHeight = FirstReductionType::height( coordinate.first() );
776 const unsigned int secondHeight = SecondReductionType::height( coordinate.second() );
778 return ((firstHeight < secondHeight) ? firstHeight : secondHeight);
781 unsigned int height ()
783 return height( dofCoordinate_ );
786 template<
class FieldType >
787 void local ( FieldVector< FieldType, dimension > &coordinate )
const
789 const FieldType factor = FieldType( 1 ) / FieldType( polynomialOrder );
790 for(
unsigned int i = 0; i < dimension; ++i )
791 coordinate[ i ] = factor * dofCoordinate_[ i ];
800 static unsigned int maxDofs (
unsigned int codim )
802 unsigned int max = 0;
803 for(
unsigned int i = 0; i <= codim; ++i ) {
805 = FirstReductionType :: maxDofs( codim - i )
806 * SecondReductionType :: maxDofs( i );
819 static unsigned int numDofs (
unsigned int codim,
unsigned int subEntity )
821 unsigned int firstCodim = codim;
822 unsigned int secondCodim = 0;
823 for( ; secondCodim <= codim; --firstCodim, ++secondCodim )
825 const unsigned int numSubEntities
826 = FirstGeometry::numSubEntities( firstCodim ) * SecondGeometry :: numSubEntities( secondCodim );
828 if( subEntity < numSubEntities )
830 subEntity -= numSubEntities;
833 const unsigned int n = FirstGeometry::numSubEntities( firstCodim );
834 const unsigned int firstSubEntity = subEntity % n;
835 const unsigned int secondSubEntity = subEntity / n;
837 return FirstReductionType :: numDofs( firstCodim, firstSubEntity )
838 * SecondReductionType :: numDofs( secondCodim, secondSubEntity );
847 static unsigned int numDofs (
unsigned int codim )
849 unsigned int count = 0;
851 unsigned int firstCodim = codim;
852 unsigned int secondCodim = 0;
853 for( ; secondCodim <= codim; --firstCodim, ++secondCodim )
854 count += FirstReductionType :: numDofs( firstCodim )
855 * SecondReductionType :: numDofs( secondCodim );
861 template<
class LocalCoordinateType >
862 static void dofCoordinate (
unsigned int index,
863 LocalCoordinateType &coordinate )
865 assert( index <= numLagrangePoints );
867 const unsigned int firstIndex
868 = index % FirstReductionType :: numLagrangePoints;
869 const unsigned int secondIndex
870 = index / FirstReductionType :: numLagrangePoints;
872 FirstReductionType :: dofCoordinate( firstIndex, coordinate.first() );
873 SecondReductionType :: dofCoordinate( secondIndex, coordinate.second() );
879 template<
class FirstGeometry,
class SecondGeometry,
unsigned int order,
bool bottom >
880 template<
unsigned int codim >
881 struct GenericLagrangePoint< ProductGeometry< FirstGeometry, SecondGeometry >, order, bottom >::Codim
883 static unsigned int maxDofs ()
885 return CodimIterator< codim, codim >::maxDofs();
892 template<
class FirstGeometry,
class SecondGeometry,
unsigned int order,
bool bottom >
893 template<
unsigned int codim,
unsigned int i >
894 struct GenericLagrangePoint< ProductGeometry< FirstGeometry, SecondGeometry >, order, bottom >::CodimIterator
896 static unsigned int maxDofs ()
898 typedef typename FirstReductionType::template Codim< codim-i > FirstCodim;
899 typedef typename SecondReductionType::template Codim< i > SecondCodim;
900 const unsigned int n = FirstCodim::maxDofs() * SecondCodim::maxDofs();
902 const unsigned int m = CodimIterator< codim, i-1 >::maxDofs();
903 return std::max( m, n );
907 template<
class FirstGeometry,
class SecondGeometry,
unsigned int order,
bool bottom >
908 template<
unsigned int codim>
909 struct GenericLagrangePoint< ProductGeometry< FirstGeometry, SecondGeometry >, order, bottom >::CodimIterator<codim, 0>
911 static unsigned int maxDofs ()
913 typedef typename FirstReductionType::template Codim< codim > FirstCodim;
914 typedef typename SecondReductionType::template
Codim< 0 > SecondCodim;
915 return FirstCodim::maxDofs() * SecondCodim::maxDofs();
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.
constexpr auto max
Function object that returns the greater of the given values.
Definition: hybridutilities.hh:484
Dune namespace.
Definition: alignedallocator.hh:13