Reference documentation for deal.II version 9.2.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\}}\)
polynomials_adini.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 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 
16 #include <deal.II/base/exceptions.h>
17 #include <deal.II/base/polynomials_adini.h>
18 #include <deal.II/base/std_cxx14/memory.h>
19 
20 #define ENTER_COEFFICIENTS( \
21  koefs, z, a0, a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, a11) \
22  koefs(0, z) = a0; \
23  koefs(1, z) = a1; \
24  koefs(2, z) = a2; \
25  koefs(3, z) = a3; \
26  koefs(4, z) = a4; \
27  koefs(5, z) = a5; \
28  koefs(6, z) = a6; \
29  koefs(7, z) = a7; \
30  koefs(8, z) = a8; \
31  koefs(9, z) = a9; \
32  koefs(10, z) = a10; \
33  koefs(11, z) = a11;
34 
35 
36 DEAL_II_NAMESPACE_OPEN
37 
38 
39 
40 template <int dim>
41 PolynomialsAdini<dim>::PolynomialsAdini()
42  : ScalarPolynomialsBase<dim>(3, 12)
43  , coef(12, 12)
44  , dx(12, 12)
45  , dy(12, 12)
46  , dxx(12, 12)
47  , dyy(12, 12)
48  , dxy(12, 12)
49 {
50  Assert(dim == 2, ExcNotImplemented());
51 
52  // 1 x y xx yy xy 3x 3y xyy xxy 3xy x3y
53  // 0 1 2 3 4 5 6 7 8 9 10 11
54  ENTER_COEFFICIENTS(coef, 0, 1, 0, 0, -3, -3, -1, 2, 2, 3, 3, -2, -2);
55  ENTER_COEFFICIENTS(coef, 1, 0, 1, 0, -2, 0, -1, 1, 0, 0, 2, -1, 0);
56  ENTER_COEFFICIENTS(coef, 2, 0, 0, 1, 0, -2, -1, 0, 1, 2, 0, 0, -1);
57  ENTER_COEFFICIENTS(coef, 3, 0, 0, 0, 3, 0, 1, -2, 0, -3, -3, 2, 2);
58  ENTER_COEFFICIENTS(coef, 4, 0, 0, 0, -1, 0, 0, 1, 0, 0, 1, -1, 0);
59  ENTER_COEFFICIENTS(coef, 5, 0, 0, 0, 0, 0, 1, 0, 0, -2, 0, 0, 1);
60  ENTER_COEFFICIENTS(coef, 6, 0, 0, 0, 0, 3, 1, 0, -2, -3, -3, 2, 2);
61  ENTER_COEFFICIENTS(coef, 7, 0, 0, 0, 0, 0, 1, 0, 0, 0, -2, 1, 0);
62  ENTER_COEFFICIENTS(coef, 8, 0, 0, 0, 0, -1, 0, 0, 1, 1, 0, 0, -1);
63  ENTER_COEFFICIENTS(coef, 9, 0, 0, 0, 0, 0, -1, 0, 0, 3, 3, -2, -2);
64  ENTER_COEFFICIENTS(coef, 10, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 1, 0);
65  ENTER_COEFFICIENTS(coef, 11, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 1);
66 
67  ENTER_COEFFICIENTS(dx, 0, 0, -6, -1, 6, 3, 6, 0, -2, 0, -6, 0, 0);
68  ENTER_COEFFICIENTS(dx, 1, 1, -4, -1, 3, 0, 4, 0, 0, 0, -3, 0, 0);
69  ENTER_COEFFICIENTS(dx, 2, 0, 0, -1, 0, 2, 0, 0, -1, 0, 0, 0, 0);
70  ENTER_COEFFICIENTS(dx, 3, 0, 6, 1, -6, -3, -6, 0, 2, 0, 6, 0, 0);
71  ENTER_COEFFICIENTS(dx, 4, 0, -2, 0, 3, 0, 2, 0, 0, 0, -3, 0, 0);
72  ENTER_COEFFICIENTS(dx, 5, 0, 0, 1, 0, -2, 0, 0, 1, 0, 0, 0, 0);
73  ENTER_COEFFICIENTS(dx, 6, 0, 0, 1, 0, -3, -6, 0, 2, 0, 6, 0, 0);
74  ENTER_COEFFICIENTS(dx, 7, 0, 0, 1, 0, 0, -4, 0, 0, 0, 3, 0, 0);
75  ENTER_COEFFICIENTS(dx, 8, 0, 0, 0, 0, 1, 0, 0, -1, 0, 0, 0, 0);
76  ENTER_COEFFICIENTS(dx, 9, 0, 0, -1, 0, 3, 6, 0, -2, 0, -6, 0, 0);
77  ENTER_COEFFICIENTS(dx, 10, 0, 0, 0, 0, 0, -2, 0, 0, 0, 3, 0, 0);
78  ENTER_COEFFICIENTS(dx, 11, 0, 0, 0, 0, -1, 0, 0, 1, 0, 0, 0, 0);
79 
80  ENTER_COEFFICIENTS(dy, 0, 0, -1, -6, 3, 6, 6, -2, 0, -6, 0, 0, 0);
81  ENTER_COEFFICIENTS(dy, 1, 0, -1, 0, 2, 0, 0, -1, 0, 0, 0, 0, 0);
82  ENTER_COEFFICIENTS(dy, 2, 1, -1, -4, 0, 3, 4, 0, 0, -3, 0, 0, 0);
83  ENTER_COEFFICIENTS(dy, 3, 0, 1, 0, -3, 0, -6, 2, 0, 6, 0, 0, 0);
84  ENTER_COEFFICIENTS(dy, 4, 0, 0, 0, 1, 0, 0, -1, 0, 0, 0, 0, 0);
85  ENTER_COEFFICIENTS(dy, 5, 0, 1, 0, 0, 0, -4, 0, 0, 3, 0, 0, 0);
86  ENTER_COEFFICIENTS(dy, 6, 0, 1, 6, -3, -6, -6, 2, 0, 6, 0, 0, 0);
87  ENTER_COEFFICIENTS(dy, 7, 0, 1, 0, -2, 0, 0, 1, 0, 0, 0, 0, 0);
88  ENTER_COEFFICIENTS(dy, 8, 0, 0, -2, 0, 3, 2, 0, 0, -3, 0, 0, 0);
89  ENTER_COEFFICIENTS(dy, 9, 0, -1, 0, 3, 0, 6, -2, 0, -6, 0, 0, 0);
90  ENTER_COEFFICIENTS(dy, 10, 0, 0, 0, -1, 0, 0, 1, 0, 0, 0, 0, 0);
91  ENTER_COEFFICIENTS(dy, 11, 0, 0, 0, 0, 0, -2, 0, 0, 3, 0, 0, 0);
92 
93  ENTER_COEFFICIENTS(dxx, 0, -6, 12, 6, 0, 0, -12, 0, 0, 0, 0, 0, 0);
94  ENTER_COEFFICIENTS(dxx, 1, -4, 6, 4, 0, 0, -6, 0, 0, 0, 0, 0, 0);
95  ENTER_COEFFICIENTS(dxx, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
96  ENTER_COEFFICIENTS(dxx, 3, 6, -12, -6, 0, 0, 12, 0, 0, 0, 0, 0, 0);
97  ENTER_COEFFICIENTS(dxx, 4, -2, 6, 2, 0, 0, -6, 0, 0, 0, 0, 0, 0);
98  ENTER_COEFFICIENTS(dxx, 5, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
99  ENTER_COEFFICIENTS(dxx, 6, 0, 0, -6, 0, 0, 12, 0, 0, 0, 0, 0, 0);
100  ENTER_COEFFICIENTS(dxx, 7, 0, 0, -4, 0, 0, 6, 0, 0, 0, 0, 0, 0);
101  ENTER_COEFFICIENTS(dxx, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
102  ENTER_COEFFICIENTS(dxx, 9, 0, 0, 6, 0, 0, -12, 0, 0, 0, 0, 0, 0);
103  ENTER_COEFFICIENTS(dxx, 10, 0, 0, -2, 0, 0, 6, 0, 0, 0, 0, 0, 0);
104  ENTER_COEFFICIENTS(dxx, 11, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
105 
106  ENTER_COEFFICIENTS(dyy, 0, -6, 6, 12, 0, 0, -12, 0, 0, 0, 0, 0, 0);
107  ENTER_COEFFICIENTS(dyy, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
108  ENTER_COEFFICIENTS(dyy, 2, -4, 4, 6, 0, 0, -6, 0, 0, 0, 0, 0, 0);
109  ENTER_COEFFICIENTS(dyy, 3, 0, -6, 0, 0, 0, 12, 0, 0, 0, 0, 0, 0);
110  ENTER_COEFFICIENTS(dyy, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
111  ENTER_COEFFICIENTS(dyy, 5, 0, -4, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0);
112  ENTER_COEFFICIENTS(dyy, 6, 6, -6, -12, 0, 0, 12, 0, 0, 0, 0, 0, 0);
113  ENTER_COEFFICIENTS(dyy, 7, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
114  ENTER_COEFFICIENTS(dyy, 8, -2, 2, 6, 0, 0, -6, 0, -0, 0, 0, 0, 0);
115  ENTER_COEFFICIENTS(dyy, 9, 0, 6, 0, 0, 0, -12, 0, 0, 0, 0, 0, 0);
116  ENTER_COEFFICIENTS(dyy, 10, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
117  ENTER_COEFFICIENTS(dyy, 11, 0, -2, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0);
118 
119  ENTER_COEFFICIENTS(dxy, 0, -1, 6, 6, -6, -6, 0, 0, 0, 0, 0, 0, 0);
120  ENTER_COEFFICIENTS(dxy, 1, -1, 4, 0, -3, 0, 0, 0, 0, 0, 0, 0, 0);
121  ENTER_COEFFICIENTS(dxy, 2, -1, 0, 4, 0, -3, 0, 0, 0, 0, 0, 0, 0);
122  ENTER_COEFFICIENTS(dxy, 3, 1, -6, -6, 6, 6, 0, 0, 0, 0, 0, 0, 0);
123  ENTER_COEFFICIENTS(dxy, 4, 0, 2, 0, -3, 0, 0, 0, 0, 0, 0, 0, 0);
124  ENTER_COEFFICIENTS(dxy, 5, 1, 0, -4, 0, 3, 0, 0, 0, 0, 0, 0, 0);
125  ENTER_COEFFICIENTS(dxy, 6, 1, -6, -6, 6, 6, 0, 0, 0, 0, 0, 0, 0);
126  ENTER_COEFFICIENTS(dxy, 7, 1, -4, 0, 3, 0, 0, 0, 0, 0, 0, 0, 0);
127  ENTER_COEFFICIENTS(dxy, 8, 0, 0, 2, 0, -3, 0, 0, 0, 0, 0, 0, 0);
128  ENTER_COEFFICIENTS(dxy, 9, -1, 6, 6, -6, -6, 0, 0, 0, 0, 0, 0, 0);
129  ENTER_COEFFICIENTS(dxy, 10, 0, -2, 0, 3, 0, 0, 0, 0, 0, 0, 0, 0);
130  ENTER_COEFFICIENTS(dxy, 11, 0, 0, -2, 0, 3, 0, 0, 0, 0, 0, 0, 0);
131 }
132 
133 
134 
135 template <int dim>
136 void
137 PolynomialsAdini<dim>::evaluate(
138  const Point<dim> & unit_point,
139  std::vector<double> & values,
140  std::vector<Tensor<1, dim>> &grads,
141  std::vector<Tensor<2, dim>> &grad_grads,
142  std::vector<Tensor<3, dim>> &third_derivatives,
143  std::vector<Tensor<4, dim>> &fourth_derivatives) const
144 {
145  const unsigned int n_pols = this->n();
146  (void)n_pols;
147 
148  Assert(values.size() == n_pols || values.size() == 0,
149  ExcDimensionMismatch(values.size(), n_pols));
150  Assert(grads.size() == n_pols || grads.size() == 0,
151  ExcDimensionMismatch(grads.size(), n_pols));
152  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
153  ExcDimensionMismatch(grad_grads.size(), n_pols));
154  (void)third_derivatives;
155  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
156  ExcDimensionMismatch(third_derivatives.size(), n_pols));
157  (void)fourth_derivatives;
158  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
159  ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
160 
161  if (values.empty() == false) // do not bother if empty
162  {
163  for (unsigned int i = 0; i < values.size(); ++i)
164  {
165  values[i] = compute_value(i, unit_point);
166  }
167  }
168 
169  if (grads.empty() == false) // do not bother if empty
170  {
171  for (unsigned int i = 0; i < grads.size(); ++i)
172  {
173  grads[i] = compute_grad(i, unit_point);
174  }
175  }
176 
177  if (grad_grads.empty() == false) // do not bother if empty
178  {
179  for (unsigned int i = 0; i < grad_grads.size(); ++i)
180  {
181  grad_grads[i] = compute_grad_grad(i, unit_point);
182  }
183  }
184 
185  return;
186 }
187 
188 
189 
190 template <int dim>
191 double
192 PolynomialsAdini<dim>::compute_value(const unsigned int i,
193  const Point<dim> & p) const
194 {
195  const double x = p(0);
196  const double y = p(1);
197  return coef(0, i) + coef(1, i) * x + coef(2, i) * y + coef(3, i) * x * x +
198  coef(4, i) * y * y + coef(5, i) * x * y + coef(6, i) * x * x * x +
199  coef(7, i) * y * y * y + coef(8, i) * x * y * y +
200  coef(9, i) * x * x * y + coef(10, i) * x * x * x * y +
201  coef(11, i) * x * y * y * y;
202 }
203 
204 
205 
206 template <int dim>
207 Tensor<1, dim>
208 PolynomialsAdini<dim>::compute_grad(const unsigned int i,
209  const Point<dim> & p) const
210 {
211  const double x = p(0);
212  const double y = p(1);
213  Tensor<1, dim> tensor;
214  tensor[0] = dx(0, i) + dx(1, i) * x + dx(2, i) * y + dx(3, i) * x * x +
215  dx(4, i) * y * y + dx(5, i) * x * y + dx(6, i) * x * x * x +
216  dx(7, i) * y * y * y + dx(8, i) * x * y * y +
217  dx(9, i) * x * x * y + dx(10, i) * x * x * x * y +
218  dx(11, i) * x * y * y * y;
219 
220  tensor[1] = dy(0, i) + dy(1, i) * x + dy(2, i) * y + dy(3, i) * x * x +
221  dy(4, i) * y * y + dy(5, i) * x * y + dy(6, i) * x * x * x +
222  dy(7, i) * y * y * y + dy(8, i) * x * y * y +
223  dy(9, i) * x * x * y + dy(10, i) * x * x * x * y +
224  dy(11, i) * x * y * y * y;
225  return tensor;
226 }
227 
228 
229 
230 template <int dim>
231 Tensor<2, dim>
232 PolynomialsAdini<dim>::compute_grad_grad(const unsigned int i,
233  const Point<dim> & p) const
234 {
235  const double x = p(0);
236  const double y = p(1);
237  Tensor<2, dim> tensor;
238  tensor[0][0] = dxx(0, i) + dxx(1, i) * x + dxx(2, i) * y + dxx(3, i) * x * x +
239  dxx(4, i) * y * y + dxx(5, i) * x * y + dxx(6, i) * x * x * x +
240  dxx(7, i) * y * y * y + dxx(8, i) * x * y * y +
241  dxx(9, i) * x * x * y + dxx(10, i) * x * x * x * y +
242  dxx(11, i) * x * y * y * y;
243  tensor[0][1] = dxy(0, i) + dxy(1, i) * x + dxy(2, i) * y + dxy(3, i) * x * x +
244  dxy(4, i) * y * y + dxy(5, i) * x * y + dxy(6, i) * x * x * x +
245  dxy(7, i) * y * y * y + dxy(8, i) * x * y * y +
246  dxy(9, i) * x * x * y + dxy(10, i) * x * x * x * y +
247  dxy(11, i) * x * y * y * y;
248  tensor[1][0] = tensor[0][1];
249  tensor[1][1] = dyy(0, i) + dyy(1, i) * x + dyy(2, i) * y + dyy(3, i) * x * x +
250  dyy(4, i) * y * y + dyy(5, i) * x * y + dyy(6, i) * x * x * x +
251  dyy(7, i) * y * y * y + dyy(8, i) * x * y * y +
252  dyy(9, i) * x * x * y + dyy(10, i) * x * x * x * y +
253  dyy(11, i) * x * y * y * y;
254  return tensor;
255 }
256 
257 
258 
259 template <int dim>
260 std::unique_ptr<ScalarPolynomialsBase<dim>>
261 PolynomialsAdini<dim>::clone() const
262 {
263  return std_cxx14::make_unique<PolynomialsAdini<dim>>(*this);
264 }
265 
266 
267 
268 template class PolynomialsAdini<0>;
269 template class PolynomialsAdini<1>;
270 template class PolynomialsAdini<2>;
271 template class PolynomialsAdini<3>;
272 
273 DEAL_II_NAMESPACE_CLOSE
PolynomialsAdini::dy
Table< 2, double > dy
Definition: polynomials_adini.h:130
PolynomialsAdini::compute_value
double compute_value(const unsigned int i, const Point< dim > &p) const
Definition: polynomials_adini.cc:192
PolynomialsAdini::dxx
Table< 2, double > dxx
Definition: polynomials_adini.h:136
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
PolynomialsAdini
Definition: polynomials_adini.h:46
DataOutBase::dx
@ dx
Definition: data_out_base.h:1562
PolynomialsAdini::dxy
Table< 2, double > dxy
Definition: polynomials_adini.h:148
PolynomialsAdini::evaluate
void evaluate(const Point< dim > &unit_point, std::vector< double > &values, std::vector< Tensor< 1, dim >> &grads, std::vector< Tensor< 2, dim >> &grad_grads, std::vector< Tensor< 3, dim >> &third_derivatives, std::vector< Tensor< 4, dim >> &fourth_derivatives) const override
Definition: polynomials_adini.cc:137
ENTER_COEFFICIENTS
#define ENTER_COEFFICIENTS( koefs, z, a0, a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, a11)
Definition: polynomials_adini.cc:20
PolynomialsAdini::dyy
Table< 2, double > dyy
Definition: polynomials_adini.h:142
Tensor< 1, dim >
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
PolynomialsAdini::compute_grad_grad
Tensor< 2, dim > compute_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: polynomials_adini.cc:232
exceptions.h
PolynomialsAdini::compute_grad
Tensor< 1, dim > compute_grad(const unsigned int i, const Point< dim > &p) const
Definition: polynomials_adini.cc:208
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
polynomials_adini.h
PolynomialsAdini::PolynomialsAdini
PolynomialsAdini()
Definition: polynomials_adini.cc:41
PolynomialsAdini::coef
Table< 2, double > coef
Definition: polynomials_adini.h:118
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Point< dim >
PolynomialsAdini::clone
virtual std::unique_ptr< ScalarPolynomialsBase< dim > > clone() const override
Definition: polynomials_adini.cc:261
PolynomialsAdini::dx
Table< 2, double > dx
Definition: polynomials_adini.h:124
memory.h
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
ScalarPolynomialsBase
Definition: scalar_polynomials_base.h:63
Generated by   doxygen 1.8.17