2 #ifndef DUNE_PDELAB_PKQKFEM_HH
3 #define DUNE_PDELAB_PKQKFEM_HH
7 #include <dune/geometry/type.hh>
9 #include <dune/localfunctions/common/virtualwrappers.hh>
10 #include <dune/localfunctions/common/virtualinterface.hh>
11 #include <dune/common/array.hh>
25 template<
class D,
class R,
int d,
int k>
26 struct InitPkQkLocalFiniteElementMap
30 static void init(C & c,
unsigned int order)
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()));
41 InitPkQkLocalFiniteElementMap<D,R,d,k-1>::init(c,order);
44 template<
class D,
class R,
int d>
45 struct InitPkQkLocalFiniteElementMap<D,R,d,-1>
48 static void init(C & c,
unsigned int order)
50 DUNE_THROW(Exception,
"Sorry, but we failed to initialize a QkPk FiniteElementMap of order " << order);
65 template<
class D,
class R,
int d,
int maxP=6>
69 typedef LocalFiniteElementVirtualInterface<Dune::LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d>,0> > FiniteElementType;
76 InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,maxP);
84 InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,order);
88 template<
class EntityType>
91 typename Dune::GeometryType geoType=e.type();
100 return *(finiteElements_[0]);
104 return *(finiteElements_[1]);
106 DUNE_THROW(
Exception,
"We can only handle cubes and simplices");
114 std::size_t
size(GeometryType gt)
const
116 assert(
false &&
"this method should never be called");
125 std::array< std::shared_ptr<FiniteElementType>, 2 > finiteElements_;
130 #endif //DUNE_PDELAB_PKQKFEM_HH
const Traits::FiniteElementType & getFEM(Dune::GeometryType gt) const
get local basis functions for a given geometrytype
Definition: pkqkfem.hh:96
PkQkLocalFiniteElementMap(unsigned int order)
Construct a space with a given order.
Definition: pkqkfem.hh:82
bool fixedSize() const
Definition: pkqkfem.hh:109
PkQkLocalFiniteElementMap()
Default constructor. Constructs a space of order maxP.
Definition: pkqkfem.hh:74
const E & e
Definition: interpolate.hh:172
FiniteElementMap which provides PkQkLocalFiniteElement instances, depending on the geometry type...
Definition: pkqkfem.hh:66
std::size_t maxLocalSize() const
Definition: pkqkfem.hh:119
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
Definition: adaptivity.hh:27
collect types exported by a finite element map
Definition: finiteelementmap.hh:27
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: pkqkfem.hh:89
std::size_t size(GeometryType gt) const
Definition: pkqkfem.hh:114
FiniteElementMapTraits< FiniteElementType > Traits
Definition: pkqkfem.hh:71