Reference documentation for deal.II version GIT relicensing-437-g81ec864850 2024-04-19 07:30: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
scalar_polynomials_base.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 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#ifndef dealii_scalar_polynomials_base_h
16#define dealii_scalar_polynomials_base_h
17
18
19#include <deal.II/base/config.h>
20
22#include <deal.II/base/point.h>
23#include <deal.II/base/tensor.h>
24
25#include <memory>
26#include <string>
27#include <vector>
28
30
61template <int dim>
63{
64public:
69 ScalarPolynomialsBase(const unsigned int deg,
70 const unsigned int n_polynomials);
71
76
81
86 virtual ~ScalarPolynomialsBase() = default;
87
100 virtual void
101 evaluate(const Point<dim> &unit_point,
102 std::vector<double> &values,
103 std::vector<Tensor<1, dim>> &grads,
104 std::vector<Tensor<2, dim>> &grad_grads,
105 std::vector<Tensor<3, dim>> &third_derivatives,
106 std::vector<Tensor<4, dim>> &fourth_derivatives) const = 0;
107
114 virtual double
115 compute_value(const unsigned int i, const Point<dim> &p) const = 0;
116
125 template <int order>
127 compute_derivative(const unsigned int i, const Point<dim> &p) const;
128
135 virtual Tensor<1, dim>
136 compute_1st_derivative(const unsigned int i, const Point<dim> &p) const = 0;
137
144 virtual Tensor<2, dim>
145 compute_2nd_derivative(const unsigned int i, const Point<dim> &p) const = 0;
146
153 virtual Tensor<3, dim>
154 compute_3rd_derivative(const unsigned int i, const Point<dim> &p) const = 0;
155
162 virtual Tensor<4, dim>
163 compute_4th_derivative(const unsigned int i, const Point<dim> &p) const = 0;
164
171 virtual Tensor<1, dim>
172 compute_grad(const unsigned int /*i*/, const Point<dim> & /*p*/) const = 0;
173
180 virtual Tensor<2, dim>
181 compute_grad_grad(const unsigned int /*i*/,
182 const Point<dim> & /*p*/) const = 0;
183
187 unsigned int
188 n() const;
189
195 virtual unsigned int
196 degree() const;
197
208 virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
209 clone() const = 0;
210
214 virtual std::string
215 name() const = 0;
216
220 virtual std::size_t
221 memory_consumption() const;
222
223private:
227 const unsigned int polynomial_degree;
228
232 const unsigned int n_pols;
233};
234
235
236
237template <int dim>
238inline unsigned int
240{
241 return n_pols;
242}
243
244
245
246template <int dim>
247inline unsigned int
249{
250 return polynomial_degree;
251}
252
253
254
255template <int dim>
256template <int order>
259 const Point<dim> &p) const
260{
261 if (order == 1)
262 {
263 auto derivative = compute_1st_derivative(i, p);
264 return *reinterpret_cast<Tensor<order, dim> *>(&derivative);
265 }
266 if (order == 2)
267 {
268 auto derivative = compute_2nd_derivative(i, p);
269 return *reinterpret_cast<Tensor<order, dim> *>(&derivative);
270 }
271 if (order == 3)
272 {
273 auto derivative = compute_3rd_derivative(i, p);
274 return *reinterpret_cast<Tensor<order, dim> *>(&derivative);
275 }
276 if (order == 4)
277 {
278 auto derivative = compute_4th_derivative(i, p);
279 return *reinterpret_cast<Tensor<order, dim> *>(&derivative);
280 }
282 Tensor<order, dim> empty;
283 return empty;
284}
285
287
288#endif
Definition point.h:111
ScalarPolynomialsBase(ScalarPolynomialsBase< dim > &&)=default
virtual Tensor< 2, dim > compute_grad_grad(const unsigned int, const Point< dim > &) const =0
virtual Tensor< 4, dim > compute_4th_derivative(const unsigned int i, const Point< dim > &p) const =0
virtual Tensor< 1, dim > compute_grad(const unsigned int, const Point< dim > &) const =0
virtual std::unique_ptr< ScalarPolynomialsBase< dim > > clone() const =0
virtual std::string name() const =0
virtual ~ScalarPolynomialsBase()=default
virtual 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 =0
virtual Tensor< 3, dim > compute_3rd_derivative(const unsigned int i, const Point< dim > &p) const =0
virtual Tensor< 1, dim > compute_1st_derivative(const unsigned int i, const Point< dim > &p) const =0
virtual unsigned int degree() const
const unsigned int polynomial_degree
virtual double compute_value(const unsigned int i, const Point< dim > &p) const =0
virtual Tensor< 2, dim > compute_2nd_derivative(const unsigned int i, const Point< dim > &p) const =0
virtual std::size_t memory_consumption() const
ScalarPolynomialsBase(const ScalarPolynomialsBase< dim > &)=default
Tensor< order, dim > compute_derivative(const unsigned int i, const Point< dim > &p) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define DEAL_II_NOT_IMPLEMENTED()