Reference documentation for deal.II version GIT 35969cdc9b 2023-12-09 01:10:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
polynomials_nedelec.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2013 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
20 
21 #include <iomanip>
22 #include <iostream>
23 #include <memory>
24 
26 
27 
28 template <int dim>
30  : TensorPolynomialsBase<dim>(k, n_polynomials(k))
31  , polynomial_space(create_polynomials(k))
32 {}
33 
34 template <int dim>
35 std::vector<std::vector<Polynomials::Polynomial<double>>>
37 {
38  std::vector<std::vector<Polynomials::Polynomial<double>>> pols(dim);
39 
41 
42  for (unsigned int i = 1; i < dim; ++i)
44 
45  return pols;
46 }
47 
48 
49 // Compute the values, gradients
50 // and double gradients of the
51 // polynomial at the given point.
52 template <int dim>
53 void
55  const Point<dim> &unit_point,
56  std::vector<Tensor<1, dim>> &values,
57  std::vector<Tensor<2, dim>> &grads,
58  std::vector<Tensor<3, dim>> &grad_grads,
59  std::vector<Tensor<4, dim>> &third_derivatives,
60  std::vector<Tensor<5, dim>> &fourth_derivatives) const
61 {
62  Assert(values.size() == this->n() || values.empty(),
63  ExcDimensionMismatch(values.size(), this->n()));
64  Assert(grads.size() == this->n() || grads.empty(),
65  ExcDimensionMismatch(grads.size(), this->n()));
66  Assert(grad_grads.size() == this->n() || grad_grads.empty(),
67  ExcDimensionMismatch(grad_grads.size(), this->n()));
68  Assert(third_derivatives.size() == this->n() || third_derivatives.empty(),
69  ExcDimensionMismatch(third_derivatives.size(), this->n()));
70  Assert(fourth_derivatives.size() == this->n() || fourth_derivatives.empty(),
71  ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
72 
73  // third and fourth derivatives not implemented
74  (void)third_derivatives;
75  Assert(third_derivatives.empty(), ExcNotImplemented());
76  (void)fourth_derivatives;
77  Assert(fourth_derivatives.empty(), ExcNotImplemented());
78 
79  // Declare the values, derivatives
80  // and second derivatives vectors of
81  // <tt>polynomial_space</tt> at
82  // <tt>unit_point</tt>
83  const unsigned int n_basis = polynomial_space.n();
84  const unsigned int my_degree = this->degree();
85  std::vector<double> unit_point_values((values.empty()) ? 0 : n_basis);
86  std::vector<Tensor<1, dim>> unit_point_grads((grads.empty()) ? 0 : n_basis);
87  std::vector<Tensor<2, dim>> unit_point_grad_grads(
88  (grad_grads.empty()) ? 0 : n_basis);
89  std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
90  std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
91 
92  switch (dim)
93  {
94  case 1:
95  {
96  polynomial_space.evaluate(unit_point,
97  unit_point_values,
98  unit_point_grads,
99  unit_point_grad_grads,
100  empty_vector_of_3rd_order_tensors,
101  empty_vector_of_4th_order_tensors);
102 
103  // Assign the correct values to the
104  // corresponding shape functions.
105  if (values.size() > 0)
106  for (unsigned int i = 0; i < unit_point_values.size(); ++i)
107  values[i][0] = unit_point_values[i];
108 
109  if (grads.size() > 0)
110  for (unsigned int i = 0; i < unit_point_grads.size(); ++i)
111  grads[i][0][0] = unit_point_grads[i][0];
112 
113  if (grad_grads.size() > 0)
114  for (unsigned int i = 0; i < unit_point_grad_grads.size(); ++i)
115  grad_grads[i][0][0][0] = unit_point_grad_grads[i][0][0];
116 
117  break;
118  }
119 
120  case 2:
121  {
122  polynomial_space.evaluate(unit_point,
123  unit_point_values,
124  unit_point_grads,
125  unit_point_grad_grads,
126  empty_vector_of_3rd_order_tensors,
127  empty_vector_of_4th_order_tensors);
128 
129  // Declare the values, derivatives and
130  // second derivatives vectors of
131  // <tt>polynomial_space</tt> at
132  // <tt>unit_point</tt> with coordinates
133  // shifted one step in positive direction
134  Point<dim> p;
135 
136  p(0) = unit_point(1);
137  p(1) = unit_point(0);
138 
139  std::vector<double> p_values((values.empty()) ? 0 : n_basis);
140  std::vector<Tensor<1, dim>> p_grads((grads.empty()) ? 0 : n_basis);
141  std::vector<Tensor<2, dim>> p_grad_grads(
142  (grad_grads.empty()) ? 0 : n_basis);
143 
144  polynomial_space.evaluate(p,
145  p_values,
146  p_grads,
147  p_grad_grads,
148  empty_vector_of_3rd_order_tensors,
149  empty_vector_of_4th_order_tensors);
150 
151  // Assign the correct values to the
152  // corresponding shape functions.
153  if (values.size() > 0)
154  {
155  for (unsigned int i = 0; i <= my_degree; ++i)
156  for (unsigned int j = 0; j < 2; ++j)
157  {
158  values[i + j * (my_degree + 1)][0] = 0.0;
159  values[i + j * (my_degree + 1)][1] =
160  p_values[i + j * (my_degree + 1)];
161  values[i + (j + 2) * (my_degree + 1)][0] =
162  unit_point_values[i + j * (my_degree + 1)];
163  values[i + (j + 2) * (my_degree + 1)][1] = 0.0;
164  }
165 
166  if (my_degree > 0)
167  for (unsigned int i = 0; i <= my_degree; ++i)
168  for (unsigned int j = 0; j < my_degree; ++j)
169  {
171  my_degree +
173  unit_point_values[i + (j + 2) * (my_degree + 1)];
175  my_degree +
176  j + GeometryInfo<dim>::lines_per_cell][1] = 0.0;
177  values[i + (j + my_degree +
179  (my_degree + 1)][0] = 0.0;
180  values[i + (j + my_degree +
182  (my_degree + 1)][1] =
183  p_values[i + (j + 2) * (my_degree + 1)];
184  }
185  }
186 
187  if (grads.size() > 0)
188  {
189  for (unsigned int i = 0; i <= my_degree; ++i)
190  for (unsigned int j = 0; j < 2; ++j)
191  {
192  for (unsigned int k = 0; k < dim; ++k)
193  {
194  grads[i + j * (my_degree + 1)][0][k] = 0.0;
195  grads[i + (j + 2) * (my_degree + 1)][0][k] =
196  unit_point_grads[i + j * (my_degree + 1)][k];
197  grads[i + (j + 2) * (my_degree + 1)][1][k] = 0.0;
198  }
199 
200  grads[i + j * (my_degree + 1)][1][0] =
201  p_grads[i + j * (my_degree + 1)][1];
202  grads[i + j * (my_degree + 1)][1][1] =
203  p_grads[i + j * (my_degree + 1)][0];
204  }
205 
206  if (my_degree > 0)
207  for (unsigned int i = 0; i <= my_degree; ++i)
208  for (unsigned int j = 0; j < my_degree; ++j)
209  {
210  for (unsigned int k = 0; k < dim; ++k)
211  {
213  my_degree +
215  unit_point_grads[i + (j + 2) * (my_degree + 1)][k];
217  my_degree +
219  0.0;
220  grads[i + (j + my_degree +
222  (my_degree + 1)][0][k] = 0.0;
223  }
224 
225  grads[i + (j + my_degree +
227  (my_degree + 1)][1][0] =
228  p_grads[i + (j + 2) * (my_degree + 1)][1];
229  grads[i + (j + my_degree +
231  (my_degree + 1)][1][1] =
232  p_grads[i + (j + 2) * (my_degree + 1)][0];
233  }
234  }
235 
236  if (grad_grads.size() > 0)
237  {
238  for (unsigned int i = 0; i <= my_degree; ++i)
239  for (unsigned int j = 0; j < 2; ++j)
240  {
241  for (unsigned int k = 0; k < dim; ++k)
242  for (unsigned int l = 0; l < dim; ++l)
243  {
244  grad_grads[i + j * (my_degree + 1)][0][k][l] = 0.0;
245  grad_grads[i + (j + 2) * (my_degree + 1)][0][k][l] =
246  unit_point_grad_grads[i + j * (my_degree + 1)][k]
247  [l];
248  grad_grads[i + (j + 2) * (my_degree + 1)][1][k][l] =
249  0.0;
250  }
251 
252  grad_grads[i + j * (my_degree + 1)][1][0][0] =
253  p_grad_grads[i + j * (my_degree + 1)][1][1];
254  grad_grads[i + j * (my_degree + 1)][1][0][1] =
255  p_grad_grads[i + j * (my_degree + 1)][1][0];
256  grad_grads[i + j * (my_degree + 1)][1][1][0] =
257  p_grad_grads[i + j * (my_degree + 1)][0][1];
258  grad_grads[i + j * (my_degree + 1)][1][1][1] =
259  p_grad_grads[i + j * (my_degree + 1)][0][0];
260  }
261 
262  if (my_degree > 0)
263  for (unsigned int i = 0; i <= my_degree; ++i)
264  for (unsigned int j = 0; j < my_degree; ++j)
265  {
266  for (unsigned int k = 0; k < dim; ++k)
267  for (unsigned int l = 0; l < dim; ++l)
268  {
269  grad_grads[(i + GeometryInfo<dim>::lines_per_cell) *
270  my_degree +
272  [k][l] = unit_point_grad_grads
273  [i + (j + 2) * (my_degree + 1)][k][l];
274  grad_grads[(i + GeometryInfo<dim>::lines_per_cell) *
275  my_degree +
277  [k][l] = 0.0;
278  grad_grads[i + (j + my_degree +
280  (my_degree + 1)][0][k][l] = 0.0;
281  }
282 
283  grad_grads[i + (j + my_degree +
285  (my_degree + 1)][1][0][0] =
286  p_grad_grads[i + (j + 2) * (my_degree + 1)][1][1];
287  grad_grads[i + (j + my_degree +
289  (my_degree + 1)][1][0][1] =
290  p_grad_grads[i + (j + 2) * (my_degree + 1)][1][0];
291  grad_grads[i + (j + my_degree +
293  (my_degree + 1)][1][1][0] =
294  p_grad_grads[i + (j + 2) * (my_degree + 1)][0][1];
295  grad_grads[i + (j + my_degree +
297  (my_degree + 1)][1][1][1] =
298  p_grad_grads[i + (j + 2) * (my_degree + 1)][0][0];
299  }
300  }
301 
302  break;
303  }
304 
305  case 3:
306  {
307  polynomial_space.evaluate(unit_point,
308  unit_point_values,
309  unit_point_grads,
310  unit_point_grad_grads,
311  empty_vector_of_3rd_order_tensors,
312  empty_vector_of_4th_order_tensors);
313 
314  // Declare the values, derivatives
315  // and second derivatives vectors of
316  // <tt>polynomial_space</tt> at
317  // <tt>unit_point</tt> with coordinates
318  // shifted two steps in positive
319  // direction
320  Point<dim> p1, p2;
321  std::vector<double> p1_values((values.empty()) ? 0 : n_basis);
322  std::vector<Tensor<1, dim>> p1_grads((grads.empty()) ? 0 : n_basis);
323  std::vector<Tensor<2, dim>> p1_grad_grads(
324  (grad_grads.empty()) ? 0 : n_basis);
325  std::vector<double> p2_values((values.empty()) ? 0 : n_basis);
326  std::vector<Tensor<1, dim>> p2_grads((grads.empty()) ? 0 : n_basis);
327  std::vector<Tensor<2, dim>> p2_grad_grads(
328  (grad_grads.empty()) ? 0 : n_basis);
329 
330  p1(0) = unit_point(1);
331  p1(1) = unit_point(2);
332  p1(2) = unit_point(0);
333  polynomial_space.evaluate(p1,
334  p1_values,
335  p1_grads,
336  p1_grad_grads,
337  empty_vector_of_3rd_order_tensors,
338  empty_vector_of_4th_order_tensors);
339  p2(0) = unit_point(2);
340  p2(1) = unit_point(0);
341  p2(2) = unit_point(1);
342  polynomial_space.evaluate(p2,
343  p2_values,
344  p2_grads,
345  p2_grad_grads,
346  empty_vector_of_3rd_order_tensors,
347  empty_vector_of_4th_order_tensors);
348 
349  // Assign the correct values to the
350  // corresponding shape functions.
351  if (values.size() > 0)
352  {
353  for (unsigned int i = 0; i <= my_degree; ++i)
354  {
355  for (unsigned int j = 0; j < 2; ++j)
356  {
357  for (unsigned int k = 0; k < 2; ++k)
358  {
359  for (unsigned int l = 0; l < 2; ++l)
360  {
361  values[i + (j + 4 * k) * (my_degree + 1)][2 * l] =
362  0.0;
363  values[i + (j + 4 * k + 2) * (my_degree + 1)]
364  [l + 1] = 0.0;
365  values[i + (j + 2 * (k + 4)) * (my_degree + 1)]
366  [l] = 0.0;
367  }
368 
369  values[i + (j + 4 * k + 2) * (my_degree + 1)][0] =
370  unit_point_values[i + (j + k * (my_degree + 2)) *
371  (my_degree + 1)];
372  values[i + (j + 2 * (k + 4)) * (my_degree + 1)][2] =
373  p2_values[i + (j + k * (my_degree + 2)) *
374  (my_degree + 1)];
375  }
376 
377  values[i + j * (my_degree + 1)][1] =
378  p1_values[i + j * (my_degree + 1) * (my_degree + 2)];
379  }
380 
381  values[i + 4 * (my_degree + 1)][1] =
382  p1_values[i + my_degree + 1];
383  values[i + 5 * (my_degree + 1)][1] =
384  p1_values[i + (my_degree + 1) * (my_degree + 3)];
385  }
386 
387  if (my_degree > 0)
388  for (unsigned int i = 0; i <= my_degree; ++i)
389  for (unsigned int j = 0; j < my_degree; ++j)
390  {
391  for (unsigned int k = 0; k < my_degree; ++k)
392  {
393  for (unsigned int l = 0; l < 2; ++l)
394  {
395  values[((i +
397  my_degree +
400  my_degree +
402  [l + 1] = 0.0;
403  values[(i +
404  (j +
406  my_degree) *
407  (my_degree + 1) +
409  my_degree +
411  [2 * l] = 0.0;
412  values[i +
413  (j +
414  (k +
416  my_degree)) *
417  my_degree +
419  (my_degree + 1)][l] = 0.0;
420  }
421 
423  my_degree +
426  my_degree +
428  unit_point_values[i +
429  (j + (k + 2) * (my_degree + 2) +
430  2) *
431  (my_degree + 1)];
432  values[(i +
434  my_degree) *
435  (my_degree + 1) +
437  my_degree +
439  p1_values[i + ((j + 2) * (my_degree + 2) + k + 2) *
440  (my_degree + 1)];
441  values[i +
442  (j +
444  my_degree)) *
445  my_degree +
447  (my_degree + 1)][2] =
448  p2_values[i + (j + (k + 2) * (my_degree + 2) + 2) *
449  (my_degree + 1)];
450  }
451 
452  for (unsigned int k = 0; k < 2; ++k)
453  {
454  for (unsigned int l = 0; l < 2; ++l)
455  {
456  for (unsigned int m = 0; m < 2; ++m)
457  {
458  values[i +
459  (j +
460  (2 * (k + 2 * l) + 1) * my_degree +
462  (my_degree + 1)][m + l] = 0.0;
463  values[(i +
464  2 * (k + 2 * (l + 1)) *
465  (my_degree + 1) +
467  my_degree +
469  [m + l] = 0.0;
470  }
471 
472  values[(i + 2 * k * (my_degree + 1) +
474  my_degree +
476  [2 * l] = 0.0;
477  values[i + (j + (2 * k + 9) * my_degree +
479  (my_degree + 1)][2 * l] = 0.0;
480  }
481 
482  values[(i + 2 * k * (my_degree + 1) +
484  my_degree +
486  p1_values[i + (j + k * (my_degree + 2) + 2) *
487  (my_degree + 1)];
488  values[i + (j + (2 * k + 1) * my_degree +
490  (my_degree + 1)][2] =
491  p2_values[i + ((j + 2) * (my_degree + 2) + k) *
492  (my_degree + 1)];
493  values[(i + 2 * (k + 2) * (my_degree + 1) +
495  my_degree +
497  p2_values[i + (j + k * (my_degree + 2) + 2) *
498  (my_degree + 1)];
499  values[i + (j + (2 * k + 5) * my_degree +
501  (my_degree + 1)][0] =
502  unit_point_values[i +
503  ((j + 2) * (my_degree + 2) + k) *
504  (my_degree + 1)];
505  values[(i + 2 * (k + 4) * (my_degree + 1) +
507  my_degree +
509  unit_point_values[i +
510  (j + k * (my_degree + 2) + 2) *
511  (my_degree + 1)];
512  values[i + (j + (2 * k + 9) * my_degree +
514  (my_degree + 1)][1] =
515  p1_values[i + ((j + 2) * (my_degree + 2) + k) *
516  (my_degree + 1)];
517  }
518  }
519  }
520 
521  if (grads.size() > 0)
522  {
523  for (unsigned int i = 0; i <= my_degree; ++i)
524  {
525  for (unsigned int j = 0; j < 2; ++j)
526  {
527  for (unsigned int k = 0; k < 2; ++k)
528  {
529  for (unsigned int l = 0; l < 2; ++l)
530  for (unsigned int m = 0; m < dim; ++m)
531  {
532  grads[i + (j + 4 * k) * (my_degree + 1)][2 * l]
533  [m] = 0.0;
534  grads[i + (j + 4 * k + 2) * (my_degree + 1)]
535  [l + 1][m] = 0.0;
536  grads[i + (j + 2 * (k + 4)) * (my_degree + 1)]
537  [l][m] = 0.0;
538  }
539 
540  for (unsigned int l = 0; l < dim; ++l)
541  grads[i + (j + 4 * k + 2) * (my_degree + 1)][0][l] =
542  unit_point_grads[i + (j + k * (my_degree + 2)) *
543  (my_degree + 1)][l];
544 
545  grads[i + (j + 2 * (k + 4)) * (my_degree + 1)][2][0] =
546  p2_grads[i + (j + k * (my_degree + 2)) *
547  (my_degree + 1)][1];
548  grads[i + (j + 2 * (k + 4)) * (my_degree + 1)][2][1] =
549  p2_grads[i + (j + k * (my_degree + 2)) *
550  (my_degree + 1)][2];
551  grads[i + (j + 2 * (k + 4)) * (my_degree + 1)][2][2] =
552  p2_grads[i + (j + k * (my_degree + 2)) *
553  (my_degree + 1)][0];
554  }
555 
556  grads[i + j * (my_degree + 1)][1][0] =
557  p1_grads[i + j * (my_degree + 1) * (my_degree + 2)][2];
558  grads[i + j * (my_degree + 1)][1][1] =
559  p1_grads[i + j * (my_degree + 1) * (my_degree + 2)][0];
560  grads[i + j * (my_degree + 1)][1][2] =
561  p1_grads[i + j * (my_degree + 1) * (my_degree + 2)][1];
562  }
563 
564  grads[i + 4 * (my_degree + 1)][1][0] =
565  p1_grads[i + my_degree + 1][2];
566  grads[i + 4 * (my_degree + 1)][1][1] =
567  p1_grads[i + my_degree + 1][0];
568  grads[i + 4 * (my_degree + 1)][1][2] =
569  p1_grads[i + my_degree + 1][1];
570  grads[i + 5 * (my_degree + 1)][1][0] =
571  p1_grads[i + (my_degree + 1) * (my_degree + 3)][2];
572  grads[i + 5 * (my_degree + 1)][1][1] =
573  p1_grads[i + (my_degree + 1) * (my_degree + 3)][0];
574  grads[i + 5 * (my_degree + 1)][1][2] =
575  p1_grads[i + (my_degree + 1) * (my_degree + 3)][1];
576  }
577 
578  if (my_degree > 0)
579  for (unsigned int i = 0; i <= my_degree; ++i)
580  for (unsigned int j = 0; j < my_degree; ++j)
581  {
582  for (unsigned int k = 0; k < my_degree; ++k)
583  {
584  for (unsigned int l = 0; l < dim; ++l)
585  {
586  for (unsigned int m = 0; m < 2; ++m)
587  {
588  grads
589  [((i +
591  my_degree +
594  my_degree +
596  [m + 1][l] = 0.0;
597  grads[(i +
598  (j +
599  2 *
601  my_degree) *
602  (my_degree + 1) +
604  my_degree +
606  [2 * m][l] = 0.0;
607  grads[i +
608  (j +
609  (k +
610  2 *
612  my_degree)) *
613  my_degree +
615  (my_degree + 1)][m][l] = 0.0;
616  }
617 
618  grads[((i +
620  my_degree +
623  my_degree +
625  [l] = unit_point_grads
626  [i + (j + (k + 2) * (my_degree + 2) + 2) *
627  (my_degree + 1)][l];
628  }
629 
630  grads[(i +
632  my_degree) *
633  (my_degree + 1) +
635  my_degree +
637  p1_grads[i + ((j + 2) * (my_degree + 2) + k + 2) *
638  (my_degree + 1)][2];
639  grads[(i +
641  my_degree) *
642  (my_degree + 1) +
644  my_degree +
646  p1_grads[i + ((j + 2) * (my_degree + 2) + k + 2) *
647  (my_degree + 1)][0];
648  grads[(i +
650  my_degree) *
651  (my_degree + 1) +
653  my_degree +
655  p1_grads[i + ((j + 2) * (my_degree + 2) + k + 2) *
656  (my_degree + 1)][1];
657  grads[i +
658  (j +
660  my_degree)) *
661  my_degree +
663  (my_degree + 1)][2][0] =
664  p2_grads[i + (j + (k + 2) * (my_degree + 2) + 2) *
665  (my_degree + 1)][1];
666  grads[i +
667  (j +
669  my_degree)) *
670  my_degree +
672  (my_degree + 1)][2][1] =
673  p2_grads[i + (j + (k + 2) * (my_degree + 2) + 2) *
674  (my_degree + 1)][2];
675  grads[i +
676  (j +
678  my_degree)) *
679  my_degree +
681  (my_degree + 1)][2][2] =
682  p2_grads[i + (j + (k + 2) * (my_degree + 2) + 2) *
683  (my_degree + 1)][0];
684  }
685 
686  for (unsigned int k = 0; k < 2; ++k)
687  {
688  for (unsigned int l = 0; l < 2; ++l)
689  for (unsigned int m = 0; m < dim; ++m)
690  {
691  for (unsigned int n = 0; n < 2; ++n)
692  {
693  grads[i +
694  (j +
695  (2 * (k + 2 * l) + 1) * my_degree +
697  (my_degree + 1)][n + l][m] = 0.0;
698  grads[(i +
699  2 * (k + 2 * (l + 1)) *
700  (my_degree + 1) +
702  my_degree +
704  [n + l][m] = 0.0;
705  }
706 
707  grads[(i + 2 * k * (my_degree + 1) +
709  my_degree +
711  [2 * l][m] = 0.0;
712  grads[i + (j + (2 * k + 9) * my_degree +
714  (my_degree + 1)][2 * l][m] = 0.0;
715  }
716 
717  for (unsigned int l = 0; l < dim; ++l)
718  {
719  grads[i + (j + (2 * k + 5) * my_degree +
721  (my_degree + 1)][0][l] =
722  unit_point_grads[i +
723  ((j + 2) * (my_degree + 2) +
724  k) *
725  (my_degree + 1)][l];
726  grads[(i + 2 * (k + 4) * (my_degree + 1) +
728  my_degree +
729  j +
731  unit_point_grads[i +
732  (j + k * (my_degree + 2) + 2) *
733  (my_degree + 1)][l];
734  }
735 
736  grads[(i + 2 * k * (my_degree + 1) +
738  my_degree +
740  p1_grads[i + (j + k * (my_degree + 2) + 2) *
741  (my_degree + 1)][2];
742  grads[(i + 2 * k * (my_degree + 1) +
744  my_degree +
746  p1_grads[i + (j + k * (my_degree + 2) + 2) *
747  (my_degree + 1)][0];
748  grads[(i + 2 * k * (my_degree + 1) +
750  my_degree +
752  p1_grads[i + (j + k * (my_degree + 2) + 2) *
753  (my_degree + 1)][1];
754  grads[i + (j + (2 * k + 1) * my_degree +
756  (my_degree + 1)][2][0] =
757  p2_grads[i + ((j + 2) * (my_degree + 2) + k) *
758  (my_degree + 1)][1];
759  grads[i + (j + (2 * k + 1) * my_degree +
761  (my_degree + 1)][2][1] =
762  p2_grads[i + ((j + 2) * (my_degree + 2) + k) *
763  (my_degree + 1)][2];
764  grads[i + (j + (2 * k + 1) * my_degree +
766  (my_degree + 1)][2][2] =
767  p2_grads[i + ((j + 2) * (my_degree + 2) + k) *
768  (my_degree + 1)][0];
769  grads[(i + 2 * (k + 2) * (my_degree + 1) +
771  my_degree +
773  p2_grads[i + (j + k * (my_degree + 2) + 2) *
774  (my_degree + 1)][1];
775  grads[(i + 2 * (k + 2) * (my_degree + 1) +
777  my_degree +
779  p2_grads[i + (j + k * (my_degree + 2) + 2) *
780  (my_degree + 1)][2];
781  grads[(i + 2 * (k + 2) * (my_degree + 1) +
783  my_degree +
785  p2_grads[i + (j + k * (my_degree + 2) + 2) *
786  (my_degree + 1)][0];
787  grads[i + (j + (2 * k + 9) * my_degree +
789  (my_degree + 1)][1][0] =
790  p1_grads[i + ((j + 2) * (my_degree + 2) + k) *
791  (my_degree + 1)][2];
792  grads[i + (j + (2 * k + 9) * my_degree +
794  (my_degree + 1)][1][1] =
795  p1_grads[i + ((j + 2) * (my_degree + 2) + k) *
796  (my_degree + 1)][0];
797  grads[i + (j + (2 * k + 9) * my_degree +
799  (my_degree + 1)][1][2] =
800  p1_grads[i + ((j + 2) * (my_degree + 2) + k) *
801  (my_degree + 1)][1];
802  }
803  }
804  }
805 
806  if (grad_grads.size() > 0)
807  {
808  for (unsigned int i = 0; i <= my_degree; ++i)
809  {
810  for (unsigned int j = 0; j < 2; ++j)
811  {
812  for (unsigned int k = 0; k < 2; ++k)
813  {
814  for (unsigned int l = 0; l < dim; ++l)
815  for (unsigned int m = 0; m < dim; ++m)
816  {
817  for (unsigned int n = 0; n < 2; ++n)
818  {
819  grad_grads[i +
820  (j + 4 * k) * (my_degree + 1)]
821  [2 * n][l][m] = 0.0;
822  grad_grads[i + (j + 4 * k + 2) *
823  (my_degree + 1)][n + 1][l]
824  [m] = 0.0;
825  grad_grads[i + (j + 2 * (k + 4)) *
826  (my_degree + 1)][n][l][m] =
827  0.0;
828  }
829 
830  grad_grads[i + (j + 4 * k + 2) *
831  (my_degree + 1)][0][l][m] =
832  unit_point_grad_grads
833  [i + (j + k * (my_degree + 2)) *
834  (my_degree + 1)][l][m];
835  }
836 
837  grad_grads[i + (j + 2 * (k + 4)) *
838  (my_degree + 1)][2][0][0] =
839  p2_grad_grads[i + (j + k * (my_degree + 2)) *
840  (my_degree + 1)][1][1];
841  grad_grads[i + (j + 2 * (k + 4)) *
842  (my_degree + 1)][2][0][1] =
843  p2_grad_grads[i + (j + k * (my_degree + 2)) *
844  (my_degree + 1)][1][2];
845  grad_grads[i + (j + 2 * (k + 4)) *
846  (my_degree + 1)][2][0][2] =
847  p2_grad_grads[i + (j + k * (my_degree + 2)) *
848  (my_degree + 1)][1][0];
849  grad_grads[i + (j + 2 * (k + 4)) *
850  (my_degree + 1)][2][1][0] =
851  p2_grad_grads[i + (j + k * (my_degree + 2)) *
852  (my_degree + 1)][2][1];
853  grad_grads[i + (j + 2 * (k + 4)) *
854  (my_degree + 1)][2][1][1] =
855  p2_grad_grads[i + (j + k * (my_degree + 2)) *
856  (my_degree + 1)][2][2];
857  grad_grads[i + (j + 2 * (k + 4)) *
858  (my_degree + 1)][2][1][2] =
859  p2_grad_grads[i + (j + k * (my_degree + 2)) *
860  (my_degree + 1)][2][0];
861  grad_grads[i + (j + 2 * (k + 4)) *
862  (my_degree + 1)][2][2][0] =
863  p2_grad_grads[i + (j + k * (my_degree + 2)) *
864  (my_degree + 1)][0][1];
865  grad_grads[i + (j + 2 * (k + 4)) *
866  (my_degree + 1)][2][2][1] =
867  p2_grad_grads[i + (j + k * (my_degree + 2)) *
868  (my_degree + 1)][0][2];
869  grad_grads[i + (j + 2 * (k + 4)) *
870  (my_degree + 1)][2][2][2] =
871  p2_grad_grads[i + (j + k * (my_degree + 2)) *
872  (my_degree + 1)][0][0];
873  }
874 
875  grad_grads[i + j * (my_degree + 1)][1][0][0] =
876  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
877  [2][2];
878  grad_grads[i + j * (my_degree + 1)][1][0][1] =
879  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
880  [2][0];
881  grad_grads[i + j * (my_degree + 1)][1][0][2] =
882  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
883  [2][1];
884  grad_grads[i + j * (my_degree + 1)][1][1][0] =
885  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
886  [0][2];
887  grad_grads[i + j * (my_degree + 1)][1][1][1] =
888  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
889  [0][0];
890  grad_grads[i + j * (my_degree + 1)][1][1][2] =
891  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
892  [0][1];
893  grad_grads[i + j * (my_degree + 1)][1][2][0] =
894  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
895  [1][2];
896  grad_grads[i + j * (my_degree + 1)][1][2][1] =
897  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
898  [1][0];
899  grad_grads[i + j * (my_degree + 1)][1][2][2] =
900  p1_grad_grads[i + j * (my_degree + 1) * (my_degree + 2)]
901  [1][1];
902  }
903 
904  grad_grads[i + 4 * (my_degree + 1)][1][0][0] =
905  p1_grad_grads[i + my_degree + 1][2][2];
906  grad_grads[i + 4 * (my_degree + 1)][1][0][1] =
907  p1_grad_grads[i + my_degree + 1][2][0];
908  grad_grads[i + 4 * (my_degree + 1)][1][0][2] =
909  p1_grad_grads[i + my_degree + 1][2][1];
910  grad_grads[i + 4 * (my_degree + 1)][1][1][0] =
911  p1_grad_grads[i + my_degree + 1][0][2];
912  grad_grads[i + 4 * (my_degree + 1)][1][1][1] =
913  p1_grad_grads[i + my_degree + 1][0][0];
914  grad_grads[i + 4 * (my_degree + 1)][1][1][2] =
915  p1_grad_grads[i + my_degree + 1][0][1];
916  grad_grads[i + 4 * (my_degree + 1)][1][2][0] =
917  p1_grad_grads[i + my_degree + 1][1][2];
918  grad_grads[i + 4 * (my_degree + 1)][1][2][1] =
919  p1_grad_grads[i + my_degree + 1][1][0];
920  grad_grads[i + 4 * (my_degree + 1)][1][2][2] =
921  p1_grad_grads[i + my_degree + 1][1][1];
922  grad_grads[i + 5 * (my_degree + 1)][1][0][0] =
923  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][2][2];
924  grad_grads[i + 5 * (my_degree + 1)][1][0][1] =
925  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][2][0];
926  grad_grads[i + 5 * (my_degree + 1)][1][0][2] =
927  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][2][1];
928  grad_grads[i + 5 * (my_degree + 1)][1][1][0] =
929  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][0][2];
930  grad_grads[i + 5 * (my_degree + 1)][1][1][1] =
931  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][0][0];
932  grad_grads[i + 5 * (my_degree + 1)][1][1][2] =
933  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][0][1];
934  grad_grads[i + 5 * (my_degree + 1)][1][2][0] =
935  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][1][2];
936  grad_grads[i + 5 * (my_degree + 1)][1][2][1] =
937  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][1][0];
938  grad_grads[i + 5 * (my_degree + 1)][1][2][2] =
939  p1_grad_grads[i + (my_degree + 1) * (my_degree + 3)][1][1];
940  }
941 
942  if (my_degree > 0)
943  for (unsigned int i = 0; i <= my_degree; ++i)
944  for (unsigned int j = 0; j < my_degree; ++j)
945  {
946  for (unsigned int k = 0; k < my_degree; ++k)
947  {
948  for (unsigned int l = 0; l < dim; ++l)
949  for (unsigned int m = 0; m < dim; ++m)
950  {
951  for (unsigned int n = 0; n < 2; ++n)
952  {
953  grad_grads
954  [((i +
955  2 *
957  my_degree +
960  my_degree +
962  [n + 1][l][m] = 0.0;
963  grad_grads
964  [(i +
965  (j +
967  my_degree) *
968  (my_degree + 1) +
970  my_degree +
972  [2 * n][l][m] = 0.0;
973  grad_grads
974  [i + (j +
975  (k + 2 * (GeometryInfo<
976  dim>::faces_per_cell +
977  my_degree)) *
978  my_degree +
980  (my_degree + 1)][n][l][m] = 0.0;
981  }
982 
983  grad_grads
984  [((i +
986  my_degree +
989  my_degree +
991  [m] = unit_point_grad_grads
992  [i + (j + (k + 2) * (my_degree + 2) + 2) *
993  (my_degree + 1)][l][m];
994  }
995 
996  grad_grads
997  [(i +
999  my_degree) *
1000  (my_degree + 1) +
1002  my_degree +
1003  k + GeometryInfo<dim>::lines_per_cell][1][0][0] =
1004  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1005  2) *
1006  (my_degree + 1)][2][2];
1007  grad_grads
1008  [(i +
1010  my_degree) *
1011  (my_degree + 1) +
1013  my_degree +
1014  k + GeometryInfo<dim>::lines_per_cell][1][0][1] =
1015  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1016  2) *
1017  (my_degree + 1)][2][0];
1018  grad_grads
1019  [(i +
1021  my_degree) *
1022  (my_degree + 1) +
1024  my_degree +
1025  k + GeometryInfo<dim>::lines_per_cell][1][0][2] =
1026  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1027  2) *
1028  (my_degree + 1)][2][1];
1029  grad_grads
1030  [(i +
1032  my_degree) *
1033  (my_degree + 1) +
1035  my_degree +
1036  k + GeometryInfo<dim>::lines_per_cell][1][1][0] =
1037  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1038  2) *
1039  (my_degree + 1)][0][2];
1040  grad_grads
1041  [(i +
1043  my_degree) *
1044  (my_degree + 1) +
1046  my_degree +
1047  k + GeometryInfo<dim>::lines_per_cell][1][1][1] =
1048  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1049  2) *
1050  (my_degree + 1)][0][0];
1051  grad_grads
1052  [(i +
1054  my_degree) *
1055  (my_degree + 1) +
1057  my_degree +
1058  k + GeometryInfo<dim>::lines_per_cell][1][1][2] =
1059  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1060  2) *
1061  (my_degree + 1)][0][1];
1062  grad_grads
1063  [(i +
1065  my_degree) *
1066  (my_degree + 1) +
1068  my_degree +
1069  k + GeometryInfo<dim>::lines_per_cell][1][2][0] =
1070  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1071  2) *
1072  (my_degree + 1)][1][2];
1073  grad_grads
1074  [(i +
1076  my_degree) *
1077  (my_degree + 1) +
1079  my_degree +
1080  k + GeometryInfo<dim>::lines_per_cell][1][2][1] =
1081  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1082  2) *
1083  (my_degree + 1)][1][0];
1084  grad_grads
1085  [(i +
1087  my_degree) *
1088  (my_degree + 1) +
1090  my_degree +
1091  k + GeometryInfo<dim>::lines_per_cell][1][2][2] =
1092  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k +
1093  2) *
1094  (my_degree + 1)][1][1];
1095  grad_grads[i +
1096  (j +
1097  (k +
1099  my_degree)) *
1100  my_degree +
1102  (my_degree + 1)][2][0][0] =
1103  p2_grad_grads[i +
1104  (j + (k + 2) * (my_degree + 2) + 2) *
1105  (my_degree + 1)][1][1];
1106  grad_grads[i +
1107  (j +
1108  (k +
1110  my_degree)) *
1111  my_degree +
1113  (my_degree + 1)][2][0][1] =
1114  p2_grad_grads[i +
1115  (j + (k + 2) * (my_degree + 2) + 2) *
1116  (my_degree + 1)][1][2];
1117  grad_grads[i +
1118  (j +
1119  (k +
1121  my_degree)) *
1122  my_degree +
1124  (my_degree + 1)][2][0][2] =
1125  p2_grad_grads[i +
1126  (j + (k + 2) * (my_degree + 2) + 2) *
1127  (my_degree + 1)][1][0];
1128  grad_grads[i +
1129  (j +
1130  (k +
1132  my_degree)) *
1133  my_degree +
1135  (my_degree + 1)][2][1][0] =
1136  p2_grad_grads[i +
1137  (j + (k + 2) * (my_degree + 2) + 2) *
1138  (my_degree + 1)][2][1];
1139  grad_grads[i +
1140  (j +
1141  (k +
1143  my_degree)) *
1144  my_degree +
1146  (my_degree + 1)][2][1][1] =
1147  p2_grad_grads[i +
1148  (j + (k + 2) * (my_degree + 2) + 2) *
1149  (my_degree + 1)][2][2];
1150  grad_grads[i +
1151  (j +
1152  (k +
1154  my_degree)) *
1155  my_degree +
1157  (my_degree + 1)][2][1][2] =
1158  p2_grad_grads[i +
1159  (j + (k + 2) * (my_degree + 2) + 2) *
1160  (my_degree + 1)][2][0];
1161  grad_grads[i +
1162  (j +
1163  (k +
1165  my_degree)) *
1166  my_degree +
1168  (my_degree + 1)][2][2][0] =
1169  p2_grad_grads[i +
1170  (j + (k + 2) * (my_degree + 2) + 2) *
1171  (my_degree + 1)][0][1];
1172  grad_grads[i +
1173  (j +
1174  (k +
1176  my_degree)) *
1177  my_degree +
1179  (my_degree + 1)][2][2][1] =
1180  p2_grad_grads[i +
1181  (j + (k + 2) * (my_degree + 2) + 2) *
1182  (my_degree + 1)][0][2];
1183  grad_grads[i +
1184  (j +
1185  (k +
1187  my_degree)) *
1188  my_degree +
1190  (my_degree + 1)][2][2][2] =
1191  p2_grad_grads[i +
1192  (j + (k + 2) * (my_degree + 2) + 2) *
1193  (my_degree + 1)][0][0];
1194  }
1195 
1196  for (unsigned int k = 0; k < 2; ++k)
1197  {
1198  for (unsigned int l = 0; l < dim; ++l)
1199  for (unsigned int m = 0; m < dim; ++m)
1200  {
1201  for (unsigned int n = 0; n < 2; ++n)
1202  {
1203  for (unsigned int o = 0; o < 2; ++o)
1204  {
1205  grad_grads
1206  [i +
1207  (j +
1208  (2 * (k + 2 * n) + 1) * my_degree +
1210  (my_degree + 1)][o + n][l][m] =
1211  0.0;
1212  grad_grads
1213  [(i +
1214  2 * (k + 2 * (n + 1)) *
1215  (my_degree + 1) +
1217  my_degree +
1218  j +
1220  [o + k][l][m] = 0.0;
1221  }
1222 
1223  grad_grads
1224  [(i + 2 * k * (my_degree + 1) +
1226  my_degree +
1228  [2 * n][l][m] = 0.0;
1229  grad_grads
1230  [i + (j + (2 * k + 9) * my_degree +
1232  (my_degree + 1)][2 * n][l][m] =
1233  0.0;
1234  }
1235 
1236  grad_grads[i +
1237  (j + (2 * k + 5) * my_degree +
1239  (my_degree + 1)][0][l][m] =
1240  unit_point_grad_grads
1241  [i + ((j + 2) * (my_degree + 2) + k) *
1242  (my_degree + 1)][l][m];
1243  grad_grads[(i + 2 * (k + 4) * (my_degree + 1) +
1245  my_degree +
1246  j +
1248  [l][m] = unit_point_grad_grads
1249  [i + (j + k * (my_degree + 2) + 2) *
1250  (my_degree + 1)][l][m];
1251  }
1252 
1253  grad_grads
1254  [(i + 2 * k * (my_degree + 1) +
1256  my_degree +
1257  j + GeometryInfo<dim>::lines_per_cell][1][0][0] =
1258  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1259  (my_degree + 1)][2][2];
1260  grad_grads
1261  [(i + 2 * k * (my_degree + 1) +
1263  my_degree +
1264  j + GeometryInfo<dim>::lines_per_cell][1][0][1] =
1265  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1266  (my_degree + 1)][2][0];
1267  grad_grads
1268  [(i + 2 * k * (my_degree + 1) +
1270  my_degree +
1271  j + GeometryInfo<dim>::lines_per_cell][1][0][2] =
1272  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1273  (my_degree + 1)][2][1];
1274  grad_grads
1275  [(i + 2 * k * (my_degree + 1) +
1277  my_degree +
1278  j + GeometryInfo<dim>::lines_per_cell][1][1][0] =
1279  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1280  (my_degree + 1)][0][2];
1281  grad_grads
1282  [(i + 2 * k * (my_degree + 1) +
1284  my_degree +
1285  j + GeometryInfo<dim>::lines_per_cell][1][1][1] =
1286  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1287  (my_degree + 1)][0][0];
1288  grad_grads
1289  [(i + 2 * k * (my_degree + 1) +
1291  my_degree +
1292  j + GeometryInfo<dim>::lines_per_cell][1][1][2] =
1293  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1294  (my_degree + 1)][0][1];
1295  grad_grads
1296  [(i + 2 * k * (my_degree + 1) +
1298  my_degree +
1299  j + GeometryInfo<dim>::lines_per_cell][1][2][0] =
1300  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1301  (my_degree + 1)][1][2];
1302  grad_grads
1303  [(i + 2 * k * (my_degree + 1) +
1305  my_degree +
1306  j + GeometryInfo<dim>::lines_per_cell][1][2][1] =
1307  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1308  (my_degree + 1)][1][0];
1309  grad_grads
1310  [(i + 2 * k * (my_degree + 1) +
1312  my_degree +
1313  j + GeometryInfo<dim>::lines_per_cell][1][2][2] =
1314  p1_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1315  (my_degree + 1)][1][1];
1316  grad_grads[i + (j + (2 * k + 1) * my_degree +
1318  (my_degree + 1)][2][0][0] =
1319  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1320  (my_degree + 1)][1][1];
1321  grad_grads[i + (j + (2 * k + 1) * my_degree +
1323  (my_degree + 1)][2][0][1] =
1324  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1325  (my_degree + 1)][1][2];
1326  grad_grads[i + (j + (2 * k + 1) * my_degree +
1328  (my_degree + 1)][2][0][2] =
1329  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1330  (my_degree + 1)][1][0];
1331  grad_grads[i + (j + (2 * k + 1) * my_degree +
1333  (my_degree + 1)][2][1][0] =
1334  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1335  (my_degree + 1)][2][1];
1336  grad_grads[i + (j + (2 * k + 1) * my_degree +
1338  (my_degree + 1)][2][1][1] =
1339  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1340  (my_degree + 1)][2][2];
1341  grad_grads[i + (j + (2 * k + 1) * my_degree +
1343  (my_degree + 1)][2][1][2] =
1344  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1345  (my_degree + 1)][2][0];
1346  grad_grads[i + (j + (2 * k + 1) * my_degree +
1348  (my_degree + 1)][2][2][0] =
1349  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1350  (my_degree + 1)][0][1];
1351  grad_grads[i + (j + (2 * k + 1) * my_degree +
1353  (my_degree + 1)][2][2][1] =
1354  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1355  (my_degree + 1)][0][2];
1356  grad_grads[i + (j + (2 * k + 1) * my_degree +
1358  (my_degree + 1)][2][2][2] =
1359  p2_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1360  (my_degree + 1)][0][0];
1361  grad_grads
1362  [(i + 2 * (k + 2) * (my_degree + 1) +
1364  my_degree +
1365  j + GeometryInfo<dim>::lines_per_cell][2][0][0] =
1366  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1367  (my_degree + 1)][1][1];
1368  grad_grads
1369  [(i + 2 * (k + 2) * (my_degree + 1) +
1371  my_degree +
1372  j + GeometryInfo<dim>::lines_per_cell][2][0][1] =
1373  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1374  (my_degree + 1)][1][2];
1375  grad_grads
1376  [(i + 2 * (k + 2) * (my_degree + 1) +
1378  my_degree +
1379  j + GeometryInfo<dim>::lines_per_cell][2][0][2] =
1380  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1381  (my_degree + 1)][1][0];
1382  grad_grads
1383  [(i + 2 * (k + 2) * (my_degree + 1) +
1385  my_degree +
1386  j + GeometryInfo<dim>::lines_per_cell][2][1][0] =
1387  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1388  (my_degree + 1)][2][1];
1389  grad_grads
1390  [(i + 2 * (k + 2) * (my_degree + 1) +
1392  my_degree +
1393  j + GeometryInfo<dim>::lines_per_cell][2][1][1] =
1394  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1395  (my_degree + 1)][2][2];
1396  grad_grads
1397  [(i + 2 * (k + 2) * (my_degree + 1) +
1399  my_degree +
1400  j + GeometryInfo<dim>::lines_per_cell][2][1][2] =
1401  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1402  (my_degree + 1)][2][0];
1403  grad_grads
1404  [(i + 2 * (k + 2) * (my_degree + 1) +
1406  my_degree +
1407  j + GeometryInfo<dim>::lines_per_cell][2][2][0] =
1408  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1409  (my_degree + 1)][0][1];
1410  grad_grads
1411  [(i + 2 * (k + 2) * (my_degree + 1) +
1413  my_degree +
1414  j + GeometryInfo<dim>::lines_per_cell][2][2][1] =
1415  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1416  (my_degree + 1)][0][2];
1417  grad_grads
1418  [(i + 2 * (k + 2) * (my_degree + 1) +
1420  my_degree +
1421  j + GeometryInfo<dim>::lines_per_cell][2][2][2] =
1422  p2_grad_grads[i + (j + k * (my_degree + 2) + 2) *
1423  (my_degree + 1)][0][0];
1424  grad_grads[i + (j + (2 * k + 9) * my_degree +
1426  (my_degree + 1)][1][0][0] =
1427  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1428  (my_degree + 1)][2][2];
1429  grad_grads[i + (j + (2 * k + 9) * my_degree +
1431  (my_degree + 1)][1][0][1] =
1432  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1433  (my_degree + 1)][2][0];
1434  grad_grads[i + (j + (2 * k + 9) * my_degree +
1436  (my_degree + 1)][1][0][2] =
1437  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1438  (my_degree + 1)][2][1];
1439  grad_grads[i + (j + (2 * k + 9) * my_degree +
1441  (my_degree + 1)][1][1][0] =
1442  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1443  (my_degree + 1)][0][2];
1444  grad_grads[i + (j + (2 * k + 9) * my_degree +
1446  (my_degree + 1)][1][1][1] =
1447  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1448  (my_degree + 1)][0][0];
1449  grad_grads[i + (j + (2 * k + 9) * my_degree +
1451  (my_degree + 1)][1][1][2] =
1452  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1453  (my_degree + 1)][0][1];
1454  grad_grads[i + (j + (2 * k + 9) * my_degree +
1456  (my_degree + 1)][1][2][0] =
1457  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1458  (my_degree + 1)][1][2];
1459  grad_grads[i + (j + (2 * k + 9) * my_degree +
1461  (my_degree + 1)][1][2][1] =
1462  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1463  (my_degree + 1)][1][0];
1464  grad_grads[i + (j + (2 * k + 9) * my_degree +
1466  (my_degree + 1)][1][2][2] =
1467  p1_grad_grads[i + ((j + 2) * (my_degree + 2) + k) *
1468  (my_degree + 1)][1][1];
1469  }
1470  }
1471  }
1472 
1473  break;
1474  }
1475 
1476  default:
1477  Assert(false, ExcNotImplemented());
1478  }
1479 }
1480 
1481 
1482 template <int dim>
1483 unsigned int
1485 {
1486  switch (dim)
1487  {
1488  case 1:
1489  return k + 1;
1490 
1491  case 2:
1492  return 2 * (k + 1) * (k + 2);
1493 
1494  case 3:
1495  return 3 * (k + 1) * (k + 2) * (k + 2);
1496 
1497  default:
1498  {
1499  Assert(false, ExcNotImplemented());
1500  return 0;
1501  }
1502  }
1503 }
1504 
1505 
1506 template <int dim>
1507 std::unique_ptr<TensorPolynomialsBase<dim>>
1509 {
1510  return std::make_unique<PolynomialsNedelec<dim>>(*this);
1511 }
1512 
1513 
1514 template class PolynomialsNedelec<1>;
1515 template class PolynomialsNedelec<2>;
1516 template class PolynomialsNedelec<3>;
1517 
1518 
Definition: point.h:112
virtual std::unique_ptr< TensorPolynomialsBase< dim > > clone() const override
PolynomialsNedelec(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 std::vector< std::vector< Polynomials::Polynomial< double > > > create_polynomials(const unsigned int k)
static unsigned int n_polynomials(const unsigned int degree)
static std::vector< Polynomial< double > > generate_complete_basis(const unsigned int degree)
Definition: polynomial.cc:730
static std::vector< Polynomial< double > > generate_complete_basis(const unsigned int p)
Definition: polynomial.cc:835
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcNotImplemented()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)