Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20: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\}}\)
Loading...
Searching...
No Matches
polynomials_rannacher_turek.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2015 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
16#ifndef dealii_polynomials_rannacher_turek_h
17#define dealii_polynomials_rannacher_turek_h
18
19#include <deal.II/base/config.h>
20
21#include <deal.II/base/point.h>
23#include <deal.II/base/tensor.h>
24
25#include <vector>
26
28
29
40template <int dim>
42{
43public:
47 static constexpr unsigned int dimension = dim;
48
53
57 double
58 compute_value(const unsigned int i, const Point<dim> &p) const override;
59
65 template <int order>
67 compute_derivative(const unsigned int i, const Point<dim> &p) const;
68
72 virtual Tensor<1, dim>
73 compute_1st_derivative(const unsigned int i,
74 const Point<dim> &p) const override;
75
79 virtual Tensor<2, dim>
80 compute_2nd_derivative(const unsigned int i,
81 const Point<dim> &p) const override;
82
86 virtual Tensor<3, dim>
87 compute_3rd_derivative(const unsigned int i,
88 const Point<dim> &p) const override;
89
93 virtual Tensor<4, dim>
94 compute_4th_derivative(const unsigned int i,
95 const Point<dim> &p) const override;
96
101 compute_grad(const unsigned int i, const Point<dim> &p) const override;
102
107 compute_grad_grad(const unsigned int i, const Point<dim> &p) const override;
108
115 void
116 evaluate(const Point<dim> &unit_point,
117 std::vector<double> &values,
118 std::vector<Tensor<1, dim>> &grads,
119 std::vector<Tensor<2, dim>> &grad_grads,
120 std::vector<Tensor<3, dim>> &third_derivatives,
121 std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
122
126 std::string
127 name() const override;
128
132 virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
133 clone() const override;
134};
135
136
137namespace internal
138{
139 namespace PolynomialsRannacherTurekImplementation
140 {
141 template <int order, int dim>
142 inline Tensor<order, dim>
143 compute_derivative(const unsigned int, const Point<dim> &)
144 {
145 Assert(dim == 2, ExcNotImplemented());
146 return Tensor<order, dim>();
147 }
148
149
150 template <int order>
151 inline Tensor<order, 2>
152 compute_derivative(const unsigned int i, const Point<2> &p)
153 {
154 const unsigned int dim = 2;
155
156 Tensor<order, dim> derivative;
157 switch (order)
158 {
159 case 1:
160 {
161 Tensor<1, dim> &grad =
162 *reinterpret_cast<Tensor<1, dim> *>(&derivative);
163 if (i == 0)
164 {
165 grad[0] = -2.5 + 3 * p[0];
166 grad[1] = 1.5 - 3 * p[1];
167 }
168 else if (i == 1)
169 {
170 grad[0] = -0.5 + 3.0 * p[0];
171 grad[1] = 1.5 - 3.0 * p[1];
172 }
173 else if (i == 2)
174 {
175 grad[0] = 1.5 - 3.0 * p[0];
176 grad[1] = -2.5 + 3.0 * p[1];
177 }
178 else if (i == 3)
179 {
180 grad[0] = 1.5 - 3.0 * p[0];
181 grad[1] = -0.5 + 3.0 * p[1];
182 }
183 else
184 {
186 }
187 return derivative;
188 }
189 case 2:
190 {
191 Tensor<2, dim> &grad_grad =
192 *reinterpret_cast<Tensor<2, dim> *>(&derivative);
193 if (i == 0)
194 {
195 grad_grad[0][0] = 3;
196 grad_grad[0][1] = 0;
197 grad_grad[1][0] = 0;
198 grad_grad[1][1] = -3;
199 }
200 else if (i == 1)
201 {
202 grad_grad[0][0] = 3;
203 grad_grad[0][1] = 0;
204 grad_grad[1][0] = 0;
205 grad_grad[1][1] = -3;
206 }
207 else if (i == 2)
208 {
209 grad_grad[0][0] = -3;
210 grad_grad[0][1] = 0;
211 grad_grad[1][0] = 0;
212 grad_grad[1][1] = 3;
213 }
214 else if (i == 3)
215 {
216 grad_grad[0][0] = -3;
217 grad_grad[0][1] = 0;
218 grad_grad[1][0] = 0;
219 grad_grad[1][1] = 3;
220 }
221 return derivative;
222 }
223 default:
224 {
225 // higher derivatives are all zero
226 return Tensor<order, dim>();
227 }
228 }
229 }
230 } // namespace PolynomialsRannacherTurekImplementation
231} // namespace internal
232
233
234
235// template functions
236template <int dim>
237template <int order>
245
246
247
248template <int dim>
249inline Tensor<1, dim>
251 const unsigned int i,
252 const Point<dim> &p) const
253{
254 return compute_derivative<1>(i, p);
255}
256
257
258
259template <int dim>
260inline Tensor<2, dim>
262 const unsigned int i,
263 const Point<dim> &p) const
264{
265 return compute_derivative<2>(i, p);
266}
267
268
269
270template <int dim>
271inline Tensor<3, dim>
273 const unsigned int i,
274 const Point<dim> &p) const
275{
276 return compute_derivative<3>(i, p);
277}
278
279
280
281template <int dim>
282inline Tensor<4, dim>
284 const unsigned int i,
285 const Point<dim> &p) const
286{
287 return compute_derivative<4>(i, p);
288}
289
290
291
292template <int dim>
293inline std::string
295{
296 return "RannacherTurek";
297}
298
299
301
302#endif
Definition point.h:111
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
virtual std::unique_ptr< ScalarPolynomialsBase< dim > > clone() const override
Tensor< order, dim > compute_derivative(const unsigned int i, const Point< dim > &p) const
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
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
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define DEAL_II_NOT_IMPLEMENTED()
Tensor< order, dim > compute_derivative(const unsigned int, const Point< dim > &)