|
[Date Index]
[Thread Index]
[Author Index]
a bug in Integrate
- To: mathgroup at smc.vnet.net
- Subject: [mg122247] a bug in Integrate
- From: dimitris <dimmechan at yahoo.com>
- Date: Sat, 22 Oct 2011 06:07:15 -0400 (EDT)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
The following has a connection with the recent post of mine called
Simplification.
http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/b6e8fc7f70a0f20f#
Let
int = 1/((x^2 + x + 1) Sqrt[x^2 - x + 1]);
Then Mathematica's result for the definite integral in the range
[-1,0] is
In[74]:= Integrate[int, {x, -1, 0}]
Out[74]= -((1/(4* Sqrt[3]))*(I*(Sqrt[
1 + I*Sqrt[3]]*(-2*I* ArcTan[(6 + 3*I*Sqrt[3])/(3*I - Sqrt[3]
+
3*Sqrt[1 - I*Sqrt[3]])] +
2*ArcTanh[(3*Sqrt[3])/(3*I + Sqrt[3] -
Sqrt[3 - 3*I*Sqrt[3]])] -
Log[11*I + 4*Sqrt[3] + 10*I*Sqrt[1 - I*Sqrt[3]]] +
Log[3*(13*I + 4*Sqrt[3] + 6*I*Sqrt[3 - 3*I*Sqrt[3]])]) +
Sqrt[1 - I*Sqrt[3]]*(2*ArcTanh[(3*(2*I + Sqrt[3]))/(3*I +
Sqrt[3] - 3*Sqrt[1 + I*Sqrt[3]])] +
2*ArcTanh[(3*Sqrt[3])/(3*I - Sqrt[3] + Sqrt[3 +
3*I*Sqrt[3]])] +
Log[(1/16)*(-11*I + 4*Sqrt[3] - 10*I*Sqrt[1 + I*Sqrt[3]])]
-
Log[(3/16)*(-13*I + 4*Sqrt[3] - 6*I*Sqrt[3 +
3*I*Sqrt[3]])]))))
which is incorect. Indeed:
In[75]:= N[%]//Chop
Out[75]= -1.294232744953466
In[76]:= NIntegrate[int, {x, -1, 0}]
Out[76]= 0.9272087241257184
I think the reason is that Mathematica firsts evaluates an
antiderivative:
In[78]:= intMath = Integrate[int, x]
Out[78]= (1/(4*Sqrt[3]))*(-2*Sqrt[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]))] -
2*Sqrt[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]))] +
I*((Sqrt[1 - I*Sqrt[3]] - Sqrt[1 + I*Sqrt[3]])*
Log[16*(1 + x + x^2)^2] +
Sqrt[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]))] -
Sqrt[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]))]))
and then simply does (intMath /. x -> 0) - (intMath /. x -> -1).
In[84]:= N[(intMath /. x -> 0) - (intMath /. x -> -1)]//Chop
Out[84]= -1.2942327449534659
However the evaluated antiderivative possesses a jump singularity at
near x = -0.725
and Mathematica fails unfortunately to take it into acount.
In[85]:= Plot[intMath, {x, -3, 3}]
Nevertheless, there are and some good news. The antiderivative is
correct in the whole complex plane.
In[89]:= FullSimplify[D[intMath, x] - int]
Out[89]= 0
If the integration range does not include the jump, Mathematica will
give correct result. E.g.
In[13]:= Integrate[int, {x, 0, Infinity}] (*after several minutes*)
Out[13]= (1/(4* Sqrt[3]))*(I*(Sqrt[
1 + I*Sqrt[3]]*(2*ArcTanh[(3*Sqrt[3])/(3*I + Sqrt[3] -
Sqrt[3 - 3*I*Sqrt[3]])] + Log[11*I + 4*Sqrt[3] -
8*I*Sqrt[1 - I*Sqrt[3]]] -
Log[11*I + 4*Sqrt[3] + 10*I*Sqrt[1 - I*Sqrt[3]]] -
Log[(-2 - I*Sqrt[3] + 2*Sqrt[1 - I*Sqrt[3]])/(-1 +
Sqrt[1 - I*Sqrt[3]])] +
Log[(-2 + I*Sqrt[3] + 2*Sqrt[1 - I*Sqrt[3]])/(-1 +
Sqrt[1 - I*Sqrt[3]])]) + Sqrt[1 - I*Sqrt[3]]*(2*
ArcTanh[(3*Sqrt[3])/(3*I - Sqrt[3] +
Sqrt[3 + 3*I*Sqrt[3]])] -
Log[-11*I + 4*Sqrt[3] + 8*I*Sqrt[1 + I*Sqrt[3]]] +
Log[-11*I + 4*Sqrt[3] - 10*I*Sqrt[1 + I*Sqrt[3]]] -
Log[(-2 - I*Sqrt[3] + 2*Sqrt[1 + I*Sqrt[3]])/(-1 +
Sqrt[1 + I*Sqrt[3]])] + Log[(-2 + I*Sqrt[3] + 2*Sqrt[1 +
I*Sqrt[3]])/(-1 +
Sqrt[1 + I*Sqrt[3]])])))
In[14]:= N[%]
Out[14]= 0.9358813101035703 + 0.*I
In[15]:= NIntegrate[int, {x, 0, Infinity}]
Out[15]= 0.9358813101035696
At first glance integrand 1/((x^2 + x + 1) Sqrt[x^2 - x + 1]) may look
simple.
This is far from true.
Dimitris
Prev by Date:
Re: count the number of sub-list in a list
Next by Date:
Re: count the number of sub-list in a list
Previous by thread:
Re: Clean up code to run faster
Next by thread:
a bug in Integrate (2nd message)
|