dune-localfunctions  2.4.1-rc2
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/common/deprecated.hh>
14 
15 #include <dune/geometry/type.hh>
16 
22 
23 namespace Dune
24 {
25 
26 
40  template<class D, class R, int d, int p, int diffOrder = p>
42  {
43  enum { static_size = MonomImp::Size<d,p>::val };
44 
45  public:
52  > Traits;
53 
55  MonomialLocalFiniteElement (const GeometryType &gt_)
56  : basis(), interpolation(gt_, basis), gt(gt_)
57  {}
58 
61  const typename Traits::LocalBasisType& localBasis () const
62  {
63  return basis;
64  }
65 
69  {
70  return coefficients;
71  }
72 
76  {
77  return interpolation;
78  }
79 
81  unsigned int size () const
82  {
83  return basis.size();
84  }
85 
88  GeometryType type () const
89  {
90  return gt;
91  }
92 
93  private:
94  MonomialLocalBasis<D,R,d,p, diffOrder> basis;
95  MonomialLocalCoefficients<static_size> coefficients;
96  MonomialLocalInterpolation<MonomialLocalBasis<D,R,d,p, diffOrder>,static_size> interpolation;
97  GeometryType gt;
98  };
99 
101 
113  template<class Geometry, class RF, std::size_t p>
115  typedef typename Geometry::ctype DF;
116  static const std::size_t dim = Geometry::mydimension;
117 
119 
120  std::vector<std::shared_ptr<const LocalFE> > localFEs;
121 
122  void init(const GeometryType &gt) {
123  std::size_t index = gt.id() >> 1;
124  if(localFEs.size() <= index)
125  localFEs.resize(index+1);
126  localFEs[index].reset(new LocalFE(gt));
127  }
128 
129  public:
132 
134 
138  template<class ForwardIterator>
139  MonomialFiniteElementFactory(const ForwardIterator &begin,
140  const ForwardIterator &end)
141  {
142  for(ForwardIterator it = begin; it != end; ++it)
143  init(*it);
144  }
145 
147 
150  MonomialFiniteElementFactory(const GeometryType &gt)
151  { init(gt); }
152 
154 
158  static_assert(dim <= 3, "MonomFiniteElementFactory knows the "
159  "available geometry types only up to dimension 3");
160 
161  GeometryType gt;
162  switch(dim) {
163  case 0 :
164  gt.makeVertex(); init(gt);
165  break;
166  case 1 :
167  gt.makeLine(); init(gt);
168  break;
169  case 2 :
170  gt.makeTriangle(); init(gt);
171  gt.makeQuadrilateral(); init(gt);
172  break;
173  case 3 :
174  gt.makeTetrahedron(); init(gt);
175  gt.makePyramid(); init(gt);
176  gt.makePrism(); init(gt);
177  gt.makeHexahedron(); init(gt);
178  break;
179  default :
180  // this should never happen -- it should be caught by the static
181  // assert above.
182  std::abort();
183  };
184  }
185 
187 
197  const FiniteElement make(const Geometry& geometry) {
198  std::size_t index = geometry.type().id() >> 1;
199  assert(localFEs.size() > index && localFEs[index]);
200  return FiniteElement(*localFEs[index], geometry);
201  }
202  };
203 }
204 
205 #endif // DUNE_LOCALFUNCTIONS_MONOMIAL_HH
Definition: monomiallocalinterpolation.hh:18
traits helper struct
Definition: localfiniteelementtraits.hh:10
MonomialFiniteElementFactory(const ForwardIterator &begin, const ForwardIterator &end)
construct a MonomialFiniteElementFactory from a list of GeometryType&#39;s
Definition: monomial.hh:139
MonomialFiniteElementFactory()
construct a MonomFiniteElementFactory for all applicable GeometryType&#39;s
Definition: monomial.hh:157
const Traits::LocalInterpolationType & localInterpolation() const
Definition: monomial.hh:75
LocalFiniteElementTraits< MonomialLocalBasis< D, R, d, p, diffOrder >, MonomialLocalCoefficients< static_size >, MonomialLocalInterpolation< MonomialLocalBasis< D, R, d, p, diffOrder >, static_size > > Traits
Definition: monomial.hh:52
ImplementationDefined FiniteElement
Type of the finite element.
Definition: interface.hh:115
Factory for global-valued MonomFiniteElement objects.
Definition: monomial.hh:114
ScalarLocalToGlobalFiniteElementAdaptor< LocalFE, Geometry > FiniteElement
Definition: monomial.hh:131
Monomial basis for discontinuous Galerkin methods.
Definition: monomial.hh:41
unsigned int size() const
Number of shape functions in this finite element.
Definition: monomial.hh:81
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
Definition: monomiallocalbasis.hh:20
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
Convert a simple scalar local finite element into a global finite element.
Definition: localtoglobaladaptors.hh:187
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
const FiniteElement make(const Geometry &geometry)
construct a global-valued MonomFiniteElement
Definition: monomial.hh:197
const Traits::LocalBasisType & localBasis() const
Definition: monomial.hh:61
GeometryType type() const
Definition: monomial.hh:88
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: monomial.hh:68
MonomialFiniteElementFactory(const GeometryType &gt)
construct a MonomialFiniteElementFactory from a single GeometryType
Definition: monomial.hh:150
Layout map for monomial finite elements.
Definition: monomiallocalcoefficients.hh:21
MonomialLocalFiniteElement(const GeometryType &gt_)
Construct a MonomLocalFiniteElement.
Definition: monomial.hh:55
Constant shape function.
Definition: monomiallocalbasis.hh:231