Reference documentation for deal.II version GIT f6a5d312c9 2023-10-04 08:50: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_rannacher_turek.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 2022 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 
17 #ifndef dealii_polynomials_rannacher_turek_h
18 #define dealii_polynomials_rannacher_turek_h
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/point.h>
24 #include <deal.II/base/tensor.h>
25 
26 #include <vector>
27 
29 
30 
41 template <int dim>
43 {
44 public:
48  static constexpr unsigned int dimension = dim;
49 
54 
58  double
59  compute_value(const unsigned int i, const Point<dim> &p) const override;
60 
66  template <int order>
68  compute_derivative(const unsigned int i, const Point<dim> &p) const;
69 
73  virtual Tensor<1, dim>
74  compute_1st_derivative(const unsigned int i,
75  const Point<dim> &p) const override;
76 
80  virtual Tensor<2, dim>
81  compute_2nd_derivative(const unsigned int i,
82  const Point<dim> &p) const override;
83 
87  virtual Tensor<3, dim>
88  compute_3rd_derivative(const unsigned int i,
89  const Point<dim> &p) const override;
90 
94  virtual Tensor<4, dim>
95  compute_4th_derivative(const unsigned int i,
96  const Point<dim> &p) const override;
97 
102  compute_grad(const unsigned int i, const Point<dim> &p) const override;
103 
108  compute_grad_grad(const unsigned int i, const Point<dim> &p) const override;
109 
116  void
117  evaluate(const Point<dim> &unit_point,
118  std::vector<double> &values,
119  std::vector<Tensor<1, dim>> &grads,
120  std::vector<Tensor<2, dim>> &grad_grads,
121  std::vector<Tensor<3, dim>> &third_derivatives,
122  std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
123 
127  std::string
128  name() const override;
129 
133  virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
134  clone() const override;
135 };
136 
137 
138 namespace internal
139 {
140  namespace PolynomialsRannacherTurekImplementation
141  {
142  template <int order, int dim>
143  inline Tensor<order, dim>
144  compute_derivative(const unsigned int, const Point<dim> &)
145  {
146  Assert(dim == 2, ExcNotImplemented());
147  return Tensor<order, dim>();
148  }
149 
150 
151  template <int order>
152  inline Tensor<order, 2>
153  compute_derivative(const unsigned int i, const Point<2> &p)
154  {
155  const unsigned int dim = 2;
156 
157  Tensor<order, dim> derivative;
158  switch (order)
159  {
160  case 1:
161  {
162  Tensor<1, dim> &grad =
163  *reinterpret_cast<Tensor<1, dim> *>(&derivative);
164  if (i == 0)
165  {
166  grad[0] = -2.5 + 3 * p(0);
167  grad[1] = 1.5 - 3 * p(1);
168  }
169  else if (i == 1)
170  {
171  grad[0] = -0.5 + 3.0 * p(0);
172  grad[1] = 1.5 - 3.0 * p(1);
173  }
174  else if (i == 2)
175  {
176  grad[0] = 1.5 - 3.0 * p(0);
177  grad[1] = -2.5 + 3.0 * p(1);
178  }
179  else if (i == 3)
180  {
181  grad[0] = 1.5 - 3.0 * p(0);
182  grad[1] = -0.5 + 3.0 * p(1);
183  }
184  else
185  {
186  Assert(false, ExcNotImplemented());
187  }
188  return derivative;
189  }
190  case 2:
191  {
192  Tensor<2, dim> &grad_grad =
193  *reinterpret_cast<Tensor<2, dim> *>(&derivative);
194  if (i == 0)
195  {
196  grad_grad[0][0] = 3;
197  grad_grad[0][1] = 0;
198  grad_grad[1][0] = 0;
199  grad_grad[1][1] = -3;
200  }
201  else if (i == 1)
202  {
203  grad_grad[0][0] = 3;
204  grad_grad[0][1] = 0;
205  grad_grad[1][0] = 0;
206  grad_grad[1][1] = -3;
207  }
208  else if (i == 2)
209  {
210  grad_grad[0][0] = -3;
211  grad_grad[0][1] = 0;
212  grad_grad[1][0] = 0;
213  grad_grad[1][1] = 3;
214  }
215  else if (i == 3)
216  {
217  grad_grad[0][0] = -3;
218  grad_grad[0][1] = 0;
219  grad_grad[1][0] = 0;
220  grad_grad[1][1] = 3;
221  }
222  return derivative;
223  }
224  default:
225  {
226  // higher derivatives are all zero
227  return Tensor<order, dim>();
228  }
229  }
230  }
231  } // namespace PolynomialsRannacherTurekImplementation
232 } // namespace internal
233 
234 
235 
236 // template functions
237 template <int dim>
238 template <int order>
241  const Point<dim> &p) const
242 {
244  order>(i, p);
245 }
246 
247 
248 
249 template <int dim>
250 inline Tensor<1, dim>
252  const unsigned int i,
253  const Point<dim> &p) const
254 {
255  return compute_derivative<1>(i, p);
256 }
257 
258 
259 
260 template <int dim>
261 inline Tensor<2, dim>
263  const unsigned int i,
264  const Point<dim> &p) const
265 {
266  return compute_derivative<2>(i, p);
267 }
268 
269 
270 
271 template <int dim>
272 inline Tensor<3, dim>
274  const unsigned int i,
275  const Point<dim> &p) const
276 {
277  return compute_derivative<3>(i, p);
278 }
279 
280 
281 
282 template <int dim>
283 inline Tensor<4, dim>
285  const unsigned int i,
286  const Point<dim> &p) const
287 {
288  return compute_derivative<4>(i, p);
289 }
290 
291 
292 
293 template <int dim>
294 inline std::string
296 {
297  return "RannacherTurek";
298 }
299 
300 
302 
303 #endif
Definition: point.h:112
virtual Tensor< 3, dim > compute_3rd_derivative(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 1, dim > compute_1st_derivative(const unsigned int i, const Point< dim > &p) const override
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
virtual std::unique_ptr< ScalarPolynomialsBase< dim > > clone() const override
Tensor< order, dim > compute_derivative(const unsigned int i, const Point< dim > &p) const
std::string name() const override
Tensor< 2, dim > compute_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 4, dim > compute_4th_derivative(const unsigned int i, const Point< dim > &p) const override
Tensor< 1, dim > compute_grad(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 2, dim > compute_2nd_derivative(const unsigned int i, const Point< dim > &p) const override
static constexpr unsigned int dimension
double compute_value(const unsigned int i, const Point< dim > &p) const override
Definition: tensor.h:516
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcNotImplemented()
Tensor< order, dim > compute_derivative(const unsigned int, const Point< dim > &)