Subject: integrate ( exp (x^3), x, 0, 1) still wrong
From: Edwin Woollett
Date: Fri, 11 May 2012 16:02:46 -0700
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.
Consider three cases: setting domain = complex
gives a different symbolic answer, as advertised.
The remaining problem is getting accurate float
values from gamma_incomplete.
Also, there is the non-trivial problem: how is the
user supposed to know that he/she is supposed to
set domain to complex??
(%i1) (display2d:false,fpprintprec:8)$
(%i2) cfloat(zzz) := expand(float(rectform(zzz)))$
(%i3) domain;
(%o3) real
case exp(x^3) over (0,1)
domain:complex cures answer since we can
neglect small fp roundoff errors:
(%i4) quad_qags(exp(x^3),x,0,1);
(%o4) [1.3419044,1.48981318E-14,21,0]
(%i5) integrate(exp(x^3),x,0,1);
(%o5) (gamma_incomplete(1/3,-1)-gamma(1/3))/3
(%i6) cfloat(%);
(%o6) -1.1621233*%i-0.670952
(%i7) domain:complex;
(%o7) complex
(%i8) integrate(exp(x^3),x,0,1);
(%o8) ((-1)^(2/3)*gamma_incomplete(1/3,-1)-(-1)^(2/3)*gamma(1/3))/3
(%i9) cfloat(%);
(%o9) 1.3419044-1.48029737E-16*%i
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
--------------------------------------------------
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.
(%i17) quad_qags(exp(x^5),x,1,2);
(%o17) [1.01323949E+12,0.661551,147,0]
(%i18) domain:real;
(%o18) real
(%i19) integrate(exp(x^5),x,1,2);
(%o19) (gamma_incomplete(1/5,-32)-gamma_incomplete(1/5,-1))/5
(%i20) cfloat(%);
(%o20) -6.21355339E+11*%i-8.55222255E+11
(%i21) domain:complex;
(%o21) complex
(%i22) integrate(exp(x^5),x,1,2);
(%o22) ((-1)^(4/5)*gamma_incomplete(1/5,-32)-(-1)^(4/5)
*gamma_incomplete(1/5,-1))/5
(%i23) cfloat(%);
(%o23) 1.05711284E+12-2.9296875E-4*%i
Ted Woollett