dune-pdelab  2.0.0
pkqkfem.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_PKQKFEM_HH
3 #define DUNE_PDELAB_PKQKFEM_HH
4 
5 #include <dune/geometry/type.hh>
6 
7 #include <dune/localfunctions/common/virtualwrappers.hh>
8 #include <dune/localfunctions/common/virtualinterface.hh>
9 #include <dune/common/array.hh>
10 #include <dune/common/shared_ptr.hh>
11 #include "finiteelementmap.hh"
12 #include "qkfem.hh"
13 #include "pkfem.hh"
14 
15 namespace Dune {
16  namespace PDELab {
17 
18  namespace {
19  template<class D, class R, int d, int k>
20  struct InitPkQkLocalFiniteElementMap
21  {
22  template<typename C>
23  static void init(C & c, unsigned int order)
24  {
25  if (order == k)
26  {
27  typedef Dune::QkLocalFiniteElement<D,R,d,k> QkLFE;
28  typedef Dune::PkLocalFiniteElement<D,R,d,k> PkLFE;
29  typedef typename C::value_type ptr;
30  c[0] = ptr(new LocalFiniteElementVirtualImp<QkLFE>(QkLFE()));
31  c[1] = ptr(new LocalFiniteElementVirtualImp<PkLFE>(PkLFE()));
32  }
33  else
34  InitPkQkLocalFiniteElementMap<D,R,d,k-1>::init(c,order);
35  }
36  };
37  template<class D, class R, int d>
38  struct InitPkQkLocalFiniteElementMap<D,R,d,-1>
39  {
40  template<typename C>
41  static void init(C & c, unsigned int order)
42  {
43  DUNE_THROW(Exception, "Sorry, but we failed to initialize a QkPk FiniteElementMap of order " << order);
44  }
45  };
46  }
47 
50  template<class D, class R, int d, int maxP=6>
52  {
54  typedef LocalFiniteElementVirtualInterface<Dune::LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d>,0> > FiniteElementType;
55  public:
57 
59  {
60  InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,maxP);
61  }
62 
63  PkQkLocalFiniteElementMap (unsigned int order)
64  {
65  InitPkQkLocalFiniteElementMap<D,R,d,maxP>::init(finiteElements_,order);
66  }
67 
69  template<class EntityType>
70  const typename Traits::FiniteElementType& find (const EntityType& e) const
71  {
72  typename Dune::GeometryType geoType=e.type();
73  return getFEM(geoType);
74  }
75 
77  const typename Traits::FiniteElementType& getFEM (Dune::GeometryType gt) const
78  {
79  if (gt.isCube())
80  {
81  return *(finiteElements_[0]);
82  }
83  if (gt.isSimplex())
84  {
85  return *(finiteElements_[1]);
86  }
87  DUNE_THROW(Exception, "We can only handle cubes and simplices");
88  }
89 
90  bool fixedSize() const
91  {
92  return false;
93  }
94 
95  std::size_t size(GeometryType gt) const
96  {
97  assert(false && "this method should never be called");
98  }
99 
100  std::size_t maxLocalSize() const
101  {
102  return (1<<d);
103  }
104 
105  private:
106  Dune::array< Dune::shared_ptr<FiniteElementType>, 2 > finiteElements_;
107  };
108  }
109 }
110 
111 #endif //DUNE_PDELAB_PKQKFEM_HH
std::size_t size(GeometryType gt) const
Definition: pkqkfem.hh:95
bool fixedSize() const
Definition: pkqkfem.hh:90
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:70
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
FiniteElementMapTraits< FiniteElementType > Traits
Definition: pkqkfem.hh:56
const Traits::FiniteElementType & getFEM(Dune::GeometryType gt) const
get local basis functions for a given geometrytype
Definition: pkqkfem.hh:77
PkQkLocalFiniteElementMap(unsigned int order)
Definition: pkqkfem.hh:63
PkQkLocalFiniteElementMap()
Definition: pkqkfem.hh:58
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
std::size_t maxLocalSize() const
Definition: pkqkfem.hh:100
const E & e
Definition: interpolate.hh:172