dune-localfunctions 3.0-git
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
21namespace Dune
22{
23
24
38 template<class D, class R, int d, int p, int diffOrder = p>
40 {
41 enum { static_size = MonomImp::Size<d,p>::val };
42
43 public:
51
53 MonomialLocalFiniteElement (const GeometryType &gt_)
54 : basis(), interpolation(gt_, basis), gt(gt_)
55 {}
56
59 const typename Traits::LocalBasisType& localBasis () const
60 {
61 return basis;
62 }
63
67 {
68 return coefficients;
69 }
70
74 {
75 return interpolation;
76 }
77
79 unsigned int size () const
80 {
81 return basis.size();
82 }
83
86 GeometryType type () const
87 {
88 return gt;
89 }
90
91 private:
95 GeometryType gt;
96 };
97
99
111 template<class Geometry, class RF, std::size_t p>
113 typedef typename Geometry::ctype DF;
114 static const std::size_t dim = Geometry::mydimension;
115
117
118 std::vector<std::shared_ptr<const LocalFE> > localFEs;
119
120 void init(const GeometryType &gt) {
121 std::size_t index = gt.id() >> 1;
122 if(localFEs.size() <= index)
123 localFEs.resize(index+1);
124 localFEs[index].reset(new LocalFE(gt));
125 }
126
127 public:
130
132
136 template<class ForwardIterator>
137 MonomialFiniteElementFactory(const ForwardIterator &begin,
138 const ForwardIterator &end)
139 {
140 for(ForwardIterator it = begin; it != end; ++it)
141 init(*it);
142 }
143
145
148 MonomialFiniteElementFactory(const GeometryType &gt)
149 { init(gt); }
150
152
156 static_assert(dim <= 3, "MonomFiniteElementFactory knows the "
157 "available geometry types only up to dimension 3");
158
159 GeometryType gt;
160 switch(dim) {
161 case 0 :
162 gt.makeVertex(); init(gt);
163 break;
164 case 1 :
165 gt.makeLine(); init(gt);
166 break;
167 case 2 :
168 gt.makeTriangle(); init(gt);
169 gt.makeQuadrilateral(); init(gt);
170 break;
171 case 3 :
172 gt.makeTetrahedron(); init(gt);
173 gt.makePyramid(); init(gt);
174 gt.makePrism(); init(gt);
175 gt.makeHexahedron(); init(gt);
176 break;
177 default :
178 // this should never happen -- it should be caught by the static
179 // assert above.
180 std::abort();
181 };
182 }
183
185
195 const FiniteElement make(const Geometry& geometry) {
196 std::size_t index = geometry.type().id() >> 1;
197 assert(localFEs.size() > index && localFEs[index]);
198 return FiniteElement(*localFEs[index], geometry);
199 }
200 };
201}
202
203#endif // DUNE_LOCALFUNCTIONS_MONOMIAL_HH
Definition brezzidouglasmarini1cube2d.hh:14
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:187
GeometryType type() const
Definition localtoglobaladaptors.hh:229
Monomial basis for discontinuous Galerkin methods.
Definition monomial.hh:40
GeometryType type() const
Definition monomial.hh:86
const Traits::LocalInterpolationType & localInterpolation() const
Definition monomial.hh:73
MonomialLocalFiniteElement(const GeometryType &gt_)
Construct a MonomLocalFiniteElement.
Definition monomial.hh:53
const Traits::LocalBasisType & localBasis() const
Definition monomial.hh:59
const Traits::LocalCoefficientsType & localCoefficients() const
Definition monomial.hh:66
LocalFiniteElementTraits< MonomialLocalBasis< D, R, d, p, diffOrder >, MonomialLocalCoefficients< static_size >, MonomialLocalInterpolation< MonomialLocalBasis< D, R, d, p, diffOrder >, static_size > > Traits
Definition monomial.hh:50
unsigned int size() const
Number of shape functions in this finite element.
Definition monomial.hh:79
Factory for global-valued MonomFiniteElement objects.
Definition monomial.hh:112
MonomialFiniteElementFactory(const ForwardIterator &begin, const ForwardIterator &end)
construct a MonomialFiniteElementFactory from a list of GeometryType's
Definition monomial.hh:137
MonomialFiniteElementFactory(const GeometryType &gt)
construct a MonomialFiniteElementFactory from a single GeometryType
Definition monomial.hh:148
const FiniteElement make(const Geometry &geometry)
construct a global-valued MonomFiniteElement
Definition monomial.hh:195
ScalarLocalToGlobalFiniteElementAdaptor< LocalFE, Geometry > FiniteElement
Definition monomial.hh:129
MonomialFiniteElementFactory()
construct a MonomFiniteElementFactory for all applicable GeometryType's
Definition monomial.hh:155
Definition monomiallocalbasis.hh:21
Constant shape function.
Definition monomiallocalbasis.hh:233
Layout map for monomial finite elements.
Definition monomiallocalcoefficients.hh:22
Definition monomiallocalinterpolation.hh:19