5#ifndef DUNE_LOCALFUNCTIONS_NEDELEC_NEDELEC1STKINDSIMPLEX_HH
6#define DUNE_LOCALFUNCTIONS_NEDELEC_NEDELEC1STKINDSIMPLEX_HH
13#include <dune/geometry/referenceelements.hh>
16#include <dune/localfunctions/common/localbasis.hh>
17#include <dune/localfunctions/common/localfiniteelementtraits.hh>
18#include <dune/localfunctions/common/localkey.hh>
36 template<
class D,
class R,
int dim,
int k>
37 class Nedelec1stKindSimplexLocalBasis
40 constexpr static std::size_t numberOfEdges = dim*(dim+1)/2;
43 using Traits = LocalBasisTraits<D,dim,FieldVector<D,dim>,
44 R,dim,FieldVector<R,dim>,
45 FieldMatrix<R,dim,dim> >;
53 Nedelec1stKindSimplexLocalBasis()
55 std::fill(edgeOrientation_.begin(), edgeOrientation_.end(), 1.0);
60 Nedelec1stKindSimplexLocalBasis(std::bitset<numberOfEdges> edgeOrientation)
61 : Nedelec1stKindSimplexLocalBasis()
63 for (std::size_t i=0; i<edgeOrientation_.size(); i++)
64 edgeOrientation_[i] *= edgeOrientation[i] ? -1.0 : 1.0;
68 static constexpr unsigned int size()
70 static_assert(dim==2 || dim==3,
"Nedelec shape functions are implemented only for 2d and 3d simplices.");
74 return k * (k+2) * (k+3) / 2;
83 std::vector<typename Traits::RangeType>& out)
const
85 static_assert(k==1,
"Evaluating Nédélec shape functions is implemented only for first order.");
93 out[0] = {D(1) - in[1], in[0]};
94 out[1] = {in[1], -in[0]+D(1)};
95 out[2] = {-in[1], in[0]};
117 out[0] = { 1 - in[1] - in[2], in[0] , in[0] };
118 out[1] = { in[1] , 1 - in[0] - in[2], in[1]};
119 out[2] = { - in[1] , in[0] , 0 };
120 out[3] = { in[2], in[2], 1 - in[0] - in[1]};
121 out[4] = { -in[2], 0 , in[0] };
122 out[5] = { 0 , -in[2], in[1]};
125 for (std::size_t i=0; i<out.size(); i++)
126 out[i] *= edgeOrientation_[i];
135 std::vector<typename Traits::JacobianType>& out)
const
140 out[0][0] = { 0, -1};
146 out[2][0] = { 0, -1};
151 out[0][0] = { 0,-1,-1};
152 out[0][1] = { 1, 0, 0};
153 out[0][2] = { 1, 0, 0};
155 out[1][0] = { 0, 1, 0};
156 out[1][1] = {-1, 0, -1};
157 out[1][2] = { 0, 1, 0};
159 out[2][0] = { 0, -1, 0};
160 out[2][1] = { 1, 0, 0};
161 out[2][2] = { 0, 0, 0};
163 out[3][0] = { 0, 0, 1};
164 out[3][1] = { 0, 0, 1};
165 out[3][2] = {-1, -1, 0};
167 out[4][0] = { 0, 0, -1};
168 out[4][1] = { 0, 0, 0};
169 out[4][2] = { 1, 0, 0};
171 out[5][0] = { 0, 0, 0};
172 out[5][1] = { 0, 0, -1};
173 out[5][2] = { 0, 1, 0};
176 for (std::size_t i=0; i<out.size(); i++)
177 out[i] *= edgeOrientation_[i];
187 void partial(
const std::array<unsigned int, dim>& order,
189 std::vector<typename Traits::RangeType>& out)
const
192 if (totalOrder == 0) {
193 evaluateFunction(in, out);
194 }
else if (totalOrder == 1) {
195 auto const direction = std::distance(order.begin(), std::find(order.begin(), order.end(), 1));
247 for (std::size_t i=0; i<out.size(); i++)
248 out[i] *= edgeOrientation_[i];
252 for (std::size_t i = 0; i <
size(); ++i)
253 for (std::size_t j = 0; j < dim; ++j)
260 unsigned int order()
const
268 std::array<R,numberOfEdges> edgeOrientation_;
276 template <
int dim,
int k>
277 class Nedelec1stKindSimplexLocalCoefficients
281 Nedelec1stKindSimplexLocalCoefficients ()
284 static_assert(k==1,
"Only first-order Nédélec local coefficients are implemented.");
287 for (std::size_t i=0; i<
size(); i++)
288 localKey_[i] = LocalKey(i,dim-1,0);
292 std::size_t
size()
const
294 static_assert(dim==2 || dim==3,
"Nédélec shape functions are implemented only for 2d and 3d simplices.");
295 return (dim==2) ? k * (k+2)
296 : k * (k+2) * (k+3) / 2;
301 const LocalKey& localKey (std::size_t i)
const
307 std::vector<LocalKey> localKey_;
315 class Nedelec1stKindSimplexLocalInterpolation
317 static constexpr auto dim = LB::Traits::dimDomain;
318 static constexpr auto size = LB::size();
321 constexpr static std::size_t numberOfEdges = dim*(dim+1)/2;
326 Nedelec1stKindSimplexLocalInterpolation (std::bitset<numberOfEdges> s = 0)
330 for (std::size_t i=0; i<numberOfEdges; i++)
331 m_[i] = refElement.position(i,dim-1);
333 for (std::size_t i=0; i<numberOfEdges; i++)
335 auto vertexIterator = refElement.subEntities(i,dim-1,dim).begin();
336 auto v0 = *vertexIterator;
337 auto v1 = *(++vertexIterator);
342 edge_[i] = refElement.position(v1,dim) - refElement.position(v0,dim);
343 edge_[i] *= (s[i]) ? -1.0 : 1.0;
352 template<
typename F,
typename C>
353 void interpolate (
const F& f, std::vector<C>& out)
const
357 for (std::size_t i=0; i<
size; i++)
361 for (
int j=0; j<dim; j++)
362 out[i] += y[j]*edge_[i][j];
368 std::array<typename LB::Traits::DomainType, numberOfEdges> m_;
370 std::array<typename LB::Traits::DomainType, numberOfEdges> edge_;
401 template<
class D,
class R,
int dim,
int k>
406 Impl::Nedelec1stKindSimplexLocalCoefficients<dim,k>,
407 Impl::Nedelec1stKindSimplexLocalInterpolation<Impl::Nedelec1stKindSimplexLocalBasis<D,R,dim,k> > >;
409 static_assert(dim==2 || dim==3,
"Nedelec elements are only implemented for 2d and 3d elements.");
410 static_assert(k==1,
"Nedelec elements of the first kind are currently only implemented for order k==1.");
433 return coefficients_;
438 return interpolation_;
441 static constexpr unsigned int size ()
443 return Traits::LocalBasisType::size();
Nédélec elements of the first kind for simplex elements.
Definition: nedelec1stkindsimplex.hh:403
Nedelec1stKindSimplexLocalFiniteElement()=default
Default constructor.
Nedelec1stKindSimplexLocalFiniteElement(std::bitset< dim *(dim+1)/2 > s)
Constructor with explicitly given edge orientations.
Definition: nedelec1stkindsimplex.hh:421
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Implements a vector constructed from a given type representing a field and a compile-time given size.
static void interpolate(const GridFunction &u, DiscreteFunction &v)
perform native interpolation of a discrete function space
Definition: interpolate.hh:54
constexpr GeometryType simplex(unsigned int dim)
Returns a GeometryType representing a simplex of dimension dim.
Definition: type.hh:453
constexpr T accumulate(Range &&range, T value, F &&f)
Accumulate values.
Definition: hybridutilities.hh:279
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
D DomainType
domain type
Definition: localbasis.hh:43
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
A unique label for each type of element that can occur in a grid.