Subject: integrate ( exp (x^3), x, 0, 1) still wrong
From: Raymond Toy
Date: Sat, 12 May 2012 07:32:56 -0700
On 5/11/12 4:02 PM, Edwin Woollett wrote:
> On May 11, 2012, I wrote
> -------------
>> This is an integral which gave a history stack
>> overflow lisp error in 5.26.0 (Jan 17, 2012).
>>
>> Using the lastest gcl 5.27.0-1, there is no longer
>> a stack problem, but a wildly wrong answer.
> ----------------------------
> The bug fix which cured the stack overflow error
> included the warning that one should use domain:complex
> to get the correct (symbolic) answer.
That seems unfortunate that you have to know, magically, to specify a
complex domain for an integral that is obviously real. If you hadn't
remembered the bug fix, would you have known? I didn't.
[snip]
> case exp(x^3) over (1,2)
> domain:complex cures answer provided we
> neglect the imaginary piece of order(E-14)??
>
> (%i10) quad_qags(exp(x^3),x,1,2);
> (%o10) [275.51098,3.23056159E-7,21,0]
>
> (%i11) domain:real;
> (%o11) real
>
> (%i12) integrate(exp(x^3),x,1,2);
> (%o12) (gamma_incomplete(1/3,-8)-gamma_incomplete(1/3,-1))/3
>
> (%i13) cfloat(%);
> (%o13) -238.59951*%i-137.75549
>
> (%i14) domain:complex;
> (%o14) complex
>
> (%i15) integrate(exp(x^3),x,1,2);
> (%o15) ((-1)^(2/3)*gamma_incomplete(1/3,-8)-(-1)^(2/3)
> *gamma_incomplete(1/3,-1))/3
>
> (%i16) cfloat(%);
> (%o16) 275.51098-1.89478063E-14*%i
I'm not sure what you're expecting. To get a real result,
gamma_incomplete obviously has to be complex to cancel out the imaginary
part from (-1)^(2/3).
>
> --------------------------------------------------
> case exp(x^5) over (1,2)
> domain:complex greatly improves answer
> but large incorrect imaginary part
> and real part only accurate to two sig. figures.
Yes, that's a problem. gamma_incomplete(1/5,-32) is not very accurate.
This can probably solved, but I'm not sure. Using bigfloats with 32
digits gives a much better answer.
Ray