Go to the source code of this file.
|
PetscErrorCode | Legendre_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim) |
| Calculate Legendre approximation basis. More...
|
|
PetscErrorCode | Jacobi_polynomials (int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim) |
| Calculate Jacobi approximation basis. More...
|
|
PetscErrorCode | IntegratedJacobi_polynomials (int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim) |
| Calculate integrated Jacobi approximation basis. More...
|
|
PetscErrorCode | Lobatto_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim) |
| Calculate Lobatto base functions [29]. More...
|
|
PetscErrorCode | LobattoKernel_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim) |
| Calculate Kernel Lobatto base functions. More...
|
|
◆ LOBATTO_PHI0
#define LOBATTO_PHI0 |
( |
|
x | ) |
(-2.0 * 1.22474487139158904909864203735) |
Definitions taken from Hermes2d code.
kernel functions for Lobatto base
Definition at line 126 of file base_functions.h.
◆ LOBATTO_PHI0X
#define LOBATTO_PHI0X |
( |
|
x | ) |
(0) |
Derivatives of kernel functions for Lobbatto base.
Definition at line 157 of file base_functions.h.
◆ LOBATTO_PHI1
#define LOBATTO_PHI1 |
( |
|
x | ) |
(-2.0 * 1.58113883008418966599944677222 * (x)) |
◆ LOBATTO_PHI1X
#define LOBATTO_PHI1X |
( |
|
x | ) |
(-2.0 * 1.58113883008418966599944677222) |
◆ LOBATTO_PHI2
#define LOBATTO_PHI2 |
( |
|
x | ) |
(-1.0 / 2.0 * 1.87082869338697069279187436616 * (5 * (x) * (x)-1)) |
◆ LOBATTO_PHI2X
#define LOBATTO_PHI2X |
( |
|
x | ) |
(-1.0 / 2.0 * 1.87082869338697069279187436616 * (10 * (x))) |
◆ LOBATTO_PHI3
#define LOBATTO_PHI3 |
( |
|
x | ) |
(-1.0 / 2.0 * 2.12132034355964257320253308631 * (7 * (x) * (x)-3) * (x)) |
◆ LOBATTO_PHI3X
#define LOBATTO_PHI3X |
( |
|
x | ) |
(-1.0 / 2.0 * 2.12132034355964257320253308631 * (21.0 * (x) * (x)-3.0)) |
◆ LOBATTO_PHI4
#define LOBATTO_PHI4 |
( |
|
x | ) |
|
Value: (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
(21 * (x) * (x) * (x) * (x)-14 * (x) * (x) + 1))
Definition at line 132 of file base_functions.h.
◆ LOBATTO_PHI4X
#define LOBATTO_PHI4X |
( |
|
x | ) |
|
Value: (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
((84.0 * (x) * (x)-28.0) * (x)))
Definition at line 163 of file base_functions.h.
◆ LOBATTO_PHI5
#define LOBATTO_PHI5 |
( |
|
x | ) |
|
Value: (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
((33 * (x) * (x)-30) * (x) * (x) + 5) * (x))
Definition at line 135 of file base_functions.h.
◆ LOBATTO_PHI5X
#define LOBATTO_PHI5X |
( |
|
x | ) |
|
Value: (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
((165.0 * (x) * (x)-90.0) * (x) * (x) + 5.0))
Definition at line 166 of file base_functions.h.
◆ LOBATTO_PHI6
#define LOBATTO_PHI6 |
( |
|
x | ) |
|
Value: (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
(((429 * (x) * (x)-495) * (x) * (x) + 135) * (x) * (x)-5))
Definition at line 138 of file base_functions.h.
◆ LOBATTO_PHI6X
#define LOBATTO_PHI6X |
( |
|
x | ) |
|
Value: (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
(((2574.0 * (x) * (x)-1980.0) * (x) * (x) + 270.0) * (x)))
Definition at line 169 of file base_functions.h.
◆ LOBATTO_PHI7
#define LOBATTO_PHI7 |
( |
|
x | ) |
|
Value: (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
(((715 * (x) * (x)-1001) * (x) * (x) + 385) * (x) * (x)-35) * (x))
Definition at line 141 of file base_functions.h.
◆ LOBATTO_PHI7X
#define LOBATTO_PHI7X |
( |
|
x | ) |
|
Value: (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
(((5005.0 * (x) * (x)-5005.0) * (x) * (x) + 1155.0) * (x) * (x)-35.0))
Definition at line 172 of file base_functions.h.
◆ LOBATTO_PHI8
#define LOBATTO_PHI8 |
( |
|
x | ) |
|
Value: (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
((((2431 * (x) * (x)-4004) * (x) * (x) + 2002) * (x) * (x)-308) * (x) * \
(x) + \
7))
Definition at line 144 of file base_functions.h.
◆ LOBATTO_PHI8X
#define LOBATTO_PHI8X |
( |
|
x | ) |
|
Value: (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
((((19448.0 * (x) * (x)-24024.0) * (x) * (x) + 8008.0) * (x) * (x)-616.0) * \
(x)))
Definition at line 175 of file base_functions.h.
◆ LOBATTO_PHI9
#define LOBATTO_PHI9 |
( |
|
x | ) |
|
Value: (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
((((4199 * (x) * (x)-7956) * (x) * (x) + 4914) * (x) * (x)-1092) * (x) * \
(x) + \
63) * \
(x))
Definition at line 149 of file base_functions.h.
◆ LOBATTO_PHI9X
#define LOBATTO_PHI9X |
( |
|
x | ) |
|
Value: (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
((((37791.0 * (x) * (x)-55692.0) * (x) * (x) + 24570.0) * (x) * \
(x)-3276.0) * \
(x) * (x)-63.0))
Definition at line 179 of file base_functions.h.
◆ IntegratedJacobi_polynomials()
Calculate integrated Jacobi approximation basis.
For more details see [30]
- Parameters
-
p | is approximation order |
alpha | polynomial parameter |
x | is position \(s\in[0,t]\) |
t | range of polynomial |
diff_x | derivatives of shape functions, i.e. \(\frac{\partial x}{\partial \xi_i}\) |
diff_t | derivatives of shape functions, i.e. \(\frac{\partial t}{\partial \xi_i}\) |
- Return values
-
L | approximation functions |
diffL | derivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\) |
- Parameters
-
- Returns
- error code
Definition at line 151 of file base_functions.c.
167 for (;
d != dim; ++
d) {
168 diffL[
d * p + 0] = diff_x[
d];
173 double jacobi[(p + 1)];
174 double diff_jacobi[(p + 1) * dim];
181 const double a = (
i + alpha) / ((2 *
i + alpha - 1) * (2 *
i + alpha));
182 const double b = alpha / ((2 *
i + alpha - 2) * (2 *
i + alpha));
183 const double c = (
i - 1) / ((2 *
i + alpha - 2) * (2 *
i + alpha - 1));
184 L[
l] =
a * jacobi[
i] + b *
t * jacobi[
i - 1] -
c *
t *
t * jacobi[
i - 2];
187 for (;
dd != dim; ++
dd) {
188 diffL[
dd * p +
l] =
a * diff_jacobi[
dd * (p + 1) +
i] +
189 b * (
t * diff_jacobi[
dd * (p + 1) +
i - 1] +
190 diff_t[
dd] * jacobi[
i - 1]) -
191 c * (
t *
t * diff_jacobi[
dd * (p + 1) +
i - 2] +
192 2 *
t * diff_t[
dd] * jacobi[
i - 2]);
◆ Jacobi_polynomials()
Calculate Jacobi approximation basis.
For more details see [30]
- Parameters
-
p | is approximation order |
alpha | polynomial parameter |
x | is position \(s\in[0,t]\) |
t | range of polynomial |
diff_x | derivatives of shape functions, i.e. \(\frac{\partial x}{\partial \xi_i}\) |
diff_t | derivatives of shape functions, i.e. \(\frac{\partial t}{\partial \xi_i}\) |
- Return values
-
L | approximation functions |
diffL | derivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\) |
- Parameters
-
- Returns
- error code
Definition at line 67 of file base_functions.c.
81 diffL[0 * (p + 1) + 0] = 0;
83 diffL[1 * (p + 1) + 0] = 0;
85 diffL[2 * (p + 1) + 0] = 0;
91 L[1] = 2 * x -
t + alpha * x;
98 double d_t = (diff_t) ? diff_t[0] : 0;
99 diffL[0 * (p + 1) + 1] = (2 + alpha) * diff_x[0] - d_t;
101 double d_t = (diff_t) ? diff_t[1] : 0;
102 diffL[1 * (p + 1) + 1] = (2 + alpha) * diff_x[1] - d_t;
104 double d_t = (diff_t) ? diff_t[2] : 0;
105 diffL[2 * (p + 1) + 1] = (2 + alpha) * diff_x[2] - d_t;
114 double a = 2 * lp1 * (lp1 + alpha) * (2 * lp1 + alpha - 2);
115 double b = 2 * lp1 + alpha - 1;
116 double c = (2 * lp1 + alpha) * (2 * lp1 + alpha - 2);
117 double d = 2 * (lp1 + alpha - 1) * (lp1 - 1) * (2 * lp1 + alpha);
118 double A = b * (
c * (2 * x -
t) + alpha * alpha *
t) /
a;
119 double B =
d *
t *
t /
a;
120 L[lp1] =
A *
L[
l] - B *
L[
l - 1];
122 double d_t = (diff_t) ? diff_t[0] : 0;
123 double diffA = b * (
c * (2 * diff_x[0] - d_t) + alpha * alpha * d_t) /
a;
124 double diffB =
d * 2 *
t * d_t /
a;
125 diffL[0 * (p + 1) + lp1] =
A * diffL[0 * (p + 1) +
l] -
126 B * diffL[0 * (p + 1) +
l - 1] + diffA *
L[
l] -
129 double d_t = (diff_t) ? diff_t[1] : 0;
131 b * (
c * (2 * diff_x[1] - d_t) + alpha * alpha * d_t) /
a;
132 double diffB =
d * 2 *
t * d_t /
a;
133 diffL[1 * (p + 1) + lp1] =
A * diffL[1 * (p + 1) +
l] -
134 B * diffL[1 * (p + 1) +
l - 1] +
135 diffA *
L[
l] - diffB *
L[
l - 1];
137 double d_t = (diff_t) ? diff_t[2] : 0;
139 b * (
c * (2 * diff_x[2] - d_t) + alpha * alpha * d_t) /
a;
140 double diffB =
d * 2 *
t * d_t /
a;
141 diffL[2 * (p + 1) + lp1] =
A * diffL[2 * (p + 1) +
l] -
142 B * diffL[2 * (p + 1) +
l - 1] +
143 diffA *
L[
l] - diffB *
L[
l - 1];
PetscErrorCode Jacobi_polynomials(int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
Calculate Jacobi approximation basis.
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)