dune-pdelab 2.7-git
Loading...
Searching...
No Matches
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 <array>
6#include <memory>
7
8#include <dune/geometry/type.hh>
9
10#include <dune/localfunctions/common/virtualwrappers.hh>
11#include <dune/localfunctions/common/virtualinterface.hh>
12#include "finiteelementmap.hh"
13#include "qkfem.hh"
14#include "pkfem.hh"
15
16namespace 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> > > FiniteElementType;
70 public:
72
74 static constexpr int dimension = d;
75
78 : order_(maxP)
79 {
80 InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,maxP);
81 }
82
86 PkQkLocalFiniteElementMap (unsigned int order)
87 : order_(order)
88 {
89 InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,order);
90 }
91
93 template<class EntityType>
94 const typename Traits::FiniteElementType& find (const EntityType& e) const
95 {
96 typename Dune::GeometryType geoType=e.type();
97 return getFEM(geoType);
98 }
99
101 const typename Traits::FiniteElementType& getFEM (Dune::GeometryType gt) const
102 {
103 if (gt.isCube())
104 {
105 return *(finiteElements_[0]);
106 }
107 if (gt.isSimplex())
108 {
109 return *(finiteElements_[1]);
110 }
111 DUNE_THROW(Exception, "We can only handle cubes and simplices");
112 }
113
114 static constexpr bool fixedSize()
115 {
116 return false;
117 }
118
119 bool hasDOFs(int codim) const
120 {
121 switch (codim)
122 {
123 case 0:
124 return order_ == 0 || order_ > 1;
125 case d:
126 return order_ > 0;
127 default:
128 return order_ > 1;
129 }
130 }
131
132 std::size_t size(GeometryType gt) const
133 {
134 DUNE_THROW(NotImplemented, "PkQkLocalFiniteElement is not fixed-size!");
135 }
136
137 static constexpr std::size_t maxLocalSize()
138 {
139 return (1<<d);
140 }
141
142 private:
143 std::array< std::shared_ptr<FiniteElementType>, 2 > finiteElements_;
144 const std::size_t order_;
145 };
146 }
147}
148
149#endif // DUNE_PDELAB_FINITEELEMENTMAP_PKQKFEM_HH
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Base class for all PDELab exceptions.
Definition: exceptions.hh:19
collect types exported by a finite element map
Definition: finiteelementmap.hh:28
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
FiniteElementMap which provides PkQkLocalFiniteElement instances, depending on the geometry type.
Definition: pkqkfem.hh:67
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: pkqkfem.hh:94
static constexpr std::size_t maxLocalSize()
Definition: pkqkfem.hh:137
PkQkLocalFiniteElementMap()
Default constructor. Constructs a space of order maxP.
Definition: pkqkfem.hh:77
const Traits::FiniteElementType & getFEM(Dune::GeometryType gt) const
get local basis functions for a given geometrytype
Definition: pkqkfem.hh:101
static constexpr int dimension
The dimension of the finite elements returned by this map.
Definition: pkqkfem.hh:74
bool hasDOFs(int codim) const
Definition: pkqkfem.hh:119
FiniteElementMapTraits< FiniteElementType > Traits
Definition: pkqkfem.hh:71
std::size_t size(GeometryType gt) const
Definition: pkqkfem.hh:132
PkQkLocalFiniteElementMap(unsigned int order)
Construct a space with a given order.
Definition: pkqkfem.hh:86
static constexpr bool fixedSize()
Definition: pkqkfem.hh:114