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