can't plot a Bessel function

[Dan Drake]

Certainly I must be doing something dumb, but I can't figure out what. I do

t = var('t')
plot(bessel_J(1, t), (t, 1, 10))

and I get:

Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/Users/dan/.sage/sage_notebook/worksheets/admin/1/code/27.py",
line 8, in <module>
plot(bessel_J(Integer(1), t), (t, Integer(1), Integer(10)))
File "/Applications/sage/local/lib/python2.5/site-packages/SQLAlchemy-0.4.6-py2.5.egg/",
line 1, in <module>

File "/Applications/sage/local/lib/python2.5/site-packages/sage/functions/special.py",
line 527, in bessel_J
b = K(nu.besselj(z))
File "real_mpfr.pyx", line 352, in
sage.rings.real_mpfr.RealField.__call__ (sage/rings/real_mpfr.c:4046)
File "real_mpfr.pyx", line 819, in
sage.rings.real_mpfr.RealNumber._set (sage/rings/real_mpfr.c:6834)
TypeError: Unable to convert x
(='1-1/8*t^2+1/192*t^4-1/9216*t^6+1/737280*t^8-1/88473600*t^10+1/14863564800*t^12-1/3329438515200*t^14+1/958878292377600*t^16+O(t^17)')
to real number.

This is with Sage 3.1.1 on OS X. I have no trouble evaluating
bessel_J, so why can't it plot the function?

Dan

--
--- Dan Drake <dr...@mathsci.kaist.ac.kr>
----- KAIST Department of Mathematical Sciences
------- http://math.kaist.ac.kr/~drake

[David Joyner]

I don't know what the problem is. Here is a workaround:

sage: list_plot([(t,bessel_J(1,1+9*t/100)) for t in range(100)],
plotjoined=True)

[John Cremona]

This works:

plot(lambda t:bessel_J(1, t), (1, 10))

so (1) a one-variable function is reequired, and the lambda
construction creates such a function from the2-variable bessel_J, and
(2) the range is a tuple (xmin,xmax) .

John Cremona

2008/9/11 Dan Drake <dr...@mathsci.kaist.ac.kr>:

[Jason Merrill]

On Sep 11, 9:07 am, "John Cremona" <john.crem...@gmail.com> wrote:
> This works:
>
> plot(lambda t:bessel_J(1, t), (1, 10))
>
> so (1) a one-variable function is reequired, and the lambda
> construction creates such a function from the2-variable bessel_J, and
> (2) the range is a tuple (xmin,xmax) .
>
> John Cremona

Is there a fundamental difference between bessel_J and, say, polylog,
for which I can do

sage: plot(polylog(1,x),(x,.1,.9))

and I get a perfectly nice looking plot?

I see that the implementation difference is that polylog is defined in
calculus.py, and is a class that derives from PrimitiveFunction,
whereas bessel_J is defined in special.py and is not a class, but just
a direct function definition. And the reason you can't do

sage: plot(bessel_J(1, t), (t, 1, 10))

is that you can't partially evaluate a bessel_J function

sage: bessel_J(1,t)
Traceback (click to the left for traceback)
...

TypeError: Unable to convert x

(='1-1/8*t^2+1/192*t^4-1/9216*t^6+1/737280*t^8-1/88473600*t^10+1/1486356\
4800*t^12-1/3329438515200*t^14+1/958878292377600*t^16+O(t^17)') to
real
number.

The lambda trick is certainly a sensible thing to do, but I can see
why it's difficult to guess that you should do that if you're just
thinking about things mathematically, rather than pythonically.

We should think hard about making things easy to partially evaluate.
Why not have all the special functions behave like polylog?
Furthermore, it would ideally be easy for a user to define their own
functions that behave like polylog.

Regards,

JM

[John Cremona]

I agree entirely. Users should not have to have to think exactly how
about polylog() and bessel_J() are defined.

John

2008/9/11 Jason Merrill <jwme...@gmail.com>:

[Dan Drake]

On Thu, Sep 11, 2008 at 10:07 PM, John Cremona <john.c...@gmail.com> wrote:
>
> This works:
>
> plot(lambda t:bessel_J(1, t), (1, 10))

That's very helpful, thanks!

[Jason Merrill]

On Sep 11, 10:32 am, Jason Merrill <jwmerr...@gmail.com> wrote:
>
> We should think hard about making things easy to partially evaluate.
> Why not have all the special functions behave like polylog?
> Furthermore, it would ideally be easy for a user to define their own
> functions that behave like polylog.

This is now http://trac.sagemath.org/sage_trac/ticket/4102

JM