Re: Bessel series

[Bill Allombert <Bill.Allombert@math.u-bordeaux.fr>]

On Tue, Oct 17, 2017 at 12:17:17AM +0000, Jacques Gélinas wrote:
> How can I compute the entire function  H_p(x) = Gamma(p+1) J_p(x) / (x/2)^p 
> where J_p denotes the Bessel-J function, and p >  -1 ?
> For example, H_{-1/2} = cos, H_{1/2} = sinc.
> 
> # besselj(-1/2,0)
>   *** besselj: domain error in gpow(0,n): n <= 0

If you need the value at 0, you need to use the power series expanssion
taking into account the documentation:

  If x converts to a power series,  the initial factor
  (x/2)^nu/Gamma(nu+1)  is omitted  (since it cannot be represented in
  PARI when nu is not integral).

So 
? besselj(-1/2,O(x))
%11 = 1+O(x^2)

Maybe you want this:
H(p,x) = if(x,gamma(p+1)*besselj(p,x)/(x/2)^p,polcoeff(besselj(p,O('x)),0))

Cheers,
Bill.