32 , raviart_thomas_space(k - 1)
37 for (
unsigned int i = 0; i <
monomials.size(); ++i)
55 Assert(grads.size() == this->n() || grads.size() == 0,
57 Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
59 Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
61 Assert(fourth_derivatives.size() == this->n() ||
62 fourth_derivatives.size() == 0,
66 (void)third_derivatives;
68 (void)fourth_derivatives;
71 const unsigned int n_sub = raviart_thomas_space.n();
72 const unsigned int my_degree = this->degree();
78 static std::mutex mutex;
79 std::lock_guard<std::mutex> lock(mutex);
81 static std::vector<Tensor<1, dim>> p_values;
82 static std::vector<Tensor<2, dim>> p_grads;
83 static std::vector<Tensor<3, dim>> p_grad_grads;
84 static std::vector<Tensor<4, dim>> p_third_derivatives;
85 static std::vector<Tensor<5, dim>> p_fourth_derivatives;
87 p_values.resize((
values.size() == 0) ? 0 : n_sub);
88 p_grads.resize((grads.size() == 0) ? 0 : n_sub);
89 p_grad_grads.resize((grad_grads.size() == 0) ? 0 : n_sub);
92 raviart_thomas_space.evaluate(unit_point,
97 p_fourth_derivatives);
98 for (
unsigned int i = 0; i < p_values.size(); ++i)
100 for (
unsigned int i = 0; i < p_grads.size(); ++i)
101 grads[i] = p_grads[i];
102 for (
unsigned int i = 0; i < p_grad_grads.size(); ++i)
103 grad_grads[i] = p_grad_grads[i];
108 const unsigned int n_derivatives = 3;
109 double monoval_plus[dim][n_derivatives + 1];
110 double monoval_i[dim][n_derivatives + 1];
119 unsigned int start = n_sub;
130 for (
unsigned int d = 0;
d < dim; ++
d)
131 monomials[my_degree + 1].value(unit_point(
d),
135 for (
unsigned int i = 0; i <= my_degree; ++i, ++start)
137 for (
unsigned int d = 0;
d < dim; ++
d)
138 monomials[i].value(unit_point(
d), n_derivatives, monoval_i[
d]);
142 values[start][0] = monoval_i[0][0] * monoval_plus[1][1];
143 values[start][1] = -monoval_i[0][1] * monoval_plus[1][0];
145 values[start + my_degree + 1][0] =
146 -monoval_plus[0][0] * monoval_i[1][1];
147 values[start + my_degree + 1][1] =
148 monoval_plus[0][1] * monoval_i[1][0];
151 if (grads.size() != 0)
153 grads[start][0][0] = monoval_i[0][1] * monoval_plus[1][1];
154 grads[start][0][1] = monoval_i[0][0] * monoval_plus[1][2];
155 grads[start][1][0] = -monoval_i[0][2] * monoval_plus[1][0];
156 grads[start][1][1] = -monoval_i[0][1] * monoval_plus[1][1];
158 grads[start + my_degree + 1][0][0] =
159 -monoval_plus[0][1] * monoval_i[1][1];
160 grads[start + my_degree + 1][0][1] =
161 -monoval_plus[0][0] * monoval_i[1][2];
162 grads[start + my_degree + 1][1][0] =
163 monoval_plus[0][2] * monoval_i[1][0];
164 grads[start + my_degree + 1][1][1] =
165 monoval_plus[0][1] * monoval_i[1][1];
168 if (grad_grads.size() != 0)
170 grad_grads[start][0][0][0] = monoval_i[0][2] * monoval_plus[1][1];
171 grad_grads[start][0][0][1] = monoval_i[0][1] * monoval_plus[1][2];
172 grad_grads[start][0][1][0] = monoval_i[0][1] * monoval_plus[1][2];
173 grad_grads[start][0][1][1] = monoval_i[0][0] * monoval_plus[1][3];
174 grad_grads[start][1][0][0] =
175 -monoval_i[0][3] * monoval_plus[1][0];
176 grad_grads[start][1][0][1] =
177 -monoval_i[0][2] * monoval_plus[1][1];
178 grad_grads[start][1][1][0] =
179 -monoval_i[0][2] * monoval_plus[1][1];
180 grad_grads[start][1][1][1] =
181 -monoval_i[0][1] * monoval_plus[1][2];
183 grad_grads[start + my_degree + 1][0][0][0] =
184 -monoval_plus[0][2] * monoval_i[1][1];
185 grad_grads[start + my_degree + 1][0][0][1] =
186 -monoval_plus[0][1] * monoval_i[1][2];
187 grad_grads[start + my_degree + 1][0][1][0] =
188 -monoval_plus[0][1] * monoval_i[1][2];
189 grad_grads[start + my_degree + 1][0][1][1] =
190 -monoval_plus[0][0] * monoval_i[1][3];
191 grad_grads[start + my_degree + 1][1][0][0] =
192 monoval_plus[0][3] * monoval_i[1][0];
193 grad_grads[start + my_degree + 1][1][0][1] =
194 monoval_plus[0][2] * monoval_i[1][1];
195 grad_grads[start + my_degree + 1][1][1][0] =
196 monoval_plus[0][2] * monoval_i[1][1];
197 grad_grads[start + my_degree + 1][1][1][1] =
198 monoval_plus[0][1] * monoval_i[1][2];
205 double monoval[dim][n_derivatives + 1];
206 double monoval_j[dim][n_derivatives + 1];
207 double monoval_jplus[dim][n_derivatives + 1];
219 for (
unsigned int d = 0;
d < dim; ++
d)
221 monomials[my_degree + 1].value(unit_point(
d),
224 monomials[my_degree].value(unit_point(
d), n_derivatives, monoval[
d]);
227 const unsigned int n_curls = (my_degree + 1) * (2 * my_degree + 1);
229 for (
unsigned int i = 0; i <= my_degree; ++i)
231 for (
unsigned int d = 0;
d < dim; ++
d)
232 monomials[i].value(unit_point(
d), n_derivatives, monoval_i[
d]);
234 for (
unsigned int j = 0; j <= my_degree; ++j)
236 for (
unsigned int d = 0;
d < dim; ++
d)
238 monomials[j].value(unit_point(
d),
241 monomials[j + 1].value(unit_point(
d),
248 values[start][0] = monoval_i[0][0] * monoval_j[1][0] *
250 static_cast<double>(j + my_degree + 2);
252 -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][0];
254 -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][0];
256 values[start + n_curls][0] =
257 -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][0];
258 values[start + n_curls][1] =
259 monoval_j[0][0] * monoval_i[1][0] * monoval[2][0] *
260 static_cast<double>(j + my_degree + 2);
261 values[start + n_curls][2] =
262 -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][0];
264 values[start + 2 * n_curls][0] =
265 -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][1];
266 values[start + 2 * n_curls][1] =
267 -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][1];
268 values[start + 2 * n_curls][2] =
269 monoval_j[0][0] * monoval[1][0] * monoval_i[2][0] *
270 static_cast<double>(j + my_degree + 2);
277 monoval_i[0][0] * monoval[1][0] * monoval_j[2][0] *
278 static_cast<double>(j + my_degree + 2);
280 -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][0];
282 -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][0];
284 values[start + n_curls + 1][0] =
285 -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][0];
286 values[start + n_curls + 1][1] =
287 monoval[0][0] * monoval_i[1][0] * monoval_j[2][0] *
288 static_cast<double>(j + my_degree + 2);
289 values[start + n_curls + 1][2] =
290 -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][0];
292 values[start + 2 * n_curls + 1][0] =
293 -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][1];
294 values[start + 2 * n_curls + 1][1] =
295 -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][1];
296 values[start + 2 * n_curls + 1][2] =
297 monoval[0][0] * monoval_j[1][0] * monoval_i[2][0] *
298 static_cast<double>(j + my_degree + 2);
302 if (grads.size() != 0)
304 grads[start][0][0] = monoval_i[0][1] * monoval_j[1][0] *
306 static_cast<double>(j + my_degree + 2);
307 grads[start][0][1] = monoval_i[0][0] * monoval_j[1][1] *
309 static_cast<double>(j + my_degree + 2);
310 grads[start][0][2] = monoval_i[0][0] * monoval_j[1][0] *
312 static_cast<double>(j + my_degree + 2);
314 -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][0];
316 -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][0];
318 -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][1];
320 -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][0];
322 -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][0];
324 -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][1];
326 grads[start + n_curls][0][0] =
327 -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][0];
328 grads[start + n_curls][0][1] =
329 -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][0];
330 grads[start + n_curls][0][2] =
331 -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][1];
332 grads[start + n_curls][1][0] =
333 monoval_j[0][1] * monoval_i[1][0] * monoval[2][0] *
334 static_cast<double>(j + my_degree + 2);
335 grads[start + n_curls][1][1] =
336 monoval_j[0][0] * monoval_i[1][1] * monoval[2][0] *
337 static_cast<double>(j + my_degree + 2);
338 grads[start + n_curls][1][2] =
339 monoval_j[0][0] * monoval_i[1][0] * monoval[2][1] *
340 static_cast<double>(j + my_degree + 2);
341 grads[start + n_curls][2][0] =
342 -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][0];
343 grads[start + n_curls][2][1] =
344 -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][0];
345 grads[start + n_curls][2][2] =
346 -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][1];
348 grads[start + 2 * n_curls][0][0] =
349 -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][1];
350 grads[start + 2 * n_curls][0][1] =
351 -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][1];
352 grads[start + 2 * n_curls][0][2] =
353 -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][2];
354 grads[start + 2 * n_curls][1][0] =
355 -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][1];
356 grads[start + 2 * n_curls][1][1] =
357 -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][1];
358 grads[start + 2 * n_curls][1][2] =
359 -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][2];
360 grads[start + 2 * n_curls][2][0] =
361 monoval_j[0][1] * monoval[1][0] * monoval_i[2][0] *
362 static_cast<double>(j + my_degree + 2);
363 grads[start + 2 * n_curls][2][1] =
364 monoval_j[0][0] * monoval[1][1] * monoval_i[2][0] *
365 static_cast<double>(j + my_degree + 2);
366 grads[start + 2 * n_curls][2][2] =
367 monoval_j[0][0] * monoval[1][0] * monoval_i[2][1] *
368 static_cast<double>(j + my_degree + 2);
372 grads[start + 1][0][0] =
373 monoval_i[0][1] * monoval[1][0] * monoval_j[2][0] *
374 static_cast<double>(j + my_degree + 2);
375 grads[start + 1][0][1] =
376 monoval_i[0][0] * monoval[1][1] * monoval_j[2][0] *
377 static_cast<double>(j + my_degree + 2);
378 grads[start + 1][0][2] =
379 monoval_i[0][0] * monoval[1][0] * monoval_j[2][1] *
380 static_cast<double>(j + my_degree + 2);
381 grads[start + 1][1][0] =
382 -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][0];
383 grads[start + 1][1][1] =
384 -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][0];
385 grads[start + 1][1][2] =
386 -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][1];
387 grads[start + 1][2][0] =
388 -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][0];
389 grads[start + 1][2][1] =
390 -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][0];
391 grads[start + 1][2][2] =
392 -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][1];
394 grads[start + n_curls + 1][0][0] =
395 -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][0];
396 grads[start + n_curls + 1][0][1] =
397 -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][0];
398 grads[start + n_curls + 1][0][2] =
399 -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][1];
400 grads[start + n_curls + 1][1][0] =
401 monoval[0][1] * monoval_i[1][0] * monoval_j[2][0] *
402 static_cast<double>(j + my_degree + 2);
403 grads[start + n_curls + 1][1][1] =
404 monoval[0][0] * monoval_i[1][1] * monoval_j[2][0] *
405 static_cast<double>(j + my_degree + 2);
406 grads[start + n_curls + 1][1][2] =
407 monoval[0][0] * monoval_i[1][0] * monoval_j[2][1] *
408 static_cast<double>(j + my_degree + 2);
409 grads[start + n_curls + 1][2][0] =
410 -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][0];
411 grads[start + n_curls + 1][2][1] =
412 -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][0];
413 grads[start + n_curls + 1][2][2] =
414 -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][1];
416 grads[start + 2 * n_curls + 1][0][0] =
417 -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][1];
418 grads[start + 2 * n_curls + 1][0][1] =
419 -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][1];
420 grads[start + 2 * n_curls + 1][0][2] =
421 -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][2];
422 grads[start + 2 * n_curls + 1][1][0] =
423 -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][1];
424 grads[start + 2 * n_curls + 1][1][1] =
425 -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][1];
426 grads[start + 2 * n_curls + 1][1][2] =
427 -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][2];
428 grads[start + 2 * n_curls + 1][2][0] =
429 monoval[0][1] * monoval_j[1][0] * monoval_i[2][0] *
430 static_cast<double>(j + my_degree + 2);
431 grads[start + 2 * n_curls + 1][2][1] =
432 monoval[0][0] * monoval_j[1][1] * monoval_i[2][0] *
433 static_cast<double>(j + my_degree + 2);
434 grads[start + 2 * n_curls + 1][2][2] =
435 monoval[0][0] * monoval_j[1][0] * monoval_i[2][1] *
436 static_cast<double>(j + my_degree + 2);
440 if (grad_grads.size() != 0)
442 grad_grads[start][0][0][0] =
443 monoval_i[0][2] * monoval_j[1][0] * monoval[2][0] *
444 static_cast<double>(j + my_degree + 2);
445 grad_grads[start][0][0][1] =
446 monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
447 static_cast<double>(j + my_degree + 2);
448 grad_grads[start][0][0][2] =
449 monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
450 static_cast<double>(j + my_degree + 2);
451 grad_grads[start][0][1][0] =
452 monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
453 static_cast<double>(j + my_degree + 2);
454 grad_grads[start][0][1][1] =
455 monoval_i[0][0] * monoval_j[1][2] * monoval[2][0] *
456 static_cast<double>(j + my_degree + 2);
457 grad_grads[start][0][1][2] =
458 monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
459 static_cast<double>(j + my_degree + 2);
460 grad_grads[start][0][2][0] =
461 monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
462 static_cast<double>(j + my_degree + 2);
463 grad_grads[start][0][2][1] =
464 monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
465 static_cast<double>(j + my_degree + 2);
466 grad_grads[start][0][2][2] =
467 monoval_i[0][0] * monoval_j[1][0] * monoval[2][2] *
468 static_cast<double>(j + my_degree + 2);
469 grad_grads[start][1][0][0] =
470 -monoval_i[0][3] * monoval_jplus[1][0] * monoval[2][0];
471 grad_grads[start][1][0][1] =
472 -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
473 grad_grads[start][1][0][2] =
474 -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
475 grad_grads[start][1][1][0] =
476 -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
477 grad_grads[start][1][1][1] =
478 -monoval_i[0][1] * monoval_jplus[1][2] * monoval[2][0];
479 grad_grads[start][1][1][2] =
480 -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
481 grad_grads[start][1][2][0] =
482 -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
483 grad_grads[start][1][2][1] =
484 -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
485 grad_grads[start][1][2][2] =
486 -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][2];
487 grad_grads[start][2][0][0] =
488 -monoval_i[0][3] * monoval_j[1][0] * monoval_plus[2][0];
489 grad_grads[start][2][0][1] =
490 -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
491 grad_grads[start][2][0][2] =
492 -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
493 grad_grads[start][2][1][0] =
494 -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
495 grad_grads[start][2][1][1] =
496 -monoval_i[0][1] * monoval_j[1][2] * monoval_plus[2][0];
497 grad_grads[start][2][1][2] =
498 -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
499 grad_grads[start][2][2][0] =
500 -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
501 grad_grads[start][2][2][1] =
502 -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
503 grad_grads[start][2][2][2] =
504 -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][2];
506 grad_grads[start + n_curls][0][0][0] =
507 -monoval_jplus[0][2] * monoval_i[1][1] * monoval[2][0];
508 grad_grads[start + n_curls][0][0][1] =
509 -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
510 grad_grads[start + n_curls][0][0][2] =
511 -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
512 grad_grads[start + n_curls][0][1][0] =
513 -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
514 grad_grads[start + n_curls][0][1][1] =
515 -monoval_jplus[0][0] * monoval_i[1][3] * monoval[2][0];
516 grad_grads[start + n_curls][0][1][2] =
517 -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
518 grad_grads[start + n_curls][0][2][0] =
519 -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
520 grad_grads[start + n_curls][0][2][1] =
521 -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
522 grad_grads[start + n_curls][0][2][2] =
523 -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][2];
524 grad_grads[start + n_curls][1][0][0] =
525 monoval_j[0][2] * monoval_i[1][0] * monoval[2][0] *
526 static_cast<double>(j + my_degree + 2);
527 grad_grads[start + n_curls][1][0][1] =
528 monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
529 static_cast<double>(j + my_degree + 2);
530 grad_grads[start + n_curls][1][0][2] =
531 monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
532 static_cast<double>(j + my_degree + 2);
533 grad_grads[start + n_curls][1][1][0] =
534 monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
535 static_cast<double>(j + my_degree + 2);
536 grad_grads[start + n_curls][1][1][1] =
537 monoval_j[0][0] * monoval_i[1][2] * monoval[2][0] *
538 static_cast<double>(j + my_degree + 2);
539 grad_grads[start + n_curls][1][1][2] =
540 monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
541 static_cast<double>(j + my_degree + 2);
542 grad_grads[start + n_curls][1][2][0] =
543 monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
544 static_cast<double>(j + my_degree + 2);
545 grad_grads[start + n_curls][1][2][1] =
546 monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
547 static_cast<double>(j + my_degree + 2);
548 grad_grads[start + n_curls][1][2][2] =
549 monoval_j[0][0] * monoval_i[1][0] * monoval[2][2] *
550 static_cast<double>(j + my_degree + 2);
551 grad_grads[start + n_curls][2][0][0] =
552 -monoval_j[0][2] * monoval_i[1][1] * monoval_plus[2][0];
553 grad_grads[start + n_curls][2][0][1] =
554 -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
555 grad_grads[start + n_curls][2][0][2] =
556 -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
557 grad_grads[start + n_curls][2][1][0] =
558 -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
559 grad_grads[start + n_curls][2][1][1] =
560 -monoval_j[0][0] * monoval_i[1][3] * monoval_plus[2][0];
561 grad_grads[start + n_curls][2][1][2] =
562 -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
563 grad_grads[start + n_curls][2][2][0] =
564 -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
565 grad_grads[start + n_curls][2][2][1] =
566 -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
567 grad_grads[start + n_curls][2][2][2] =
568 -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][2];
570 grad_grads[start + 2 * n_curls][0][0][0] =
571 -monoval_jplus[0][2] * monoval[1][0] * monoval_i[2][1];
572 grad_grads[start + 2 * n_curls][0][0][1] =
573 -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
574 grad_grads[start + 2 * n_curls][0][0][2] =
575 -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
576 grad_grads[start + 2 * n_curls][0][1][0] =
577 -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
578 grad_grads[start + 2 * n_curls][0][1][1] =
579 -monoval_jplus[0][0] * monoval[1][2] * monoval_i[2][1];
580 grad_grads[start + 2 * n_curls][0][1][2] =
581 -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
582 grad_grads[start + 2 * n_curls][0][2][0] =
583 -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
584 grad_grads[start + 2 * n_curls][0][2][1] =
585 -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
586 grad_grads[start + 2 * n_curls][0][2][2] =
587 -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][3];
588 grad_grads[start + 2 * n_curls][1][0][0] =
589 -monoval_j[0][2] * monoval_plus[1][0] * monoval_i[2][1];
590 grad_grads[start + 2 * n_curls][1][0][1] =
591 -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
592 grad_grads[start + 2 * n_curls][1][0][2] =
593 -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
594 grad_grads[start + 2 * n_curls][1][1][0] =
595 -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
596 grad_grads[start + 2 * n_curls][1][1][1] =
597 -monoval_j[0][0] * monoval_plus[1][2] * monoval_i[2][1];
598 grad_grads[start + 2 * n_curls][1][1][2] =
599 -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
600 grad_grads[start + 2 * n_curls][1][2][0] =
601 -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
602 grad_grads[start + 2 * n_curls][1][2][1] =
603 -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
604 grad_grads[start + 2 * n_curls][1][2][2] =
605 -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][3];
606 grad_grads[start + 2 * n_curls][2][0][0] =
607 monoval_j[0][2] * monoval[1][0] * monoval_i[2][0] *
608 static_cast<double>(j + my_degree + 2);
609 grad_grads[start + 2 * n_curls][2][0][1] =
610 monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
611 static_cast<double>(j + my_degree + 2);
612 grad_grads[start + 2 * n_curls][2][0][2] =
613 monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
614 static_cast<double>(j + my_degree + 2);
615 grad_grads[start + 2 * n_curls][2][1][0] =
616 monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
617 static_cast<double>(j + my_degree + 2);
618 grad_grads[start + 2 * n_curls][2][1][1] =
619 monoval_j[0][0] * monoval[1][2] * monoval_i[2][0] *
620 static_cast<double>(j + my_degree + 2);
621 grad_grads[start + 2 * n_curls][2][1][2] =
622 monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
623 static_cast<double>(j + my_degree + 2);
624 grad_grads[start + 2 * n_curls][2][2][0] =
625 monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
626 static_cast<double>(j + my_degree + 2);
627 grad_grads[start + 2 * n_curls][2][2][1] =
628 monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
629 static_cast<double>(j + my_degree + 2);
630 grad_grads[start + 2 * n_curls][2][2][2] =
631 monoval_j[0][0] * monoval[1][0] * monoval_i[2][2] *
632 static_cast<double>(j + my_degree + 2);
636 grad_grads[start + 1][0][0][0] =
637 monoval_i[0][2] * monoval[1][0] * monoval_j[2][0] *
638 static_cast<double>(j + my_degree + 2);
639 grad_grads[start + 1][0][0][1] =
640 monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
641 static_cast<double>(j + my_degree + 2);
642 grad_grads[start + 1][0][0][2] =
643 monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
644 static_cast<double>(j + my_degree + 2);
645 grad_grads[start + 1][0][1][0] =
646 monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
647 static_cast<double>(j + my_degree + 2);
648 grad_grads[start + 1][0][1][1] =
649 monoval_i[0][0] * monoval[1][2] * monoval_j[2][0] *
650 static_cast<double>(j + my_degree + 2);
651 grad_grads[start + 1][0][1][2] =
652 monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
653 static_cast<double>(j + my_degree + 2);
654 grad_grads[start + 1][0][2][0] =
655 monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
656 static_cast<double>(j + my_degree + 2);
657 grad_grads[start + 1][0][2][1] =
658 monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
659 static_cast<double>(j + my_degree + 2);
660 grad_grads[start + 1][0][2][2] =
661 monoval_i[0][0] * monoval[1][0] * monoval_j[2][2] *
662 static_cast<double>(j + my_degree + 2);
663 grad_grads[start + 1][1][0][0] =
664 -monoval_i[0][3] * monoval_plus[1][0] * monoval_j[2][0];
665 grad_grads[start + 1][1][0][1] =
666 -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
667 grad_grads[start + 1][1][0][2] =
668 -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
669 grad_grads[start + 1][1][1][0] =
670 -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
671 grad_grads[start + 1][1][1][1] =
672 -monoval_i[0][1] * monoval_plus[1][2] * monoval_j[2][0];
673 grad_grads[start + 1][1][1][2] =
674 -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
675 grad_grads[start + 1][1][2][0] =
676 -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
677 grad_grads[start + 1][1][2][1] =
678 -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
679 grad_grads[start + 1][1][2][2] =
680 -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][2];
681 grad_grads[start + 1][2][0][0] =
682 -monoval_i[0][3] * monoval[1][0] * monoval_jplus[2][0];
683 grad_grads[start + 1][2][0][1] =
684 -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
685 grad_grads[start + 1][2][0][2] =
686 -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
687 grad_grads[start + 1][2][1][0] =
688 -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
689 grad_grads[start + 1][2][1][1] =
690 -monoval_i[0][1] * monoval[1][2] * monoval_jplus[2][0];
691 grad_grads[start + 1][2][1][2] =
692 -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
693 grad_grads[start + 1][2][2][0] =
694 -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
695 grad_grads[start + 1][2][2][1] =
696 -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
697 grad_grads[start + 1][2][2][2] =
698 -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][2];
700 grad_grads[start + n_curls + 1][0][0][0] =
701 -monoval_plus[0][2] * monoval_i[1][1] * monoval_j[2][0];
702 grad_grads[start + n_curls + 1][0][0][1] =
703 -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
704 grad_grads[start + n_curls + 1][0][0][2] =
705 -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
706 grad_grads[start + n_curls + 1][0][1][0] =
707 -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
708 grad_grads[start + n_curls + 1][0][1][1] =
709 -monoval_plus[0][0] * monoval_i[1][3] * monoval_j[2][0];
710 grad_grads[start + n_curls + 1][0][1][2] =
711 -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
712 grad_grads[start + n_curls + 1][0][2][0] =
713 -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
714 grad_grads[start + n_curls + 1][0][2][1] =
715 -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
716 grad_grads[start + n_curls + 1][0][2][2] =
717 -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][2];
718 grad_grads[start + n_curls + 1][1][0][0] =
719 monoval[0][2] * monoval_i[1][0] * monoval_j[2][0] *
720 static_cast<double>(j + my_degree + 2);
721 grad_grads[start + n_curls + 1][1][0][1] =
722 monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
723 static_cast<double>(j + my_degree + 2);
724 grad_grads[start + n_curls + 1][1][0][2] =
725 monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
726 static_cast<double>(j + my_degree + 2);
727 grad_grads[start + n_curls + 1][1][1][0] =
728 monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
729 static_cast<double>(j + my_degree + 2);
730 grad_grads[start + n_curls + 1][1][1][1] =
731 monoval[0][0] * monoval_i[1][2] * monoval_j[2][0] *
732 static_cast<double>(j + my_degree + 2);
733 grad_grads[start + n_curls + 1][1][1][2] =
734 monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
735 static_cast<double>(j + my_degree + 2);
736 grad_grads[start + n_curls + 1][1][2][0] =
737 monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
738 static_cast<double>(j + my_degree + 2);
739 grad_grads[start + n_curls + 1][1][2][1] =
740 monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
741 static_cast<double>(j + my_degree + 2);
742 grad_grads[start + n_curls + 1][1][2][2] =
743 monoval[0][0] * monoval_i[1][0] * monoval_j[2][2] *
744 static_cast<double>(j + my_degree + 2);
745 grad_grads[start + n_curls + 1][2][0][0] =
746 -monoval[0][2] * monoval_i[1][1] * monoval_jplus[2][0];
747 grad_grads[start + n_curls + 1][2][0][1] =
748 -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
749 grad_grads[start + n_curls + 1][2][0][2] =
750 -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
751 grad_grads[start + n_curls + 1][2][1][0] =
752 -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
753 grad_grads[start + n_curls + 1][2][1][1] =
754 -monoval[0][0] * monoval_i[1][3] * monoval_jplus[2][0];
755 grad_grads[start + n_curls + 1][2][1][2] =
756 -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
757 grad_grads[start + n_curls + 1][2][2][0] =
758 -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
759 grad_grads[start + n_curls + 1][2][2][1] =
760 -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
761 grad_grads[start + n_curls + 1][2][2][2] =
762 -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][2];
764 grad_grads[start + 2 * n_curls + 1][0][0][0] =
765 -monoval_plus[0][2] * monoval_j[1][0] * monoval_i[2][1];
766 grad_grads[start + 2 * n_curls + 1][0][0][1] =
767 -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
768 grad_grads[start + 2 * n_curls + 1][0][0][2] =
769 -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
770 grad_grads[start + 2 * n_curls + 1][0][1][0] =
771 -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
772 grad_grads[start + 2 * n_curls + 1][0][1][1] =
773 -monoval_plus[0][0] * monoval_j[1][2] * monoval_i[2][1];
774 grad_grads[start + 2 * n_curls + 1][0][1][2] =
775 -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
776 grad_grads[start + 2 * n_curls + 1][0][2][0] =
777 -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
778 grad_grads[start + 2 * n_curls + 1][0][2][1] =
779 -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
780 grad_grads[start + 2 * n_curls + 1][0][2][2] =
781 -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][3];
782 grad_grads[start + 2 * n_curls + 1][1][0][0] =
783 -monoval[0][2] * monoval_jplus[1][0] * monoval_i[2][1];
784 grad_grads[start + 2 * n_curls + 1][1][0][1] =
785 -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
786 grad_grads[start + 2 * n_curls + 1][1][0][2] =
787 -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
788 grad_grads[start + 2 * n_curls + 1][1][1][0] =
789 -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
790 grad_grads[start + 2 * n_curls + 1][1][1][1] =
791 -monoval[0][0] * monoval_jplus[1][2] * monoval_i[2][1];
792 grad_grads[start + 2 * n_curls + 1][1][1][2] =
793 -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
794 grad_grads[start + 2 * n_curls + 1][1][2][0] =
795 -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
796 grad_grads[start + 2 * n_curls + 1][1][2][1] =
797 -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
798 grad_grads[start + 2 * n_curls + 1][1][2][2] =
799 -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][3];
800 grad_grads[start + 2 * n_curls + 1][2][0][0] =
801 monoval[0][2] * monoval_j[1][0] * monoval_i[2][0] *
802 static_cast<double>(j + my_degree + 2);
803 grad_grads[start + 2 * n_curls + 1][2][0][1] =
804 monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
805 static_cast<double>(j + my_degree + 2);
806 grad_grads[start + 2 * n_curls + 1][2][0][2] =
807 monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
808 static_cast<double>(j + my_degree + 2);
809 grad_grads[start + 2 * n_curls + 1][2][1][0] =
810 monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
811 static_cast<double>(j + my_degree + 2);
812 grad_grads[start + 2 * n_curls + 1][2][1][1] =
813 monoval[0][0] * monoval_j[1][2] * monoval_i[2][0] *
814 static_cast<double>(j + my_degree + 2);
815 grad_grads[start + 2 * n_curls + 1][2][1][2] =
816 monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
817 static_cast<double>(j + my_degree + 2);
818 grad_grads[start + 2 * n_curls + 1][2][2][0] =
819 monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
820 static_cast<double>(j + my_degree + 2);
821 grad_grads[start + 2 * n_curls + 1][2][2][1] =
822 monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
823 static_cast<double>(j + my_degree + 2);
824 grad_grads[start + 2 * n_curls + 1][2][2][2] =
825 monoval[0][0] * monoval_j[1][0] * monoval_i[2][2] *
826 static_cast<double>(j + my_degree + 2);
846 if (dim == 1 || dim == 2 || dim == 3)
847 return dim * Utilities::fixed_power<dim>(k + 1);
855 std::unique_ptr<TensorPolynomialsBase<dim>>
858 return std::make_unique<PolynomialsRT_Bubbles<dim>>(*this);
virtual std::unique_ptr< TensorPolynomialsBase< dim > > clone() const override
PolynomialsRT_Bubbles(const unsigned int k)
void evaluate(const Point< dim > &unit_point, std::vector< Tensor< 1, dim >> &values, std::vector< Tensor< 2, dim >> &grads, std::vector< Tensor< 3, dim >> &grad_grads, std::vector< Tensor< 4, dim >> &third_derivatives, std::vector< Tensor< 5, dim >> &fourth_derivatives) const override
static unsigned int n_polynomials(const unsigned int degree)
std::vector< Polynomials::Polynomial< double > > monomials
#define DEAL_II_CONSTEXPR_IN_CONDITIONAL
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)