|
| v0.14.0
|
Go to the documentation of this file.
5 #ifndef __BASE_FUNCTIONS_H__
6 #define __BASE_FUNCTIONS_H__
40 double *diffL,
const int dim);
59 double *diff_x,
double *diff_t,
double *
L,
60 double *diffL,
const int dim);
79 double t,
double *diff_x,
80 double *diff_t,
double *
L,
81 double *diffL,
const int dim);
100 double *diffL,
const int dim);
120 double *
L,
double *diffL,
126 #define LOBATTO_PHI0(x) (-2.0 * 1.22474487139158904909864203735)
127 #define LOBATTO_PHI1(x) (-2.0 * 1.58113883008418966599944677222 * (x))
128 #define LOBATTO_PHI2(x) \
129 (-1.0 / 2.0 * 1.87082869338697069279187436616 * (5 * (x) * (x)-1))
130 #define LOBATTO_PHI3(x) \
131 (-1.0 / 2.0 * 2.12132034355964257320253308631 * (7 * (x) * (x)-3) * (x))
132 #define LOBATTO_PHI4(x) \
133 (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
134 (21 * (x) * (x) * (x) * (x)-14 * (x) * (x) + 1))
135 #define LOBATTO_PHI5(x) \
136 (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
137 ((33 * (x) * (x)-30) * (x) * (x) + 5) * (x))
138 #define LOBATTO_PHI6(x) \
139 (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
140 (((429 * (x) * (x)-495) * (x) * (x) + 135) * (x) * (x)-5))
141 #define LOBATTO_PHI7(x) \
142 (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
143 (((715 * (x) * (x)-1001) * (x) * (x) + 385) * (x) * (x)-35) * (x))
144 #define LOBATTO_PHI8(x) \
145 (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
146 ((((2431 * (x) * (x)-4004) * (x) * (x) + 2002) * (x) * (x)-308) * (x) * \
149 #define LOBATTO_PHI9(x) \
150 (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
151 ((((4199 * (x) * (x)-7956) * (x) * (x) + 4914) * (x) * (x)-1092) * (x) * \
157 #define LOBATTO_PHI0X(x) (0)
158 #define LOBATTO_PHI1X(x) (-2.0 * 1.58113883008418966599944677222)
159 #define LOBATTO_PHI2X(x) \
160 (-1.0 / 2.0 * 1.87082869338697069279187436616 * (10 * (x)))
161 #define LOBATTO_PHI3X(x) \
162 (-1.0 / 2.0 * 2.12132034355964257320253308631 * (21.0 * (x) * (x)-3.0))
163 #define LOBATTO_PHI4X(x) \
164 (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
165 ((84.0 * (x) * (x)-28.0) * (x)))
166 #define LOBATTO_PHI5X(x) \
167 (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
168 ((165.0 * (x) * (x)-90.0) * (x) * (x) + 5.0))
169 #define LOBATTO_PHI6X(x) \
170 (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
171 (((2574.0 * (x) * (x)-1980.0) * (x) * (x) + 270.0) * (x)))
172 #define LOBATTO_PHI7X(x) \
173 (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
174 (((5005.0 * (x) * (x)-5005.0) * (x) * (x) + 1155.0) * (x) * (x)-35.0))
175 #define LOBATTO_PHI8X(x) \
176 (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
177 ((((19448.0 * (x) * (x)-24024.0) * (x) * (x) + 8008.0) * (x) * (x)-616.0) * \
179 #define LOBATTO_PHI9X(x) \
180 (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
181 ((((37791.0 * (x) * (x)-55692.0) * (x) * (x) + 24570.0) * (x) * \
189 #endif //__BASE_FUNCTIONS_H__
PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Lobatto base functions .
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.
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
constexpr double t
plate stiffness
PetscErrorCode LobattoKernel_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Kernel Lobatto base functions.
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.