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_wedge.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2021 - 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_base_polynomials_wedge_h
18#define dealii_base_polynomials_wedge_h
19
20#include <deal.II/base/config.h>
21
25
27
28
29namespace internal
30{
36 static const constexpr ::ndarray<unsigned int, 6, 2> wedge_table_1{
37 {{{0, 0}}, {{1, 0}}, {{2, 0}}, {{0, 1}}, {{1, 1}}, {{2, 1}}}};
38
44 static const constexpr ::ndarray<unsigned int, 18, 2> wedge_table_2{
45 {{{0, 0}},
46 {{1, 0}},
47 {{2, 0}},
48 {{0, 1}},
49 {{1, 1}},
50 {{2, 1}},
51 {{3, 0}},
52 {{4, 0}},
53 {{5, 0}},
54 {{3, 1}},
55 {{4, 1}},
56 {{5, 1}},
57 {{0, 2}},
58 {{1, 2}},
59 {{2, 2}},
60 {{3, 2}},
61 {{4, 2}},
62 {{5, 2}}}};
63} // namespace internal
64
65
75template <int dim>
77{
78public:
82 static constexpr unsigned int dimension = dim;
83
84 /*
85 * Constructor taking the polynomial @p degree as input.
86 *
87 * @note Currently, only linear (degree=1) and quadratic polynomials
88 * (degree=2) are implemented.
89 */
90 ScalarLagrangePolynomialWedge(const unsigned int degree);
91
97 void
98 evaluate(const Point<dim> & unit_point,
99 std::vector<double> & values,
100 std::vector<Tensor<1, dim>> &grads,
101 std::vector<Tensor<2, dim>> &grad_grads,
102 std::vector<Tensor<3, dim>> &third_derivatives,
103 std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
104
105 double
106 compute_value(const unsigned int i, const Point<dim> &p) const override;
107
113 template <int order>
115 compute_derivative(const unsigned int i, const Point<dim> &p) const;
116
118 compute_1st_derivative(const unsigned int i,
119 const Point<dim> & p) const override;
120
127 compute_2nd_derivative(const unsigned int i,
128 const Point<dim> & p) const override;
129
136 compute_3rd_derivative(const unsigned int i,
137 const Point<dim> & p) const override;
138
145 compute_4th_derivative(const unsigned int i,
146 const Point<dim> & p) const override;
147
154 compute_grad(const unsigned int i, const Point<dim> &p) const override;
155
162 compute_grad_grad(const unsigned int i, const Point<dim> &p) const override;
163
164 std::string
165 name() const override;
166
167 virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
168 clone() const override;
169
170private:
175
180};
181
182
183
184template <int dim>
185template <int order>
188 const unsigned int i,
189 const Point<dim> & p) const
190{
192
193 AssertDimension(order, 1);
194 const auto grad = compute_grad(i, p);
195
196 for (unsigned int i = 0; i < dim; ++i)
197 der[i] = grad[i];
198
199 return der;
200}
201
203
204#endif
Definition point.h:112
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
Tensor< 1, dim > compute_grad(const unsigned int i, const Point< dim > &p) const override
Tensor< 1, dim > compute_1st_derivative(const unsigned int i, const Point< dim > &p) const override
Tensor< 4, dim > compute_4th_derivative(const unsigned int i, const Point< dim > &p) const override
std::string name() const override
Tensor< 3, dim > compute_3rd_derivative(const unsigned int i, const Point< dim > &p) const override
Tensor< 2, dim > compute_2nd_derivative(const unsigned int i, const Point< dim > &p) const override
Tensor< order, dim > compute_derivative(const unsigned int i, const Point< dim > &p) const
virtual std::unique_ptr< ScalarPolynomialsBase< dim > > clone() const override
static constexpr unsigned int dimension
double compute_value(const unsigned int i, const Point< dim > &p) const override
const BarycentricPolynomials< 1 > poly_line
const BarycentricPolynomials< 2 > poly_tri
Tensor< 2, dim > compute_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual unsigned int degree() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define AssertDimension(dim1, dim2)
static const constexpr ::ndarray< unsigned int, 18, 2 > wedge_table_2
static const constexpr ::ndarray< unsigned int, 6, 2 > wedge_table_1