Group Abstract Group Abstract

Message Boards Message Boards

1
|
10.1K Views
|
24 Replies
|
4 Total Likes
View groups...
Share
Share this post:

Problem with Integrate

Posted 11 years ago

Hello.

I was trying to do the following integral:

Integrate[1/Sqrt[2 (1 + x)^3 + 1], {x, 0, 1}] // N

The output in Mathematica 9 is a complex number that is not correct

-2.96267 - 1.92762 I

When I try to solve by using NIntegrate, the output is the correct one

NIntegrate[1/Sqrt[2 (1 + x)^3 + 1], {x, 0, 1}]

0.376065

Is this a bug in Integrate or I am missing something? Can someone test in Mathematica 10?

Thanks

POSTED BY: otsocol
24 Replies
POSTED BY: Daniel Lichtblau

A side observation...

If one computes

theIntegral = Integrate[1/Sqrt[2 (1 + x)^3 + 1], x]

and then computes

theDerivative = D[theIntegral, x] // FullSimplify

which yields

-(((-1)^(5/6) Sqrt[3] Sqrt[
     1 + 2 (1 + x)^3])/(Sqrt[-(-1)^(5/6) (1 + 2^(1/3) + 2^(1/3) x)]
      Sqrt[(-1)^(5/6) (-1 + (-2)^(1/3) (1 + x))] Sqrt[
     1 + (-2)^(1/3) (1 + x) + (-2)^(2/3) (1 + x)^2] Sqrt[
     3 - (-1)^(1/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 1]] Sqrt[
     3 + (-1)^(5/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 2]])) 

and then attempts to Plot this using

Plot[-(((-1)^(5/6) Sqrt[3] Sqrt[
      1 + 2 (1 + x)^3])/(Sqrt[-(-1)^(5/6) (1 + 2^(1/3) + 2^(1/3) x)]
       Sqrt[(-1)^(5/6) (-1 + (-2)^(1/3) (1 + x))] Sqrt[
      1 + (-2)^(1/3) (1 + x) + (-2)^(2/3) (1 + x)^2] Sqrt[
      3 - (-1)^(1/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 1]] Sqrt[
      3 + (-1)^(5/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 2]])) , {x, 0, 1}]

the result is that it crashes the Mathematica Kernel. (I will report this.)

In[1]:= SystemInformation["Small"]

Out[1]= {"Kernel" -> {"SystemID" -> "MacOSX-x86-64", 
   "ReleaseID" -> "10.0.0.0 (5098698, 5098537)", 
   "CreationDate" -> 
    DateObject[{2014, 6, 29}, TimeObject[{20, 38, 32}]]}, 
 "FrontEnd" -> {"OperatingSystem" -> "MacOSX", 
   "ReleaseID" -> "10.0.0.0 (5098698, 2014062704)", 
   "CreationDate" -> 
    DateObject[{2014, 6, 27}, TimeObject[{23, 17, 40.}]]}}
POSTED BY: David Reiss

Let me add two remarks

1) the exact value of the jump position seems to be xj = 2^(2/3) - 1

N[2^(2/3) - 1, 20]
(* 0.58740105196819947475 *)

2) an antiderivative without a jump can be calculated as follows

ad2F1[x_] = Integrate[1/Sqrt[1 + 2 x^3], {x, 0, t}, Assumptions -> t > 0] /. t -> 1 + x

(* Out[116]= (1 + x) Hypergeometric2F1[1/3, 1/2, 4/3, -2 (1 + x)^3] *)

Hence we find a hypergeometric function instead of the "bizarre" elliptic integral.

Proof

D[ad2F1[x], x] // Simplify

(* Out[109]= 1/Sqrt[1 + 2 (1 + x)^3] *)

Did somebody notice already that the simple substitution x+1->y, dx->dy, {x,0,1} -> {y,1,2} leads to the correct result?

In fact:

sol2 = Integrate[1/Sqrt[2 y^3 + 1], {y, 1, 2}]

(* Out[27]= Sqrt[2] Hypergeometric2F1[1/6, 1/2, 7/6, -(1/2)] -  Hypergeometric2F1[1/6, 1/2, 7/6, -(1/16)] *)

sol2 // N

(* Out[28]= 0.376065 *)

solN = NIntegrate[1/Sqrt[2 (1 + x)^3 + 1], {x, 0, 1}]

(* Out[39]= 0.376065 *)


$Version

(* Out[40]= "8.0 for Microsoft Windows (64-bit) (October 7, 2011)" *)

Best regards, Wolfgang

Yes, for the indefinite integral you're right; it is still in Version 10.0.0 the case. The difference of valued indefinite integrals does not convert into the Hypergeometric2F1[] form and the bug reappears

