On 20/10/2011 7:31 AM, Edwin Woollett wrote:
> I am using integrate to evaluate one the the
> Mathematica list of NIntegrate examples:
>
> integrate(bessel_y(2,x),x,1,10^4);
>
> which produces
>
> -2*bessel_y(2,10000)+(4*bessel_y(2,1)-%pi*struve_h(0,1)*bessel_y(1,1)
>
> -%pi*struve_h(-1,1)*bessel_y(0,1))
>
> /2+5000*%pi*struve_h(0,10000)*bessel_y(1,10000)
> +5000*%pi*struve_h(-1,10000)*bessel_y(0,10000)
>
> The floats of the bessel functions are no problem, but
> the struve_h float chokes:
>
> (%i66) float(struve_h(0,10));
> (%o66) 0.11874368368746
> (%i67) float(struve_h(0,100));
> (%o67) -0.070878751689647
> (%i68) float(struve_h(0,1000));
> Exceeded maximum allowed fpprec.Exceeded maximum allowed
> fpprec.Exceeded maximum allowed fpprec.Exceeded maximum allowed
> fpprec.Exceeded maximum allowed fpprec.Exceeded maximum allowed
> fpprec.Exceeded maximum allowed fpprec.Exceeded maximum allowed fpprec.
> (%o68) 636.6197723675813*hypergeometric([1.0],[1.5,1.5],-250000.0)
>
> How can I get a numerical value for the output
> of integrate for this case?
>
> Ted Woollett
Sorry for the glib answer earlier. When I look at the code I see that
you need to increase the variable max_fpprec, defined in
share/hypergeometric/hypergeometric.lisp. It defaults to 1000
In this case 10000 is sufficient:
(%i1) load('hypergeometric);
(%o1)
C:/msys/1.0/local/src/maxima-20111002/share/hypergeometric/hypergeometric.lisp
(%i2) max_fpprec;
(%o2) 1000
(%i3) struve_h(0,100.0);
(%o3) - 0.070878751689647
(%i4) struve_h(0,1000.0);
Exceeded maximum allowed fpprec.Exceeded maximum allowed fpprec.Exceeded
maximum allowed fpprec.Exceeded maximum allowed fpprec.Exceeded maximum
allowed fpprec.Exceeded maximum allowed fpprec.Exceeded maximum allowed
fpprec.Exceeded maximum allowed fpprec.(%o4) 636.6197723675813
hypergeometric([1.0], [1.5, 1.5], - 250000.0)
(%i5) max_fpprec:10000;
(%o5) 10000
(%i6) struve_h(0,1000.0);
(%o6) 0.0053525371133764