1 /* specfunc/gsl_sf_bessel.h
2  *
3  * Copyright (C) 1996, 1997, 1998, 1999, 2000 Gerard Jungman
4  *
5  * This program is free software; you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation; either version 3 of the License, or (at
8  * your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful, but
11  * WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
13  * General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
18  */
19 
20 /* Author:  G. Jungman */
21 
22 #ifndef __GSL_SF_BESSEL_H__
23 #define __GSL_SF_BESSEL_H__
24 
25 #include <stdlib.h>
26 #include <gsl/gsl_mode.h>
27 #include <gsl/gsl_precision.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 /* Regular Bessel Function J_0(x)
44  *
45  * exceptions: none
46  */
47 int gsl_sf_bessel_J0_e(const double x,  gsl_sf_result * result);
48 double gsl_sf_bessel_J0(const double x);
49 
50 
51 /* Regular Bessel Function J_1(x)
52  *
53  * exceptions: GSL_EUNDRFLW
54  */
55 int gsl_sf_bessel_J1_e(const double x, gsl_sf_result * result);
56 double gsl_sf_bessel_J1(const double x);
57 
58 
59 /* Regular Bessel Function J_n(x)
60  *
61  * exceptions: GSL_EUNDRFLW
62  */
63 int gsl_sf_bessel_Jn_e(int n, double x, gsl_sf_result * result);
64 double gsl_sf_bessel_Jn(const int n, const double x);
65 
66 
67 /* Regular Bessel Function J_n(x),  nmin <= n <= nmax
68  *
69  * exceptions: GSL_EDOM, GSL_EUNDRFLW
70  */
71 int gsl_sf_bessel_Jn_array(int nmin, int nmax, double x, double * result_array);
72 
73 
74 /* Irregular Bessel function Y_0(x)
75  *
76  * x > 0.0
77  * exceptions: GSL_EDOM, GSL_EUNDRFLW
78  */
79 int gsl_sf_bessel_Y0_e(const double x, gsl_sf_result * result);
80 double gsl_sf_bessel_Y0(const double x);
81 
82 
83 /* Irregular Bessel function Y_1(x)
84  *
85  * x > 0.0
86  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
87  */
88 int gsl_sf_bessel_Y1_e(const double x, gsl_sf_result * result);
89 double gsl_sf_bessel_Y1(const double x);
90 
91 
92 /* Irregular Bessel function Y_n(x)
93  *
94  * x > 0.0
95  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
96  */
97 int gsl_sf_bessel_Yn_e(int n,const double x, gsl_sf_result * result);
98 double gsl_sf_bessel_Yn(const int n,const double x);
99 
100 
101 /* Irregular Bessel function Y_n(x), nmin <= n <= nmax
102  *
103  * x > 0.0
104  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
105  */
106 int gsl_sf_bessel_Yn_array(const int nmin, const int nmax, const double x, double * result_array);
107 
108 
109 /* Regular modified Bessel function I_0(x)
110  *
111  * exceptions: GSL_EOVRFLW
112  */
113 int gsl_sf_bessel_I0_e(const double x, gsl_sf_result * result);
114 double gsl_sf_bessel_I0(const double x);
115 
116 
117 /* Regular modified Bessel function I_1(x)
118  *
119  * exceptions: GSL_EOVRFLW, GSL_EUNDRFLW
120  */
121 int gsl_sf_bessel_I1_e(const double x, gsl_sf_result * result);
122 double gsl_sf_bessel_I1(const double x);
123 
124 
125 /* Regular modified Bessel function I_n(x)
126  *
127  * exceptions: GSL_EOVRFLW, GSL_EUNDRFLW
128  */
129 int gsl_sf_bessel_In_e(const int n, const double x, gsl_sf_result * result);
130 double gsl_sf_bessel_In(const int n, const double x);
131 
132 
133 /* Regular modified Bessel function  I_n(x) for n=nmin,...,nmax
134  *
135  * nmin >=0, nmax >= nmin
136  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
137  */
138 int gsl_sf_bessel_In_array(const int nmin, const int nmax, const double x, double * result_array);
139 
140 
141 /* Scaled regular modified Bessel function
142  *  exp(-|x|) I_0(x)
143  *
144  * exceptions: none
145  */
146 int gsl_sf_bessel_I0_scaled_e(const double x, gsl_sf_result * result);
147 double gsl_sf_bessel_I0_scaled(const double x);
148 
149 
150 /* Scaled regular modified Bessel function
151  *  exp(-|x|) I_1(x)
152  *
153  * exceptions: GSL_EUNDRFLW
154  */
155 int gsl_sf_bessel_I1_scaled_e(const double x, gsl_sf_result * result);
156 double gsl_sf_bessel_I1_scaled(const double x);
157 
158 
159 /* Scaled regular modified Bessel function
160  *  exp(-|x|) I_n(x)
161  *
162  * exceptions: GSL_EUNDRFLW
163  */
164 int gsl_sf_bessel_In_scaled_e(int n, const double x, gsl_sf_result * result);
165 double gsl_sf_bessel_In_scaled(const int n, const double x);
166 
167 
168 /* Scaled regular modified Bessel function
169  *  exp(-|x|) I_n(x)  for n=nmin,...,nmax
170  *
171  * nmin >=0, nmax >= nmin
172  * exceptions: GSL_EUNDRFLW
173  */
174 int gsl_sf_bessel_In_scaled_array(const int nmin, const int nmax, const double x, double * result_array);
175 
176 
177 /* Irregular modified Bessel function K_0(x)
178  *
179  * x > 0.0
180  * exceptions: GSL_EDOM, GSL_EUNDRFLW
181  */
182 int gsl_sf_bessel_K0_e(const double x, gsl_sf_result * result);
183 double gsl_sf_bessel_K0(const double x);
184 
185 
186 /* Irregular modified Bessel function K_1(x)
187  *
188  * x > 0.0
189  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
190  */
191 int gsl_sf_bessel_K1_e(const double x, gsl_sf_result * result);
192 double gsl_sf_bessel_K1(const double x);
193 
194 
195 /* Irregular modified Bessel function K_n(x)
196  *
197  * x > 0.0
198  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
199  */
200 int gsl_sf_bessel_Kn_e(const int n, const double x, gsl_sf_result * result);
201 double gsl_sf_bessel_Kn(const int n, const double x);
202 
203 
204 /* Irregular modified Bessel function  K_n(x)  for n=nmin,...,nmax
205  *
206  * x > 0.0, nmin >=0, nmax >= nmin
207  * exceptions: GSL_EDOM, GSL_EOVRFLW, GSL_EUNDRFLW
208  */
209 int gsl_sf_bessel_Kn_array(const int nmin, const int nmax, const double x, double * result_array);
210 
211 
212 /* Scaled irregular modified Bessel function
213  *  exp(x) K_0(x)
214  *
215  * x > 0.0
216  * exceptions: GSL_EDOM
217  */
218 int gsl_sf_bessel_K0_scaled_e(const double x, gsl_sf_result * result);
219 double gsl_sf_bessel_K0_scaled(const double x);
220 
221 
222 /* Scaled irregular modified Bessel function
223  *  exp(x) K_1(x)
224  *
225  * x > 0.0
226  * exceptions: GSL_EDOM, GSL_EUNDRFLW
227  */
228 int gsl_sf_bessel_K1_scaled_e(const double x, gsl_sf_result * result);
229 double gsl_sf_bessel_K1_scaled(const double x);
230 
231 
232 /* Scaled irregular modified Bessel function
233  *  exp(x) K_n(x)
234  *
235  * x > 0.0
236  * exceptions: GSL_EDOM, GSL_EUNDRFLW
237  */
238 int gsl_sf_bessel_Kn_scaled_e(int n, const double x, gsl_sf_result * result);
239 double gsl_sf_bessel_Kn_scaled(const int n, const double x);
240 
241 
242 /* Scaled irregular modified Bessel function  exp(x) K_n(x)  for n=nmin,...,nmax
243  *
244  * x > 0.0, nmin >=0, nmax >= nmin
245  * exceptions: GSL_EDOM, GSL_EUNDRFLW
246  */
247 int gsl_sf_bessel_Kn_scaled_array(const int nmin, const int nmax, const double x, double * result_array);
248 
249 
250 /* Regular spherical Bessel function j_0(x) = sin(x)/x
251  *
252  * exceptions: none
253  */
254 int gsl_sf_bessel_j0_e(const double x, gsl_sf_result * result);
255 double gsl_sf_bessel_j0(const double x);
256 
257 
258 /* Regular spherical Bessel function j_1(x) = (sin(x)/x - cos(x))/x
259  *
260  * exceptions: GSL_EUNDRFLW
261  */
262 int gsl_sf_bessel_j1_e(const double x, gsl_sf_result * result);
263 double gsl_sf_bessel_j1(const double x);
264 
265 
266 /* Regular spherical Bessel function j_2(x) = ((3/x^2 - 1)sin(x) - 3cos(x)/x)/x
267  *
268  * exceptions: GSL_EUNDRFLW
269  */
270 int gsl_sf_bessel_j2_e(const double x, gsl_sf_result * result);
271 double gsl_sf_bessel_j2(const double x);
272 
273 
274 /* Regular spherical Bessel function j_l(x)
275  *
276  * l >= 0, x >= 0.0
277  * exceptions: GSL_EDOM, GSL_EUNDRFLW
278  */
279 int gsl_sf_bessel_jl_e(const int l, const double x, gsl_sf_result * result);
280 double gsl_sf_bessel_jl(const int l, const double x);
281 
282 
283 /* Regular spherical Bessel function j_l(x) for l=0,1,...,lmax
284  *
285  * exceptions: GSL_EDOM, GSL_EUNDRFLW
286  */
287 int gsl_sf_bessel_jl_array(const int lmax, const double x, double * result_array);
288 
289 
290 /* Regular spherical Bessel function j_l(x) for l=0,1,...,lmax
291  * Uses Steed's method.
292  *
293  * exceptions: GSL_EDOM, GSL_EUNDRFLW
294  */
295 int gsl_sf_bessel_jl_steed_array(const int lmax, const double x, double * jl_x_array);
296 
297 
298 /* Irregular spherical Bessel function y_0(x)
299  *
300  * exceptions: none
301  */
302 int gsl_sf_bessel_y0_e(const double x, gsl_sf_result * result);
303 double gsl_sf_bessel_y0(const double x);
304 
305 
306 /* Irregular spherical Bessel function y_1(x)
307  *
308  * exceptions: GSL_EUNDRFLW
309  */
310 int gsl_sf_bessel_y1_e(const double x, gsl_sf_result * result);
311 double gsl_sf_bessel_y1(const double x);
312 
313 
314 /* Irregular spherical Bessel function y_2(x)
315  *
316  * exceptions: GSL_EUNDRFLW
317  */
318 int gsl_sf_bessel_y2_e(const double x, gsl_sf_result * result);
319 double gsl_sf_bessel_y2(const double x);
320 
321 
322 /* Irregular spherical Bessel function y_l(x)
323  *
324  * exceptions: GSL_EUNDRFLW
325  */
326 int gsl_sf_bessel_yl_e(int l, const double x, gsl_sf_result * result);
327 double gsl_sf_bessel_yl(const int l, const double x);
328 
329 
330 /* Irregular spherical Bessel function y_l(x) for l=0,1,...,lmax
331  *
332  * exceptions: GSL_EUNDRFLW
333  */
334 int gsl_sf_bessel_yl_array(const int lmax, const double x, double * result_array);
335 
336 
337 /* Regular scaled modified spherical Bessel function
338  *
339  * Exp[-|x|] i_0(x)
340  *
341  * exceptions: none
342  */
343 int gsl_sf_bessel_i0_scaled_e(const double x, gsl_sf_result * result);
344 double gsl_sf_bessel_i0_scaled(const double x);
345 
346 
347 /* Regular scaled modified spherical Bessel function
348  *
349  * Exp[-|x|] i_1(x)
350  *
351  * exceptions: GSL_EUNDRFLW
352  */
353 int gsl_sf_bessel_i1_scaled_e(const double x, gsl_sf_result * result);
354 double gsl_sf_bessel_i1_scaled(const double x);
355 
356 
357 /* Regular scaled modified spherical Bessel function
358  *
359  * Exp[-|x|] i_2(x)
360  *
361  * exceptions: GSL_EUNDRFLW
362  */
363 int gsl_sf_bessel_i2_scaled_e(const double x, gsl_sf_result * result);
364 double gsl_sf_bessel_i2_scaled(const double x);
365 
366 
367 /* Regular scaled modified spherical Bessel functions
368  *
369  * Exp[-|x|] i_l(x)
370  *
371  * i_l(x) = Sqrt[Pi/(2x)] BesselI[l+1/2,x]
372  *
373  * l >= 0
374  * exceptions: GSL_EDOM, GSL_EUNDRFLW
375  */
376 int gsl_sf_bessel_il_scaled_e(const int l, double x, gsl_sf_result * result);
377 double gsl_sf_bessel_il_scaled(const int l, const double x);
378 
379 
380 /* Regular scaled modified spherical Bessel functions
381  *
382  * Exp[-|x|] i_l(x)
383  * for l=0,1,...,lmax
384  *
385  * exceptions: GSL_EUNDRFLW
386  */
387 int gsl_sf_bessel_il_scaled_array(const int lmax, const double x, double * result_array);
388 
389 
390 /* Irregular scaled modified spherical Bessel function
391  * Exp[x] k_0(x)
392  *
393  * x > 0.0
394  * exceptions: GSL_EDOM, GSL_EUNDRFLW
395  */
396 int gsl_sf_bessel_k0_scaled_e(const double x, gsl_sf_result * result);
397 double gsl_sf_bessel_k0_scaled(const double x);
398 
399 
400 /* Irregular modified spherical Bessel function
401  * Exp[x] k_1(x)
402  *
403  * x > 0.0
404  * exceptions: GSL_EDOM, GSL_EUNDRFLW, GSL_EOVRFLW
405  */
406 int gsl_sf_bessel_k1_scaled_e(const double x, gsl_sf_result * result);
407 double gsl_sf_bessel_k1_scaled(const double x);
408 
409 
410 /* Irregular modified spherical Bessel function
411  * Exp[x] k_2(x)
412  *
413  * x > 0.0
414  * exceptions: GSL_EDOM, GSL_EUNDRFLW, GSL_EOVRFLW
415  */
416 int gsl_sf_bessel_k2_scaled_e(const double x, gsl_sf_result * result);
417 double gsl_sf_bessel_k2_scaled(const double x);
418 
419 
420 /* Irregular modified spherical Bessel function
421  * Exp[x] k_l[x]
422  *
423  * k_l(x) = Sqrt[Pi/(2x)] BesselK[l+1/2,x]
424  *
425  * exceptions: GSL_EDOM, GSL_EUNDRFLW
426  */
427 int gsl_sf_bessel_kl_scaled_e(int l, const double x, gsl_sf_result * result);
428 double gsl_sf_bessel_kl_scaled(const int l, const double x);
429 
430 
431 /* Irregular scaled modified spherical Bessel function
432  * Exp[x] k_l(x)
433  *
434  * for l=0,1,...,lmax
435  * exceptions: GSL_EDOM, GSL_EUNDRFLW
436  */
437 int gsl_sf_bessel_kl_scaled_array(const int lmax, const double x, double * result_array);
438 
439 
440 /* Regular cylindrical Bessel function J_nu(x)
441  *
442  * exceptions: GSL_EDOM, GSL_EUNDRFLW
443  */
444 int gsl_sf_bessel_Jnu_e(const double nu, const double x, gsl_sf_result * result);
445 double gsl_sf_bessel_Jnu(const double nu, const double x);
446 
447 
448 /* Irregular cylindrical Bessel function Y_nu(x)
449  *
450  * exceptions:
451  */
452 int gsl_sf_bessel_Ynu_e(double nu, double x, gsl_sf_result * result);
453 double gsl_sf_bessel_Ynu(const double nu, const double x);
454 
455 
456 /* Regular cylindrical Bessel function J_nu(x)
457  * evaluated at a series of x values. The array
458  * contains the x values. They are assumed to be
459  * strictly ordered and positive. The array is
460  * over-written with the values of J_nu(x_i).
461  *
462  * exceptions: GSL_EDOM, GSL_EINVAL
463  */
464 int gsl_sf_bessel_sequence_Jnu_e(double nu, gsl_mode_t mode, size_t size, double * v);
465 
466 
467 /* Scaled modified cylindrical Bessel functions
468  *
469  * Exp[-|x|] BesselI[nu, x]
470  * x >= 0, nu >= 0
471  *
472  * exceptions: GSL_EDOM
473  */
474 int gsl_sf_bessel_Inu_scaled_e(double nu, double x, gsl_sf_result * result);
475 double gsl_sf_bessel_Inu_scaled(double nu, double x);
476 
477 
478 /* Modified cylindrical Bessel functions
479  *
480  * BesselI[nu, x]
481  * x >= 0, nu >= 0
482  *
483  * exceptions: GSL_EDOM, GSL_EOVRFLW
484  */
485 int gsl_sf_bessel_Inu_e(double nu, double x, gsl_sf_result * result);
486 double gsl_sf_bessel_Inu(double nu, double x);
487 
488 
489 /* Scaled modified cylindrical Bessel functions
490  *
491  * Exp[+|x|] BesselK[nu, x]
492  * x > 0, nu >= 0
493  *
494  * exceptions: GSL_EDOM
495  */
496 int gsl_sf_bessel_Knu_scaled_e(const double nu, const double x, gsl_sf_result * result);
497 double gsl_sf_bessel_Knu_scaled(const double nu, const double x);
498 
499 int gsl_sf_bessel_Knu_scaled_e10_e(const double nu, const double x, gsl_sf_result_e10 * result);
500 
501 /* Modified cylindrical Bessel functions
502  *
503  * BesselK[nu, x]
504  * x > 0, nu >= 0
505  *
506  * exceptions: GSL_EDOM, GSL_EUNDRFLW
507  */
508 int gsl_sf_bessel_Knu_e(const double nu, const double x, gsl_sf_result * result);
509 double gsl_sf_bessel_Knu(const double nu, const double x);
510 
511 
512 /* Logarithm of modified cylindrical Bessel functions.
513  *
514  * Log[BesselK[nu, x]]
515  * x > 0, nu >= 0
516  *
517  * exceptions: GSL_EDOM
518  */
519 int gsl_sf_bessel_lnKnu_e(const double nu, const double x, gsl_sf_result * result);
520 double gsl_sf_bessel_lnKnu(const double nu, const double x);
521 
522 
523 /* s'th positive zero of the Bessel function J_0(x).
524  *
525  * exceptions:
526  */
527 int gsl_sf_bessel_zero_J0_e(unsigned int s, gsl_sf_result * result);
528 double gsl_sf_bessel_zero_J0(unsigned int s);
529 
530 
531 /* s'th positive zero of the Bessel function J_1(x).
532  *
533  * exceptions:
534  */
535 int gsl_sf_bessel_zero_J1_e(unsigned int s, gsl_sf_result * result);
536 double gsl_sf_bessel_zero_J1(unsigned int s);
537 
538 
539 /* s'th positive zero of the Bessel function J_nu(x).
540  *
541  * exceptions:
542  */
543 int gsl_sf_bessel_zero_Jnu_e(double nu, unsigned int s, gsl_sf_result * result);
544 double gsl_sf_bessel_zero_Jnu(double nu, unsigned int s);
545 
546 
547 __END_DECLS
548 
549 #endif /* __GSL_SF_BESSEL_H__ */
550