GNU Octave  4.4.1
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
besselj.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 1997-2018 John W. Eaton
4 
5 This file is part of Octave.
6 
7 Octave is free software: you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11 
12 Octave is distributed in the hope that it will be useful, but
13 WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with Octave; see the file COPYING. If not, see
19 <https://www.gnu.org/licenses/>.
20 
21 */
22 
23 #if defined (HAVE_CONFIG_H)
24 # include "config.h"
25 #endif
26 
27 #include "lo-specfun.h"
28 #include "quit.h"
29 
30 #include "defun.h"
31 #include "error.h"
32 #include "ovl.h"
33 #include "utils.h"
34 
35 enum bessel_type
36 {
37  BESSEL_J,
38  BESSEL_Y,
39  BESSEL_I,
40  BESSEL_K,
41  BESSEL_H1,
42  BESSEL_H2
43 };
44 
45 #define DO_BESSEL(type, alpha, x, scaled, ierr, result) \
46  do \
47  { \
48  switch (type) \
49  { \
50  case BESSEL_J: \
51  result = octave::math::besselj (alpha, x, scaled, ierr); \
52  break; \
53  \
54  case BESSEL_Y: \
55  result = octave::math::bessely (alpha, x, scaled, ierr); \
56  break; \
57  \
58  case BESSEL_I: \
59  result = octave::math::besseli (alpha, x, scaled, ierr); \
60  break; \
61  \
62  case BESSEL_K: \
63  result = octave::math::besselk (alpha, x, scaled, ierr); \
64  break; \
65  \
66  case BESSEL_H1: \
67  result = octave::math::besselh1 (alpha, x, scaled, ierr); \
68  break; \
69  \
70  case BESSEL_H2: \
71  result = octave::math::besselh2 (alpha, x, scaled, ierr); \
72  break; \
73  \
74  default: \
75  break; \
76  } \
77  } \
78  while (0)
79 
80 octave_value_list
81 do_bessel (enum bessel_type type, const char *fn,
82  const octave_value_list& args, int nargout)
83 {
84  int nargin = args.length ();
85 
86  if (nargin < 2 || nargin > 3)
87  print_usage ();
88 
89  octave_value_list retval (nargout > 1 ? 2 : 1);
90 
91  bool scaled = false;
92  if (nargin == 3)
93  {
94  octave_value opt_arg = args(2);
95  bool rpt_error = false;
96 
97  if (! opt_arg.is_scalar_type ())
98  rpt_error = true;
99  else if (opt_arg.isnumeric ())
100  {
101  double opt_val = opt_arg.double_value ();
102  if (opt_val != 0.0 && opt_val != 1.0)
103  rpt_error = true;
104  scaled = (opt_val == 1.0);
105  }
106  else if (opt_arg.islogical ())
107  scaled = opt_arg.bool_value ();
108 
109  if (rpt_error)
110  error ("%s: OPT must be 0 (or false) or 1 (or true)", fn);
111  }
112 
113  octave_value alpha_arg = args(0);
114  octave_value x_arg = args(1);
115 
116  if (alpha_arg.is_single_type () || x_arg.is_single_type ())
117  {
118  if (alpha_arg.is_scalar_type ())
119  {
120  float alpha = args(0).xfloat_value ("%s: ALPHA must be a scalar or matrix", fn);
121 
122  if (x_arg.is_scalar_type ())
123  {
124  FloatComplex x = x_arg.xfloat_complex_value ("%s: X must be a scalar or matrix", fn);
125 
126  octave_idx_type ierr;
127  octave_value result;
128 
129  DO_BESSEL (type, alpha, x, scaled, ierr, result);
130 
131  retval(0) = result;
132  if (nargout > 1)
133  retval(1) = static_cast<float> (ierr);
134  }
135  else
136  {
137  FloatComplexNDArray x = x_arg.xfloat_complex_array_value ("%s: X must be a scalar or matrix", fn);
138 
139  Array<octave_idx_type> ierr;
140  octave_value result;
141 
142  DO_BESSEL (type, alpha, x, scaled, ierr, result);
143 
144  retval(0) = result;
145  if (nargout > 1)
146  retval(1) = NDArray (ierr);
147  }
148  }
149  else
150  {
151  dim_vector dv0 = args(0).dims ();
152  dim_vector dv1 = args(1).dims ();
153 
154  bool args0_is_row_vector = (dv0(1) == dv0.numel ());
155  bool args1_is_col_vector = (dv1(0) == dv1.numel ());
156 
157  if (args0_is_row_vector && args1_is_col_vector)
158  {
159  FloatRowVector ralpha = args(0).xfloat_row_vector_value ("%s: ALPHA must be a scalar or matrix", fn);
160 
161  FloatComplexColumnVector cx =
162  x_arg.xfloat_complex_column_vector_value ("%s: X must be a scalar or matrix", fn);
163 
164  Array<octave_idx_type> ierr;
165  octave_value result;
166 
167  DO_BESSEL (type, ralpha, cx, scaled, ierr, result);
168 
169  retval(0) = result;
170  if (nargout > 1)
171  retval(1) = NDArray (ierr);
172  }
173  else
174  {
175  FloatNDArray alpha = args(0).xfloat_array_value ("%s: ALPHA must be a scalar or matrix", fn);
176 
177  if (x_arg.is_scalar_type ())
178  {
179  FloatComplex x = x_arg.xfloat_complex_value ("%s: X must be a scalar or matrix", fn);
180 
181  Array<octave_idx_type> ierr;
182  octave_value result;
183 
184  DO_BESSEL (type, alpha, x, scaled, ierr, result);
185 
186  retval(0) = result;
187  if (nargout > 1)
188  retval(1) = NDArray (ierr);
189  }
190  else
191  {
192  FloatComplexNDArray x = x_arg.xfloat_complex_array_value ("%s: X must be a scalar or matrix", fn);
193 
194  Array<octave_idx_type> ierr;
195  octave_value result;
196 
197  DO_BESSEL (type, alpha, x, scaled, ierr, result);
198 
199  retval(0) = result;
200  if (nargout > 1)
201  retval(1) = NDArray (ierr);
202  }
203  }
204  }
205  }
206  else
207  {
208  if (alpha_arg.is_scalar_type ())
209  {
210  double alpha = args(0).xdouble_value ("%s: ALPHA must be a scalar or matrix", fn);
211 
212  if (x_arg.is_scalar_type ())
213  {
214  Complex x = x_arg.xcomplex_value ("%s: X must be a scalar or matrix", fn);
215 
216  octave_idx_type ierr;
217  octave_value result;
218 
219  DO_BESSEL (type, alpha, x, scaled, ierr, result);
220 
221  retval(0) = result;
222  if (nargout > 1)
223  retval(1) = static_cast<double> (ierr);
224  }
225  else
226  {
227  ComplexNDArray x = x_arg.xcomplex_array_value ("%s: X must be a scalar or matrix", fn);
228 
229  Array<octave_idx_type> ierr;
230  octave_value result;
231 
232  DO_BESSEL (type, alpha, x, scaled, ierr, result);
233 
234  retval(0) = result;
235  if (nargout > 1)
236  retval(1) = NDArray (ierr);
237  }
238  }
239  else
240  {
241  dim_vector dv0 = args(0).dims ();
242  dim_vector dv1 = args(1).dims ();
243 
244  bool args0_is_row_vector = (dv0(1) == dv0.numel ());
245  bool args1_is_col_vector = (dv1(0) == dv1.numel ());
246 
247  if (args0_is_row_vector && args1_is_col_vector)
248  {
249  RowVector ralpha = args(0).xrow_vector_value ("%s: ALPHA must be a scalar or matrix", fn);
250 
251  ComplexColumnVector cx =
252  x_arg.xcomplex_column_vector_value ("%s: X must be a scalar or matrix", fn);
253 
254  Array<octave_idx_type> ierr;
255  octave_value result;
256 
257  DO_BESSEL (type, ralpha, cx, scaled, ierr, result);
258 
259  retval(0) = result;
260  if (nargout > 1)
261  retval(1) = NDArray (ierr);
262  }
263  else
264  {
265  NDArray alpha = args(0).xarray_value ("%s: ALPHA must be a scalar or matrix", fn);
266 
267  if (x_arg.is_scalar_type ())
268  {
269  Complex x = x_arg.xcomplex_value ("%s: X must be a scalar or matrix", fn);
270 
271  Array<octave_idx_type> ierr;
272  octave_value result;
273 
274  DO_BESSEL (type, alpha, x, scaled, ierr, result);
275 
276  retval(0) = result;
277  if (nargout > 1)
278  retval(1) = NDArray (ierr);
279  }
280  else
281  {
282  ComplexNDArray x = x_arg.xcomplex_array_value ("%s: X must be a scalar or matrix", fn);
283 
284  Array<octave_idx_type> ierr;
285  octave_value result;
286 
287  DO_BESSEL (type, alpha, x, scaled, ierr, result);
288 
289  retval(0) = result;
290  if (nargout > 1)
291  retval(1) = NDArray (ierr);
292  }
293  }
294  }
295  }
296 
297  return retval;
298 }
299 
300 DEFUN (besselj, args, nargout,
301  doc: /* -*- texinfo -*-
302 @deftypefn {} {@var{J} =} besselj (@var{alpha}, @var{x})
303 @deftypefnx {} {@var{J} =} besselj (@var{alpha}, @var{x}, @var{opt})
304 @deftypefnx {} {[@var{J}, @var{ierr}] =} besselj (@dots{})
305 Compute Bessel functions of the first kind.
306 
307 The order of the Bessel function @var{alpha} must be real. The points for
308 evaluation @var{x} may be complex.
309 
310 If the optional argument @var{opt} is 1 or true, the result @var{J} is
311 multiplied by @w{@code{exp (-abs (imag (@var{x})))}}.
312 
313 If @var{alpha} is a scalar, the result is the same size as @var{x}. If @var{x}
314 is a scalar, the result is the same size as @var{alpha}. If @var{alpha} is a
315 row vector and @var{x} is a column vector, the result is a matrix with
316 @code{length (@var{x})} rows and @code{length (@var{alpha})} columns.
317 Otherwise, @var{alpha} and @var{x} must conform and the result will be the same
318 size.
319 
320 If requested, @var{ierr} contains the following status information and is the
321 same size as the result.
322 
323 @enumerate 0
324 @item
325 Normal return.
326 
327 @item
328 Input error, return @code{NaN}.
329 
330 @item
331 Overflow, return @code{Inf}.
332 
333 @item
334 Loss of significance by argument reduction results in less than half of machine
335 accuracy.
336 
337 @item
338 Loss of significance by argument reduction, output may be inaccurate.
339 
340 @item
341 Error---no computation, algorithm termination condition not met, return
342 @code{NaN}.
343 @end enumerate
344 
345 @seealso{bessely, besseli, besselk, besselh}
346 @end deftypefn */)
347 {
348  return do_bessel (BESSEL_J, "besselj", args, nargout);
349 }
350 
351 /*
352 %!# Function besselj is tested along with other bessels at the end of this file
353 */
354 
355 DEFUN (bessely, args, nargout,
356  doc: /* -*- texinfo -*-
357 @deftypefn {} {@var{Y} =} bessely (@var{alpha}, @var{x})
358 @deftypefnx {} {@var{Y} =} bessely (@var{alpha}, @var{x}, @var{opt})
359 @deftypefnx {} {[@var{Y}, @var{ierr}] =} bessely (@dots{})
360 Compute Bessel functions of the second kind.
361 
362 The order of the Bessel function @var{alpha} must be real. The points for
363 evaluation @var{x} may be complex.
364 
365 If the optional argument @var{opt} is 1 or true, the result @var{Y} is
366 multiplied by @w{@code{exp (-abs (imag (@var{x})))}}.
367 
368 If @var{alpha} is a scalar, the result is the same size as @var{x}. If @var{x}
369 is a scalar, the result is the same size as @var{alpha}. If @var{alpha} is a
370 row vector and @var{x} is a column vector, the result is a matrix with
371 @code{length (@var{x})} rows and @code{length (@var{alpha})} columns.
372 Otherwise, @var{alpha} and @var{x} must conform and the result will be the same
373 size.
374 
375 If requested, @var{ierr} contains the following status information and is the
376 same size as the result.
377 
378 @enumerate 0
379 @item
380 Normal return.
381 
382 @item
383 Input error, return @code{NaN}.
384 
385 @item
386 Overflow, return @code{Inf}.
387 
388 @item
389 Loss of significance by argument reduction results in less than half of machine
390 accuracy.
391 
392 @item
393 Complete loss of significance by argument reduction, return @code{NaN}.
394 
395 @item
396 Error---no computation, algorithm termination condition not met, return
397 @code{NaN}.
398 @end enumerate
399 
400 @seealso{besselj, besseli, besselk, besselh}
401 @end deftypefn */)
402 {
403  return do_bessel (BESSEL_Y, "bessely", args, nargout);
404 }
405 
406 /*
407 %!# Function bessely is tested along with other bessels at the end of this file
408 */
409 
410 DEFUN (besseli, args, nargout,
411  doc: /* -*- texinfo -*-
412 @deftypefn {} {@var{I} =} besseli (@var{alpha}, @var{x})
413 @deftypefnx {} {@var{I} =} besseli (@var{alpha}, @var{x}, @var{opt})
414 @deftypefnx {} {[@var{I}, @var{ierr}] =} besseli (@dots{})
415 Compute modified Bessel functions of the first kind.
416 
417 The order of the Bessel function @var{alpha} must be real. The points for
418 evaluation @var{x} may be complex.
419 
420 If the optional argument @var{opt} is 1 or true, the result @var{I} is
421 multiplied by @w{@code{exp (-abs (real (@var{x})))}}.
422 
423 If @var{alpha} is a scalar, the result is the same size as @var{x}. If @var{x}
424 is a scalar, the result is the same size as @var{alpha}. If @var{alpha} is a
425 row vector and @var{x} is a column vector, the result is a matrix with
426 @code{length (@var{x})} rows and @code{length (@var{alpha})} columns.
427 Otherwise, @var{alpha} and @var{x} must conform and the result will be the same
428 size.
429 
430 If requested, @var{ierr} contains the following status information and is the
431 same size as the result.
432 
433 @enumerate 0
434 @item
435 Normal return.
436 
437 @item
438 Input error, return @code{NaN}.
439 
440 @item
441 Overflow, return @code{Inf}.
442 
443 @item
444 Loss of significance by argument reduction results in less than half of machine
445 accuracy.
446 
447 @item
448 Complete loss of significance by argument reduction, return @code{NaN}.
449 
450 @item
451 Error---no computation, algorithm termination condition not met, return
452 @code{NaN}.
453 @end enumerate
454 
455 @seealso{besselk, besselj, bessely, besselh}
456 @end deftypefn */)
457 
458 {
459  return do_bessel (BESSEL_I, "besseli", args, nargout);
460 }
461 
462 /*
463 %!# Function besseli is tested along with other bessels at the end of this file
464 */
465 
466 DEFUN (besselk, args, nargout,
467  doc: /* -*- texinfo -*-
468 @deftypefn {} {@var{K} =} besselk (@var{alpha}, @var{x})
469 @deftypefnx {} {@var{K} =} besselk (@var{alpha}, @var{x}, @var{opt})
470 @deftypefnx {} {[@var{K}, @var{ierr}] =} besselk (@dots{})
471 
472 Compute modified Bessel functions of the second kind.
473 
474 The order of the Bessel function @var{alpha} must be real. The points for
475 evaluation @var{x} may be complex.
476 
477 If the optional argument @var{opt} is 1 or true, the result @var{K} is
478 multiplied by @w{@code{exp (@var{x})}}.
479 
480 If @var{alpha} is a scalar, the result is the same size as @var{x}. If @var{x}
481 is a scalar, the result is the same size as @var{alpha}. If @var{alpha} is a
482 row vector and @var{x} is a column vector, the result is a matrix with
483 @code{length (@var{x})} rows and @code{length (@var{alpha})} columns.
484 Otherwise, @var{alpha} and @var{x} must conform and the result will be the same
485 size.
486 
487 If requested, @var{ierr} contains the following status information and is the
488 same size as the result.
489 
490 @enumerate 0
491 @item
492 Normal return.
493 
494 @item
495 Input error, return @code{NaN}.
496 
497 @item
498 Overflow, return @code{Inf}.
499 
500 @item
501 Loss of significance by argument reduction results in less than half of machine
502 accuracy.
503 
504 @item
505 Complete loss of significance by argument reduction, return @code{NaN}.
506 
507 @item
508 Error---no computation, algorithm termination condition not met, return
509 @code{NaN}.
510 @end enumerate
511 
512 @seealso{besseli, besselj, bessely, besselh}
513 @end deftypefn */)
514 {
515  return do_bessel (BESSEL_K, "besselk", args, nargout);
516 }
517 
518 /*
519 %!# Function besselk is tested along with other bessels at the end of this file
520 */
521 
522 DEFUN (besselh, args, nargout,
523  doc: /* -*- texinfo -*-
524 @deftypefn {} {@var{H} =} besselh (@var{alpha}, @var{x})
525 @deftypefnx {} {@var{H} =} besselh (@var{alpha}, @var{k}, @var{x})
526 @deftypefnx {} {@var{H} =} besselh (@var{alpha}, @var{k}, @var{x}, @var{opt})
527 @deftypefnx {} {[@var{H}, @var{ierr}] =} besselh (@dots{})
528 Compute Bessel functions of the third kind (Hankel functions).
529 
530 The order of the Bessel function @var{alpha} must be real. The kind of Hankel
531 function is specified by @var{k} and may be either first (@var{k} = 1) or
532 second (@var{k} = 2). The default is Hankel functions of the first kind. The
533 points for evaluation @var{x} may be complex.
534 
535 If the optional argument @var{opt} is 1 or true, the result is multiplied
536 by @code{exp (-I*@var{x})} for @var{k} = 1 or @code{exp (I*@var{x})} for
537 @var{k} = 2.
538 
539 If @var{alpha} is a scalar, the result is the same size as @var{x}. If @var{x}
540 is a scalar, the result is the same size as @var{alpha}. If @var{alpha} is a
541 row vector and @var{x} is a column vector, the result is a matrix with
542 @code{length (@var{x})} rows and @code{length (@var{alpha})} columns.
543 Otherwise, @var{alpha} and @var{x} must conform and the result will be the same
544 size.
545 
546 If requested, @var{ierr} contains the following status information and is the
547 same size as the result.
548 
549 @enumerate 0
550 @item
551 Normal return.
552 
553 @item
554 Input error, return @code{NaN}.
555 
556 @item
557 Overflow, return @code{Inf}.
558 
559 @item
560 Loss of significance by argument reduction results in less than half of machine
561 accuracy.
562 
563 @item
564 Complete loss of significance by argument reduction, return @code{NaN}.
565 
566 @item
567 Error---no computation, algorithm termination condition not met, return
568 @code{NaN}.
569 @end enumerate
570 
571 @seealso{besselj, bessely, besseli, besselk}
572 @end deftypefn */)
573 {
574  int nargin = args.length ();
575 
576  if (nargin < 2 || nargin > 4)
577  print_usage ();
578 
579  octave_value_list retval;
580 
581  if (nargin == 2)
582  {
583  retval = do_bessel (BESSEL_H1, "besselh", args, nargout);
584  }
585  else
586  {
587  octave_idx_type kind = args(1).xint_value ("besselh: invalid value of K");
588 
589  octave_value_list tmp_args;
590 
591  if (nargin == 4)
592  tmp_args(2) = args(3);
593 
594  tmp_args(1) = args(2);
595  tmp_args(0) = args(0);
596 
597  if (kind == 1)
598  retval = do_bessel (BESSEL_H1, "besselh", tmp_args, nargout);
599  else if (kind == 2)
600  retval = do_bessel (BESSEL_H2, "besselh", tmp_args, nargout);
601  else
602  error ("besselh: K must be 1 or 2");
603  }
604 
605  return retval;
606 }
607 
608 /*
609 %!# Function besselh is tested along with other bessels at the end of this file
610 */
611 
612 DEFUN (airy, args, nargout,
613  doc: /* -*- texinfo -*-
614 @deftypefn {} {[@var{a}, @var{ierr}] =} airy (@var{k}, @var{z}, @var{opt})
615 Compute Airy functions of the first and second kind, and their derivatives.
616 
617 @example
618 @group
619  K Function Scale factor (if "opt" is supplied)
620 --- -------- ---------------------------------------
621  0 Ai (Z) exp ((2/3) * Z * sqrt (Z))
622  1 dAi(Z)/dZ exp ((2/3) * Z * sqrt (Z))
623  2 Bi (Z) exp (-abs (real ((2/3) * Z * sqrt (Z))))
624  3 dBi(Z)/dZ exp (-abs (real ((2/3) * Z * sqrt (Z))))
625 @end group
626 @end example
627 
628 The function call @code{airy (@var{z})} is equivalent to
629 @code{airy (0, @var{z})}.
630 
631 The result is the same size as @var{z}.
632 
633 If requested, @var{ierr} contains the following status information and
634 is the same size as the result.
635 
636 @enumerate 0
637 @item
638 Normal return.
639 
640 @item
641 Input error, return @code{NaN}.
642 
643 @item
644 Overflow, return @code{Inf}.
645 
646 @item
647 Loss of significance by argument reduction results in less than half
648  of machine accuracy.
649 
650 @item
651 Loss of significance by argument reduction, output may be inaccurate.
652 
653 @item
654 Error---no computation, algorithm termination condition not met,
655 return @code{NaN}.
656 @end enumerate
657 @end deftypefn */)
658 {
659  int nargin = args.length ();
660 
661  if (nargin < 1 || nargin > 3)
662  print_usage ();
663 
664  octave_value_list retval (nargout > 1 ? 2 : 1);
665 
666  int kind = 0;
667  if (nargin > 1)
668  {
669  kind = args(0).xint_value ("airy: K must be an integer value");
670 
671  if (kind < 0 || kind > 3)
672  error ("airy: K must be 0, 1, 2, or 3");
673  }
674 
675  bool scale = (nargin == 3);
676 
677  int idx = (nargin == 1 ? 0 : 1);
678 
679  if (args(idx).is_single_type ())
680  {
681  FloatComplexNDArray z = args(idx).xfloat_complex_array_value ("airy: Z must be a complex matrix");
682 
683  Array<octave_idx_type> ierr;
684  octave_value result;
685 
686  if (kind > 1)
687  result = octave::math::biry (z, kind == 3, scale, ierr);
688  else
689  result = octave::math::airy (z, kind == 1, scale, ierr);
690 
691  retval(0) = result;
692  if (nargout > 1)
693  retval(1) = NDArray (ierr);
694  }
695  else
696  {
697  ComplexNDArray z = args(idx).xcomplex_array_value ("airy: Z must be a complex matrix");
698 
699  Array<octave_idx_type> ierr;
700  octave_value result;
701 
702  if (kind > 1)
703  result = octave::math::biry (z, kind == 3, scale, ierr);
704  else
705  result = octave::math::airy (z, kind == 1, scale, ierr);
706 
707  retval(0) = result;
708  if (nargout > 1)
709  retval(1) = NDArray (ierr);
710  }
711 
712  return retval;
713 }
714 
715 /*
716 FIXME: Function airy does not yet have BIST tests
717 */
718 
719 /*
720 ## Test values computed with GP/PARI version 2.3.3
721 %!shared alpha, x, jx, yx, ix, kx, nix
722 %!
723 %! ## Bessel functions, even order, positive and negative x
724 %! alpha = 2; x = 1.25;
725 %! jx = 0.1710911312405234823613091417;
726 %! yx = -1.193199310178553861283790424;
727 %! ix = 0.2220184483766341752692212604;
728 %! kx = 0.9410016167388185767085460540;
729 %!
730 %!assert (besselj (alpha,x), jx, 100*eps)
731 %!assert (bessely (alpha,x), yx, 100*eps)
732 %!assert (besseli (alpha,x), ix, 100*eps)
733 %!assert (besselk (alpha,x), kx, 100*eps)
734 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
735 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
736 %!
737 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
738 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
739 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
740 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
741 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
742 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
743 %!
744 %!assert (besselj (-alpha,x), jx, 100*eps)
745 %!assert (bessely (-alpha,x), yx, 100*eps)
746 %!assert (besseli (-alpha,x), ix, 100*eps)
747 %!assert (besselk (-alpha,x), kx, 100*eps)
748 %!assert (besselh (-alpha,1,x), jx + I*yx, 100*eps)
749 %!assert (besselh (-alpha,2,x), jx - I*yx, 100*eps)
750 %!
751 %!assert (besselj (-alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
752 %!assert (bessely (-alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
753 %!assert (besseli (-alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
754 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
755 %!assert (besselh (-alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
756 %!assert (besselh (-alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
757 %!
758 %! x *= -1;
759 %! yx = -1.193199310178553861283790424 + 0.3421822624810469647226182835*I;
760 %! kx = 0.9410016167388185767085460540 - 0.6974915263814386815610060884*I;
761 %!
762 %!assert (besselj (alpha,x), jx, 100*eps)
763 %!assert (bessely (alpha,x), yx, 100*eps)
764 %!assert (besseli (alpha,x), ix, 100*eps)
765 %!assert (besselk (alpha,x), kx, 100*eps)
766 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
767 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
768 %!
769 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
770 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
771 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
772 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
773 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
774 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
775 %!
776 %! ## Bessel functions, odd order, positive and negative x
777 %! alpha = 3; x = 2.5;
778 %! jx = 0.2166003910391135247666890035;
779 %! yx = -0.7560554967536709968379029772;
780 %! ix = 0.4743704087780355895548240179;
781 %! kx = 0.2682271463934492027663765197;
782 %!
783 %!assert (besselj (alpha,x), jx, 100*eps)
784 %!assert (bessely (alpha,x), yx, 100*eps)
785 %!assert (besseli (alpha,x), ix, 100*eps)
786 %!assert (besselk (alpha,x), kx, 100*eps)
787 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
788 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
789 %!
790 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
791 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
792 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
793 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
794 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
795 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
796 %!
797 %!assert (besselj (-alpha,x), -jx, 100*eps)
798 %!assert (bessely (-alpha,x), -yx, 100*eps)
799 %!assert (besseli (-alpha,x), ix, 100*eps)
800 %!assert (besselk (-alpha,x), kx, 100*eps)
801 %!assert (besselh (-alpha,1,x), -(jx + I*yx), 100*eps)
802 %!assert (besselh (-alpha,2,x), -(jx - I*yx), 100*eps)
803 %!
804 %!assert (besselj (-alpha,x,1), -jx*exp(-abs(imag(x))), 100*eps)
805 %!assert (bessely (-alpha,x,1), -yx*exp(-abs(imag(x))), 100*eps)
806 %!assert (besseli (-alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
807 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
808 %!assert (besselh (-alpha,1,x,1), -(jx + I*yx)*exp(-I*x), 100*eps)
809 %!assert (besselh (-alpha,2,x,1), -(jx - I*yx)*exp(I*x), 100*eps)
810 %!
811 %! x *= -1;
812 %! jx = -jx;
813 %! yx = 0.7560554967536709968379029772 - 0.4332007820782270495333780070*I;
814 %! ix = -ix;
815 %! kx = -0.2682271463934492027663765197 - 1.490278591297463775542004240*I;
816 %!
817 %!assert (besselj (alpha,x), jx, 100*eps)
818 %!assert (bessely (alpha,x), yx, 100*eps)
819 %!assert (besseli (alpha,x), ix, 100*eps)
820 %!assert (besselk (alpha,x), kx, 100*eps)
821 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
822 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
823 %!
824 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
825 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
826 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
827 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
828 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
829 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
830 %!
831 %! ## Bessel functions, fractional order, positive and negative x
832 %!
833 %! alpha = 3.5; x = 2.75;
834 %! jx = 0.1691636439842384154644784389;
835 %! yx = -0.8301381935499356070267953387;
836 %! ix = 0.3930540878794826310979363668;
837 %! kx = 0.2844099013460621170288192503;
838 %!
839 %!assert (besselj (alpha,x), jx, 100*eps)
840 %!assert (bessely (alpha,x), yx, 100*eps)
841 %!assert (besseli (alpha,x), ix, 100*eps)
842 %!assert (besselk (alpha,x), kx, 100*eps)
843 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
844 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
845 %!
846 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
847 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
848 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
849 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
850 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
851 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
852 %!
853 %! nix = 0.2119931212254662995364461998;
854 %!
855 %!assert (besselj (-alpha,x), yx, 100*eps)
856 %!assert (bessely (-alpha,x), -jx, 100*eps)
857 %!assert (besseli (-alpha,x), nix, 100*eps)
858 %!assert (besselk (-alpha,x), kx, 100*eps)
859 %!assert (besselh (-alpha,1,x), -I*(jx + I*yx), 100*eps)
860 %!assert (besselh (-alpha,2,x), I*(jx - I*yx), 100*eps)
861 %!
862 %!assert (besselj (-alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
863 %!assert (bessely (-alpha,x,1), -jx*exp(-abs(imag(x))), 100*eps)
864 %!assert (besseli (-alpha,x,1), nix*exp(-abs(real(x))), 100*eps)
865 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
866 %!assert (besselh (-alpha,1,x,1), -I*(jx + I*yx)*exp(-I*x), 100*eps)
867 %!assert (besselh (-alpha,2,x,1), I*(jx - I*yx)*exp(I*x), 100*eps)
868 %!
869 %! x *= -1;
870 %! jx *= -I;
871 %! yx = -0.8301381935499356070267953387*I;
872 %! ix *= -I;
873 %! kx = -0.9504059335995575096509874508*I;
874 %!
875 %!assert (besselj (alpha,x), jx, 100*eps)
876 %!assert (bessely (alpha,x), yx, 100*eps)
877 %!assert (besseli (alpha,x), ix, 100*eps)
878 %!assert (besselk (alpha,x), kx, 100*eps)
879 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
880 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
881 %!
882 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
883 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
884 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
885 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
886 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
887 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
888 %!
889 %! ## Bessel functions, even order, complex x
890 %!
891 %! alpha = 2; x = 1.25 + 3.625 * I;
892 %! jx = -1.299533366810794494030065917 + 4.370833116012278943267479589*I;
893 %! yx = -4.370357232383223896393056727 - 1.283083391453582032688834041*I;
894 %! ix = -0.6717801680341515541002273932 - 0.2314623443930774099910228553*I;
895 %! kx = -0.01108009888623253515463783379 + 0.2245218229358191588208084197*I;
896 %!
897 %!assert (besselj (alpha,x), jx, 100*eps)
898 %!assert (bessely (alpha,x), yx, 100*eps)
899 %!assert (besseli (alpha,x), ix, 100*eps)
900 %!assert (besselk (alpha,x), kx, 100*eps)
901 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
902 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
903 %!
904 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
905 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
906 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
907 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
908 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
909 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
910 %!
911 %!assert (besselj (-alpha,x), jx, 100*eps)
912 %!assert (bessely (-alpha,x), yx, 100*eps)
913 %!assert (besseli (-alpha,x), ix, 100*eps)
914 %!assert (besselk (-alpha,x), kx, 100*eps)
915 %!assert (besselh (-alpha,1,x), jx + I*yx, 100*eps)
916 %!assert (besselh (-alpha,2,x), jx - I*yx, 100*eps)
917 %!
918 %!assert (besselj (-alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
919 %!assert (bessely (-alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
920 %!assert (besseli (-alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
921 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
922 %!assert (besselh (-alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
923 %!assert (besselh (-alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
924 %!
925 %! ## Bessel functions, odd order, complex x
926 %!
927 %! alpha = 3; x = 2.5 + 1.875 * I;
928 %! jx = 0.1330721523048277493333458596 + 0.5386295217249660078754395597*I;
929 %! yx = -0.6485072392105829901122401551 + 0.2608129289785456797046996987*I;
930 %! ix = -0.6182064685486998097516365709 + 0.4677561094683470065767989920*I;
931 %! kx = -0.1568585587733540007867882337 - 0.05185853709490846050505141321*I;
932 %!
933 %!assert (besselj (alpha,x), jx, 100*eps)
934 %!assert (bessely (alpha,x), yx, 100*eps)
935 %!assert (besseli (alpha,x), ix, 100*eps)
936 %!assert (besselk (alpha,x), kx, 100*eps)
937 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
938 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
939 %!
940 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
941 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
942 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
943 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
944 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
945 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
946 %!
947 %!assert (besselj (-alpha,x), -jx, 100*eps)
948 %!assert (bessely (-alpha,x), -yx, 100*eps)
949 %!assert (besseli (-alpha,x), ix, 100*eps)
950 %!assert (besselk (-alpha,x), kx, 100*eps)
951 %!assert (besselh (-alpha,1,x), -(jx + I*yx), 100*eps)
952 %!assert (besselh (-alpha,2,x), -(jx - I*yx), 100*eps)
953 %!
954 %!assert (besselj (-alpha,x,1), -jx*exp(-abs(imag(x))), 100*eps)
955 %!assert (bessely (-alpha,x,1), -yx*exp(-abs(imag(x))), 100*eps)
956 %!assert (besseli (-alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
957 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
958 %!assert (besselh (-alpha,1,x,1), -(jx + I*yx)*exp(-I*x), 100*eps)
959 %!assert (besselh (-alpha,2,x,1), -(jx - I*yx)*exp(I*x), 100*eps)
960 %!
961 %! ## Bessel functions, fractional order, complex x
962 %!
963 %! alpha = 3.5; x = 1.75 + 4.125 * I;
964 %! jx = -3.018566131370455929707009100 - 0.7585648436793900607704057611*I;
965 %! yx = 0.7772278839106298215614791107 - 3.018518722313849782683792010*I;
966 %! ix = 0.2100873577220057189038160913 - 0.6551765604618246531254970926*I;
967 %! kx = 0.1757147290513239935341488069 + 0.08772348296883849205562558311*I;
968 %!
969 %!assert (besselj (alpha,x), jx, 100*eps)
970 %!assert (bessely (alpha,x), yx, 100*eps)
971 %!assert (besseli (alpha,x), ix, 100*eps)
972 %!assert (besselk (alpha,x), kx, 100*eps)
973 %!assert (besselh (alpha,1,x), jx + I*yx, 100*eps)
974 %!assert (besselh (alpha,2,x), jx - I*yx, 100*eps)
975 %!
976 %!assert (besselj (alpha,x,1), jx*exp(-abs(imag(x))), 100*eps)
977 %!assert (bessely (alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
978 %!assert (besseli (alpha,x,1), ix*exp(-abs(real(x))), 100*eps)
979 %!assert (besselk (alpha,x,1), kx*exp(x), 100*eps)
980 %!assert (besselh (alpha,1,x,1), (jx + I*yx)*exp(-I*x), 100*eps)
981 %!assert (besselh (alpha,2,x,1), (jx - I*yx)*exp(I*x), 100*eps)
982 %!
983 %! nix = 0.09822388691172060573913739253 - 0.7110230642207380127317227407*I;
984 %!
985 %!assert (besselj (-alpha,x), yx, 100*eps)
986 %!assert (bessely (-alpha,x), -jx, 100*eps)
987 %!assert (besseli (-alpha,x), nix, 100*eps)
988 %!assert (besselk (-alpha,x), kx, 100*eps)
989 %!assert (besselh (-alpha,1,x), -I*(jx + I*yx), 100*eps)
990 %!assert (besselh (-alpha,2,x), I*(jx - I*yx), 100*eps)
991 %!
992 %!assert (besselj (-alpha,x,1), yx*exp(-abs(imag(x))), 100*eps)
993 %!assert (bessely (-alpha,x,1), -jx*exp(-abs(imag(x))), 100*eps)
994 %!assert (besseli (-alpha,x,1), nix*exp(-abs(real(x))), 100*eps)
995 %!assert (besselk (-alpha,x,1), kx*exp(x), 100*eps)
996 %!assert (besselh (-alpha,1,x,1), -I*(jx + I*yx)*exp(-I*x), 100*eps)
997 %!assert (besselh (-alpha,2,x,1), I*(jx - I*yx)*exp(I*x), 100*eps)
998 
999 Tests contributed by Robert T. Short.
1000 Tests are based on the properties and tables in A&S:
1001  Abramowitz and Stegun, "Handbook of Mathematical Functions",
1002  1972.
1003 
1004 For regular Bessel functions, there are 3 tests. These compare octave
1005 results against Tables 9.1, 9.2, and 9.4 in A&S. Tables 9.1 and 9.2
1006 are good to only a few decimal places, so any failures should be
1007 considered a broken implementation. Table 9.4 is an extended table
1008 for larger orders and arguments. There are some differences between
1009 Octave and Table 9.4, mostly in the last decimal place but in a very
1010 few instances the errors are in the last two places. The comparison
1011 tolerance has been changed to reflect this.
1012 
1013 Similarly for modified Bessel functions, there are 3 tests. These
1014 compare octave results against Tables 9.8, 9.9, and 9.11 in A&S.
1015 Tables 9.8 and 9.9 are good to only a few decimal places, so any
1016 failures should be considered a broken implementation. Table 9.11 is
1017 an extended table for larger orders and arguments. There are some
1018 differences between octave and Table 9.11, mostly in the last decimal
1019 place but in a very few instances the errors are in the last two
1020 places. The comparison tolerance has been changed to reflect this.
1021 
1022 For spherical Bessel functions, there are also three tests, comparing
1023 octave results to Tables 10.1, 10.2, and 10.4 in A&S. Very similar
1024 comments may be made here as in the previous lines. At this time,
1025 modified spherical Bessel function tests are not included.
1026 
1027 % Table 9.1 - J and Y for integer orders 0, 1, 2.
1028 % Compare against excerpts of Table 9.1, Abramowitz and Stegun.
1029 %!test
1030 %! n = 0:2;
1031 %! z = (0:2.5:17.5)';
1032 %!
1033 %! Jt = [[ 1.000000000000000, 0.0000000000, 0.0000000000];
1034 %! [-0.048383776468198, 0.4970941025, 0.4460590584];
1035 %! [-0.177596771314338, -0.3275791376, 0.0465651163];
1036 %! [ 0.266339657880378, 0.1352484276, -0.2302734105];
1037 %! [-0.245935764451348, 0.0434727462, 0.2546303137];
1038 %! [ 0.146884054700421, -0.1654838046, -0.1733614634];
1039 %! [-0.014224472826781, 0.2051040386, 0.0415716780];
1040 %! [-0.103110398228686, -0.1634199694, 0.0844338303]];
1041 %!
1042 %! Yt = [[-Inf, -Inf, -Inf ];
1043 %! [ 0.4980703596, 0.1459181380, -0.38133585 ];
1044 %! [-0.3085176252, 0.1478631434, 0.36766288 ];
1045 %! [ 0.1173132861, -0.2591285105, -0.18641422 ];
1046 %! [ 0.0556711673, 0.2490154242, -0.00586808 ];
1047 %! [-0.1712143068, -0.1538382565, 0.14660019 ];
1048 %! [ 0.2054642960, 0.0210736280, -0.20265448 ];
1049 %! [-0.1604111925, 0.0985727987, 0.17167666 ]];
1050 %!
1051 %! J = besselj (n,z);
1052 %! Y = bessely (n,z);
1053 %! assert (Jt(:,1), J(:,1), 0.5e-10);
1054 %! assert (Yt(:,1), Y(:,1), 0.5e-10);
1055 %! assert (Jt(:,2:3), J(:,2:3), 0.5e-10);
1056 
1057 Table 9.2 - J and Y for integer orders 3-9.
1058 
1059 %!test
1060 %! n = (3:9);
1061 %! z = (0:2:20).';
1062 %!
1063 %! Jt = [[ 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00, 0.0000e+00];
1064 %! [ 1.2894e-01, 3.3996e-02, 7.0396e-03, 1.2024e-03, 1.7494e-04, 2.2180e-05, 2.4923e-06];
1065 %! [ 4.3017e-01, 2.8113e-01, 1.3209e-01, 4.9088e-02, 1.5176e-02, 4.0287e-03, 9.3860e-04];
1066 %! [ 1.1477e-01, 3.5764e-01, 3.6209e-01, 2.4584e-01, 1.2959e-01, 5.6532e-02, 2.1165e-02];
1067 %! [-2.9113e-01,-1.0536e-01, 1.8577e-01, 3.3758e-01, 3.2059e-01, 2.2345e-01, 1.2632e-01];
1068 %! [ 5.8379e-02,-2.1960e-01,-2.3406e-01,-1.4459e-02, 2.1671e-01, 3.1785e-01, 2.9186e-01];
1069 %! [ 1.9514e-01, 1.8250e-01,-7.3471e-02,-2.4372e-01,-1.7025e-01, 4.5095e-02, 2.3038e-01];
1070 %! [-1.7681e-01, 7.6244e-02, 2.2038e-01, 8.1168e-02,-1.5080e-01,-2.3197e-01,-1.1431e-01];
1071 %! [-4.3847e-02,-2.0264e-01,-5.7473e-02, 1.6672e-01, 1.8251e-01,-7.0211e-03,-1.8953e-01];
1072 %! [ 1.8632e-01, 6.9640e-02,-1.5537e-01,-1.5596e-01, 5.1399e-02, 1.9593e-01, 1.2276e-01];
1073 %! [-9.8901e-02, 1.3067e-01, 1.5117e-01,-5.5086e-02,-1.8422e-01,-7.3869e-02, 1.2513e-01]];
1074 %!
1075 %! Yt = [[ -Inf, -Inf, -Inf, -Inf, -Inf, -Inf, -Inf];
1076 %! [-1.1278e+00,-2.7659e+00,-9.9360e+00,-4.6914e+01,-2.7155e+02,-1.8539e+03,-1.4560e+04];
1077 %! [-1.8202e-01,-4.8894e-01,-7.9585e-01,-1.5007e+00,-3.7062e+00,-1.1471e+01,-4.2178e+01];
1078 %! [ 3.2825e-01, 9.8391e-02,-1.9706e-01,-4.2683e-01,-6.5659e-01,-1.1052e+00,-2.2907e+00];
1079 %! [ 2.6542e-02, 2.8294e-01, 2.5640e-01, 3.7558e-02,-2.0006e-01,-3.8767e-01,-5.7528e-01];
1080 %! [-2.5136e-01,-1.4495e-01, 1.3540e-01, 2.8035e-01, 2.0102e-01, 1.0755e-03,-1.9930e-01];
1081 %! [ 1.2901e-01,-1.5122e-01,-2.2982e-01,-4.0297e-02, 1.8952e-01, 2.6140e-01, 1.5902e-01];
1082 %! [ 1.2350e-01, 2.0393e-01,-6.9717e-03,-2.0891e-01,-1.7209e-01, 3.6816e-02, 2.1417e-01];
1083 %! [-1.9637e-01,-7.3222e-05, 1.9633e-01, 1.2278e-01,-1.0425e-01,-2.1399e-01,-1.0975e-01];
1084 %! [ 3.3724e-02,-1.7722e-01,-1.1249e-01, 1.1472e-01, 1.8897e-01, 3.2253e-02,-1.6030e-01];
1085 %! [ 1.4967e-01, 1.2409e-01,-1.0004e-01,-1.7411e-01,-4.4312e-03, 1.7101e-01, 1.4124e-01]];
1086 %!
1087 %! n = (3:9);
1088 %! z = (0:2:20).';
1089 %! J = besselj (n,z);
1090 %! Y = bessely (n,z);
1091 %!
1092 %! assert (J(1,:), zeros (1, columns (J)));
1093 %! assert (J(2:end,:), Jt(2:end,:), -5e-5);
1094 %! assert (Yt(1,:), Y(1,:));
1095 %! assert (Y(2:end,:), Yt(2:end,:), -5e-5);
1096 
1097 Table 9.4 - J and Y for various integer orders and arguments.
1098 
1099 %!test
1100 %! Jt = [[ 7.651976866e-01, 2.238907791e-01, -1.775967713e-01, -2.459357645e-01, 5.581232767e-02, 1.998585030e-02];
1101 %! [ 2.497577302e-04, 7.039629756e-03, 2.611405461e-01, -2.340615282e-01, -8.140024770e-02, -7.419573696e-02];
1102 %! [ 2.630615124e-10, 2.515386283e-07, 1.467802647e-03, 2.074861066e-01, -1.138478491e-01, -5.473217694e-02];
1103 %! [ 2.297531532e-17, 7.183016356e-13, 4.796743278e-07, 4.507973144e-03, -1.082255990e-01, 1.519812122e-02];
1104 %! [ 3.873503009e-25, 3.918972805e-19, 2.770330052e-11, 1.151336925e-05, -1.167043528e-01, 6.221745850e-02];
1105 %! [ 3.482869794e-42, 3.650256266e-33, 2.671177278e-21, 1.551096078e-12, 4.843425725e-02, 8.146012958e-02];
1106 %! [ 1.107915851e-60, 1.196077458e-48, 8.702241617e-33, 6.030895312e-21, -1.381762812e-01, 7.270175482e-02];
1107 %! [ 2.906004948e-80, 3.224095839e-65, 2.294247616e-45, 1.784513608e-30, 1.214090219e-01, -3.869833973e-02];
1108 %! [ 8.431828790e-189, 1.060953112e-158, 6.267789396e-119, 6.597316064e-89, 1.115927368e-21, 9.636667330e-02]];
1109 %!
1110 %! Yt = [[ 8.825696420e-02, 5.103756726e-01, -3.085176252e-01, 5.567116730e-02, -9.806499547e-02, -7.724431337e-02]
1111 %! [-2.604058666e+02, -9.935989128e+00, -4.536948225e-01, 1.354030477e-01, -7.854841391e-02, -2.948019628e-02]
1112 %! [-1.216180143e+08, -1.291845422e+05, -2.512911010e+01, -3.598141522e-01, 5.723897182e-03, 5.833157424e-02]
1113 %! [-9.256973276e+14, -2.981023646e+10, -4.694049564e+04, -6.364745877e+00, 4.041280205e-02, 7.879068695e-02]
1114 %! [-4.113970315e+22, -4.081651389e+16, -5.933965297e+08, -1.597483848e+03, 1.644263395e-02, 5.124797308e-02]
1115 %! [-3.048128783e+39, -2.913223848e+30, -4.028568418e+18, -7.256142316e+09, -1.164572349e-01, 6.138839212e-03]
1116 %! [-7.184874797e+57, -6.661541235e+45, -9.216816571e+29, -1.362803297e+18, -4.530801120e-02, 4.074685217e-02]
1117 %! [-2.191142813e+77, -1.976150576e+62, -2.788837017e+42, -3.641066502e+27, -2.103165546e-01, 7.650526394e-02]
1118 %! [-3.775287810e+185, -3.000826049e+155, -5.084863915e+115, -4.849148271e+85, -3.293800188e+18, -1.669214114e-01]];
1119 %!
1120 %! n = [(0:5:20).';30;40;50;100];
1121 %! z = [1,2,5,10,50,100];
1122 %! J = besselj (n.', z.').';
1123 %! Y = bessely (n.', z.').';
1124 %! assert (J, Jt, -1e-9);
1125 %! assert (Y, Yt, -1e-9);
1126 
1127 Table 9.8 - I and K for integer orders 0, 1, 2.
1128 
1129 %!test
1130 %! n = 0:2;
1131 %! z1 = [0.1;2.5;5.0];
1132 %! z2 = [7.5;10.0;15.0;20.0];
1133 %! rtbl = [[ 0.9071009258 0.0452984468 0.1251041992 2.6823261023 10.890182683 1.995039646 ];
1134 %! [ 0.2700464416 0.2065846495 0.2042345837 0.7595486903 0.9001744239 0.759126289 ];
1135 %! [ 0.1835408126 0.1639722669 0.7002245988 0.5478075643 0.6002738588 0.132723593 ];
1136 %! [ 0.1483158301 0.1380412115 0.111504840 0.4505236991 0.4796689336 0.57843541 ];
1137 %! [ 0.1278333372 0.1212626814 0.103580801 0.3916319344 0.4107665704 0.47378525 ];
1138 %! [ 0.1038995314 0.1003741751 0.090516308 0.3210023535 0.3315348950 0.36520701 ];
1139 %! [ 0.0897803119 0.0875062222 0.081029690 0.2785448768 0.2854254970 0.30708743 ]];
1140 %!
1141 %! tbl = [besseli(n,z1,1), besselk(n,z1,1)];
1142 %! tbl(:,3) = tbl(:,3) .* (exp (z1) .* z1.^(-2));
1143 %! tbl(:,6) = tbl(:,6) .* (exp (-z1) .* z1.^(2));
1144 %! tbl = [tbl;[besseli(n,z2,1),besselk(n,z2,1)]];
1145 %!
1146 %! assert (tbl, rtbl, -2e-8);
1147 
1148 Table 9.9 - I and K for orders 3-9.
1149 
1150 %!test
1151 %! It = [[ 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00 0.0000e+00];
1152 %! [ 2.8791e-02 6.8654e-03 1.3298e-03 2.1656e-04 3.0402e-05 3.7487e-06 4.1199e-07];
1153 %! [ 6.1124e-02 2.5940e-02 9.2443e-03 2.8291e-03 7.5698e-04 1.7968e-04 3.8284e-05];
1154 %! [ 7.4736e-02 4.1238e-02 1.9752e-02 8.3181e-03 3.1156e-03 1.0484e-03 3.1978e-04];
1155 %! [ 7.9194e-02 5.0500e-02 2.8694e-02 1.4633e-02 6.7449e-03 2.8292e-03 1.0866e-03];
1156 %! [ 7.9830e-02 5.5683e-02 3.5284e-02 2.0398e-02 1.0806e-02 5.2694e-03 2.3753e-03];
1157 %! [ 7.8848e-02 5.8425e-02 3.9898e-02 2.5176e-02 1.4722e-02 8.0010e-03 4.0537e-03];
1158 %! [ 7.7183e-02 5.9723e-02 4.3056e-02 2.8969e-02 1.8225e-02 1.0744e-02 5.9469e-03];
1159 %! [ 7.5256e-02 6.0155e-02 4.5179e-02 3.1918e-02 2.1240e-02 1.3333e-02 7.9071e-03];
1160 %! [ 7.3263e-02 6.0059e-02 4.6571e-02 3.4186e-02 2.3780e-02 1.5691e-02 9.8324e-03];
1161 %! [ 7.1300e-02 5.9640e-02 4.7444e-02 3.5917e-02 2.5894e-02 1.7792e-02 1.1661e-02]];
1162 %!
1163 %! Kt = [[ Inf Inf Inf Inf Inf Inf Inf];
1164 %! [ 4.7836e+00 1.6226e+01 6.9687e+01 3.6466e+02 2.2576e+03 1.6168e+04 1.3160e+05];
1165 %! [ 1.6317e+00 3.3976e+00 8.4268e+00 2.4465e+01 8.1821e+01 3.1084e+02 1.3252e+03];
1166 %! [ 9.9723e-01 1.6798e+00 3.2370e+00 7.0748e+00 1.7387e+01 4.7644e+01 1.4444e+02];
1167 %! [ 7.3935e-01 1.1069e+00 1.8463e+00 3.4148e+00 6.9684e+00 1.5610e+01 3.8188e+01];
1168 %! [ 6.0028e-01 8.3395e-01 1.2674e+00 2.1014e+00 3.7891e+00 7.4062e+00 1.5639e+01];
1169 %! [ 5.1294e-01 6.7680e-01 9.6415e-01 1.4803e+00 2.4444e+00 4.3321e+00 8.2205e+00];
1170 %! [ 4.5266e-01 5.7519e-01 7.8133e-01 1.1333e+00 1.7527e+00 2.8860e+00 5.0510e+00];
1171 %! [ 4.0829e-01 5.0414e-01 6.6036e-01 9.1686e-01 1.3480e+00 2.0964e+00 3.4444e+00];
1172 %! [ 3.7411e-01 4.5162e-01 5.7483e-01 7.7097e-01 1.0888e+00 1.6178e+00 2.5269e+00];
1173 %! [ 3.4684e-01 4.1114e-01 5.1130e-01 6.6679e-01 9.1137e-01 1.3048e+00 1.9552e+00]];
1174 %!
1175 %! n = (3:9);
1176 %! z = (0:2:20).';
1177 %! I = besseli (n,z,1);
1178 %! K = besselk (n,z,1);
1179 %!
1180 %! assert (abs (I(1,:)), zeros (1, columns (I)));
1181 %! assert (I(2:end,:), It(2:end,:), -5e-5);
1182 %! assert (Kt(1,:), K(1,:));
1183 %! assert (K(2:end,:), Kt(2:end,:), -5e-5);
1184 
1185 Table 9.11 - I and K for various integer orders and arguments.
1186 
1187 %!test
1188 %! It = [[ 1.266065878e+00 2.279585302e+00 2.723987182e+01 2.815716628e+03 2.93255378e+20 1.07375171e+42 ];
1189 %! [ 2.714631560e-04 9.825679323e-03 2.157974547e+00 7.771882864e+02 2.27854831e+20 9.47009387e+41 ];
1190 %! [ 2.752948040e-10 3.016963879e-07 4.580044419e-03 2.189170616e+01 1.07159716e+20 6.49897552e+41 ];
1191 %! [ 2.370463051e-17 8.139432531e-13 1.047977675e-06 1.043714907e-01 3.07376455e+19 3.47368638e+41 ];
1192 %! [ 3.966835986e-25 4.310560576e-19 5.024239358e-11 1.250799736e-04 5.44200840e+18 1.44834613e+41 ];
1193 %! [ 3.539500588e-42 3.893519664e-33 3.997844971e-21 7.787569783e-12 4.27499365e+16 1.20615487e+40 ];
1194 %! [ 1.121509741e-60 1.255869192e-48 1.180426980e-32 2.042123274e-20 6.00717897e+13 3.84170550e+38 ];
1195 %! [ 2.934635309e-80 3.353042830e-65 2.931469647e-45 4.756894561e-30 1.76508024e+10 4.82195809e+36 ];
1196 %! [ 8.473674008e-189 1.082171475e-158 7.093551489e-119 1.082344202e-88 2.72788795e-16 4.64153494e+21 ]];
1197 %!
1198 %! Kt = [[ 4.210244382e-01 1.138938727e-01 3.691098334e-03 1.778006232e-05 3.41016774e-23 4.65662823e-45 ];
1199 %! [ 3.609605896e+02 9.431049101e+00 3.270627371e-02 5.754184999e-05 4.36718224e-23 5.27325611e-45 ];
1200 %! [ 1.807132899e+08 1.624824040e+05 9.758562829e+00 1.614255300e-03 9.15098819e-23 7.65542797e-45 ];
1201 %! [ 1.403066801e+15 4.059213332e+10 3.016976630e+04 2.656563849e-01 3.11621117e-22 1.42348325e-44 ];
1202 %! [ 6.294369360e+22 5.770856853e+16 4.827000521e+08 1.787442782e+02 1.70614838e-21 3.38520541e-44 ];
1203 %! [ 4.706145527e+39 4.271125755e+30 4.112132063e+18 2.030247813e+09 2.00581681e-19 3.97060205e-43 ];
1204 %! [ 1.114220651e+58 9.940839886e+45 1.050756722e+30 5.938224681e+17 1.29986971e-16 1.20842080e-41 ];
1205 %! [ 3.406896854e+77 2.979981740e+62 3.394322243e+42 2.061373775e+27 4.00601349e-13 9.27452265e-40 ];
1206 %! [ 5.900333184e+185 4.619415978e+155 7.039860193e+115 4.596674084e+85 1.63940352e+13 7.61712963e-25 ]];
1207 %!
1208 %! n = [(0:5:20).';30;40;50;100];
1209 %! z = [1,2,5,10,50,100];
1210 %! I = besseli (n.', z.').';
1211 %! K = besselk (n.', z.').';
1212 %! assert (I, It, -5e-9);
1213 %! assert (K, Kt, -5e-9);
1214 
1215 The next section checks that negative integer orders and positive
1216 integer orders are appropriately related.
1217 
1218 %!test
1219 %! n = (0:2:20);
1220 %! assert (besselj (n,1), besselj (-n,1), 1e-8);
1221 %! assert (-besselj (n+1,1), besselj (-n-1,1), 1e-8);
1222 
1223 besseli (n,z) = besseli (-n,z);
1224 
1225 %!test
1226 %! n = (0:2:20);
1227 %! assert (besseli (n,1), besseli (-n,1), 1e-8);
1228 
1229 Table 10.1 - j and y for integer orders 0, 1, 2.
1230 Compare against excerpts of Table 10.1, Abramowitz and Stegun.
1231 
1232 %!test
1233 %! n = (0:2);
1234 %! z = [0.1;(2.5:2.5:10.0).'];
1235 %!
1236 %! jt = [[ 9.9833417e-01 3.33000119e-02 6.6619061e-04 ];
1237 %! [ 2.3938886e-01 4.16212989e-01 2.6006673e-01 ];
1238 %! [-1.9178485e-01 -9.50894081e-02 1.3473121e-01 ];
1239 %! [ 1.2507e-01 -2.9542e-02 -1.3688e-01 ];
1240 %! [ -5.4402e-02 7.8467e-02 7.7942e-02 ]];
1241 %!
1242 %! yt = [[-9.9500417e+00 -1.0049875e+02 -3.0050125e+03 ];
1243 %! [ 3.2045745e-01 -1.1120588e-01 -4.5390450e-01 ];
1244 %! [-5.6732437e-02 1.8043837e-01 1.6499546e-01 ];
1245 %! [ -4.6218e-02 -1.3123e-01 -6.2736e-03 ];
1246 %! [ 8.3907e-02 6.2793e-02 -6.5069e-02 ]];
1247 %!
1248 %! j = sqrt ((pi/2)./z) .* besselj (n+1/2,z);
1249 %! y = sqrt ((pi/2)./z) .* bessely (n+1/2,z);
1250 %! assert (jt, j, -5e-5);
1251 %! assert (yt, y, -5e-5);
1252 
1253 Table 10.2 - j and y for orders 3-8.
1254 Compare against excerpts of Table 10.2, Abramowitzh and Stegun.
1255 
1256  Important note: In A&S, y_4(0.1) = -1.0507e+7, but Octave returns
1257  y_4(0.1) = -1.0508e+07 (-10507503.75). If I compute the same term using
1258  a series, the difference is in the eighth significant digit so I left
1259  the Octave results in place.
1260 
1261 %!test
1262 %! n = (3:8);
1263 %! z = (0:2.5:10).'; z(1) = 0.1;
1264 %!
1265 %! jt = [[ 9.5185e-06 1.0577e-07 9.6163e-10 7.3975e-12 4.9319e-14 2.9012e-16];
1266 %! [ 1.0392e-01 3.0911e-02 7.3576e-03 1.4630e-03 2.5009e-04 3.7516e-05];
1267 %! [ 2.2982e-01 1.8702e-01 1.0681e-01 4.7967e-02 1.7903e-02 5.7414e-03];
1268 %! [-6.1713e-02 7.9285e-02 1.5685e-01 1.5077e-01 1.0448e-01 5.8188e-02];
1269 %! [-3.9496e-02 -1.0559e-01 -5.5535e-02 4.4501e-02 1.1339e-01 1.2558e-01]];
1270 %!
1271 %! yt = [[-1.5015e+05 -1.0508e+07 -9.4553e+08 -1.0400e+11 -1.3519e+13 -2.0277e+15];
1272 %! [-7.9660e-01 -1.7766e+00 -5.5991e+00 -2.2859e+01 -1.1327e+02 -6.5676e+02];
1273 %! [-1.5443e-02 -1.8662e-01 -3.2047e-01 -5.1841e-01 -1.0274e+00 -2.5638e+00];
1274 %! [ 1.2705e-01 1.2485e-01 2.2774e-02 -9.1449e-02 -1.8129e-01 -2.7112e-01];
1275 %! [-9.5327e-02 -1.6599e-03 9.3834e-02 1.0488e-01 4.2506e-02 -4.1117e-02]];
1276 %!
1277 %! j = sqrt ((pi/2)./z) .* besselj (n+1/2,z);
1278 %! y = sqrt ((pi/2)./z) .* bessely (n+1/2,z);
1279 %!
1280 %! assert (jt, j, -5e-5);
1281 %! assert (yt, y, -5e-5);
1282 
1283 Table 10.4 - j and y for various integer orders and arguments.
1284 
1285 %!test
1286 %! jt = [[ 8.414709848e-01 4.546487134e-01 -1.917848549e-01 -5.440211109e-02 -5.247497074e-03 -5.063656411e-03];
1287 %! [ 9.256115861e-05 2.635169770e-03 1.068111615e-01 -5.553451162e-02 -2.004830056e-02 -9.290148935e-03];
1288 %! [ 7.116552640e-11 6.825300865e-08 4.073442442e-04 6.460515449e-02 -1.503922146e-02 -1.956578597e-04];
1289 %! [ 5.132686115e-18 1.606982166e-13 1.084280182e-07 1.063542715e-03 -1.129084539e-02 7.877261748e-03];
1290 %! [ 7.537795722e-26 7.632641101e-20 5.427726761e-12 2.308371961e-06 -1.578502990e-02 1.010767128e-02];
1291 %! [ 5.566831267e-43 5.836617888e-34 4.282730217e-22 2.512057385e-13 -1.494673454e-03 8.700628514e-03];
1292 %! [ 1.538210374e-61 1.660978779e-49 1.210347583e-33 8.435671634e-22 -2.606336952e-02 1.043410851e-02];
1293 %! [ 3.615274717e-81 4.011575290e-66 2.857479350e-46 2.230696023e-31 1.882910737e-02 5.797140882e-04];
1294 %! [7.444727742e-190 9.367832591e-160 5.535650303e-120 5.832040182e-90 1.019012263e-22 1.088047701e-02]];
1295 %!
1296 %! yt = [[ -5.403023059e-01 2.080734183e-01 -5.673243709e-02 8.390715291e-02 -1.929932057e-02 -8.623188723e-03]
1297 %! [ -9.994403434e+02 -1.859144531e+01 -3.204650467e-01 9.383354168e-02 -6.971131965e-04 3.720678486e-03]
1298 %! [ -6.722150083e+08 -3.554147201e+05 -2.665611441e+01 -1.724536721e-01 1.352468751e-02 1.002577737e-02]
1299 %! [ -6.298007233e+15 -1.012182944e+11 -6.288146513e+04 -3.992071745e+00 1.712319725e-02 6.258641510e-03]
1300 %! [ -3.239592219e+23 -1.605436493e+17 -9.267951403e+08 -1.211210605e+03 1.375953130e-02 5.631729379e-05]
1301 %! [ -2.946428547e+40 -1.407393871e+31 -7.760717570e+18 -6.908318646e+09 -2.241226812e-02 -5.412929349e-03]
1302 %! [ -8.028450851e+58 -3.720929322e+46 -2.055758716e+30 -1.510304919e+18 4.978797221e-05 -7.048420407e-04]
1303 %! [ -2.739192285e+78 -1.235021944e+63 -6.964109188e+42 -4.528227272e+27 -4.190000150e-02 1.074782297e-02]
1304 %! [-6.683079463e+186 -2.655955830e+156 -1.799713983e+116 -8.573226309e+85 -1.125692891e+18 -2.298385049e-02]];
1305 %!
1306 %! n = [(0:5:20).';30;40;50;100];
1307 %! z = [1,2,5,10,50,100];
1308 %! j = sqrt ((pi/2)./z) .* besselj ((n+1/2).', z.').';
1309 %! y = sqrt ((pi/2)./z) .* bessely ((n+1/2).', z.').';
1310 %! assert (j, jt, -1e-9);
1311 %! assert (y, yt, -1e-9);
1312 */
octave_value::xcomplex_column_vector_value
ComplexColumnVector xcomplex_column_vector_value(const char *fmt,...) const
octave_value
Definition: ov.h:75
octave_value::xfloat_complex_column_vector_value
FloatComplexColumnVector xfloat_complex_column_vector_value(const char *fmt,...) const
octave_value::islogical
bool islogical(void) const
Definition: ov.h:696
print_usage
OCTINTERP_API void print_usage(void)
Definition: defun.cc:54
DO_BESSEL
#define DO_BESSEL(type, alpha, x, scaled, ierr, result)
Definition: besselj.cc:45
octave::math::besselk
Complex besselk(double alpha, const Complex &x, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:818
DEFUN
#define DEFUN(name, args_name, nargout_name, doc)
Macro to define a builtin function.
Definition: defun.h:53
error
void error(const char *fmt,...)
Definition: error.cc:578
octave_value_list
Definition: ovl.h:39
octave_value::bool_value
bool bool_value(bool warn=false) const
Definition: ov.h:866
nargout
OCTAVE_EXPORT octave_value_list return the number of command line arguments passed to Octave If called with the optional argument the function xample nargout(@histc)
Definition: ov-usr-fcn.cc:997
octave::math::besseli
Complex besseli(double alpha, const Complex &x, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:817
do_bessel
octave_value_list do_bessel(enum bessel_type type, const char *fn, const octave_value_list &args, int nargout)
Definition: besselj.cc:81
octave_value::is_single_type
bool is_single_type(void) const
Definition: ov.h:651
octave_value::xcomplex_value
Complex xcomplex_value(const char *fmt,...) const
retval
octave_value retval
Definition: data.cc:6246
octave_value::xfloat_complex_value
FloatComplex xfloat_complex_value(const char *fmt,...) const
octave_value::xcomplex_array_value
ComplexNDArray xcomplex_array_value(const char *fmt,...) const
type
idx type
Definition: ov.cc:3114
result
With real return the complex result
Definition: data.cc:3260
octave_value::xfloat_complex_array_value
FloatComplexNDArray xfloat_complex_array_value(const char *fmt,...) const
dim_vector::numel
octave_idx_type numel(int n=0) const
Number of elements that a matrix with this dimensions would have.
Definition: dim-vector.h:362
octave_value::double_value
double double_value(bool frc_str_conv=false) const
Definition: ov.h:822
octave_value_list::length
octave_idx_type length(void) const
Definition: ovl.h:96
scale
void scale(Matrix &m, double x, double y, double z)
Definition: graphics.cc:5442
octave_idx_type
nargin
args.length() nargin
Definition: file-io.cc:589
FloatComplex
std::complex< float > FloatComplex
Definition: oct-cmplx.h:32
octave::math::airy
Complex airy(const Complex &z, bool deriv, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:132
Complex
std::complex< double > Complex
Definition: oct-cmplx.h:31
bessel_type
bessel_type
Definition: besselj.cc:35
dim_vector
Vector representing the dimensions (size) of an Array.
Definition: dim-vector.h:87
ierr
F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE const F77_INT F77_INT & ierr
Definition: lo-slatec-proto.h:55
octave_value::is_scalar_type
bool is_scalar_type(void) const
Definition: ov.h:717
octave_value::isnumeric
bool isnumeric(void) const
Definition: ov.h:723
x
F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T const F77_DBLE * x
Definition: lo-slatec-proto.h:55
octave::math::besselj
Complex besselj(double alpha, const Complex &x, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:815
octave::math::bessely
Complex bessely(double alpha, const Complex &x, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:816
octave::math::biry
Complex biry(const Complex &z, bool deriv, bool scaled, octave_idx_type &ierr)
Definition: lo-specfun.cc:1378