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