4 #ifndef DUNE_LOCALFUNCTIONS_QK_LOCALCOEFFICIENTS_HH 5 #define DUNE_LOCALFUNCTIONS_QK_LOCALCOEFFICIENTS_HH 11 #include <dune/common/exceptions.hh> 12 #include <dune/common/power.hh> 23 template<
int k,
int d>
26 static const unsigned unsignedK = k;
29 static std::array<unsigned int,d> multiindex (
unsigned int i)
31 std::array<unsigned int,d> alpha;
32 for (
int j=0; j<d; j++)
41 void setup1d(std::vector<unsigned int>& subEntity)
57 subEntity[lastIndex++] = 0;
58 for (
unsigned i = 0; i < unsignedK - 1; ++i)
59 subEntity[lastIndex++] = 0;
61 subEntity[lastIndex++] = 1;
63 assert((StaticPower<k+1,d>::power==lastIndex));
66 void setup2d(std::vector<unsigned int>& subEntity)
89 subEntity[lastIndex++] = 0;
90 for (
unsigned i = 0; i < unsignedK - 1; ++i)
91 subEntity[lastIndex++] = 2;
93 subEntity[lastIndex++] = 1;
96 for (
unsigned e = 0; e < unsignedK - 1; ++e) {
97 subEntity[lastIndex++] = 0;
98 for (
unsigned i = 0; i < unsignedK - 1; ++i)
99 subEntity[lastIndex++] = 0;
100 subEntity[lastIndex++] = 1;
104 subEntity[lastIndex++] = 2;
105 for (
unsigned i = 0; i < unsignedK - 1; ++i)
106 subEntity[lastIndex++] = 3;
108 subEntity[lastIndex++] = 3;
110 assert((StaticPower<k+1,d>::power==lastIndex));
115 void setup3d(std::vector<unsigned int>& subEntity)
124 unsigned lastIndex=0;
126 const unsigned numIndices = StaticPower<k+1,d>::power;
127 const unsigned numFaceIndices = StaticPower<k+1,d-1>::power;
129 const unsigned numInnerEdgeDofs = k-1;
150 subEntity[lastIndex++] = 0;
151 for (
unsigned i = 0; i < numInnerEdgeDofs; ++i)
152 subEntity[lastIndex++] = 6;
154 subEntity[lastIndex++] = 1;
157 for (
unsigned e = 0; e < numInnerEdgeDofs; ++e) {
158 subEntity[lastIndex++] = 4;
159 for (
unsigned i = 0; i < numInnerEdgeDofs; ++i)
160 subEntity[lastIndex++] = 4;
161 subEntity[lastIndex++] = 5;
165 subEntity[lastIndex++] = 2;
166 for (
unsigned i = 0; i < unsignedK - 1; ++i)
167 subEntity[lastIndex++] = 7;
168 subEntity[lastIndex++] = 3;
170 assert(numFaceIndices==lastIndex);
174 for(
unsigned f = 0; f < numInnerEdgeDofs; ++f) {
177 subEntity[lastIndex++] = 0;
178 for (
unsigned i = 0; i < numInnerEdgeDofs; ++i)
179 subEntity[lastIndex++] = 2;
180 subEntity[lastIndex++] = 1;
183 for (
unsigned e = 0; e < numInnerEdgeDofs; ++e) {
184 subEntity[lastIndex++] = 0;
185 for (
unsigned i = 0; i < numInnerEdgeDofs; ++i)
186 subEntity[lastIndex++] = 0;
187 subEntity[lastIndex++] = 1;
191 subEntity[lastIndex++] = 2;
192 for (
unsigned i = 0; i < numInnerEdgeDofs; ++i)
193 subEntity[lastIndex++] = 3;
194 subEntity[lastIndex++] = 3;
196 assert(lastIndex==(f+1+1)*numFaceIndices);
201 subEntity[lastIndex++] = 4;
202 for (
unsigned i = 0; i < unsignedK - 1; ++i)
203 subEntity[lastIndex++] = 10;
204 subEntity[lastIndex++] = 5;
207 for (
unsigned e = 0; e < unsignedK - 1; ++e) {
208 subEntity[lastIndex++] = 8;
209 for (
unsigned i = 0; i < unsignedK - 1; ++i)
210 subEntity[lastIndex++] = 5;
211 subEntity[lastIndex++] = 9;
215 subEntity[lastIndex++] = 6;
216 for (
unsigned i = 0; i < unsignedK - 1; ++i)
217 subEntity[lastIndex++] = 11;
218 subEntity[lastIndex++] = 7;
220 assert(numIndices==lastIndex);
228 std::vector<unsigned int> codim(li.size());
230 for (std::size_t i=0; i<codim.size(); i++) {
236 std::array<unsigned int,d> mIdx = multiindex(i);
237 for (
int j=0; j<d; j++)
238 if (mIdx[j]==0 or mIdx[j]==k)
247 std::vector<unsigned int> index(
size());
249 for (std::size_t i=0; i<
size(); i++) {
253 std::array<unsigned int,d> mIdx = multiindex(i);
255 for (
int j=d-1; j>=0; j--)
256 if (mIdx[j]>0 and mIdx[j]<unsignedK)
257 index[i] = (k-1)*index[i] + (mIdx[j]-1);
262 std::vector<unsigned int> subEntity(li.size());
266 for (std::size_t i=0; i<
size(); i++)
282 DUNE_THROW(Dune::NotImplemented,
"QkLocalCoefficients for k==" << k <<
" and d==" << d);
284 for (
size_t i=0; i<li.size(); i++)
285 li[i] =
LocalKey(subEntity[i], codim[i], index[i]);
291 return StaticPower<k+1,d>::power;
301 std::vector<LocalKey> li;
QkLocalCoefficients()
Default constructor.
Definition: qklocalcoefficients.hh:225
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
std::size_t size() const
number of coefficients
Definition: qklocalcoefficients.hh:289
Attaches a shape function to an entity.
Definition: qklocalcoefficients.hh:24
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qklocalcoefficients.hh:295
Describe position of one degree of freedom.
Definition: localkey.hh:20