dune-localfunctions 3.0-git
pk1d.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#ifndef DUNE_PK1DLOCALFINITEELEMENT_HH
4#define DUNE_PK1DLOCALFINITEELEMENT_HH
5
6#include <cstddef>
7
8#include <dune/geometry/type.hh>
9
15
16namespace Dune
17{
18
21 template<class D, class R, unsigned int k>
23 {
24 public:
30
34 {
35 gt.makeLine();
36 }
37
40 Pk1DLocalFiniteElement (int variant) : coefficients(variant)
41 {
42 gt.makeLine();
43 }
44
49 Pk1DLocalFiniteElement (const unsigned int vertexmap[3]) : coefficients(vertexmap)
50 {
51 gt.makeLine();
52 }
53
56 const typename Traits::LocalBasisType& localBasis () const
57 {
58 return basis;
59 }
60
64 {
65 return coefficients;
66 }
67
71 {
72 return interpolation;
73 }
74
76 unsigned int size () const
77 {
78 return basis.size();
79 }
80
83 GeometryType type () const
84 {
85 return gt;
86 }
87
89 {
90 return new Pk1DLocalFiniteElement(*this);
91 }
92
93 private:
95 Pk1DLocalCoefficients<k> coefficients;
97 GeometryType gt;
98 };
99
101
108 template<class Geometry, class RF, std::size_t k>
110 typedef typename Geometry::ctype DF;
113
114 public:
126
127 private:
128 static const GeometryType gt;
129 static const LocalBasis localBasis;
130 static const LocalInterpolation localInterpolation;
131
132 typename Traits::Basis basis_;
133 typename Traits::Interpolation interpolation_;
134 typename Traits::Coefficients coefficients_;
135
136 public:
138
151 template<class VertexOrder>
152 Pk1DFiniteElement(const Geometry &geometry,
153 const VertexOrder& vertexOrder) :
154 basis_(localBasis, geometry), interpolation_(localInterpolation),
155 coefficients_(vertexOrder.begin(0, 0))
156 { }
157
158 const typename Traits::Basis& basis() const { return basis_; }
159 const typename Traits::Interpolation& interpolation() const
160 { return interpolation_; }
161 const typename Traits::Coefficients& coefficients() const
162 { return coefficients_; }
163 const GeometryType &type() const { return gt; }
164 };
165
166 template<class Geometry, class RF, std::size_t k>
167 const GeometryType
168 Pk1DFiniteElement<Geometry, RF, k>::gt(GeometryType::simplex, 2);
169
170 template<class Geometry, class RF, std::size_t k>
171 const typename Pk1DFiniteElement<Geometry, RF, k>::LocalBasis
172 Pk1DFiniteElement<Geometry, RF, k>::localBasis = LocalBasis();
173
174 template<class Geometry, class RF, std::size_t k>
175 const typename Pk1DFiniteElement<Geometry, RF, k>::LocalInterpolation
176 Pk1DFiniteElement<Geometry, RF, k>::localInterpolation =
177 LocalInterpolation();
178
180
190 template<class Geometry, class RF, std::size_t k>
193
195
209 template<class VertexOrder>
210 const FiniteElement make(const Geometry& geometry,
211 const VertexOrder& vertexOrder)
212 { return FiniteElement(geometry, vertexOrder); }
213 };
214}
215
216#endif
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
Traits class for local-to-global basis adaptors.
Definition localtoglobaladaptors.hh:28
Convert a simple scalar local basis into a global basis.
Definition localtoglobaladaptors.hh:65
Convert a local interpolation into a global interpolation.
Definition localtoglobaladaptors.hh:149
Definition pk1d.hh:23
Pk1DLocalFiniteElement(const unsigned int vertexmap[3])
Definition pk1d.hh:49
Pk1DLocalFiniteElement()
Definition pk1d.hh:33
LocalFiniteElementTraits< Pk1DLocalBasis< D, R, k >, Pk1DLocalCoefficients< k >, Pk1DLocalInterpolation< Pk1DLocalBasis< D, R, k > > > Traits
Definition pk1d.hh:29
unsigned int size() const
Number of shape functions in this finite element.
Definition pk1d.hh:76
const Traits::LocalCoefficientsType & localCoefficients() const
Definition pk1d.hh:63
Pk1DLocalFiniteElement(int variant)
Definition pk1d.hh:40
Pk1DLocalFiniteElement * clone() const
Definition pk1d.hh:88
GeometryType type() const
Definition pk1d.hh:83
const Traits::LocalBasisType & localBasis() const
Definition pk1d.hh:56
const Traits::LocalInterpolationType & localInterpolation() const
Definition pk1d.hh:70
Langrange finite element of arbitrary order on triangles.
Definition pk1d.hh:109
Pk1DFiniteElement(const Geometry &geometry, const VertexOrder &vertexOrder)
construct a Pk1DFiniteElement
Definition pk1d.hh:152
const Traits::Basis & basis() const
Definition pk1d.hh:158
const GeometryType & type() const
Definition pk1d.hh:163
const Traits::Coefficients & coefficients() const
Definition pk1d.hh:161
const Traits::Interpolation & interpolation() const
Definition pk1d.hh:159
Definition pk1d.hh:118
ScalarLocalToGlobalBasisAdaptor< LocalBasis, Geometry > Basis
Definition pk1d.hh:119
Pk1DLocalCoefficients< k > Coefficients
Definition pk1d.hh:124
LocalToGlobalInterpolationAdaptor< LocalInterpolation, typename Basis::Traits > Interpolation
Definition pk1d.hh:123
Factory for Pk1DFiniteElement objects.
Definition pk1d.hh:191
const FiniteElement make(const Geometry &geometry, const VertexOrder &vertexOrder)
construct Pk1DFiniteElementFactory
Definition pk1d.hh:210
Pk1DFiniteElement< Geometry, RF, k > FiniteElement
Definition pk1d.hh:192
Lagrange shape functions of arbitrary order on the 1D reference triangle.
Definition pk1dlocalbasis.hh:26
Layout map for Pk elements.
Definition pk1dlocalcoefficients.hh:23
Definition pk1dlocalinterpolation.hh:12