#2191 error in integrating exp(-x)*sinh(sqrt(x))

closed
5
2011-07-22
2011-04-23
Anonymous
No

Maxima 5.24.0 http://maxima.sourceforge.net
using Lisp SBCL 1.0.24
(%i3) integrate(exp(-x)*sinh(sqrt(x)),x,0,inf);
(%o3) 0
(%i4) quad_qagi(exp(-x)*sinh(sqrt(x)),x,0,inf);
(%o4) [1.137937897234377, 5.171862937913829e-11, 345, 0]

Discussion

  • Aleksas

    Aleksas - 2011-05-18

    solving with Maxima 5.24.0:
    Dwight Tables of integrals 862.04
    (%i1) declare(integrate,linear)$ assume(y>0)$
    (%i3) S:'integrate(exp(-x)*sinh(sqrt(x)),x,0,inf)$
    (%i4) exponentialize(%);
    (%o4) integrate((%e^sqrt(x)-%e^(-sqrt(x)))*%e^(-x),x,0,inf)/2
    (%i5) expand(%);
    (%o5) integrate(%e^(sqrt(x)-x),x,0,inf)/2-integrate(%e^(-x-sqrt(x)),x,0,inf)/2
    (%i6) S1:changevar(%, y=sqrt(x), y, x);
    (%o6) integrate(y*%e^(y-y^2),y,0,inf)-integrate(y*%e^(-y^2-y),y,0,inf)
    (%i10) changevar(part(S1,1),2*y-1=z,z,y)+
    changevar(part(S1,2),2*y+1=z,z,y);
    (%o10) integrate((z+1)*%e^(-(z^2-1)/4),z,-1,inf)/4-integrate((z-1)*%e^(-(z^2-1)/4),z,1,inf)/4
    (%i8) ev(%, nouns),ratsimp;
    (%o8) (%e^(1/4)*sqrt(%pi))/2
    (%i9) float(%), numer;
    (%o9) 1.137937897234373

    Where is the error?
    wrong:
    (%i1) R:integrate(y*%e^(y-y^2),y,0,inf);
    (%o1) (%e^(1/4)*gamma_incomplete(1,1/4))/2-(%e^(1/4)*gamma_incomplete(1/2,1/4))/4
    (%i2) float(R), numer;
    (%o2) 0.22717931961748
    (%i3) quad_qagi(y*%e^(y-y^2), y, 0, inf);
    (%o3) [1.365117216851849,7.6857545173630158*10^-10,165,0]

    correct:
    (%i4) R1:integrate(y*%e^(-y-y^2),y,0,inf);
    (%o4) (%e^(1/4)*gamma_incomplete(1,1/4))/2-(%e^(1/4)*gamma_incomplete(1/2,1/4))/4
    (%i5) float(R1), numer;
    (%o5) 0.22717931961748
    (%i6) quad_qagi(y*%e^(-y-y^2), y, 0, inf);
    (%o6) [0.22717931961748,1.4118104693665376*10^-9,135,0]

    Aleksas D

     
  • Dan Gildea

    Dan Gildea - 2011-05-24

    A simpler test case:
    (%i16) integrate(exp(sqrt(x)-x), x, 0, inf);
    (%o16) %e^(1/4)*gamma_incomplete(1,1/4)-%e^(1/4)*gamma_incomplete(1/2,1/4)/2
    +2*%e^(1/4)*sqrt(%pi)
    (%i17) float(%);
    (%o17) 5.006110228172447
    (%i18) quad_qagi(exp(sqrt(x)-x), x, 0, inf);
    (%o18) [2.730234433703704,9.207745677031198e-11,345,0]

    The problem seems to be caused by the fact that the indefinite integral:
    (%i15) integrate(exp(sqrt(x)-x), x);
    (%o15) %e^(1/4)*%i
    *(%i*gamma_incomplete(1,(1-2*sqrt(x))^2/4)*(1-2*sqrt(x))^2
    /(2*sqrt(x)-1)^2
    -%i*gamma_incomplete(1/2,(1-2*sqrt(x))^2/4)*(1-2*sqrt(x))
    /(2*abs(2*sqrt(x)-1)))
    has a discontinuity at x=1/4. This integral is computed
    by case m2-exp-type-5 in integrate-exp-special.

    Mathematica gives the indefinite integral as:
    -exp(sqrt(x)-x)+exp(1/4)*sqrt(%pi)*erf((-1+2*sqrt(x))/2)/2;
    which results in the correct value for the definite integral in this case.

     
  • Dan Gildea

    Dan Gildea - 2011-07-22
    • assigned_to: nobody --> dgildea
    • status: open --> closed
     
  • Dan Gildea

    Dan Gildea - 2011-07-22

    Fixed in defint.lisp.

    Look for discontinuities in antiderivative.

     

Log in to post a comment.

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:

JavaScript is required for this form.





No, thanks