[SciPy-user] Possible bug in bessel special function h2vp
Robert Kern
rkern at ucsd.edu
Mon Oct 3 08:33:41 CDT 2005
Julian Swartz wrote:
> Hi,
>
> I've encountered what appears to be a bug in the function h2vp, which
> calculates the nth order derivative of a hankel function of the second kind.
> Here is an example with a second derivative being requested. The
> function returns two different values on successive calls. Any
> subsequent calls will yield the same result as the second call shown
> below. This problem is also evident in function h1vp, calculating the
> derivatives of a hankel function of the first kind.
>
>
>>>> from scipy import *
>>>> from scipy.special import *
>>>> n = 1.0
>>>> x = 2*0.3*pi
>>>> h2vp(n, x, 2)
> (2.9279764804171628e+165+0.22981683442664608j)
>>>> h2vp(n, x, 2)
> (-0.23515420288453626+0.54623524288750835j)
>
> Even more interesting is that given that the hankel function of the
> second kind is defined as:
>
> (2)
> H (x) = J (x) - jY (x)
> n n n
>
> Where J and Y are the standard cylindrical bessel functions of the first
> and second kind. One would therefore expect that the second derivative
> of the Hankel function should be the sum of the second derivatives of
> the ordinary bessel functions. Using the scipy functions jvp and yvp,
>
>>>> jvp(n, x, 2)
> -0.40831349982163045
>>>> yvp(n, x, 2)
> -0.18651604740953789
>
> Which is quite different to either of the results above.
>
> Any ideas on what the problem could be. (or am I just missing something?)
Not a clue. On OS X with the last CVS of scipy before the migration to
SVN, I don't see differences when repeating the call, but I do see a
difference from the stated definition for derivatives of order > 1:
In [1]: n = 1.0
In [2]: x = 2*0.3*pi
In [3]: special.h2vp(n,x,2)
Out[3]: (-0.26294530063194937+0.22981683442664608j)
In [4]: special.h2vp(n,x,2)
Out[4]: (-0.26294530063194937+0.22981683442664608j)
In [8]: special.jvp(n,x,2) - 1j*special.yvp(n,x,2)
Out[8]: (-0.40831349982163045+0.18651604740953792j)
In [23]: special.h2vp(n,x,1)
Out[23]: (-0.017916685502942564-0.58616758520739864j)
In [24]: special.jvp(n,x,1) - 1j*special.yvp(n,x,1)
Out[24]: (-0.017916685502942648-0.58616758520739876j)
In [21]: special.h2vp(n,x,3)
Out[21]: (-0.2382065153405013+0.02078698962754455j)
In [22]: special.jvp(n,x,3) - 1j*special.yvp(n,x,3)
Out[22]: (0.050806004596450884+0.10554382826770611j)
In [25]: special.hankel2(n,x)
Out[25]: (0.58147279675872476+0.1732031480684324j)
In [26]: special.jv(n,x) - 1j*special.yv(n,x)
Out[26]: (0.58147279675872487+0.17320314806843251j)
--
Robert Kern
rkern at ucsd.edu
"In the fields of hell where the grass grows high
Are the graves of dreams allowed to die."
-- Richard Harter
More information about the SciPy-user
mailing list