dune-pdelab  2.5-dev
pkqkfem.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_FINITEELEMENTMAP_PKQKFEM_HH
3 #define DUNE_PDELAB_FINITEELEMENTMAP_PKQKFEM_HH
4 
5 #include <memory>
6 
7 #include <dune/geometry/type.hh>
8 
9 #include <dune/localfunctions/common/virtualwrappers.hh>
10 #include <dune/localfunctions/common/virtualinterface.hh>
11 #include <dune/common/array.hh>
12 #include "finiteelementmap.hh"
13 #include "qkfem.hh"
14 #include "pkfem.hh"
15 
16 namespace Dune {
17  namespace PDELab {
18 
19  namespace {
20 
25  template<class D, class R, int d, int k>
26  struct InitPkQkLocalFiniteElementMap
27  {
29  template<typename C>
30  static void init(C & c, unsigned int order)
31  {
32  if (order == k)
33  {
34  typedef Dune::QkLocalFiniteElement<D,R,d,k> QkLFE;
35  typedef Dune::PkLocalFiniteElement<D,R,d,k> PkLFE;
36  typedef typename C::value_type ptr;
37  c[0] = ptr(new LocalFiniteElementVirtualImp<QkLFE>(QkLFE()));
38  c[1] = ptr(new LocalFiniteElementVirtualImp<PkLFE>(PkLFE()));
39  }
40  else
41  InitPkQkLocalFiniteElementMap<D,R,d,k-1>::init(c,order);
42  }
43  };
44  template<class D, class R, int d>
45  struct InitPkQkLocalFiniteElementMap<D,R,d,-1>
46  {
47  template<typename C>
48  static void init(C & c, unsigned int order)
49  {
50  DUNE_THROW(Exception, "Sorry, but we failed to initialize a QkPk FiniteElementMap of order " << order);
51  }
52  };
53  }
54 
65  template<class D, class R, int d, int maxP=6>
67  {
69  typedef LocalFiniteElementVirtualInterface<Dune::LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d>,0> > FiniteElementType;
70  public:
72 
75  : order_(maxP)
76  {
77  InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,maxP);
78  }
79 
83  PkQkLocalFiniteElementMap (unsigned int order)
84  : order_(order)
85  {
86  InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,order);
87  }
88 
90  template<class EntityType>
91  const typename Traits::FiniteElementType& find (const EntityType& e) const
92  {
93  typename Dune::GeometryType geoType=e.type();
94  return getFEM(geoType);
95  }
96 
98  const typename Traits::FiniteElementType& getFEM (Dune::GeometryType gt) const
99  {
100  if (gt.isCube())
101  {
102  return *(finiteElements_[0]);
103  }
104  if (gt.isSimplex())
105  {
106  return *(finiteElements_[1]);
107  }
108  DUNE_THROW(Exception, "We can only handle cubes and simplices");
109  }
110 
111  bool fixedSize() const
112  {
113  return false;
114  }
115 
116  bool hasDOFs(int codim) const
117  {
118  switch (codim)
119  {
120  case 0:
121  return order_ == 0 || order_ > 1;
122  case d:
123  return order_ > 0;
124  default:
125  return order_ > 1;
126  }
127  }
128 
129  std::size_t size(GeometryType gt) const
130  {
131  DUNE_THROW(NotImplemented, "PkQkLocalFiniteElement is not fixed-size!");
132  }
133 
134  std::size_t maxLocalSize() const
135  {
136  return (1<<d);
137  }
138 
139  private:
140  std::array< std::shared_ptr<FiniteElementType>, 2 > finiteElements_;
141  const std::size_t order_;
142  };
143  }
144 }
145 
146 #endif // DUNE_PDELAB_FINITEELEMENTMAP_PKQKFEM_HH
FiniteElementMap which provides PkQkLocalFiniteElement instances, depending on the geometry type...
Definition: pkqkfem.hh:66
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
bool hasDOFs(int codim) const
Definition: pkqkfem.hh:116
PkQkLocalFiniteElementMap()
Default constructor. Constructs a space of order maxP.
Definition: pkqkfem.hh:74
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Entity & e
Definition: localfunctionspace.hh:111
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: pkqkfem.hh:91
std::size_t size(GeometryType gt) const
Definition: pkqkfem.hh:129
bool fixedSize() const
Definition: pkqkfem.hh:111
PkQkLocalFiniteElementMap(unsigned int order)
Construct a space with a given order.
Definition: pkqkfem.hh:83
collect types exported by a finite element map
Definition: finiteelementmap.hh:27
std::size_t maxLocalSize() const
Definition: pkqkfem.hh:134
const Traits::FiniteElementType & getFEM(Dune::GeometryType gt) const
get local basis functions for a given geometrytype
Definition: pkqkfem.hh:98
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
FiniteElementMapTraits< FiniteElementType > Traits
Definition: pkqkfem.hh:71