\n
\n
\n
\n-
5#ifndef DUNE_LOCALFUNCTIONS_NEDELEC_NEDELEC1STKINDCUBE_HH \n-
6#define DUNE_LOCALFUNCTIONS_NEDELEC_NEDELEC1STKINDCUBE_HH \n+
5#ifndef DUNE_RANNACHER_TUREK_3D_LOCALBASIS_HH \n+
6#define DUNE_RANNACHER_TUREK_3D_LOCALBASIS_HH \n
\n
\n-
\n-
10#include <dune/common/fmatrix.hh> \n+
\n+
\n
11#include <dune/common/fvector.hh> \n-
12#include <dune/common/math.hh> \n+
12#include <dune/common/fmatrix.hh> \n
\n-
14#include <dune/geometry/referenceelements.hh> \n-
15#include <dune/geometry/type.hh> \n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
36 template<
class D,
class R,
int dim,
int k>
\n-
37 class Nedelec1stKindCubeLocalBasis
\n-
\n-
\n-
40 constexpr static std::size_t numberOfEdges = power(2,dim-1)*dim;
\n-
\n-
\n-
43 using Traits = LocalBasisTraits<D,dim,FieldVector<D,dim>,
\n-
44 R,dim,FieldVector<R,dim>,
\n-
45 FieldMatrix<R,dim,dim> >;
\n-
\n-
53 Nedelec1stKindCubeLocalBasis()
\n-
\n-
55 std::fill(edgeOrientation_.begin(), edgeOrientation_.end(), 1.0);
\n-
\n-
\n-
60 Nedelec1stKindCubeLocalBasis(std::bitset<numberOfEdges> edgeOrientation)
\n-
61 : Nedelec1stKindCubeLocalBasis()
\n-
\n-
63 for (std::size_t i=0; i<edgeOrientation_.size(); i++)
\n-
64 edgeOrientation_[i] *= edgeOrientation[i] ? -1.0 : 1.0;
\n-
\n-
\n-
68 static constexpr unsigned int size()
\n-
\n-
70 static_assert(dim==2 || dim==3,
"Nedelec shape functions are implemented only for 2d and 3d cubes.");
\n-
\n-
\n-
\n-
74 return 3*k * (k+1) * (k+1);
\n-
\n-
\n-
\n-
83 std::vector<typename Traits::RangeType>& out)
const \n-
\n-
85 static_assert(k==1,
"Evaluating N\u00e9d\u00e9lec shape functions is implemented only for first order.");
\n-
\n+
\n+
\n+
\n+
\n+
\n+
19 template<
class D,
class R >
\n+
\n+
\n+
\n+
22 static const int coefficients[ 6 ][ 6 ];
\n+
\n+
\n+
\n+
26 R, 1, FieldVector< R, 1 >,
\n+
\n+
\n+
\n+
\n+
\n+
\n+
37 std::vector< typename Traits::RangeType > &out )
const \n+
\n+
\n+
40 RangeFieldType y[ 6 ] = { 1, in[ 0 ], in[ 1 ], in[ 2 ],
\n+
41 in[ 0 ]*in[ 0 ] - in[ 1 ]*in[ 1 ],
\n+
42 in[ 1 ]*in[ 1 ] - in[ 2 ]*in[ 2 ] };
\n+
\n+
44 for(
unsigned int i = 0; i <
size(); ++i )
\n+
\n+
46 out[ i ] = RangeFieldType( 0 );
\n+
47 for(
unsigned int j = 0; j < 6; ++j )
\n+
48 out[ i ] += coefficients[ i ][ j ]*y[ j ];
\n+
49 out[ i ] /= RangeFieldType( 3 );
\n+
\n+
\n+
\n+
\n+
\n+
\n+
55 std::vector< typename Traits::JacobianType > &out )
const \n+
\n+
\n+
58 RangeFieldType y0[ 5 ] = { 1, 0, 0, 2*in[ 0 ], 0 };
\n+
59 RangeFieldType y1[ 5 ] = { 0, 1, 0, -2*in[ 1 ], 2*in[ 1 ] };
\n+
60 RangeFieldType y2[ 5 ] = { 0, 0, 1, 0, -2*in[ 2 ] };
\n+
\n+
\n+
63 for(
unsigned int i = 0; i <
size(); ++i )
\n+
\n+
65 out[ i ] = RangeFieldType( 0 );
\n+
66 for(
unsigned int j = 0; j < 5; ++j )
\n+
\n+
68 out[ i ][ 0 ][ 0 ] += coefficients[ i ][ j+1 ]*y0[ j ];
\n+
69 out[ i ][ 0 ][ 1 ] += coefficients[ i ][ j+1 ]*y1[ j ];
\n+
70 out[ i ][ 0 ][ 2 ] += coefficients[ i ][ j+1 ]*y2[ j ];
\n+
\n+
72 out[ i ] /= RangeFieldType( 3 );
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
79 std::vector<typename Traits::RangeType>& out)
const \n+
\n+
81 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
\n+
82 if (totalOrder == 0) {
\n+
\n+
84 }
else if (totalOrder == 1) {
\n+
\n+
86 auto const direction = std::distance(
order.begin(), std::find(
order.begin(),
order.end(), 1));
\n
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
102 out[0] = { 0, D(1) - in[0]};
\n-
103 out[1] = { 0, in[0]};
\n-
104 out[2] = { D(1) - in[1], 0};
\n-
105 out[3] = { in[1], 0};
\n-
\n-
\n-
108 if constexpr (dim==3)
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n+
\n+
89 RangeFieldType y[3][5] = { { 1.0, 0.0, 0.0, 2*in[0], 0.0 },
\n+
90 { 0.0, 1.0, 0.0, -2*in[1], 2*in[1] },
\n+
91 { 0.0, 0.0, 1.0, 0.0, -2*in[2] } };
\n+
\n+
93 for (std::size_t i = 0; i <
size(); ++i) {
\n+
94 out[i] = RangeFieldType{0};
\n+
95 for (std::size_t j = 0; j < 5; ++j)
\n+
96 out[i] += coefficients[i][j+1] * y[direction][j];
\n+
97 out[i] /= RangeFieldType{3};
\n+
\n+
\n+
100 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
\n+
116 template<
class D,
class R >
\n+
117 const int RannacherTurek3DLocalBasis< D, R >
\n+
118 ::coefficients[ 6 ][ 6 ] = {{ 2, -7, 2, 2, 4, 2 },
\n+
119 { -1, -1, 2, 2, 4, 2 },
\n+
120 { 2, 2, -7, 2, -2, 2 },
\n+
121 { -1, 2, -1, 2, -2, 2 },
\n+
122 { 2, 2, 2, -7, -2, -4 },
\n+
123 { -1, 2, 2, -1, -2, -4 }};
\n+
\n+
\n
\n-
127 for (std::size_t i=0; i<out.size(); i++)
\n-
\n-
\n-
130 out[0][2] = { 1 - in[0] - in[1] + in[0]*in[1]};
\n-
131 out[1][2] = { in[0] - in[0]*in[1]};
\n-
132 out[2][2] = { in[1] - in[0]*in[1]};
\n-
133 out[3][2] = { in[0]*in[1]};
\n-
\n-
135 out[4][1] = { 1 - in[0] - in[2] + in[0]*in[2]};
\n-
136 out[5][1] = { in[0] - in[0]*in[2]};
\n-
137 out[8][1] = { in[2] - in[0]*in[2]};
\n-
138 out[9][1] = { in[0]*in[2]};
\n-
\n-
140 out[6][0] = { 1 - in[1] - in[2] + in[1]*in[2]};
\n-
141 out[7][0] = { in[1] - in[1]*in[2]};
\n-
142 out[10][0] = { in[2] - in[1]*in[2]};
\n-
143 out[11][0] = { in[1]*in[2]};
\n-
\n-
\n-
146 for (std::size_t i=0; i<out.size(); i++)
\n-
147 out[i] *= edgeOrientation_[i];
\n-
\n-
\n-
\n-
156 std::vector<typename Traits::JacobianType>& out)
const \n-
\n-
\n-
\n-
\n-
161 for (std::size_t i=0; i<out.size(); i++)
\n-
162 for (std::size_t j=0; j<dim; j++)
\n-
\n-
\n-
165 out[0][1] = { -1, 0};
\n-
\n-
\n-
168 out[2][0] = { 0, -1};
\n-
\n-
\n-
\n-
\n-
173 for (std::size_t i=0; i<out.size(); i++)
\n-
174 for(std::size_t j=0;j<dim; j++)
\n-
\n-
\n-
\n-
178 out[0][2] = {-1 +in[1], -1 + in[0], 0};
\n-
179 out[1][2] = { 1 -in[1], - in[0], 0};
\n-
180 out[2][2] = { -in[1], 1 - in[0], 0};
\n-
181 out[3][2] = { in[1], in[0], 0};
\n-
\n-
183 out[4][1] = {-1 +in[2], 0, -1 + in[0]};
\n-
184 out[5][1] = { 1 -in[2], 0, - in[0]};
\n-
185 out[8][1] = { -in[2], 0, 1 - in[0]};
\n-
186 out[9][1] = { in[2], 0, in[0]};
\n-
\n-
188 out[6][0] = { 0, -1 + in[2], -1 + in[1]};
\n-
189 out[7][0] = { 0, 1 - in[2], - in[1]};
\n-
190 out[10][0] = { 0, - in[2], 1 - in[1]};
\n-
191 out[11][0] = { 0, in[2], in[1]};
\n-
\n-
\n-
\n-
195 for (std::size_t i=0; i<out.size(); i++)
\n-
196 out[i] *= edgeOrientation_[i];
\n-
\n-
\n-
\n-
206 void partial(
const std::array<unsigned int, dim>& order,
\n-
\n-
208 std::vector<typename Traits::RangeType>& out)
const \n-
\n-
210 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
\n-
211 if (totalOrder == 0) {
\n-
212 evaluateFunction(in, out);
\n-
213 }
else if (totalOrder == 1) {
\n-
214 auto const direction = std::distance(order.begin(), std::find(order.begin(), order.end(), 1));
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
240 out[0] = { 0, 0, -1 +in[1]};
\n-
241 out[1] = { 0, 0, 1 -in[1]};
\n-
242 out[2] = { 0, 0, -in[1]};
\n-
243 out[3] = { 0, 0, in[1]};
\n-
\n-
245 out[4] = { 0, -1 +in[2], 0};
\n-
246 out[5] = { 0, 1 -in[2], 0};
\n-
247 out[8] = { 0, -in[2], 0};
\n-
248 out[9] = { 0, in[2], 0};
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
257 out[0] = { 0, 0, -1 + in[0]};
\n-
258 out[1] = { 0, 0, - in[0]};
\n-
259 out[2] = { 0, 0, 1 - in[0]};
\n-
260 out[3] = { 0, 0, in[0]};
\n-
\n-
\n-
\n-
\n-
\n-
\n-
267 out[6] = { -1 + in[2], 0, 0};
\n-
268 out[7] = { 1 - in[2], 0, 0};
\n-
269 out[10] = { - in[2], 0, 0};
\n-
270 out[11] = { in[2], 0, 0};
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
279 out[4] = { 0, -1 + in[0], 0};
\n-
280 out[5] = { 0, - in[0], 0};
\n-
281 out[8] = { 0, 1 - in[0], 0};
\n-
282 out[9] = { 0, in[0], 0};
\n-
\n-
284 out[6] = { -1 + in[1], 0, 0};
\n-
285 out[7] = { - in[1], 0, 0};
\n-
286 out[10] = { 1 - in[1], 0, 0};
\n-
287 out[11] = { in[1], 0, 0};
\n-
\n-
\n-
\n-
\n-
292 for (std::size_t i=0; i<out.size(); i++)
\n-
293 out[i] *= edgeOrientation_[i];
\n-
\n-
295 }
else if (totalOrder == 2) {
\n-
\n-
\n-
\n-
299 for (std::size_t i = 0; i < size(); ++i)
\n-
300 for (std::size_t j = 0; j < dim; ++j)
\n-
\n-
\n-
\n-
\n-
305 for(
size_t i=0; i<out.size(); i++)
\n-
\n-
\n-
\n-
309 if( order[0] == 1 and order[1]==1)
\n-
\n-
\n-
312 out[1] = { 0, 0, -1};
\n-
313 out[2] = { 0, 0, -1};
\n-
\n-
\n-
\n-
\n-
318 if( order[0] == 1 and order[2]==1)
\n-
\n-
\n-
321 out[5] = { 0, -1, 0};
\n-
322 out[8] = { 0, -1, 0};
\n-
\n-
\n-
\n-
\n-
327 if( order[1] == 1 and order[2]==1)
\n-
\n-
\n-
330 out[7] = { -1, 0, 0};
\n-
331 out[10] = { -1, 0, 0};
\n-
332 out[11] = { 1, 0, 0};
\n-
\n-
\n-
335 for (std::size_t i=0; i<out.size(); i++)
\n-
336 out[i] *= edgeOrientation_[i];
\n-
\n-
\n-
\n-
\n-
\n-
342 for (std::size_t i = 0; i < size(); ++i)
\n-
343 for (std::size_t j = 0; j < dim; ++j)
\n-
\n-
\n-
\n-
\n-
\n-
350 unsigned int order()
const \n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
361 std::array<R,numberOfEdges> edgeOrientation_;
\n-
\n-
\n-
\n-
369 template <
int dim,
int k>
\n-
370 class Nedelec1stKindCubeLocalCoefficients
\n-
\n-
\n-
374 Nedelec1stKindCubeLocalCoefficients ()
\n-
\n-
\n-
377 static_assert(k==1,
"Only first-order N\u00e9d\u00e9lec local coefficients are implemented.");
\n-
\n-
\n-
380 for (std::size_t i=0; i<size(); i++)
\n-
381 localKey_[i] = LocalKey(i,dim-1,0);
\n-
\n-
\n-
385 std::size_t size()
const \n-
\n-
387 static_assert(dim==2 || dim==3,
"N\u00e9d\u00e9lec shape functions are implemented only for 2d and 3d cubes.");
\n-
388 return (dim==2) ? 2*k * (k+1)
\n-
389 : 3*k * (k+1) * (k+1);
\n-
\n-
\n-
394 const LocalKey& localKey (std::size_t i)
const \n-
\n-
\n-
\n-
\n-
\n-
400 std::vector<LocalKey> localKey_;
\n-
\n-
\n-
\n-
408 class Nedelec1stKindCubeLocalInterpolation
\n-
\n-
410 static constexpr auto dim = LB::Traits::dimDomain;
\n-
411 static constexpr auto size = LB::size();
\n-
\n-
\n-
414 constexpr static std::size_t numberOfEdges = power(2,dim-1)*dim;
\n-
\n-
\n-
\n-
419 Nedelec1stKindCubeLocalInterpolation (std::bitset<numberOfEdges> s = 0)
\n-
\n-
421 auto refElement = Dune::referenceElement<double,dim>(GeometryTypes::cube(dim));
\n-
\n-
423 for (std::size_t i=0; i<numberOfEdges; i++)
\n-
424 m_[i] = refElement.position(i,dim-1);
\n-
\n-
426 for (std::size_t i=0; i<numberOfEdges; i++)
\n-
\n-
428 auto vertexIterator = refElement.subEntities(i,dim-1,dim).begin();
\n-
429 auto v0 = *vertexIterator;
\n-
430 auto v1 = *(++vertexIterator);
\n-
\n-
\n-
\n-
\n-
\n-
436 edge_[i] = refElement.position(v1,dim) - refElement.position(v0,dim);
\n-
437 edge_[i] *= (s[i]) ? -1.0 : 1.0;
\n-
\n-
\n-
\n-
446 template<
typename F,
typename C>
\n-
447 void interpolate (
const F& ff, std::vector<C>& out)
const \n-
\n-
\n-
450 auto&& f = Impl::makeFunctionWithCallOperator<typename LB::Traits::DomainType>(ff);
\n-
\n-
452 for (std::size_t i=0; i<size; i++)
\n-
\n-
\n-
\n-
456 for (
int j=0; j<dim; j++)
\n-
457 out[i] += y[j]*edge_[i][j];
\n-
\n-
\n-
\n-
\n-
\n-
463 std::array<typename LB::Traits::DomainType, numberOfEdges> m_;
\n-
\n-
465 std::array<typename LB::Traits::DomainType, numberOfEdges> edge_;
\n-
\n-
\n-
\n-
\n-
\n-
494 template<
class D,
class R,
int dim,
int k>
\n-
\n-
\n-
\n-
\n-
\n-
499 Impl::Nedelec1stKindCubeLocalCoefficients<dim,k>,
\n-
500 Impl::Nedelec1stKindCubeLocalInterpolation<Impl::Nedelec1stKindCubeLocalBasis<D,R,dim,k> > >;
\n-
\n-
502 static_assert(dim==2 || dim==3,
"Nedelec elements are only implemented for 2d and 3d elements.");
\n-
503 static_assert(k==1,
"Nedelec elements of the first kind are currently only implemented for order k==1.");
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
526 return coefficients_;
\n-
\n-
\n-
\n-
\n-
\n-
\n-
531 return interpolation_;
\n-
\n-
\n-
\n-
\n-
534 static constexpr unsigned int size ()
\n-
\n-
536 return Traits::LocalBasisType::size();
\n-
\n-
\n-
\n-
\n-
539 static constexpr GeometryType
type ()
\n-
\n-
541 return GeometryTypes::cube(dim);
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n-
\n+
\n
Definition bdfmcube.hh:18
\n+
Type traits for LocalBasisVirtualInterface.
Definition common/localbasis.hh:34
\n
D DomainType
domain type
Definition common/localbasis.hh:42
\n-
traits helper struct
Definition localfiniteelementtraits.hh:13
\n-
LB LocalBasisType
Definition localfiniteelementtraits.hh:16
\n-
LC LocalCoefficientsType
Definition localfiniteelementtraits.hh:20
\n-
LI LocalInterpolationType
Definition localfiniteelementtraits.hh:24
\n-
N\u00e9d\u00e9lec elements of the first kind for cube elements.
Definition nedelec1stkindcube.hh:496
\n-
const Traits::LocalInterpolationType & localInterpolation() const
Definition nedelec1stkindcube.hh:529
\n-
static constexpr unsigned int size()
Definition nedelec1stkindcube.hh:534
\n-
Nedelec1stKindCubeLocalFiniteElement()=default
Default constructor.
\n-
static constexpr GeometryType type()
Definition nedelec1stkindcube.hh:539
\n-
const Traits::LocalCoefficientsType & localCoefficients() const
Definition nedelec1stkindcube.hh:524
\n-
const Traits::LocalBasisType & localBasis() const
Definition nedelec1stkindcube.hh:519
\n-
Nedelec1stKindCubeLocalFiniteElement(std::bitset< power(2, dim-1) *dim > s)
Constructor with explicitly given edge orientations.
Definition nedelec1stkindcube.hh:514
\n+
RF RangeFieldType
Export type for range field.
Definition common/localbasis.hh:45
\n+
Definition rannacherturek3dlocalbasis.hh:21
\n+
LocalBasisTraits< D, 3, FieldVector< D, 3 >, R, 1, FieldVector< R, 1 >, FieldMatrix< R, 1, 3 > > Traits
Definition rannacherturek3dlocalbasis.hh:27
\n+
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
evaluate all shape functions
Definition rannacherturek3dlocalbasis.hh:36
\n+
unsigned int size() const
number of shape functions
Definition rannacherturek3dlocalbasis.hh:30
\n+
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
evaluate jacobian of all shape functions
Definition rannacherturek3dlocalbasis.hh:54
\n+
unsigned int order() const
polynomial order of the shape functions
Definition rannacherturek3dlocalbasis.hh:105
\n+
void partial(const std::array< unsigned int, 3 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition rannacherturek3dlocalbasis.hh:77
\n
\n-
\n
![\"doxygen\"/](\"doxygen.svg\")
1.9.8\n \n \n \n", "details": [{"source1": "html2text {}", "source2": "html2text {}", "unified_diff": "@@ -1,26 +1,12 @@\n dune-localfunctions\u00a02.9.0\n Loading...\n Searching...\n No Matches\n * _\bd_\bu_\bn_\be\n * _\bl_\bo_\bc_\ba_\bl_\bf_\bu_\bn_\bc_\bt_\bi_\bo_\bn_\bs\n- * _\br_\be_\bf_\bi_\bn_\be_\bd\n- * _\br_\be_\bf_\bi_\bn_\be_\bd_\bp_\b0\n-_\bC_\bl_\ba_\bs_\bs_\be_\bs | _\bN_\ba_\bm_\be_\bs_\bp_\ba_\bc_\be_\bs\n-refinedp0localbasis.hh File Reference\n-#include