dune-fem  2.8-git
basisfunctionset/hpdg/orthogonal.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_HPDG_SPACE_BASISFUNCTIONSETS_ORTHOGONAL_HH
2 #define DUNE_FEM_HPDG_SPACE_BASISFUNCTIONSETS_ORTHOGONAL_HH
3 
4 #include <cassert>
5 #include <cstddef>
6 
7 #include <algorithm>
8 #include <array>
9 #include <memory>
10 #include <type_traits>
11 #include <vector>
12 
13 #include <dune/geometry/type.hh>
14 
15 #include <dune/grid/common/capabilities.hh>
16 
22 
24 
25 #include "typeindexset.hh"
26 #include "typemap.hh"
27 #include "basisfunctionsets.hh"
28 
29 namespace Dune
30 {
31 
32  namespace Fem
33  {
34 
35  namespace hpDG
36  {
37 
38  // Internal forward declarations
39  // -----------------------------
40 
41  template< class FunctionSpace, class GridPart, int maxOrder, class Storage >
42  class OrthogonalBasisFunctionSets;
43 
44 
45 
46 #ifndef DOXYGEN
47 
48  // OrthonormalShapeFunctionSets
49  // ----------------------------
50 
51  template< class FunctionSpace, int order, class Storage >
52  class OrthonormalShapeFunctionSets
53  {
54  using ThisType = OrthonormalShapeFunctionSets< FunctionSpace, order, Storage >;
55 
56  public:
58  using ShapeFunctionSetType =
60 
61  private:
62  static const int dimension = ShapeFunctionSetType::DomainType::dimension;
63 
64  protected:
65  OrthonormalShapeFunctionSets ()
66  {
67  using GeometryTypeIndexSet = LocalGeometryTypeIndexSet< dimension, true >;
68  const std::size_t types = GeometryTypeIndexSet::size();
69  for( std::size_t i = 0u; i < types; ++i )
70  {
71  GeometryType type = GeometryTypeIndexSet::type( i );
72  assert( !type.isNone() );
73  types_.push_back( type );
74  for( int p = 0; p <= order; ++p )
75  shapeFunctionSets_[ type ][ p ].reset( new ShapeFunctionSetType( type, typename ShapeFunctionSetType::ImplementationType( type, p ) ) );
76  }
77  }
78 
79  static const ThisType &instance ()
80  {
81  static ThisType instance;
82  return instance;
83  }
84 
85  public:
87  static const std::vector< GeometryType > &types ()
88  {
89  return instance().types();
90  }
91 
93  static std::size_t maxSize ()
94  {
95  return OrthonormalShapeFunctions< dimension >::size( order );
96  }
97 
99  static std::size_t size ( GeometryType type, int p )
100  {
101  return OrthonormalShapeFunctions< dimension >::size( p );
102  }
103 
105  static const ShapeFunctionSetType &get ( GeometryType type, int p )
106  {
107  return *instance().shapeFunctionSets_[ type ][ p ];
108  }
109 
110  private:
111  std::vector< GeometryType > types_;
112  LocalGeometryTypeMap< std::array< std::unique_ptr< ShapeFunctionSetType >, order+1 >, dimension > shapeFunctionSets_;
113  };
114 
115 
116 
117  // OrthogonalBasisFunctionSetsTraits
118  // ---------------------------------
119 
120  template< class FunctionSpace, class GridPart, int maxOrder, class Storage >
121  class OrthogonalBasisFunctionSetsTraits
122  {
123  public:
124  using ImplementationType = OrthogonalBasisFunctionSets< FunctionSpace, GridPart, maxOrder, Storage >;
125 
126  using GridPartType = GridPart;
127  using Types = const std::vector< GeometryType > &;
128 
129  using KeyType = int;
130  using DataType = KeyType;
131 
132  using EntityType = typename GridPartType::template Codim< 0 >::EntityType;
133 
134  using ShapeFunctionSetsType = OrthonormalShapeFunctionSets< Dune::Fem::FunctionSpace< typename FunctionSpace::DomainFieldType, typename FunctionSpace::RangeFieldType, EntityType::mydimension, 1 >, maxOrder, Storage >;
136 
137  using BasisFunctionSetType = DefaultBasisFunctionSet< EntityType, ShapeFunctionSetType >;
138 
139  static const int localBlockSize = BasisFunctionSetType::RangeType::dimension;
140  };
141 
142 #endif // #ifndef DOXYGEN
143 
144 
145 
146  // OrthogonalBasisFunctionSets
147  // ---------------------------
148 
158  template< class FunctionSpace, class GridPart, int maxOrder, class Storage >
160  : public BasisFunctionSets< OrthogonalBasisFunctionSetsTraits< FunctionSpace, GridPart, maxOrder, Storage > >
161  {
164 
165  public:
170 
171  private:
172  using ShapeFunctionSetsType = typename BaseType::Traits::ShapeFunctionSetsType;
173  using ShapeFunctionSetType = typename BaseType::Traits::ShapeFunctionSetType;
174 
175  public:
178 
180  using KeyType = typename BaseType::KeyType;
182  using DataType = typename BaseType::DataType;
183 
186 
189 
192 
199 
201  OrthogonalBasisFunctionSets ( const ThisType & ) = default;
202 
205 
213  ThisType &operator= ( const ThisType & ) = default;
214 
216  ThisType &operator= ( ThisType && ) = default;
217 
221  static const std::vector< GeometryType > &types ()
222  {
223  return ShapeFunctionSetsType::types();
224  }
225 
227  static std::size_t maxBlocks ()
228  {
229  return ShapeFunctionSetsType::maxSize();
230  }
231 
233  static std::size_t maxBlocks ( GeometryType type )
234  {
235  return contains( type ) ? maxBlocks() : 0u;
236  }
237 
239  static std::size_t blocks ( GeometryType type, KeyType key )
240  {
241  return contains( type ) ? ShapeFunctionSetsType::size( type, key ) : 0u;
242  }
243 
245  static DataType encode ( const KeyType &key ) noexcept { return key; }
246 
248  static KeyType decode ( const DataType &data ) noexcept { return data; }
249 
251  static constexpr bool orthogonal () noexcept
252  {
253  //using GridType = typename GridPartType::GridType;
254  //return Dune::Capabilities::isCartesian< GridType >::v;
255  return true;
256  }
257 
259  static constexpr int order () noexcept { return maxOrder; }
260 
262  static int order ( GeometryType type ) noexcept { return order(); }
263 
265  static int order ( GeometryType type, KeyType key )
266  {
267  assert( contains( type ) );
268  return static_cast< int >( key );
269  }
270 
273  {
274  assert( contains( entity.type() ) );
275  return BasisFunctionSetType( entity, shapeFunctionSet( entity.type(), key ) );
276  }
277 
278  private:
279  static bool contains ( GeometryType type )
280  {
281  return (type.dim() <= 3 && !type.isNone());
282  }
283 
284  static ShapeFunctionSetType shapeFunctionSet ( GeometryType type, KeyType key )
285  {
286  return ShapeFunctionSetType( &ShapeFunctionSetsType::get( type, static_cast< int >( key ) ) );
287  }
288  };
289 
290  } // namespace hpDG
291 
292  } // namespace Fem
293 
294 } // namespace Dune
295 
296 #endif // #ifndef DUNE_FEM_HPDG_SPACE_BASISFUNCTIONSETS_ORTHOGONAL_HH
Provides a proxy class for pointers to a shape function set.
Definition: bindguard.hh:11
std::tuple_element< i, Tuple >::type & get(Dune::TypeIndexedTuple< Tuple, Types > &tuple)
Definition: typeindexedtuple.hh:122
abstract interface class for a family of local basis function sets
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:30
typename Traits::DataType DataType
data type
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:53
typename Traits::BasisFunctionSetType BasisFunctionSetType
basis function set type
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:42
typename BasisFunctionSetType::EntityType EntityType
entity type
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:44
typename Traits::KeyType KeyType
key type
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:39
typename Traits::GridPartType GridPartType
grid part type
Definition: basisfunctionset/hpdg/basisfunctionsets.hh:36
A family of orthogonal local basis function sets.
Definition: basisfunctionset/hpdg/orthogonal.hh:161
static std::size_t maxBlocks(GeometryType type)
Definition: basisfunctionset/hpdg/orthogonal.hh:233
OrthogonalBasisFunctionSets(const ThisType &)=default
copy constructor
typename BaseType::GridPartType GridPartType
Definition: basisfunctionset/hpdg/orthogonal.hh:167
typename BaseType::DataType DataType
Definition: basisfunctionset/hpdg/orthogonal.hh:182
static int order(GeometryType type, KeyType key)
return maximum order
Definition: basisfunctionset/hpdg/orthogonal.hh:265
static const std::vector< GeometryType > & types()
Definition: basisfunctionset/hpdg/orthogonal.hh:221
ThisType & operator=(const ThisType &)=default
assignment operator
static BasisFunctionSetType basisFunctionSet(const EntityType &entity, KeyType key)
return basis function set for given entity
Definition: basisfunctionset/hpdg/orthogonal.hh:272
static std::size_t blocks(GeometryType type, KeyType key)
Definition: basisfunctionset/hpdg/orthogonal.hh:239
static KeyType decode(const DataType &data) noexcept
Definition: basisfunctionset/hpdg/orthogonal.hh:248
static std::size_t maxBlocks()
Definition: basisfunctionset/hpdg/orthogonal.hh:227
typename BaseType::EntityType EntityType
entity type
Definition: basisfunctionset/hpdg/orthogonal.hh:169
typename FunctionSpaceType ::DomainType DomainType
range type of basis functions
Definition: basisfunctionset/hpdg/orthogonal.hh:188
static constexpr bool orthogonal() noexcept
Definition: basisfunctionset/hpdg/orthogonal.hh:251
typename BaseType::BasisFunctionSetType BasisFunctionSetType
basis function set
Definition: basisfunctionset/hpdg/orthogonal.hh:177
typename BaseType::KeyType KeyType
Definition: basisfunctionset/hpdg/orthogonal.hh:180
static constexpr int order() noexcept
return maximum order
Definition: basisfunctionset/hpdg/orthogonal.hh:259
static int order(GeometryType type) noexcept
return maximum order
Definition: basisfunctionset/hpdg/orthogonal.hh:262
OrthogonalBasisFunctionSets(ThisType &&)=default
move constructor
static DataType encode(const KeyType &key) noexcept
Definition: basisfunctionset/hpdg/orthogonal.hh:245
typename FunctionSpaceType ::RangeType RangeType
range type of basis functions
Definition: basisfunctionset/hpdg/orthogonal.hh:191
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
Definition: selectcaching.hh:26
Definition: shapefunctionset/vectorial.hh:447
Please doc me.