dune-localfunctions  2.9.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 // SPDX-FileCopyrightInfo: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
4 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5 
6 #ifndef DUNE_LOCALFUNCTIONS_MONOMIAL_HH
7 #define DUNE_LOCALFUNCTIONS_MONOMIAL_HH
8 
9 #include <cassert>
10 #include <cstddef>
11 #include <cstdlib>
12 #include <memory>
13 #include <vector>
14 
15 #include <dune/geometry/type.hh>
16 
22 
23 namespace Dune
24 {
25 
26 
39  template<class D, class R, int d, int p>
41  {
42  constexpr static int static_size = MonomialLocalBasis<D,R,d,p>::size();
43 
44  public:
51  > Traits;
52 
54  MonomialLocalFiniteElement (const GeometryType &gt_)
55  : basis(), interpolation(gt_, basis), gt(gt_)
56  {}
57 
60  const typename Traits::LocalBasisType& localBasis () const
61  {
62  return basis;
63  }
64 
68  {
69  return coefficients;
70  }
71 
75  {
76  return interpolation;
77  }
78 
80  unsigned int size () const
81  {
82  return basis.size();
83  }
84 
87  GeometryType type () const
88  {
89  return gt;
90  }
91 
92  private:
96  GeometryType gt;
97  };
98 
100 
112  template<class Geometry, class RF, std::size_t p>
114  typedef typename Geometry::ctype DF;
115  static const std::size_t dim = Geometry::mydimension;
116 
118 
119  std::vector<std::shared_ptr<const LocalFE> > localFEs;
120 
121  void init(const GeometryType &gt) {
122  std::size_t index = gt.id() >> 1;
123  if(localFEs.size() <= index)
124  localFEs.resize(index+1);
125  localFEs[index].reset(new LocalFE(gt));
126  }
127 
128  public:
131 
133 
137  template<class ForwardIterator>
138  MonomialFiniteElementFactory(const ForwardIterator &begin,
139  const ForwardIterator &end)
140  {
141  for(ForwardIterator it = begin; it != end; ++it)
142  init(*it);
143  }
144 
146 
149  MonomialFiniteElementFactory(const GeometryType &gt)
150  { init(gt); }
151 
153 
157  static_assert(dim <= 3, "MonomFiniteElementFactory knows the "
158  "available geometry types only up to dimension 3");
159 
160  GeometryType gt;
161  switch(dim) {
162  case 0 :
163  gt = Dune::GeometryTypes::vertex; init(gt);
164  break;
165  case 1 :
166  gt = Dune::GeometryTypes::line; init(gt);
167  break;
168  case 2 :
169  gt = Dune::GeometryTypes::triangle; init(gt);
170  gt = Dune::GeometryTypes::quadrilateral; init(gt);
171  break;
172  case 3 :
173  gt = Dune::GeometryTypes::tetrahedron; init(gt);
174  gt = Dune::GeometryTypes::pyramid; init(gt);
175  gt = Dune::GeometryTypes::prism; init(gt);
176  gt = Dune::GeometryTypes::hexahedron; init(gt);
177  break;
178  default :
179  // this should never happen -- it should be caught by the static
180  // assert above.
181  std::abort();
182  };
183  }
184 
186 
196  const FiniteElement make(const Geometry& geometry) {
197  std::size_t index = geometry.type().id() >> 1;
198  assert(localFEs.size() > index && localFEs[index]);
199  return FiniteElement(*localFEs[index], geometry);
200  }
201  };
202 }
203 
204 #endif // DUNE_LOCALFUNCTIONS_MONOMIAL_HH
Definition: bdfmcube.hh:18
ImplementationDefined FiniteElement
Type of the finite element.
Definition: interface.hh:117
traits helper struct
Definition: localfiniteelementtraits.hh:13
LB LocalBasisType
Definition: localfiniteelementtraits.hh:16
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:20
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:24
Convert a simple scalar local finite element into a global finite element.
Definition: localtoglobaladaptors.hh:187
GeometryType type() const
Definition: localtoglobaladaptors.hh:229
Constant shape function.
Definition: monomiallocalbasis.hh:201
static constexpr unsigned int size()
Number of shape functions.
Definition: monomiallocalbasis.hh:216
Layout map for monomial finite elements.
Definition: monomiallocalcoefficients.hh:24
Definition: monomiallocalinterpolation.hh:22
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:80
const Traits::LocalBasisType & localBasis() const
Definition: monomial.hh:60
GeometryType type() const
Definition: monomial.hh:87
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: monomial.hh:67
const Traits::LocalInterpolationType & localInterpolation() const
Definition: monomial.hh:74
LocalFiniteElementTraits< MonomialLocalBasis< D, R, d, p >, MonomialLocalCoefficients< static_size >, MonomialLocalInterpolation< MonomialLocalBasis< D, R, d, p >, static_size > > Traits
Definition: monomial.hh:51
MonomialLocalFiniteElement(const GeometryType &gt_)
Construct a MonomLocalFiniteElement.
Definition: monomial.hh:54
Factory for global-valued MonomFiniteElement objects.
Definition: monomial.hh:113
MonomialFiniteElementFactory(const ForwardIterator &begin, const ForwardIterator &end)
construct a MonomialFiniteElementFactory from a list of GeometryType's
Definition: monomial.hh:138
MonomialFiniteElementFactory(const GeometryType &gt)
construct a MonomialFiniteElementFactory from a single GeometryType
Definition: monomial.hh:149
const FiniteElement make(const Geometry &geometry)
construct a global-valued MonomFiniteElement
Definition: monomial.hh:196
ScalarLocalToGlobalFiniteElementAdaptor< LocalFE, Geometry > FiniteElement
Definition: monomial.hh:130
MonomialFiniteElementFactory()
construct a MonomFiniteElementFactory for all applicable GeometryType's
Definition: monomial.hh:156