In[11]:= N[(ReplaceAll[#, y -> 2] - ReplaceAll[#, y -> 1]) &[Integrate[1/Sqrt[2 y^3 + 1], y]]]
Out[11]= -2.96267 - 1.92762 I

In[15]:= N[FullSimplify[(ReplaceAll[#, y -> 2] - ReplaceAll[#, y -> 1])] &[Integrate[1/Sqrt[2 y^3 + 1], y]]]
Out[15]= -2.96267 - 1.92762 I

In[17]:= (Integrate[1/Sqrt[2 y^3 + 1], y] /. y -> x + 1) === Integrate[1/Sqrt[2 (1 + x)^3 + 1], x]
Out[17]= True

that means this is a bug not about a path singularity, because the indefinite integral does not go along a path. It's something about turning back to the right leaf. For example, look at this, the derivative does not turn back to the original integrand, but presents an equation of six-th degree instead:

otsocol6Degree

POSTED BY: Udo Krause
Posted 11 years ago
POSTED BY: otsocol
POSTED BY: Udo Krause
sln = DSolve[{y'[x] == 1/Sqrt[2 (1 + x)^3 + 1], y[0] == 0}, y]

sln represents a syntax error as such. If one is concerned about newbies and rookies not sticking to the syntax, one should stick to it at least oneself.

POSTED BY: Udo Krause

Much simpler example that causes the crash:

Plot[Root[108 + #1^6 &, 2], {x, 0, 1}]
POSTED BY: David Reiss

Interesting. I've reported the crash and it's in the bug system. And I just added your example to the report.

POSTED BY: David Reiss

Slipping in Evaluate prevents the crash. fyi.

Plot[ Evaluate[-(((-1)^(5/6) Sqrt[3] Sqrt[
        1 + 2 (1 + x)^3])/(Sqrt[-(-1)^(5/6) (1 + 2^(1/3) + 
           2^(1/3) x)] Sqrt[(-1)^(5/
            6) (-1 + (-2)^(1/3) (1 + x))] Sqrt[
        1 + (-2)^(1/3) (1 + x) + (-2)^(2/3) (1 + x)^2] Sqrt[
        3 - (-1)^(1/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 1]] Sqrt[
        3 + (-1)^(5/6) Sqrt[3] + (1 + x) Root[108 + #1^6 &, 
           2]]))], {x, 0, 1}]

Ditto pre-processing the expression with FullSimplify.

expr2 = FullSimplify[  expr , x >= 0] 

I don't want to discourage your reporting the crash.

POSTED BY: Bruce Miller
POSTED BY: David Reiss
POSTED BY: W. Craig Carter
Posted 11 years ago
POSTED BY: otsocol

I filed a bug report.

POSTED BY: David Reiss

Yes, that's what happens on my computer also. Weird.

POSTED BY: Frank Kampas

Plotting theDerivative using

Plot[theDerivative,{x,0,1}] 

does not crash my Kernel, but plotting the full expression for theDerivative as written in my post does crash it reproducibly.

POSTED BY: David Reiss
POSTED BY: Frank Kampas
POSTED BY: Daniel Lichtblau
Posted 11 years ago

Hi Daniel. Thanks for your answer.

When you say that the antiderivative has a jump discontinuity near x=0.59, I think you are looking at the Mathematica answer

In[6]:= Integrate[1/Sqrt[2 (1 + x)^3 + 1], x]

Out[6]= ((-1)^(1/6) 2^(2/3) Sqrt[(-1)^(5/6) (-1 + (-2)^(1/3) (1 + x))] Sqrt[1 + (-2)^(1/3) (1 + x) + (-2)^(2/3) (1 + x)^2] EllipticF[ ArcSin[Sqrt[-(-1)^(5/6) - I (-2)^(1/3) (1 + x)]/3^(1/4)], (-1)^( 1/3)])/(3^(1/4) Sqrt[1 + 2 (1 + x)^3])

Is that correct? But is it possible to trust in this result in this case?

POSTED BY: otsocol

I see what you mean

sln = DSolve[{y'[x] == 1/Sqrt[2 (1 + x)^3 + 1], y[0] == 0}, y]

Plot [Im[y[x] /. sln], {x, 0, 1}, PlotRange -> All ] shows that y[x] develops an imaginary part at 0.59

POSTED BY: Frank Kampas
POSTED BY: Frank Kampas

There is a path singularity around .59 that is not found. Will investigate.

POSTED BY: Daniel Lichtblau

Using Mathematica 10

In[2]:= Integrate[1/Sqrt[2 (1 + x)^3 + 1], {x, 0, 1}]

Out[2]= (1/(3^(1/4)))(-1)^(1/12) 2^( 2/3) (EllipticF[ ArcSin[((-1)^(11/12) Sqrt[1 + 2^(1/3)])/3^(1/4)], (-1)^(1/3)] - EllipticF[ ArcSin[((-1)^(11/12) Sqrt[1 + 2 2^(1/3)])/3^(1/4)], (-1)^(1/3)])

In[3]:= N[%]

Out[3]= -2.96267 - 1.92762 I

Seeing -1 raised to fractional powers makes a branch cut problem a possibility.

POSTED BY: Frank Kampas
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard