On Tue, Dec 20, 2011 at 2:26 PM, Barton Willis <willisb at unk.edu> wrote:
> Numerically evaluating this expression results in subtractive
> cancellation; try:
>
> (%i9) e : integrate(integrate(cos(x*y*log(2)),y,1,10),x,1,10)$
>
> (%i10) makelist(block([fpprec : 15], bfloat(block([fpprec : 10 * k],
> expand(bfloat(e))))), k,1,10);
>
> (%o10)
> [-2.238107472896b14,6.40685300762984b-16*%i+6.58088693446743b6,-9.59754643393143b-1,
>
> -9.5975928150589b-1,-9.59759281505811b-1,-9.59759281505811b-1,-9.59759281505811b-1,
> -9.59759281505811b-1,-9.59759281505811b-1,-9.59759281505811b-1]
>
> The nfloat function (defined in hypergeometric) knows nothing about
> gamma_incomplete :(
>
>
I think this is a pretty good example where, if the goal is a numerical
answer, it's better to use a numerical algorithm instead of taking the
analytic solution and numerically evaluating that. :-)
Curiously, Wolfram Alpha says gamma_incomplete(0, 100*%i*log(2)) is
infinity. I didn't check to see if that is true or not.
Ray