v0.9.0
base_functions.h
Go to the documentation of this file.
1 /** \file base_functions.c
2 
3 */
4 
5 /**
6  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
7  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
8  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
9  * License for more details.
10  *
11  * You should have received a copy of the GNU Lesser General Public
12  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
13  */
14 
15 #ifndef __BASE_FUNCTIONS_H__
16 #define __BASE_FUNCTIONS_H__
17 
18 #ifdef __cplusplus
19 extern "C" {
20 #endif
21 
22 /**
23 \brief Calculate Legendre approximation basis
24 
25 \ingroup mofem_base_functions
26 
27 Lagrange polynomial is given by
28 \f[
29 L_0(s)=1;\quad L_1(s) = s
30 \f]
31 and following terms are generated inductively
32 \f[
33 L_{l+1}=\frac{2l+1}{l+1}sL_l(s)-\frac{l}{l+1}L_{l-1}(s)
34 \f]
35 
36 Note that:
37 \f[
38 s\in[-1,1] \quad \textrm{and}\; s=s(\xi_0,\xi_1,\xi_2)
39 \f]
40 where \f$\xi_i\f$ are barycentric coordinates of element.
41 
42 \param p is approximation order
43 \param s is position \f$s\in[-1,1]\f$
44 \param diff_s derivatives of shape functions, i.e. \f$\frac{\partial s}{\partial
45 \xi_i}\f$ \retval L approximation functions \retval diffL derivatives, i.e.
46 \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension \return error code
47 
48 */
49 PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L,
50  double *diffL, const int dim);
51 
52 /**
53 \brief Calculate Jacobi approximation basis
54 
55 For more details see \cite fuentes2015orientation
56 
57 \param p is approximation order
58 \param alpha polynomial parameter
59 \param x is position \f$s\in[0,t]\f$
60 \param t range of polynomial
61 \param diff_x derivatives of shape functions, i.e. \f$\frac{\partial x}{\partial
62 \xi_i}\f$ \param diff_t derivatives of shape functions, i.e. \f$\frac{\partial
63 t}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
64 derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
65 \return error code
66 
67 */
68 PetscErrorCode Jacobi_polynomials(int p, double alpha, double x, double t,
69  double *diff_x, double *diff_t, double *L,
70  double *diffL, const int dim);
71 
72 /**
73 \brief Calculate integrated Jacobi approximation basis
74 
75 For more details see \cite fuentes2015orientation
76 
77 \param p is approximation order
78 \param alpha polynomial parameter
79 \param x is position \f$s\in[0,t]\f$
80 \param t range of polynomial
81 \param diff_x derivatives of shape functions, i.e. \f$\frac{\partial x}{\partial
82 \xi_i}\f$ \param diff_t derivatives of shape functions, i.e. \f$\frac{\partial
83 t}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
84 derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
85 \return error code
86 
87 */
88 PetscErrorCode IntegratedJacobi_polynomials(int p, double alpha, double x,
89  double t, double *diff_x,
90  double *diff_t, double *L,
91  double *diffL, const int dim);
92 
93 /**
94  \brief Calculate Lobatto base functions.
95 
96  \ingroup mofem_base_functions
97 
98  Order of first function is 2.
99 
100  \param p is approximation order
101  \param s is position \f$s\in[-1,1]\f$
102  \param diff_s derivatives of shape functions, i.e. \f$\frac{\partial
103  s}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
104  derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
105  \return error code
106 
107 */
108 PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L,
109  double *diffL, const int dim);
110 
111 /**
112  \brief Calculate Kernel Lobatto base functions.
113 
114  \ingroup mofem_base_functions
115 
116  This is implemented using definitions from Hermes2d
117  <https://github.com/hpfem/hermes> following book by Pavel Solin et al \cite
118  solin2003higher.
119 
120  \param p is approximation order
121  \param s is position \f$s\in[-1,1]\f$
122  \param diff_s derivatives of shape functions, i.e. \f$\frac{\partial
123  s}{\partial \xi_i}\f$ \retval L approximation functions \retval diffL
124  derivatives, i.e. \f$\frac{\partial L}{\partial \xi_i}\f$ \param dim dimension
125  \return error code
126 
127 */
128 PetscErrorCode LobattoKernel_polynomials(int p, double s, double *diff_s,
129  double *L, double *diffL,
130  const int dim);
131 
132 /// Definitions taken from Hermes2d code
133 
134 /// kernel functions for Lobatto base
135 #define LOBATTO_PHI0(x) (-2.0 * 1.22474487139158904909864203735)
136 #define LOBATTO_PHI1(x) (-2.0 * 1.58113883008418966599944677222 * (x))
137 #define LOBATTO_PHI2(x) \
138  (-1.0 / 2.0 * 1.87082869338697069279187436616 * (5 * (x) * (x)-1))
139 #define LOBATTO_PHI3(x) \
140  (-1.0 / 2.0 * 2.12132034355964257320253308631 * (7 * (x) * (x)-3) * (x))
141 #define LOBATTO_PHI4(x) \
142  (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
143  (21 * (x) * (x) * (x) * (x)-14 * (x) * (x) + 1))
144 #define LOBATTO_PHI5(x) \
145  (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
146  ((33 * (x) * (x)-30) * (x) * (x) + 5) * (x))
147 #define LOBATTO_PHI6(x) \
148  (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
149  (((429 * (x) * (x)-495) * (x) * (x) + 135) * (x) * (x)-5))
150 #define LOBATTO_PHI7(x) \
151  (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
152  (((715 * (x) * (x)-1001) * (x) * (x) + 385) * (x) * (x)-35) * (x))
153 #define LOBATTO_PHI8(x) \
154  (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
155  ((((2431 * (x) * (x)-4004) * (x) * (x) + 2002) * (x) * (x)-308) * (x) * \
156  (x) + \
157  7))
158 #define LOBATTO_PHI9(x) \
159  (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
160  ((((4199 * (x) * (x)-7956) * (x) * (x) + 4914) * (x) * (x)-1092) * (x) * \
161  (x) + \
162  63) * \
163  (x))
164 
165 /// Derivatives of kernel functions for Lobbatto base
166 #define LOBATTO_PHI0X(x) (0)
167 #define LOBATTO_PHI1X(x) (-2.0 * 1.58113883008418966599944677222)
168 #define LOBATTO_PHI2X(x) \
169  (-1.0 / 2.0 * 1.87082869338697069279187436616 * (10 * (x)))
170 #define LOBATTO_PHI3X(x) \
171  (-1.0 / 2.0 * 2.12132034355964257320253308631 * (21.0 * (x) * (x)-3.0))
172 #define LOBATTO_PHI4X(x) \
173  (-1.0 / 4.0 * 2.34520787991171477728281505677 * \
174  ((84.0 * (x) * (x)-28.0) * (x)))
175 #define LOBATTO_PHI5X(x) \
176  (-1.0 / 4.0 * 2.54950975679639241501411205451 * \
177  ((165.0 * (x) * (x)-90.0) * (x) * (x) + 5.0))
178 #define LOBATTO_PHI6X(x) \
179  (-1.0 / 32.0 * 2.73861278752583056728484891400 * \
180  (((2574.0 * (x) * (x)-1980.0) * (x) * (x) + 270.0) * (x)))
181 #define LOBATTO_PHI7X(x) \
182  (-1.0 / 32.0 * 2.91547594742265023543707643877 * \
183  (((5005.0 * (x) * (x)-5005.0) * (x) * (x) + 1155.0) * (x) * (x)-35.0))
184 #define LOBATTO_PHI8X(x) \
185  (-1.0 / 64.0 * 3.08220700148448822512509619073 * \
186  ((((19448.0 * (x) * (x)-24024.0) * (x) * (x) + 8008.0) * (x) * (x)-616.0) * \
187  (x)))
188 #define LOBATTO_PHI9X(x) \
189  (-1.0 / 128.0 * 6.4807406984078603784382721642 * \
190  ((((37791.0 * (x) * (x)-55692.0) * (x) * (x) + 24570.0) * (x) * \
191  (x)-3276.0) * \
192  (x) * (x)-63.0))
193 
194 #ifdef __cplusplus
195 }
196 #endif
197 
198 #endif //__BASE_FUNCTIONS_H__
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
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.
PetscErrorCode LobattoKernel_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Kernel Lobatto base functions.
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.