Comparison of quad_qag with bessel functions using integer or float orders
Subject: Comparison of quad_qag with bessel functions using integer or float orders
From: Edwin Woollett
Date: Sat, 5 Nov 2011 16:14:59 -0700
We compare using integer order with using
float order for the realpart (bessel_u (n,x)) and
the imagpart (bessel_u (n,x)) for u = j, i, k, and y,
when using quag_qag.
-------------------------
case bessel_j(0,z):
can use either integer or float for
realpart(bessel_j(0,z)), for z = x or %i*x
imagpart(bessel_j(0,z)), for z = x or %i*x
--------------------------------
case bessel_j(1,z):
can use either integer or float for
realpart(bessel_j(1,z)) for z = x or %i*x
but not for imagpart:
arg x: (yes,yes)
(%i110) quad_qag(imagpart(bessel_j(1,x)),x,0,3,3,limit=700);
(%o110) [0.0,0.0,31,0]
(%i111) quad_qag(imagpart(bessel_j(1.0,x)),x,0,3,3,limit=700);
(%o111) [0.0,0.0,31,0]
arg %i*x: (no,yes)
(%i112) quad_qag(imagpart(bessel_j(1,%i*x)),x,0,3,3,limit=700);
(%o112) quad_qag(bessel_j(1,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i113) quad_qag(imagpart(bessel_j(1.0,%i*x)),x,0,3,3,limit=700);
(%o113) [3.880792585865024,4.3085452826219503E-14,31,0]
-------------------------------
case bessel_j(2,z):
can use either integer or float for
realpart(bessel_j(2,z)), for z = x or %i*x
or for
imagpart(bessel_j(2,z)), for z = x or %i*x
------------------------------
case bessel_j(3,z)
can use either integer or float for
realpart(bessel_j(3,z)) for z = x or %i*x
but not for imagpart:
arg x: (yes,yes)
(%i126) quad_qag(imagpart(bessel_j(3,x)),x,0,3,3,limit=700);
(%o126) [0.0,0.0,31,0]
(%i127) quad_qag(imagpart(bessel_j(3.0,x)),x,0,3,3,limit=700);
(%o127) [0.0,0.0,31,0]
arg %i*x: (no,yes)
(%i128) quad_qag(imagpart(bessel_j(3,%i*x)),x,0,3,3,limit=700);
(%o128) quad_qag(bessel_j(3,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i129) quad_qag(imagpart(bessel_j(3.0,%i*x)),x,0,3,3,limit=700);
(%o129) [-0.60963229599488,6.7682781156861259E-15,31,0]
---------------------------------------------
================================
==============================
bessel_i(n,z) comparison of integer or float
-------------------------------------------
case bessel_i(0,z) (yes,yes), (yes,no), (no,yes)
realpart (yes,yes) and (no,yes)
arg x: (yes,yes)
(%i45) quad_qag(realpart(bessel_i(0,x)),x,0,3,3,limit=700);
(%o45) [6.160961491502393,6.8400413016949048E-14,31,0]
(%i48) quad_qag(realpart(bessel_i(0.0,x)),x,0,3,3,limit=700);
(%o48) [6.160961491502393,6.8400413016949048E-14,31,0]
arg %i*x: (no,yes)
(%i46) quad_qag(realpart(bessel_i(0,%i*x)),x,0,3,3,limit=700);
(%o46) quad_qag(bessel_i(0,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i47) quad_qag(realpart(bessel_i(0.0,%i*x)),x,0,3,3,limit=700);
(%o47) [1.387567252009865,1.7248014181203231E-14,31,0]
imagpart (yes,yes) and (yes,no)
arg x: (yes,yes)
(%i49) quad_qag(imagpart(bessel_i(0,x)),x,0,3,3,limit=700);
(%o49) [0.0,0.0,31,0]
(%i50) quad_qag(imagpart(bessel_i(0.0,x)),x,0,3,3,limit=700);
(%o50) [0.0,0.0,31,0]
arg %i*x (yes, no and hangs)
(%i51) quad_qag(imagpart(bessel_i(0,%i*x)),x,0,3,3,limit=700);
(%o51) [0.0,0.0,31,0]
(%i52) quad_qag(imagpart(bessel_i(0.0,%i*x)),x,0,3,3,limit=700);
(%o52) quad_qag('imagpart(bessel_i(0.0,%i*x)),x,0,3,3,epsrel = 1.0E-8,
epsabs = 0.0,limit = 700)
(had to use ctrl-g to end hang on last one)
-----------------------------------------------
case bessel_i(1,z) (yes,yes), (yes,no) and (no,yes)
realpart (yes,yes) and (yes,no)
arg x: (yes,yes)
(%i53) quad_qag(realpart(bessel_i(1,x)),x,0,3,3,limit=700);
(%o53) [3.880792585865024,4.3085452826219503E-14,31,0]
(%i54) quad_qag(realpart(bessel_i(1.0,x)),x,0,3,3,limit=700);
(%o54) [3.880792585865024,4.3085452826219503E-14,31,0]
arg %i*x (yes, no and hangs)
(%i55) quad_qag(realpart(bessel_i(1,%i*x)),x,0,3,3,limit=700);
(%o55) [0.0,0.0,31,0]
(%i56) quad_qag(realpart(bessel_i(1.0,%i*x)),x,0,3,3,limit=700);
(%o56) quad_qag('realpart(bessel_i(1.0,%i*x)),x,0,3,3,epsrel = 1.0E-8,
epsabs = 0.0,limit = 700)
(had to use ctrl-g to escape last one)
imagpart (yes,yes) and (no,yes)
arg x :(yes,yes)
(%i57) quad_qag(imagpart(bessel_i(1,x)),x,0,3,3,limit=700);
(%o57) [0.0,0.0,31,0]
(%i58) quad_qag(imagpart(bessel_i(1.0,x)),x,0,3,3,limit=700);
(%o58) [0.0,0.0,31,0]
arg %i*x: (no,yes)
(%i59) quad_qag(imagpart(bessel_i(1,%i*x)),x,0,3,3,limit=700);
(%o59) quad_qag(bessel_i(1,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i60) quad_qag(imagpart(bessel_i(1.0,%i*x)),x,0,3,3,limit=700);
(%o60) [1.260051954901933,1.3989386925560662E-14,31,0]
-------------------------------------------
case bessel_i(2,z) (yes,yes) and (no,yes)
can use either integer or float for
imagpart(bessel_j(1,z)) for z = x or %i*x
but not for realpart:
arg x: (yes,yes)
(%i61) quad_qag(realpart(bessel_i(2,x)),x,0,3,3,limit=700);
(%o61) [1.745778943302825,1.9382039787605724E-14,31,0]
(%i62) quad_qag(realpart(bessel_i(2.0,x)),x,0,3,3,limit=700);
(%o62) [1.745778943302825,1.9382039787605724E-14,31,0]
arg %i*x: (no, yes)
(%i63) quad_qag(realpart(bessel_i(2,%i*x)),x,0,3,3,limit=700);
(%o63) quad_qag(bessel_i(2,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i64) quad_qag(realpart(bessel_i(2.0,%i*x)),x,0,3,3,limit=700);
(%o64) [-0.70944933495799,7.8764698647536805E-15,31,0]
-----------------------------------
case bessel_i(3,z)
can use either integer or float for
realpart(bessel_j(1,z)) for z = x or %i*x
but not for imagpart:
arg x: (yes,yes)
(%i73) quad_qag(imagpart(bessel_i(3,x)),x,0,3,3,limit=700);
(%o73) [0.0,0.0,31,0]
(%i74) quad_qag(imagpart(bessel_i(3.0,x)),x,0,3,3,limit=700);
(%o74) [0.0,0.0,31,0]
arg %i*x (no,yes)
(%i75) quad_qag(imagpart(bessel_i(3,%i*x)),x,0,3,3,limit=700);
(%o75) quad_qag(bessel_i(3,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i76) quad_qag(imagpart(bessel_i(3.0,%i*x)),x,0,3,3,limit=700);
(%o76) [-0.28786943373015,3.1959927341301972E-15,31,0]
-----------------------------------------------
case bessel_i(4,z)
can use either integer or float for
imagpart(bessel_j(1,z)) for z = x or %i*x
but not for realpart:
arg x: (yes,yes)
(%i130) quad_qag(realpart(bessel_i(4,x)),x,0,3,3,limit=700);
(%o130) [0.17372831568919,1.9287717610748791E-15,31,0]
(%i131) quad_qag(realpart(bessel_i(4.0,x)),x,0,3,3,limit=700);
(%o131) [0.17372831568919,1.9287717610748791E-15,31,0]
arg %i*x: (no,yes)
(%i132) quad_qag(realpart(bessel_i(4,%i*x)),x,0,3,3,limit=700);
(%o132) quad_qag(bessel_i(4,%i*x),x,0,3,3,epsrel = 1.0E-8,epsabs = 0.0,
limit = 700)
(%i133) quad_qag(realpart(bessel_i(4.0,%i*x)),x,0,3,3,limit=700);
(%o133) [0.091323890447489,1.0138988587314746E-15,31,0]
================================
==================================
bessel_k and bessel_y (for examples tested) work
equally well with either integer or float orders.
=====================
Ted Woollett