dune-localfunctions  2.8.0
monomial.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_LOCALFUNCTIONS_MONOMIAL_HH
5 #define DUNE_LOCALFUNCTIONS_MONOMIAL_HH
6 
7 #include <cassert>
8 #include <cstddef>
9 #include <cstdlib>
10 #include <memory>
11 #include <vector>
12 
13 #include <dune/geometry/type.hh>
14 
20 
21 namespace Dune
22 {
23 
24 
37  template<class D, class R, int d, int p>
39  {
40  enum { static_size = MonomialLocalBasis<D,R,d,p>::size() };
41 
42  public:
49  > Traits;
50 
52  MonomialLocalFiniteElement (const GeometryType &gt_)
53  : basis(), interpolation(gt_, basis), gt(gt_)
54  {}
55 
58  const typename Traits::LocalBasisType& localBasis () const
59  {
60  return basis;
61  }
62 
66  {
67  return coefficients;
68  }
69 
73  {
74  return interpolation;
75  }
76 
78  unsigned int size () const
79  {
80  return basis.size();
81  }
82 
85  GeometryType type () const
86  {
87  return gt;
88  }
89 
90  private:
94  GeometryType gt;
95  };
96 
98 
110  template<class Geometry, class RF, std::size_t p>
112  typedef typename Geometry::ctype DF;
113  static const std::size_t dim = Geometry::mydimension;
114 
116 
117  std::vector<std::shared_ptr<const LocalFE> > localFEs;
118 
119  void init(const GeometryType &gt) {
120  std::size_t index = gt.id() >> 1;
121  if(localFEs.size() <= index)
122  localFEs.resize(index+1);
123  localFEs[index].reset(new LocalFE(gt));
124  }
125 
126  public:
129 
131 
135  template<class ForwardIterator>
136  MonomialFiniteElementFactory(const ForwardIterator &begin,
137  const ForwardIterator &end)
138  {
139  for(ForwardIterator it = begin; it != end; ++it)
140  init(*it);
141  }
142 
144 
147  MonomialFiniteElementFactory(const GeometryType &gt)
148  { init(gt); }
149 
151 
155  static_assert(dim <= 3, "MonomFiniteElementFactory knows the "
156  "available geometry types only up to dimension 3");
157 
158  GeometryType gt;
159  switch(dim) {
160  case 0 :
161  gt = Dune::GeometryTypes::vertex; init(gt);
162  break;
163  case 1 :
164  gt = Dune::GeometryTypes::line; init(gt);
165  break;
166  case 2 :
167  gt = Dune::GeometryTypes::triangle; init(gt);
168  gt = Dune::GeometryTypes::quadrilateral; init(gt);
169  break;
170  case 3 :
171  gt = Dune::GeometryTypes::tetrahedron; init(gt);
172  gt = Dune::GeometryTypes::pyramid; init(gt);
173  gt = Dune::GeometryTypes::prism; init(gt);
174  gt = Dune::GeometryTypes::hexahedron; init(gt);
175  break;
176  default :
177  // this should never happen -- it should be caught by the static
178  // assert above.
179  std::abort();
180  };
181  }
182 
184 
194  const FiniteElement make(const Geometry& geometry) {
195  std::size_t index = geometry.type().id() >> 1;
196  assert(localFEs.size() > index && localFEs[index]);
197  return FiniteElement(*localFEs[index], geometry);
198  }
199  };
200 }
201 
202 #endif // DUNE_LOCALFUNCTIONS_MONOMIAL_HH
Definition: bdfmcube.hh:16
ImplementationDefined FiniteElement
Type of the finite element.
Definition: interface.hh:115
traits helper struct
Definition: localfiniteelementtraits.hh:11
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
Convert a simple scalar local finite element into a global finite element.
Definition: localtoglobaladaptors.hh:185
GeometryType type() const
Definition: localtoglobaladaptors.hh:227
Constant shape function.
Definition: monomiallocalbasis.hh:200
static constexpr unsigned int size()
Number of shape functions.
Definition: monomiallocalbasis.hh:215
Layout map for monomial finite elements.
Definition: monomiallocalcoefficients.hh:22
Definition: monomiallocalinterpolation.hh:20
Monomial basis for discontinuous Galerkin methods.
Definition: monomial.hh:39
unsigned int size() const
Number of shape functions in this finite element.
Definition: monomial.hh:78
const Traits::LocalBasisType & localBasis() const
Definition: monomial.hh:58
GeometryType type() const
Definition: monomial.hh:85
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: monomial.hh:65
const Traits::LocalInterpolationType & localInterpolation() const
Definition: monomial.hh:72
LocalFiniteElementTraits< MonomialLocalBasis< D, R, d, p >, MonomialLocalCoefficients< static_size >, MonomialLocalInterpolation< MonomialLocalBasis< D, R, d, p >, static_size > > Traits
Definition: monomial.hh:49
MonomialLocalFiniteElement(const GeometryType &gt_)
Construct a MonomLocalFiniteElement.
Definition: monomial.hh:52
Factory for global-valued MonomFiniteElement objects.
Definition: monomial.hh:111
MonomialFiniteElementFactory(const ForwardIterator &begin, const ForwardIterator &end)
construct a MonomialFiniteElementFactory from a list of GeometryType's
Definition: monomial.hh:136
MonomialFiniteElementFactory(const GeometryType &gt)
construct a MonomialFiniteElementFactory from a single GeometryType
Definition: monomial.hh:147
const FiniteElement make(const Geometry &geometry)
construct a global-valued MonomFiniteElement
Definition: monomial.hh:194
ScalarLocalToGlobalFiniteElementAdaptor< LocalFE, Geometry > FiniteElement
Definition: monomial.hh:128
MonomialFiniteElementFactory()
construct a MonomFiniteElementFactory for all applicable GeometryType's
Definition: monomial.hh:154