v0.14.0
base_functions.h
Go to the documentation of this file.
1 /** \file base_functions.c
2 
3 */
4 
5 #ifndef __BASE_FUNCTIONS_H__
6 #define __BASE_FUNCTIONS_H__
7 
8 #ifdef __cplusplus
9 extern "C" {
10 #endif
11 
12 /**
13 \brief Calculate Legendre approximation basis
14 
15 \ingroup mofem_base_functions
16 
17 Lagrange polynomial is given by
18 \f[
19 L_0(s)=1;\quad L_1(s) = s
20 \f]
21 and following terms are generated inductively
22 \f[
23 L_{l+1}=\frac{2l+1}{l+1}sL_l(s)-\frac{l}{l+1}L_{l-1}(s)
24 \f]
25 
26 Note that:
27 \f[
28 s\in[-1,1] \quad \textrm{and}\; s=s(\xi_0,\xi_1,\xi_2)
29 \f]
30 where \f$\xi_i\f$ are barycentric coordinates of element.
31 
32 \param p is approximation order
33 \param s is position \f$s\in[-1,1]\f$
34 \param diff_s derivatives of shape functions, i.e. \f$\frac{\partial s}{\partial
35 \xi_i}\f$ \retval L approximation functions \retval diffL derivatives, i.e.
36 \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension \return error code
37 
38 */
39 PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L,
40  double *diffL, const int dim);
41 
42 /**
43 \brief Calculate Jacobi approximation basis
44 
45 For more details see \cite fuentes2015orientation
46 
47 \param p is approximation order
48 \param alpha polynomial parameter
49 \param x is position \f$s\in[0,t]\f$
50 \param t range of polynomial
51 \param diff_x derivatives of shape functions, i.e. \f$\frac{\partial x}{\partial
52 \xi_i}\f$ \param diff_t derivatives of shape functions, i.e. \f$\frac{\partial
53 t}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
54 derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
55 \return error code
56 
57 */
58 PetscErrorCode Jacobi_polynomials(int p, double alpha, double x, double t,
59  double *diff_x, double *diff_t, double *L,
60  double *diffL, const int dim);
61 
62 /**
63 \brief Calculate integrated Jacobi approximation basis
64 
65 For more details see \cite fuentes2015orientation
66 
67 \param p is approximation order
68 \param alpha polynomial parameter
69 \param x is position \f$s\in[0,t]\f$
70 \param t range of polynomial
71 \param diff_x derivatives of shape functions, i.e. \f$\frac{\partial x}{\partial
72 \xi_i}\f$ \param diff_t derivatives of shape functions, i.e. \f$\frac{\partial
73 t}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
74 derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
75 \return error code
76 
77 */
78 PetscErrorCode IntegratedJacobi_polynomials(int p, double alpha, double x,
79  double t, double *diff_x,
80  double *diff_t, double *L,
81  double *diffL, const int dim);
82 
83 /**
84  \brief Calculate Lobatto base functions \cite FUENTES2015353.
85 
86  \ingroup mofem_base_functions
87 
88  Order of first function is 2 and goes to p.
89 
90  \param p is approximation order
91  \param s is a mapping of coordinates of edge to \f$[-1, 1]\f$, i.e., \f$s(\xi_1,\cdot,\xi_{dim})\in[-1,1]\f$
92  \param diff_s jacobian of the transformation, i.e. \f$\frac{\partial
93  s}{\partial \xi_i}\f$
94  - output
95  \retval L values basis functions at s
96  \retval diffL derivatives of basis functions at s, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
97  \return error code
98 */
99 PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L,
100  double *diffL, const int dim);
101 
102 /**
103  \brief Calculate Kernel Lobatto base functions.
104 
105  \ingroup mofem_base_functions
106 
107  This is implemented using definitions from Hermes2d
108  <https://github.com/hpfem/hermes> following book by Pavel Solin et al \cite
109  solin2003higher.
110 
111  \param p is approximation order
112  \param s is position \f$s\in[-1,1]\f$
113  \param diff_s derivatives of shape functions, i.e. \f$\frac{\partial
114  s}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
115  derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
116  \return error code
117 
118 */
119 PetscErrorCode LobattoKernel_polynomials(int p, double s, double *diff_s,
120  double *L, double *diffL,
121  const int dim);
122 
123 /// Definitions taken from Hermes2d code
124 
125 /// kernel functions for Lobatto base
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) * \
147  (x) + \
148  7))
149 #define LOBATTO_PHI9(x) \
150  (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
151  ((((4199 * (x) * (x)-7956) * (x) * (x) + 4914) * (x) * (x)-1092) * (x) * \
152  (x) + \
153  63) * \
154  (x))
155 
156 /// Derivatives of kernel functions for Lobbatto base
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) * \
178  (x)))
179 #define LOBATTO_PHI9X(x) \
180  (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
181  ((((37791.0 * (x) * (x)-55692.0) * (x) * (x) + 24570.0) * (x) * \
182  (x)-3276.0) * \
183  (x) * (x)-63.0))
184 
185 #ifdef __cplusplus
186 }
187 #endif
188 
189 #endif //__BASE_FUNCTIONS_H__
Lobatto_polynomials
PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Lobatto base functions .
Definition: base_functions.c:182
IntegratedJacobi_polynomials
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.
Definition: base_functions.c:134
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
Legendre_polynomials
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
Definition: base_functions.c:15
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
LobattoKernel_polynomials
PetscErrorCode LobattoKernel_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Kernel Lobatto base functions.
Definition: base_functions.c:328
Jacobi_polynomials
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.
Definition: base_functions.c:67