1/* specfunc/gsl_sf_legendre.h
2 *
3 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001, 2002, 2004 Gerard Jungman
4 * Copyright (C) 2019 Patrick Alken
5 *
6 * This program is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 3 of the License, or (at
9 * your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful, but
12 * WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 * General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
19 */
20
21/* Author: G. Jungman */
22
23#ifndef __GSL_SF_LEGENDRE_H__
24#define __GSL_SF_LEGENDRE_H__
25
26#include <stdlib.h>
27#include <gsl/gsl_inline.h>
28#include <gsl/gsl_sf_result.h>
29
30#undef __BEGIN_DECLS
31#undef __END_DECLS
32#ifdef __cplusplus
33# define __BEGIN_DECLS extern "C" {
34# define __END_DECLS }
35#else
36# define __BEGIN_DECLS /* empty */
37# define __END_DECLS /* empty */
38#endif
39
40__BEGIN_DECLS
41
42
43/* P_l(x) l >= 0; |x| <= 1
44 *
45 * exceptions: GSL_EDOM
46 */
47int gsl_sf_legendre_Pl_e(const int l, const double x, gsl_sf_result * result);
48double gsl_sf_legendre_Pl(const int l, const double x);
49
50
51/* P_l(x) for l=0,...,lmax; |x| <= 1
52 *
53 * exceptions: GSL_EDOM
54 */
55int gsl_sf_legendre_Pl_array(
56 const int lmax, const double x,
57 double * result_array
58 );
59
60
61/* P_l(x) and P_l'(x) for l=0,...,lmax; |x| <= 1
62 *
63 * exceptions: GSL_EDOM
64 */
65int gsl_sf_legendre_Pl_deriv_array(
66 const int lmax, const double x,
67 double * result_array,
68 double * result_deriv_array
69 );
70
71
72/* P_l(x), l=1,2,3
73 *
74 * exceptions: none
75 */
76int gsl_sf_legendre_P1_e(double x, gsl_sf_result * result);
77int gsl_sf_legendre_P2_e(double x, gsl_sf_result * result);
78int gsl_sf_legendre_P3_e(double x, gsl_sf_result * result);
79double gsl_sf_legendre_P1(const double x);
80double gsl_sf_legendre_P2(const double x);
81double gsl_sf_legendre_P3(const double x);
82
83
84/* Q_0(x), x > -1, x != 1
85 *
86 * exceptions: GSL_EDOM
87 */
88int gsl_sf_legendre_Q0_e(const double x, gsl_sf_result * result);
89double gsl_sf_legendre_Q0(const double x);
90
91
92/* Q_1(x), x > -1, x != 1
93 *
94 * exceptions: GSL_EDOM
95 */
96int gsl_sf_legendre_Q1_e(const double x, gsl_sf_result * result);
97double gsl_sf_legendre_Q1(const double x);
98
99
100/* Q_l(x), x > -1, x != 1, l >= 0
101 *
102 * exceptions: GSL_EDOM
103 */
104int gsl_sf_legendre_Ql_e(const int l, const double x, gsl_sf_result * result);
105double gsl_sf_legendre_Ql(const int l, const double x);
106
107
108/* P_l^m(x) m >= 0; l >= m; |x| <= 1.0
109 *
110 * Note that this function grows combinatorially with l.
111 * Therefore we can easily generate an overflow for l larger
112 * than about 150.
113 *
114 * There is no trouble for small m, but when m and l are both large,
115 * then there will be trouble. Rather than allow overflows, these
116 * functions refuse to calculate when they can sense that l and m are
117 * too big.
118 *
119 * If you really want to calculate a spherical harmonic, then DO NOT
120 * use this. Instead use legendre_sphPlm() below, which uses a similar
121 * recursion, but with the normalized functions.
122 *
123 * exceptions: GSL_EDOM, GSL_EOVRFLW
124 */
125int gsl_sf_legendre_Plm_e(const int l, const int m, const double x, gsl_sf_result * result);
126double gsl_sf_legendre_Plm(const int l, const int m, const double x);
127
128
129/* P_l^m(x) m >= 0; l >= m; |x| <= 1.0
130 * l=|m|,...,lmax
131 *
132 * exceptions: GSL_EDOM, GSL_EOVRFLW
133 */
134int gsl_sf_legendre_Plm_array(
135 const int lmax, const int m, const double x,
136 double * result_array
137 );
138
139
140/* P_l^m(x) and d(P_l^m(x))/dx; m >= 0; lmax >= m; |x| <= 1.0
141 * l=|m|,...,lmax
142 *
143 * exceptions: GSL_EDOM, GSL_EOVRFLW
144 */
145int gsl_sf_legendre_Plm_deriv_array(
146 const int lmax, const int m, const double x,
147 double * result_array,
148 double * result_deriv_array
149 );
150
151
152/* P_l^m(x), normalized properly for use in spherical harmonics
153 * m >= 0; l >= m; |x| <= 1.0
154 *
155 * There is no overflow problem, as there is for the
156 * standard normalization of P_l^m(x).
157 *
158 * Specifically, it returns:
159 *
160 * sqrt((2l+1)/(4pi)) sqrt((l-m)!/(l+m)!) P_l^m(x)
161 *
162 * exceptions: GSL_EDOM
163 */
164int gsl_sf_legendre_sphPlm_e(const int l, int m, const double x, gsl_sf_result * result);
165double gsl_sf_legendre_sphPlm(const int l, const int m, const double x);
166
167
168/* sphPlm(l,m,x) values
169 * m >= 0; l >= m; |x| <= 1.0
170 * l=|m|,...,lmax
171 *
172 * exceptions: GSL_EDOM
173 */
174int gsl_sf_legendre_sphPlm_array(
175 const int lmax, int m, const double x,
176 double * result_array
177 );
178
179
180/* sphPlm(l,m,x) and d(sphPlm(l,m,x))/dx values
181 * m >= 0; l >= m; |x| <= 1.0
182 * l=|m|,...,lmax
183 *
184 * exceptions: GSL_EDOM
185 */
186int gsl_sf_legendre_sphPlm_deriv_array(
187 const int lmax, const int m, const double x,
188 double * result_array,
189 double * result_deriv_array
190 );
191
192
193
194/* size of result_array[] needed for the array versions of Plm
195 * (lmax - m + 1)
196 */
197int gsl_sf_legendre_array_size(const int lmax, const int m);
198
199/* Irregular Spherical Conical Function
200 * P^{1/2}_{-1/2 + I lambda}(x)
201 *
202 * x > -1.0
203 * exceptions: GSL_EDOM
204 */
205int gsl_sf_conicalP_half_e(const double lambda, const double x, gsl_sf_result * result);
206double gsl_sf_conicalP_half(const double lambda, const double x);
207
208
209/* Regular Spherical Conical Function
210 * P^{-1/2}_{-1/2 + I lambda}(x)
211 *
212 * x > -1.0
213 * exceptions: GSL_EDOM
214 */
215int gsl_sf_conicalP_mhalf_e(const double lambda, const double x, gsl_sf_result * result);
216double gsl_sf_conicalP_mhalf(const double lambda, const double x);
217
218
219/* Conical Function
220 * P^{0}_{-1/2 + I lambda}(x)
221 *
222 * x > -1.0
223 * exceptions: GSL_EDOM
224 */
225int gsl_sf_conicalP_0_e(const double lambda, const double x, gsl_sf_result * result);
226double gsl_sf_conicalP_0(const double lambda, const double x);
227
228
229/* Conical Function
230 * P^{1}_{-1/2 + I lambda}(x)
231 *
232 * x > -1.0
233 * exceptions: GSL_EDOM
234 */
235int gsl_sf_conicalP_1_e(const double lambda, const double x, gsl_sf_result * result);
236double gsl_sf_conicalP_1(const double lambda, const double x);
237
238
239/* Regular Spherical Conical Function
240 * P^{-1/2-l}_{-1/2 + I lambda}(x)
241 *
242 * x > -1.0, l >= -1
243 * exceptions: GSL_EDOM
244 */
245int gsl_sf_conicalP_sph_reg_e(const int l, const double lambda, const double x, gsl_sf_result * result);
246double gsl_sf_conicalP_sph_reg(const int l, const double lambda, const double x);
247
248
249/* Regular Cylindrical Conical Function
250 * P^{-m}_{-1/2 + I lambda}(x)
251 *
252 * x > -1.0, m >= -1
253 * exceptions: GSL_EDOM
254 */
255int gsl_sf_conicalP_cyl_reg_e(const int m, const double lambda, const double x, gsl_sf_result * result);
256double gsl_sf_conicalP_cyl_reg(const int m, const double lambda, const double x);
257
258
259/* The following spherical functions are specializations
260 * of Legendre functions which give the regular eigenfunctions
261 * of the Laplacian on a 3-dimensional hyperbolic space.
262 * Of particular interest is the flat limit, which is
263 * Flat-Lim := {lambda->Inf, eta->0, lambda*eta fixed}.
264 */
265
266/* Zeroth radial eigenfunction of the Laplacian on the
267 * 3-dimensional hyperbolic space.
268 *
269 * legendre_H3d_0(lambda,eta) := sin(lambda*eta)/(lambda*sinh(eta))
270 *
271 * Normalization:
272 * Flat-Lim legendre_H3d_0(lambda,eta) = j_0(lambda*eta)
273 *
274 * eta >= 0.0
275 * exceptions: GSL_EDOM
276 */
277int gsl_sf_legendre_H3d_0_e(const double lambda, const double eta, gsl_sf_result * result);
278double gsl_sf_legendre_H3d_0(const double lambda, const double eta);
279
280
281/* First radial eigenfunction of the Laplacian on the
282 * 3-dimensional hyperbolic space.
283 *
284 * legendre_H3d_1(lambda,eta) :=
285 * 1/sqrt(lambda^2 + 1) sin(lam eta)/(lam sinh(eta))
286 * (coth(eta) - lambda cot(lambda*eta))
287 *
288 * Normalization:
289 * Flat-Lim legendre_H3d_1(lambda,eta) = j_1(lambda*eta)
290 *
291 * eta >= 0.0
292 * exceptions: GSL_EDOM
293 */
294int gsl_sf_legendre_H3d_1_e(const double lambda, const double eta, gsl_sf_result * result);
295double gsl_sf_legendre_H3d_1(const double lambda, const double eta);
296
297
298/* l'th radial eigenfunction of the Laplacian on the
299 * 3-dimensional hyperbolic space.
300 *
301 * Normalization:
302 * Flat-Lim legendre_H3d_l(l,lambda,eta) = j_l(lambda*eta)
303 *
304 * eta >= 0.0, l >= 0
305 * exceptions: GSL_EDOM
306 */
307int gsl_sf_legendre_H3d_e(const int l, const double lambda, const double eta, gsl_sf_result * result);
308double gsl_sf_legendre_H3d(const int l, const double lambda, const double eta);
309
310
311/* Array of H3d(ell), 0 <= ell <= lmax
312 */
313int gsl_sf_legendre_H3d_array(const int lmax, const double lambda, const double eta, double * result_array);
314
315/* associated legendre P_{lm} routines */
316
317typedef enum
318{
319 GSL_SF_LEGENDRE_SCHMIDT,
320 GSL_SF_LEGENDRE_SPHARM,
321 GSL_SF_LEGENDRE_FULL,
322 GSL_SF_LEGENDRE_NONE
323} gsl_sf_legendre_t;
324
325int gsl_sf_legendre_array(const gsl_sf_legendre_t norm,
326 const size_t lmax, const double x,
327 double result_array[]);
328int gsl_sf_legendre_array_e(const gsl_sf_legendre_t norm,
329 const size_t lmax, const double x,
330 const double csphase,
331 double result_array[]);
332int gsl_sf_legendre_deriv_array(const gsl_sf_legendre_t norm,
333 const size_t lmax, const double x,
334 double result_array[],
335 double result_deriv_array[]);
336int gsl_sf_legendre_deriv_array_e(const gsl_sf_legendre_t norm,
337 const size_t lmax, const double x,
338 const double csphase,
339 double result_array[],
340 double result_deriv_array[]);
341int gsl_sf_legendre_deriv_alt_array(const gsl_sf_legendre_t norm,
342 const size_t lmax, const double x,
343 double result_array[],
344 double result_deriv_array[]);
345int gsl_sf_legendre_deriv_alt_array_e(const gsl_sf_legendre_t norm,
346 const size_t lmax, const double x,
347 const double csphase,
348 double result_array[],
349 double result_deriv_array[]);
350int gsl_sf_legendre_deriv2_array(const gsl_sf_legendre_t norm,
351 const size_t lmax, const double x,
352 double result_array[],
353 double result_deriv_array[],
354 double result_deriv2_array[]);
355int gsl_sf_legendre_deriv2_array_e(const gsl_sf_legendre_t norm,
356 const size_t lmax, const double x,
357 const double csphase,
358 double result_array[],
359 double result_deriv_array[],
360 double result_deriv2_array[]);
361int gsl_sf_legendre_deriv2_alt_array(const gsl_sf_legendre_t norm,
362 const size_t lmax, const double x,
363 double result_array[],
364 double result_deriv_array[],
365 double result_deriv2_array[]);
366int gsl_sf_legendre_deriv2_alt_array_e(const gsl_sf_legendre_t norm,
367 const size_t lmax, const double x,
368 const double csphase,
369 double result_array[],
370 double result_deriv_array[],
371 double result_deriv2_array[]);
372size_t gsl_sf_legendre_array_n(const size_t lmax);
373size_t gsl_sf_legendre_nlm(const size_t lmax);
374
375INLINE_DECL size_t gsl_sf_legendre_array_index(const size_t l, const size_t m);
376
377#ifdef HAVE_INLINE
378
379/*
380gsl_sf_legendre_array_index()
381This routine computes the index into a result_array[] corresponding
382to a given (l,m)
383*/
384INLINE_FUN
385size_t
386gsl_sf_legendre_array_index(const size_t l, const size_t m)
387{
388 return (((l * (l + 1)) >> 1) + m);
389}
390
391#endif /* HAVE_INLINE */
392
393__END_DECLS
394
395#endif /* __GSL_SF_LEGENDRE_H__ */
396