Reference documentation for deal.II version GIT e22cb6a53e 2023-09-21 14:00: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_rt_bubbles.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 - 2023 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 
16 
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  , raviart_thomas_space(k - 1)
32  , monomials(k + 2)
33 {
34  Assert(dim >= 2, ExcImpossibleInDim(dim));
35 
36  for (unsigned int i = 0; i < monomials.size(); ++i)
38 }
39 
40 
41 
42 template <int dim>
43 void
45  const Point<dim> &unit_point,
46  std::vector<Tensor<1, dim>> &values,
47  std::vector<Tensor<2, dim>> &grads,
48  std::vector<Tensor<3, dim>> &grad_grads,
49  std::vector<Tensor<4, dim>> &third_derivatives,
50  std::vector<Tensor<5, dim>> &fourth_derivatives) const
51 {
52  Assert(values.size() == this->n() || values.empty(),
53  ExcDimensionMismatch(values.size(), this->n()));
54  Assert(grads.size() == this->n() || grads.empty(),
55  ExcDimensionMismatch(grads.size(), this->n()));
56  Assert(grad_grads.size() == this->n() || grad_grads.empty(),
57  ExcDimensionMismatch(grad_grads.size(), this->n()));
58  Assert(third_derivatives.size() == this->n() || third_derivatives.empty(),
59  ExcDimensionMismatch(third_derivatives.size(), this->n()));
60  Assert(fourth_derivatives.size() == this->n() || fourth_derivatives.empty(),
61  ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
62 
63  // Third and fourth derivatives are not implemented
64  (void)third_derivatives;
65  Assert(third_derivatives.empty(), ExcNotImplemented());
66  (void)fourth_derivatives;
67  Assert(fourth_derivatives.empty(), ExcNotImplemented());
68 
69  const unsigned int n_sub = raviart_thomas_space.n();
70  const unsigned int my_degree = this->degree();
71 
72  // Guard access to the scratch arrays in the following block
73  // using a mutex to make sure they are not used by multiple threads
74  // at once
75  {
76  static std::mutex mutex;
77  std::lock_guard<std::mutex> lock(mutex);
78 
79  static std::vector<Tensor<1, dim>> p_values;
80  static std::vector<Tensor<2, dim>> p_grads;
81  static std::vector<Tensor<3, dim>> p_grad_grads;
82  static std::vector<Tensor<4, dim>> p_third_derivatives;
83  static std::vector<Tensor<5, dim>> p_fourth_derivatives;
84 
85  p_values.resize((values.empty()) ? 0 : n_sub);
86  p_grads.resize((grads.empty()) ? 0 : n_sub);
87  p_grad_grads.resize((grad_grads.empty()) ? 0 : n_sub);
88 
89  // This is the Raviart-Thomas part of the space
90  raviart_thomas_space.evaluate(unit_point,
91  p_values,
92  p_grads,
93  p_grad_grads,
94  p_third_derivatives,
95  p_fourth_derivatives);
96  for (unsigned int i = 0; i < p_values.size(); ++i)
97  values[i] = p_values[i];
98  for (unsigned int i = 0; i < p_grads.size(); ++i)
99  grads[i] = p_grads[i];
100  for (unsigned int i = 0; i < p_grad_grads.size(); ++i)
101  grad_grads[i] = p_grad_grads[i];
102  }
103 
104  // Next we compute the polynomials and derivatives
105  // of the curl part of the space
106  const unsigned int n_derivatives = 3;
107  double monoval_plus[dim][n_derivatives + 1];
108  double monoval_i[dim][n_derivatives + 1];
109 
110 
111  if constexpr (dim <= 1)
112  {
113  (void)monoval_plus;
114  (void)monoval_i;
115  }
116 
117  unsigned int start = n_sub;
118  if constexpr (dim == 2)
119  {
120  // In 2d the curl part of the space is spanned by the vectors
121  // of two types. The first one is
122  // [ x^i * [y^(k+1)]' ]
123  // [ -[x^i]' * y^(k+1) ]
124  // The second one can be obtained from the first by a cyclic
125  // rotation of the coordinates.
126  // monoval_i = x^i,
127  // monoval_plus = x^(k+1)
128  for (unsigned int d = 0; d < dim; ++d)
129  monomials[my_degree + 1].value(unit_point(d),
130  n_derivatives,
131  monoval_plus[d]);
132 
133  for (unsigned int i = 0; i <= my_degree; ++i, ++start)
134  {
135  for (unsigned int d = 0; d < dim; ++d)
136  monomials[i].value(unit_point(d), n_derivatives, monoval_i[d]);
137 
138  if (values.size() != 0)
139  {
140  values[start][0] = monoval_i[0][0] * monoval_plus[1][1];
141  values[start][1] = -monoval_i[0][1] * monoval_plus[1][0];
142 
143  values[start + my_degree + 1][0] =
144  -monoval_plus[0][0] * monoval_i[1][1];
145  values[start + my_degree + 1][1] =
146  monoval_plus[0][1] * monoval_i[1][0];
147  }
148 
149  if (grads.size() != 0)
150  {
151  grads[start][0][0] = monoval_i[0][1] * monoval_plus[1][1];
152  grads[start][0][1] = monoval_i[0][0] * monoval_plus[1][2];
153  grads[start][1][0] = -monoval_i[0][2] * monoval_plus[1][0];
154  grads[start][1][1] = -monoval_i[0][1] * monoval_plus[1][1];
155 
156  grads[start + my_degree + 1][0][0] =
157  -monoval_plus[0][1] * monoval_i[1][1];
158  grads[start + my_degree + 1][0][1] =
159  -monoval_plus[0][0] * monoval_i[1][2];
160  grads[start + my_degree + 1][1][0] =
161  monoval_plus[0][2] * monoval_i[1][0];
162  grads[start + my_degree + 1][1][1] =
163  monoval_plus[0][1] * monoval_i[1][1];
164  }
165 
166  if (grad_grads.size() != 0)
167  {
168  grad_grads[start][0][0][0] = monoval_i[0][2] * monoval_plus[1][1];
169  grad_grads[start][0][0][1] = monoval_i[0][1] * monoval_plus[1][2];
170  grad_grads[start][0][1][0] = monoval_i[0][1] * monoval_plus[1][2];
171  grad_grads[start][0][1][1] = monoval_i[0][0] * monoval_plus[1][3];
172  grad_grads[start][1][0][0] =
173  -monoval_i[0][3] * monoval_plus[1][0];
174  grad_grads[start][1][0][1] =
175  -monoval_i[0][2] * monoval_plus[1][1];
176  grad_grads[start][1][1][0] =
177  -monoval_i[0][2] * monoval_plus[1][1];
178  grad_grads[start][1][1][1] =
179  -monoval_i[0][1] * monoval_plus[1][2];
180 
181  grad_grads[start + my_degree + 1][0][0][0] =
182  -monoval_plus[0][2] * monoval_i[1][1];
183  grad_grads[start + my_degree + 1][0][0][1] =
184  -monoval_plus[0][1] * monoval_i[1][2];
185  grad_grads[start + my_degree + 1][0][1][0] =
186  -monoval_plus[0][1] * monoval_i[1][2];
187  grad_grads[start + my_degree + 1][0][1][1] =
188  -monoval_plus[0][0] * monoval_i[1][3];
189  grad_grads[start + my_degree + 1][1][0][0] =
190  monoval_plus[0][3] * monoval_i[1][0];
191  grad_grads[start + my_degree + 1][1][0][1] =
192  monoval_plus[0][2] * monoval_i[1][1];
193  grad_grads[start + my_degree + 1][1][1][0] =
194  monoval_plus[0][2] * monoval_i[1][1];
195  grad_grads[start + my_degree + 1][1][1][1] =
196  monoval_plus[0][1] * monoval_i[1][2];
197  }
198  }
199  Assert(start == this->n() - my_degree - 1, ExcInternalError());
200  }
201  else if constexpr (dim == 3)
202  {
203  double monoval[dim][n_derivatives + 1];
204  double monoval_j[dim][n_derivatives + 1];
205  double monoval_jplus[dim][n_derivatives + 1];
206 
207  // In 3d the first type of basis vector is
208  // [ x^i * y^j * z^k * (j+k+2) ]
209  // [ -[x^i]' * y^(j+1) * z^k ]
210  // [ -[x^i]' * y^j * z^(k+1) ],
211  // For the second type of basis vector y and z
212  // are swapped. Then for each of these,
213  // two more are obtained by the cyclic rotation
214  // of the coordinates.
215  // monoval = x^k, monoval_plus = x^(k+1)
216  // monoval_* = x^*, monoval_jplus = x^(j+1)
217  for (unsigned int d = 0; d < dim; ++d)
218  {
219  monomials[my_degree + 1].value(unit_point(d),
220  n_derivatives,
221  monoval_plus[d]);
222  monomials[my_degree].value(unit_point(d), n_derivatives, monoval[d]);
223  }
224 
225  const unsigned int n_curls = (my_degree + 1) * (2 * my_degree + 1);
226  // Span of @f$\tilde{B}@f$
227  for (unsigned int i = 0; i <= my_degree; ++i)
228  {
229  for (unsigned int d = 0; d < dim; ++d)
230  monomials[i].value(unit_point(d), n_derivatives, monoval_i[d]);
231 
232  for (unsigned int j = 0; j <= my_degree; ++j)
233  {
234  for (unsigned int d = 0; d < dim; ++d)
235  {
236  monomials[j].value(unit_point(d),
237  n_derivatives,
238  monoval_j[d]);
239  monomials[j + 1].value(unit_point(d),
240  n_derivatives,
241  monoval_jplus[d]);
242  }
243 
244  if (values.size() != 0)
245  {
246  values[start][0] = monoval_i[0][0] * monoval_j[1][0] *
247  monoval[2][0] *
248  static_cast<double>(j + my_degree + 2);
249  values[start][1] =
250  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][0];
251  values[start][2] =
252  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][0];
253 
254  values[start + n_curls][0] =
255  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][0];
256  values[start + n_curls][1] =
257  monoval_j[0][0] * monoval_i[1][0] * monoval[2][0] *
258  static_cast<double>(j + my_degree + 2);
259  values[start + n_curls][2] =
260  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][0];
261 
262  values[start + 2 * n_curls][0] =
263  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][1];
264  values[start + 2 * n_curls][1] =
265  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][1];
266  values[start + 2 * n_curls][2] =
267  monoval_j[0][0] * monoval[1][0] * monoval_i[2][0] *
268  static_cast<double>(j + my_degree + 2);
269 
270  // Only unique triples of powers (i j k)
271  // and (i k j) are allowed, 0 <= i,j <= k
272  if (j != my_degree)
273  {
274  values[start + 1][0] =
275  monoval_i[0][0] * monoval[1][0] * monoval_j[2][0] *
276  static_cast<double>(j + my_degree + 2);
277  values[start + 1][1] =
278  -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][0];
279  values[start + 1][2] =
280  -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][0];
281 
282  values[start + n_curls + 1][0] =
283  -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][0];
284  values[start + n_curls + 1][1] =
285  monoval[0][0] * monoval_i[1][0] * monoval_j[2][0] *
286  static_cast<double>(j + my_degree + 2);
287  values[start + n_curls + 1][2] =
288  -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][0];
289 
290  values[start + 2 * n_curls + 1][0] =
291  -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][1];
292  values[start + 2 * n_curls + 1][1] =
293  -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][1];
294  values[start + 2 * n_curls + 1][2] =
295  monoval[0][0] * monoval_j[1][0] * monoval_i[2][0] *
296  static_cast<double>(j + my_degree + 2);
297  }
298  }
299 
300  if (grads.size() != 0)
301  {
302  grads[start][0][0] = monoval_i[0][1] * monoval_j[1][0] *
303  monoval[2][0] *
304  static_cast<double>(j + my_degree + 2);
305  grads[start][0][1] = monoval_i[0][0] * monoval_j[1][1] *
306  monoval[2][0] *
307  static_cast<double>(j + my_degree + 2);
308  grads[start][0][2] = monoval_i[0][0] * monoval_j[1][0] *
309  monoval[2][1] *
310  static_cast<double>(j + my_degree + 2);
311  grads[start][1][0] =
312  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][0];
313  grads[start][1][1] =
314  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][0];
315  grads[start][1][2] =
316  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][1];
317  grads[start][2][0] =
318  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][0];
319  grads[start][2][1] =
320  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][0];
321  grads[start][2][2] =
322  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][1];
323 
324  grads[start + n_curls][0][0] =
325  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][0];
326  grads[start + n_curls][0][1] =
327  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][0];
328  grads[start + n_curls][0][2] =
329  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][1];
330  grads[start + n_curls][1][0] =
331  monoval_j[0][1] * monoval_i[1][0] * monoval[2][0] *
332  static_cast<double>(j + my_degree + 2);
333  grads[start + n_curls][1][1] =
334  monoval_j[0][0] * monoval_i[1][1] * monoval[2][0] *
335  static_cast<double>(j + my_degree + 2);
336  grads[start + n_curls][1][2] =
337  monoval_j[0][0] * monoval_i[1][0] * monoval[2][1] *
338  static_cast<double>(j + my_degree + 2);
339  grads[start + n_curls][2][0] =
340  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][0];
341  grads[start + n_curls][2][1] =
342  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][0];
343  grads[start + n_curls][2][2] =
344  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][1];
345 
346  grads[start + 2 * n_curls][0][0] =
347  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][1];
348  grads[start + 2 * n_curls][0][1] =
349  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][1];
350  grads[start + 2 * n_curls][0][2] =
351  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][2];
352  grads[start + 2 * n_curls][1][0] =
353  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][1];
354  grads[start + 2 * n_curls][1][1] =
355  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][1];
356  grads[start + 2 * n_curls][1][2] =
357  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][2];
358  grads[start + 2 * n_curls][2][0] =
359  monoval_j[0][1] * monoval[1][0] * monoval_i[2][0] *
360  static_cast<double>(j + my_degree + 2);
361  grads[start + 2 * n_curls][2][1] =
362  monoval_j[0][0] * monoval[1][1] * monoval_i[2][0] *
363  static_cast<double>(j + my_degree + 2);
364  grads[start + 2 * n_curls][2][2] =
365  monoval_j[0][0] * monoval[1][0] * monoval_i[2][1] *
366  static_cast<double>(j + my_degree + 2);
367 
368  if (j != my_degree)
369  {
370  grads[start + 1][0][0] =
371  monoval_i[0][1] * monoval[1][0] * monoval_j[2][0] *
372  static_cast<double>(j + my_degree + 2);
373  grads[start + 1][0][1] =
374  monoval_i[0][0] * monoval[1][1] * monoval_j[2][0] *
375  static_cast<double>(j + my_degree + 2);
376  grads[start + 1][0][2] =
377  monoval_i[0][0] * monoval[1][0] * monoval_j[2][1] *
378  static_cast<double>(j + my_degree + 2);
379  grads[start + 1][1][0] =
380  -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][0];
381  grads[start + 1][1][1] =
382  -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][0];
383  grads[start + 1][1][2] =
384  -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][1];
385  grads[start + 1][2][0] =
386  -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][0];
387  grads[start + 1][2][1] =
388  -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][0];
389  grads[start + 1][2][2] =
390  -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][1];
391 
392  grads[start + n_curls + 1][0][0] =
393  -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][0];
394  grads[start + n_curls + 1][0][1] =
395  -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][0];
396  grads[start + n_curls + 1][0][2] =
397  -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][1];
398  grads[start + n_curls + 1][1][0] =
399  monoval[0][1] * monoval_i[1][0] * monoval_j[2][0] *
400  static_cast<double>(j + my_degree + 2);
401  grads[start + n_curls + 1][1][1] =
402  monoval[0][0] * monoval_i[1][1] * monoval_j[2][0] *
403  static_cast<double>(j + my_degree + 2);
404  grads[start + n_curls + 1][1][2] =
405  monoval[0][0] * monoval_i[1][0] * monoval_j[2][1] *
406  static_cast<double>(j + my_degree + 2);
407  grads[start + n_curls + 1][2][0] =
408  -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][0];
409  grads[start + n_curls + 1][2][1] =
410  -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][0];
411  grads[start + n_curls + 1][2][2] =
412  -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][1];
413 
414  grads[start + 2 * n_curls + 1][0][0] =
415  -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][1];
416  grads[start + 2 * n_curls + 1][0][1] =
417  -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][1];
418  grads[start + 2 * n_curls + 1][0][2] =
419  -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][2];
420  grads[start + 2 * n_curls + 1][1][0] =
421  -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][1];
422  grads[start + 2 * n_curls + 1][1][1] =
423  -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][1];
424  grads[start + 2 * n_curls + 1][1][2] =
425  -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][2];
426  grads[start + 2 * n_curls + 1][2][0] =
427  monoval[0][1] * monoval_j[1][0] * monoval_i[2][0] *
428  static_cast<double>(j + my_degree + 2);
429  grads[start + 2 * n_curls + 1][2][1] =
430  monoval[0][0] * monoval_j[1][1] * monoval_i[2][0] *
431  static_cast<double>(j + my_degree + 2);
432  grads[start + 2 * n_curls + 1][2][2] =
433  monoval[0][0] * monoval_j[1][0] * monoval_i[2][1] *
434  static_cast<double>(j + my_degree + 2);
435  }
436  }
437 
438  if (grad_grads.size() != 0)
439  {
440  grad_grads[start][0][0][0] =
441  monoval_i[0][2] * monoval_j[1][0] * monoval[2][0] *
442  static_cast<double>(j + my_degree + 2);
443  grad_grads[start][0][0][1] =
444  monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
445  static_cast<double>(j + my_degree + 2);
446  grad_grads[start][0][0][2] =
447  monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
448  static_cast<double>(j + my_degree + 2);
449  grad_grads[start][0][1][0] =
450  monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
451  static_cast<double>(j + my_degree + 2);
452  grad_grads[start][0][1][1] =
453  monoval_i[0][0] * monoval_j[1][2] * monoval[2][0] *
454  static_cast<double>(j + my_degree + 2);
455  grad_grads[start][0][1][2] =
456  monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
457  static_cast<double>(j + my_degree + 2);
458  grad_grads[start][0][2][0] =
459  monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
460  static_cast<double>(j + my_degree + 2);
461  grad_grads[start][0][2][1] =
462  monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
463  static_cast<double>(j + my_degree + 2);
464  grad_grads[start][0][2][2] =
465  monoval_i[0][0] * monoval_j[1][0] * monoval[2][2] *
466  static_cast<double>(j + my_degree + 2);
467  grad_grads[start][1][0][0] =
468  -monoval_i[0][3] * monoval_jplus[1][0] * monoval[2][0];
469  grad_grads[start][1][0][1] =
470  -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
471  grad_grads[start][1][0][2] =
472  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
473  grad_grads[start][1][1][0] =
474  -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
475  grad_grads[start][1][1][1] =
476  -monoval_i[0][1] * monoval_jplus[1][2] * monoval[2][0];
477  grad_grads[start][1][1][2] =
478  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
479  grad_grads[start][1][2][0] =
480  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
481  grad_grads[start][1][2][1] =
482  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
483  grad_grads[start][1][2][2] =
484  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][2];
485  grad_grads[start][2][0][0] =
486  -monoval_i[0][3] * monoval_j[1][0] * monoval_plus[2][0];
487  grad_grads[start][2][0][1] =
488  -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
489  grad_grads[start][2][0][2] =
490  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
491  grad_grads[start][2][1][0] =
492  -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
493  grad_grads[start][2][1][1] =
494  -monoval_i[0][1] * monoval_j[1][2] * monoval_plus[2][0];
495  grad_grads[start][2][1][2] =
496  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
497  grad_grads[start][2][2][0] =
498  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
499  grad_grads[start][2][2][1] =
500  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
501  grad_grads[start][2][2][2] =
502  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][2];
503 
504  grad_grads[start + n_curls][0][0][0] =
505  -monoval_jplus[0][2] * monoval_i[1][1] * monoval[2][0];
506  grad_grads[start + n_curls][0][0][1] =
507  -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
508  grad_grads[start + n_curls][0][0][2] =
509  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
510  grad_grads[start + n_curls][0][1][0] =
511  -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
512  grad_grads[start + n_curls][0][1][1] =
513  -monoval_jplus[0][0] * monoval_i[1][3] * monoval[2][0];
514  grad_grads[start + n_curls][0][1][2] =
515  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
516  grad_grads[start + n_curls][0][2][0] =
517  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
518  grad_grads[start + n_curls][0][2][1] =
519  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
520  grad_grads[start + n_curls][0][2][2] =
521  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][2];
522  grad_grads[start + n_curls][1][0][0] =
523  monoval_j[0][2] * monoval_i[1][0] * monoval[2][0] *
524  static_cast<double>(j + my_degree + 2);
525  grad_grads[start + n_curls][1][0][1] =
526  monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
527  static_cast<double>(j + my_degree + 2);
528  grad_grads[start + n_curls][1][0][2] =
529  monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
530  static_cast<double>(j + my_degree + 2);
531  grad_grads[start + n_curls][1][1][0] =
532  monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
533  static_cast<double>(j + my_degree + 2);
534  grad_grads[start + n_curls][1][1][1] =
535  monoval_j[0][0] * monoval_i[1][2] * monoval[2][0] *
536  static_cast<double>(j + my_degree + 2);
537  grad_grads[start + n_curls][1][1][2] =
538  monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
539  static_cast<double>(j + my_degree + 2);
540  grad_grads[start + n_curls][1][2][0] =
541  monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
542  static_cast<double>(j + my_degree + 2);
543  grad_grads[start + n_curls][1][2][1] =
544  monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
545  static_cast<double>(j + my_degree + 2);
546  grad_grads[start + n_curls][1][2][2] =
547  monoval_j[0][0] * monoval_i[1][0] * monoval[2][2] *
548  static_cast<double>(j + my_degree + 2);
549  grad_grads[start + n_curls][2][0][0] =
550  -monoval_j[0][2] * monoval_i[1][1] * monoval_plus[2][0];
551  grad_grads[start + n_curls][2][0][1] =
552  -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
553  grad_grads[start + n_curls][2][0][2] =
554  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
555  grad_grads[start + n_curls][2][1][0] =
556  -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
557  grad_grads[start + n_curls][2][1][1] =
558  -monoval_j[0][0] * monoval_i[1][3] * monoval_plus[2][0];
559  grad_grads[start + n_curls][2][1][2] =
560  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
561  grad_grads[start + n_curls][2][2][0] =
562  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
563  grad_grads[start + n_curls][2][2][1] =
564  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
565  grad_grads[start + n_curls][2][2][2] =
566  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][2];
567 
568  grad_grads[start + 2 * n_curls][0][0][0] =
569  -monoval_jplus[0][2] * monoval[1][0] * monoval_i[2][1];
570  grad_grads[start + 2 * n_curls][0][0][1] =
571  -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
572  grad_grads[start + 2 * n_curls][0][0][2] =
573  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
574  grad_grads[start + 2 * n_curls][0][1][0] =
575  -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
576  grad_grads[start + 2 * n_curls][0][1][1] =
577  -monoval_jplus[0][0] * monoval[1][2] * monoval_i[2][1];
578  grad_grads[start + 2 * n_curls][0][1][2] =
579  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
580  grad_grads[start + 2 * n_curls][0][2][0] =
581  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
582  grad_grads[start + 2 * n_curls][0][2][1] =
583  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
584  grad_grads[start + 2 * n_curls][0][2][2] =
585  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][3];
586  grad_grads[start + 2 * n_curls][1][0][0] =
587  -monoval_j[0][2] * monoval_plus[1][0] * monoval_i[2][1];
588  grad_grads[start + 2 * n_curls][1][0][1] =
589  -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
590  grad_grads[start + 2 * n_curls][1][0][2] =
591  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
592  grad_grads[start + 2 * n_curls][1][1][0] =
593  -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
594  grad_grads[start + 2 * n_curls][1][1][1] =
595  -monoval_j[0][0] * monoval_plus[1][2] * monoval_i[2][1];
596  grad_grads[start + 2 * n_curls][1][1][2] =
597  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
598  grad_grads[start + 2 * n_curls][1][2][0] =
599  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
600  grad_grads[start + 2 * n_curls][1][2][1] =
601  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
602  grad_grads[start + 2 * n_curls][1][2][2] =
603  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][3];
604  grad_grads[start + 2 * n_curls][2][0][0] =
605  monoval_j[0][2] * monoval[1][0] * monoval_i[2][0] *
606  static_cast<double>(j + my_degree + 2);
607  grad_grads[start + 2 * n_curls][2][0][1] =
608  monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
609  static_cast<double>(j + my_degree + 2);
610  grad_grads[start + 2 * n_curls][2][0][2] =
611  monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
612  static_cast<double>(j + my_degree + 2);
613  grad_grads[start + 2 * n_curls][2][1][0] =
614  monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
615  static_cast<double>(j + my_degree + 2);
616  grad_grads[start + 2 * n_curls][2][1][1] =
617  monoval_j[0][0] * monoval[1][2] * monoval_i[2][0] *
618  static_cast<double>(j + my_degree + 2);
619  grad_grads[start + 2 * n_curls][2][1][2] =
620  monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
621  static_cast<double>(j + my_degree + 2);
622  grad_grads[start + 2 * n_curls][2][2][0] =
623  monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
624  static_cast<double>(j + my_degree + 2);
625  grad_grads[start + 2 * n_curls][2][2][1] =
626  monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
627  static_cast<double>(j + my_degree + 2);
628  grad_grads[start + 2 * n_curls][2][2][2] =
629  monoval_j[0][0] * monoval[1][0] * monoval_i[2][2] *
630  static_cast<double>(j + my_degree + 2);
631 
632  if (j != my_degree)
633  {
634  grad_grads[start + 1][0][0][0] =
635  monoval_i[0][2] * monoval[1][0] * monoval_j[2][0] *
636  static_cast<double>(j + my_degree + 2);
637  grad_grads[start + 1][0][0][1] =
638  monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
639  static_cast<double>(j + my_degree + 2);
640  grad_grads[start + 1][0][0][2] =
641  monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
642  static_cast<double>(j + my_degree + 2);
643  grad_grads[start + 1][0][1][0] =
644  monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
645  static_cast<double>(j + my_degree + 2);
646  grad_grads[start + 1][0][1][1] =
647  monoval_i[0][0] * monoval[1][2] * monoval_j[2][0] *
648  static_cast<double>(j + my_degree + 2);
649  grad_grads[start + 1][0][1][2] =
650  monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
651  static_cast<double>(j + my_degree + 2);
652  grad_grads[start + 1][0][2][0] =
653  monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
654  static_cast<double>(j + my_degree + 2);
655  grad_grads[start + 1][0][2][1] =
656  monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
657  static_cast<double>(j + my_degree + 2);
658  grad_grads[start + 1][0][2][2] =
659  monoval_i[0][0] * monoval[1][0] * monoval_j[2][2] *
660  static_cast<double>(j + my_degree + 2);
661  grad_grads[start + 1][1][0][0] =
662  -monoval_i[0][3] * monoval_plus[1][0] * monoval_j[2][0];
663  grad_grads[start + 1][1][0][1] =
664  -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
665  grad_grads[start + 1][1][0][2] =
666  -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
667  grad_grads[start + 1][1][1][0] =
668  -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
669  grad_grads[start + 1][1][1][1] =
670  -monoval_i[0][1] * monoval_plus[1][2] * monoval_j[2][0];
671  grad_grads[start + 1][1][1][2] =
672  -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
673  grad_grads[start + 1][1][2][0] =
674  -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
675  grad_grads[start + 1][1][2][1] =
676  -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
677  grad_grads[start + 1][1][2][2] =
678  -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][2];
679  grad_grads[start + 1][2][0][0] =
680  -monoval_i[0][3] * monoval[1][0] * monoval_jplus[2][0];
681  grad_grads[start + 1][2][0][1] =
682  -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
683  grad_grads[start + 1][2][0][2] =
684  -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
685  grad_grads[start + 1][2][1][0] =
686  -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
687  grad_grads[start + 1][2][1][1] =
688  -monoval_i[0][1] * monoval[1][2] * monoval_jplus[2][0];
689  grad_grads[start + 1][2][1][2] =
690  -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
691  grad_grads[start + 1][2][2][0] =
692  -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
693  grad_grads[start + 1][2][2][1] =
694  -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
695  grad_grads[start + 1][2][2][2] =
696  -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][2];
697 
698  grad_grads[start + n_curls + 1][0][0][0] =
699  -monoval_plus[0][2] * monoval_i[1][1] * monoval_j[2][0];
700  grad_grads[start + n_curls + 1][0][0][1] =
701  -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
702  grad_grads[start + n_curls + 1][0][0][2] =
703  -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
704  grad_grads[start + n_curls + 1][0][1][0] =
705  -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
706  grad_grads[start + n_curls + 1][0][1][1] =
707  -monoval_plus[0][0] * monoval_i[1][3] * monoval_j[2][0];
708  grad_grads[start + n_curls + 1][0][1][2] =
709  -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
710  grad_grads[start + n_curls + 1][0][2][0] =
711  -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
712  grad_grads[start + n_curls + 1][0][2][1] =
713  -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
714  grad_grads[start + n_curls + 1][0][2][2] =
715  -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][2];
716  grad_grads[start + n_curls + 1][1][0][0] =
717  monoval[0][2] * monoval_i[1][0] * monoval_j[2][0] *
718  static_cast<double>(j + my_degree + 2);
719  grad_grads[start + n_curls + 1][1][0][1] =
720  monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
721  static_cast<double>(j + my_degree + 2);
722  grad_grads[start + n_curls + 1][1][0][2] =
723  monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
724  static_cast<double>(j + my_degree + 2);
725  grad_grads[start + n_curls + 1][1][1][0] =
726  monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
727  static_cast<double>(j + my_degree + 2);
728  grad_grads[start + n_curls + 1][1][1][1] =
729  monoval[0][0] * monoval_i[1][2] * monoval_j[2][0] *
730  static_cast<double>(j + my_degree + 2);
731  grad_grads[start + n_curls + 1][1][1][2] =
732  monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
733  static_cast<double>(j + my_degree + 2);
734  grad_grads[start + n_curls + 1][1][2][0] =
735  monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
736  static_cast<double>(j + my_degree + 2);
737  grad_grads[start + n_curls + 1][1][2][1] =
738  monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
739  static_cast<double>(j + my_degree + 2);
740  grad_grads[start + n_curls + 1][1][2][2] =
741  monoval[0][0] * monoval_i[1][0] * monoval_j[2][2] *
742  static_cast<double>(j + my_degree + 2);
743  grad_grads[start + n_curls + 1][2][0][0] =
744  -monoval[0][2] * monoval_i[1][1] * monoval_jplus[2][0];
745  grad_grads[start + n_curls + 1][2][0][1] =
746  -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
747  grad_grads[start + n_curls + 1][2][0][2] =
748  -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
749  grad_grads[start + n_curls + 1][2][1][0] =
750  -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
751  grad_grads[start + n_curls + 1][2][1][1] =
752  -monoval[0][0] * monoval_i[1][3] * monoval_jplus[2][0];
753  grad_grads[start + n_curls + 1][2][1][2] =
754  -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
755  grad_grads[start + n_curls + 1][2][2][0] =
756  -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
757  grad_grads[start + n_curls + 1][2][2][1] =
758  -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
759  grad_grads[start + n_curls + 1][2][2][2] =
760  -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][2];
761 
762  grad_grads[start + 2 * n_curls + 1][0][0][0] =
763  -monoval_plus[0][2] * monoval_j[1][0] * monoval_i[2][1];
764  grad_grads[start + 2 * n_curls + 1][0][0][1] =
765  -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
766  grad_grads[start + 2 * n_curls + 1][0][0][2] =
767  -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
768  grad_grads[start + 2 * n_curls + 1][0][1][0] =
769  -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
770  grad_grads[start + 2 * n_curls + 1][0][1][1] =
771  -monoval_plus[0][0] * monoval_j[1][2] * monoval_i[2][1];
772  grad_grads[start + 2 * n_curls + 1][0][1][2] =
773  -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
774  grad_grads[start + 2 * n_curls + 1][0][2][0] =
775  -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
776  grad_grads[start + 2 * n_curls + 1][0][2][1] =
777  -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
778  grad_grads[start + 2 * n_curls + 1][0][2][2] =
779  -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][3];
780  grad_grads[start + 2 * n_curls + 1][1][0][0] =
781  -monoval[0][2] * monoval_jplus[1][0] * monoval_i[2][1];
782  grad_grads[start + 2 * n_curls + 1][1][0][1] =
783  -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
784  grad_grads[start + 2 * n_curls + 1][1][0][2] =
785  -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
786  grad_grads[start + 2 * n_curls + 1][1][1][0] =
787  -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
788  grad_grads[start + 2 * n_curls + 1][1][1][1] =
789  -monoval[0][0] * monoval_jplus[1][2] * monoval_i[2][1];
790  grad_grads[start + 2 * n_curls + 1][1][1][2] =
791  -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
792  grad_grads[start + 2 * n_curls + 1][1][2][0] =
793  -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
794  grad_grads[start + 2 * n_curls + 1][1][2][1] =
795  -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
796  grad_grads[start + 2 * n_curls + 1][1][2][2] =
797  -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][3];
798  grad_grads[start + 2 * n_curls + 1][2][0][0] =
799  monoval[0][2] * monoval_j[1][0] * monoval_i[2][0] *
800  static_cast<double>(j + my_degree + 2);
801  grad_grads[start + 2 * n_curls + 1][2][0][1] =
802  monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
803  static_cast<double>(j + my_degree + 2);
804  grad_grads[start + 2 * n_curls + 1][2][0][2] =
805  monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
806  static_cast<double>(j + my_degree + 2);
807  grad_grads[start + 2 * n_curls + 1][2][1][0] =
808  monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
809  static_cast<double>(j + my_degree + 2);
810  grad_grads[start + 2 * n_curls + 1][2][1][1] =
811  monoval[0][0] * monoval_j[1][2] * monoval_i[2][0] *
812  static_cast<double>(j + my_degree + 2);
813  grad_grads[start + 2 * n_curls + 1][2][1][2] =
814  monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
815  static_cast<double>(j + my_degree + 2);
816  grad_grads[start + 2 * n_curls + 1][2][2][0] =
817  monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
818  static_cast<double>(j + my_degree + 2);
819  grad_grads[start + 2 * n_curls + 1][2][2][1] =
820  monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
821  static_cast<double>(j + my_degree + 2);
822  grad_grads[start + 2 * n_curls + 1][2][2][2] =
823  monoval[0][0] * monoval_j[1][0] * monoval_i[2][2] *
824  static_cast<double>(j + my_degree + 2);
825  }
826  }
827 
828  if (j == my_degree)
829  start += 1;
830  else
831  start += 2;
832  }
833  }
834  Assert(start == this->n() - 2 * n_curls, ExcInternalError());
835  }
836 }
837 
838 
839 
840 template <int dim>
841 unsigned int
843 {
844  if (dim == 1 || dim == 2 || dim == 3)
845  return dim * Utilities::fixed_power<dim>(k + 1);
846 
847  Assert(false, ExcNotImplemented());
848  return 0;
849 }
850 
851 
852 template <int dim>
853 std::unique_ptr<TensorPolynomialsBase<dim>>
855 {
856  return std::make_unique<PolynomialsRT_Bubbles<dim>>(*this);
857 }
858 
859 
860 template class PolynomialsRT_Bubbles<1>;
861 template class PolynomialsRT_Bubbles<2>;
862 template class PolynomialsRT_Bubbles<3>;
863 
864 
Definition: point.h:112
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_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1616
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)