Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

Re: a bug in Integrate (2nd message)

• To: mathgroup at smc.vnet.net
• Subject: [mg122276] Re: a bug in Integrate (2nd message)
• From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
• Date: Sun, 23 Oct 2011 06:23:39 -0400 (EDT)
• Delivered-to: l-mathgroup@mail-archive0.wolfram.com
• References: <201110221007.GAA29655@smc.vnet.net>

I don' think the jump point is quite at the point where you suggest it 
is (and I don't think that your claim that Mathematica "fails to 
evaluate" these limits is correct. You can actually see that the jump 
can't be quite at the point you suggest by evaluating:

N[intMath /. x -> jump - 10^(-2), 10] // Chop

0.3647401105

Chop[N[intMath /. x -> jump - 10^(-3), 10]]

 -1.84923250024426519561423619339918941442`10.014331869778037

where jump is the value you computed. As you can see, the actual jump 
takes place slightly to the left of your "jump".

So where does the jump take place? Well, the function ArcTan has branch 
cuts on the imaginary axis running from I to I Infinity and -I to - I 
Infinity. Hence the jump will occur when the argument of ArcTan becomes 
purely Imaginary. One can see this for one of your ArcTan parts by 
looking at the graphs of the real and imaginary parts:

Plot[{Re[(3 (x^2 - x + 1) (x^2 - x +
        I Sqrt[3]))/(-(2 Sqrt[3]
          x^4) + (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] + Sqrt[3] -
         3 I) x^2 + (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] -
         3 Sqrt[3] + 3 I) x -
      Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] +
      2 (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] + 3 I) x^3 + Sqrt[
      3] + 3 I)],
  Im[(3 (x^2 - x + 1) (x^2 - x +
        I Sqrt[3]))/(-(2 Sqrt[3]
          x^4) + (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] + Sqrt[3] -
         3 I) x^2 + (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] -
         3 Sqrt[3] + 3 I) x -
      Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] +
      2 (Sqrt[3 - 3 I Sqrt[3]] Sqrt[x^2 - x + 1] + 3 I) x^3 + Sqrt[
      3] + 3 I)]}, {x, -1, 0}]

One can see from the graph that this point is around x = -0.7248, 
somewhat to the left of the value -0.7224 you computed. In principle one 
could try to compute it by solving the equation f==0, where f is the 
first of the plotted functions, but I don't think any Mathematica solver 
will return this value in a reasonable time.

Andrzej Kozlowski






On 22 Oct 2011, at 12:07, dimitris wrote:

> int = 1/((x^2 + x + 1)*Sqrt[x^2 - x + 1]);
>
> As I told in the previous message Mathematica gives wrong result for
> the following integral
>
> Integrate[int, {x, -1, 0}]
>
> the problem is (I guess from what I figured out in the previous
> message) that after applying the Risch algorithm and gets an
> antiderivative
>
> intMath = Integrate[int, x]
>
> Mathematica uses the fundamental theorem of Calculus incorrectly
> failing to take into account the jump singularity that possesses the
> antiderivative in the integration range.
>
> Let's figure out where is this jump.
>
> In[43]:= ff = Expand[intMath]
>
> Out[43]:=(-(1/2))*Sqrt[(1/3)*(1 + I*Sqrt[3])]*
>  ArcTan[(3*(1 - x + x^2)*(I*Sqrt[3] - x + x^2))/(3*I + Sqrt[3] -
>      2*Sqrt[3]*x^4 -
>            Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>      2*x^3*(3*I + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>      x*(3*I - 3*Sqrt[3] + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2])
> +
>      x^2*(-3*I + Sqrt[3] +
>         Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]))] - (1/2)*
>  Sqrt[(1/3)*(1 - I*Sqrt[3])]*
>     ArcTan[(3*(I*Sqrt[3] + x - x^2)*(1 - x + x^2))/(3*I - Sqrt[3] +
>      2*Sqrt[3]*x^4 + Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>            x^3*(6*I - 2*Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>      x*(3*I + 3*Sqrt[3] - Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2])
> -
>      x^2*(3*I + Sqrt[3] +
>         Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]))] + (1/4)*I*
>  Sqrt[(1/3)*(1 - I*Sqrt[3])]*
>     Log[16*(1 + x + x^2)^2] - (1/4)*I*Sqrt[(1/3)*(1 + I*Sqrt[3])]*
>  Log[16*(1 + x + x^2)^2] +
>   (1/4)*I*Sqrt[(1/3)*(1 + I*Sqrt[3])]*
>  Log[(1 + x + x^2)*(11*I + 4*Sqrt[3] + (11*I + 4*Sqrt[3])*x^2 +
>            10*I*Sqrt[1 - I*Sqrt[3]]*Sqrt[1 - x + x^2] -
>      x*(17*I + 4*Sqrt[3] +
>         8*I*Sqrt[1 - I*Sqrt[3]]*Sqrt[1 - x + x^2]))] -
>   (1/4)*I*Sqrt[(1/3)*(1 - I*Sqrt[3])]*
>  Log[(1 + x + x^2)*(-11*I + 4*Sqrt[3] + (-11*I + 4*Sqrt[3])*x^2 -
>            10*I*Sqrt[1 + I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>      x*(17*I - 4*Sqrt[3] +
>         8*I*Sqrt[1 + I*Sqrt[3]]*Sqrt[1 - x + x^2]))]
>
> In[44]:= Length[ff]
>
> Out[44]= 6
>
> In[45]:= ({#1, Plot[Evaluate[Re[ff[[#1]]]], {x, -3, 3}]} & ) /@
> Range[Length[ff]]
>
> The first two terms of ff (which contain ArcTan) have jump in the same
> point.
>
> In[52]:= ff[[1]]
> ff[[2]]
>
> Out[52]= (-(1/2))*Sqrt[(1/3)*(1 + I*Sqrt[3])]*
> ArcTan[(3*(1 - x + x^2)*(I*Sqrt[3] - x + x^2))/
>       (3*I + Sqrt[3] - 2*Sqrt[3]*x^4 -
>     Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>          2*x^3*(3*I + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>     x*(3*I - 3*Sqrt[3] + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>
>     x^2*(-3*I + Sqrt[3] + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]))]
>
> Out[53]= (-(1/2))*Sqrt[(1/3)*(1 - I*Sqrt[3])]*
> ArcTan[(3*(I*Sqrt[3] + x - x^2)*(1 - x + x^2))/
>       (3*I - Sqrt[3] + 2*Sqrt[3]*x^4 +
>     Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>          x^3*(6*I - 2*Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>     x*(3*I + 3*Sqrt[3] - Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) -
>
>     x^2*(3*I + Sqrt[3] + Sqrt[3 + 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]))]
>
> In particular the jump is at the point:
>
> In[100]:= (N[#1, 20] & )[Re[Select[
>    x /. Solve[3*I + Sqrt[3] - 2*Sqrt[3]*x^4 -
>        Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2] +
>        2*x^3*(3*I + Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>        x*(3*I - 3*Sqrt[3] +
>           Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) +
>        x^2*(-3*I + Sqrt[3] +
>           Sqrt[3 - 3*I*Sqrt[3]]*Sqrt[1 - x + x^2]) == 0], -1 <
>       Re[#1] < 0 && Abs[Abs[#1] - 0.7] < 0.1 & ][[
>        1]]]]
>
> Out[100]= -0.72241927849008646432307690472106948257`20.
>
> or in symbolic form
>
> jump = 1/6 + Re[(I*(1 + I*Sqrt[3])*(1 -
>         5*I*Sqrt[3]))/(12*(-71*I + 9*Sqrt[3] +
>          Sqrt[-5022 - 918*I*Sqrt[3]])^(1/3)) -
>        (1/12)*I*(1 - I*Sqrt[3])*(-71*I + 9*Sqrt[3] +
>        Sqrt[-5022 - 918*I*Sqrt[3]])^(1/3)];
>
> Unfortunately, Mathematica fails to evaluate the side limits
>
> Limit[intMath,x->jump,Direction->1]
> and
> Limit[intMath,x->jump,Direction->-1]
>
> The package NLimit fails as well.
>
> So I think the only solution is to help a bit Mathematica in order to
> get an antiderivative continuous in the whole real axis.
> The question is how?
>
> Dimitris Anagnostou
>
> P.S. In another forum a user gave the following antiderivative:
>
> intUser=ArcTan[(Sqrt[2]*(1 + x))/Sqrt[1 - x + x^2]]/Sqrt[2] +
> ArcTanh[(Sqrt[2/3]*(-1 + x))/Sqrt[1 - x + x^2]]/Sqrt[6];
>
> which is correct in the whole complex plane
>
> D[intUser, x] - 1/((x^2 + x + 1) Sqrt[x^2 - x + 1]) // FullSimplify
> 0
>
> and possesses no jump in the whole real axis.
>
> Hence the definite integral in the range [-1,0] is
>
> In[111]:= FullSimplify[int /. x -> 0 - int /. x -> -1]
>
> Out[111]= (3/13)*Sqrt[3*(4 + Sqrt[3])]
>
> In[112]:= N[%] (*check*)
>
> Out[112]= 0.956959
>
>
>

• References:
  • a bug in Integrate (2nd message)
    • From: dimitris <dimmechan@yahoo.com>