Dune Core Modules (2.8.0)

raviartthomas0cube2dall.hh
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_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE2D_ALL_HH
4#define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE2D_ALL_HH
5
6#include <cstddef>
7#include <numeric>
8#include <vector>
9
11
12#include <dune/localfunctions/common/localbasis.hh>
13#include <dune/localfunctions/common/localkey.hh>
14#include <dune/localfunctions/common/localinterpolation.hh>
15
16namespace Dune
17{
26 template<class D, class R>
28 {
29 public:
32
34 RT0Cube2DLocalBasis (std::bitset<4> s = 0)
35 {
36 for (int i=0; i<4; i++)
37 sign_[i] = s[i] ? -1.0 : 1.0;
38 }
39
41 unsigned int size () const
42 {
43 return 4;
44 }
45
47 inline void evaluateFunction (const typename Traits::DomainType& in,
48 std::vector<typename Traits::RangeType>& out) const
49 {
50 out.resize(4);
51 out[0] = {sign_[0]*(in[0]-1.0), 0.0};
52 out[1] = {sign_[1]*(in[0]), 0.0};
53 out[2] = {0.0, sign_[2]*(in[1]-1.0)};
54 out[3] = {0.0, sign_[3]*(in[1])};
55 }
56
58 inline void
59 evaluateJacobian (const typename Traits::DomainType& in, // position
60 std::vector<typename Traits::JacobianType>& out) const // return value
61 {
62 out.resize(4);
63 out[0][0] = {sign_[0], 0};
64 out[0][1] = {0, 0};
65
66 out[1][0] = {sign_[1], 0};
67 out[1][1] = {0, 0};
68
69 out[2][0] = {0, 0};
70 out[2][1] = {0, sign_[2]};
71
72 out[3][0] = {0, 0};
73 out[3][1] = {0, sign_[3]};
74 }
75
77 void partial (const std::array<unsigned int, 2>& order,
78 const typename Traits::DomainType& in, // position
79 std::vector<typename Traits::RangeType>& out) const // return value
80 {
81 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
82 if (totalOrder == 0) {
83 evaluateFunction(in, out);
84 } else if (totalOrder == 1) {
85 auto const direction = std::distance(order.begin(), std::find(order.begin(), order.end(), 1));
86 out.resize(size());
87
88 for (std::size_t i = 0; i < size(); ++i)
89 out[i] = {0, 0};
90
91 switch (direction) {
92 case 0:
93 out[0][0] = sign_[0];
94 out[1][0] = sign_[1];
95 break;
96 case 1:
97 out[2][1] = sign_[2];
98 out[3][1] = sign_[3];
99 break;
100 default:
101 DUNE_THROW(RangeError, "Component out of range.");
102 }
103 } else {
104 out.resize(size());
105 for (std::size_t i = 0; i < size(); ++i)
106 for (std::size_t j = 0; j < 2; ++j)
107 out[i][j] = 0;
108 }
109
110 }
111
113 unsigned int order () const
114 {
115 return 1;
116 }
117
118 private:
119 std::array<R,4> sign_;
120 };
121
122
130 template<class LB>
132 {
133 public:
134
136 RT0Cube2DLocalInterpolation (std::bitset<4> s = 0)
137 {
138 for (int i=0; i<4; i++)
139 sign_[i] = s[i] ? -1.0 : 1.0;
140
141 m0 = {0.0, 0.5};
142 m1 = {1.0, 0.5};
143 m2 = {0.5, 0.0};
144 m3 = {0.5, 1.0};
145
146 n0 = {-1.0, 0.0};
147 n1 = { 1.0, 0.0};
148 n2 = { 0.0, -1.0};
149 n3 = { 0.0, 1.0};
150 }
151
152 template<typename F, typename C>
153 void interpolate (const F& ff, std::vector<C>& out) const
154 {
155 // f gives v*outer normal at a point on the edge!
156 auto&& f = Impl::makeFunctionWithCallOperator<typename LB::Traits::DomainType>(ff);
157
158 out.resize(4);
159
160 // Evaluate the normal components at the edge midpoints
161 auto y = f(m0); out[0] = (y[0]*n0[0]+y[1]*n0[1])*sign_[0];
162 y = f(m1); out[1] = (y[0]*n1[0]+y[1]*n1[1])*sign_[1];
163 y = f(m2); out[2] = (y[0]*n2[0]+y[1]*n2[1])*sign_[2];
164 y = f(m3); out[3] = (y[0]*n3[0]+y[1]*n3[1])*sign_[3];
165 }
166
167 private:
168 std::array<typename LB::Traits::RangeFieldType,4> sign_;
169
170 // The four edge midpoints of the reference quadrilateral
171 typename LB::Traits::DomainType m0,m1,m2,m3;
172
173 // The four edge normals of the reference quadrilateral
174 typename LB::Traits::DomainType n0,n1,n2,n3;
175 };
176
184 {
185 public:
188 {
189 for (std::size_t i=0; i<4; i++)
190 li[i] = LocalKey(i,1,0);
191 }
192
194 std::size_t size () const
195 {
196 return 4;
197 }
198
200 const LocalKey& localKey (std::size_t i) const
201 {
202 return li[i];
203 }
204
205 private:
206 std::vector<LocalKey> li;
207 };
208
209}
210#endif // DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE2D_ALL_HH
A dense n x m matrix.
Definition: fmatrix.hh:69
vector space out of a tensor product of fields.
Definition: fvector.hh:95
Describe position of one degree of freedom.
Definition: localkey.hh:21
Lowest order Raviart-Thomas shape functions on the reference quadrilateral.
Definition: raviartthomas0cube2dall.hh:28
RT0Cube2DLocalBasis(std::bitset< 4 > s=0)
Constructor with a set of edge orientations.
Definition: raviartthomas0cube2dall.hh:34
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: raviartthomas0cube2dall.hh:47
void partial(const std::array< unsigned int, 2 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition: raviartthomas0cube2dall.hh:77
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: raviartthomas0cube2dall.hh:59
unsigned int order() const
Polynomial order of the shape functions.
Definition: raviartthomas0cube2dall.hh:113
unsigned int size() const
number of shape functions
Definition: raviartthomas0cube2dall.hh:41
Layout map for RT0 elements on quadrilaterals.
Definition: raviartthomas0cube2dall.hh:184
RT0Cube2DLocalCoefficients()
Standard constructor.
Definition: raviartthomas0cube2dall.hh:187
std::size_t size() const
number of coefficients
Definition: raviartthomas0cube2dall.hh:194
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: raviartthomas0cube2dall.hh:200
Lowest order Raviart-Thomas shape functions on the reference quadrilateral.
Definition: raviartthomas0cube2dall.hh:132
RT0Cube2DLocalInterpolation(std::bitset< 4 > s=0)
Constructor with explicitly given edge orientations.
Definition: raviartthomas0cube2dall.hh:136
Default exception class for range errors.
Definition: exceptions.hh:252
Implements a matrix constructed from a given type representing a field and compile-time given number ...
#define DUNE_THROW(E, m)
Definition: exceptions.hh:216
constexpr T accumulate(Range &&range, T value, F &&f)
Accumulate values.
Definition: hybridutilities.hh:289
Dune namespace.
Definition: alignedallocator.hh:11
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:32
D DomainType
domain type
Definition: localbasis.hh:43
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 21, 23:30, 2024)