Re: integrate
- To: mathgroup at smc.vnet.net
- Subject: [mg72453] Re: integrate
- From: "David W. Cantrell" <DWCantrell at sigmaxi.net>
- Date: Thu, 28 Dec 2006 05:36:00 -0500 (EST)
- References: <emgg3i$2ca$1@smc.vnet.net>
"dimitris" <dimmechan at yahoo.com> wrote:
For your third question:
> Integrate[Sqrt[x^2 + 1]/Sqrt[1 + x^6], {x, 0, 10}]
> (-(-1)^(1/6))*EllipticF[I*ArcSinh[10*(-1)^(1/3)], (-1)^(2/3)]
>
> N[%]
> -0.10016641038463325 - 1.6857503548125956*I
>
> NIntegrate[Sqrt[x^2 + 1]/Sqrt[1 + x^6], {x, 0, 10}]
> 2.056349237110889
>
> --------> Any ideas to "help" Mathematica to give a correct answer?
I have an answer, but getting there was not pleasant. I hope others
will have better suggestions.
Part 1:
I wanted to know why the symbolic answer was wrong, so I looked at.
In[7]:= Integrate[Sqrt[(1 + x^2)/(1 + x^6)], x]
Out[7]=
(-(-1)^(1/6))*Sqrt[1 - (-1)^(1/3)*x^2]*Sqrt[1 + (-1)^(2/3)*x^2]*
Sqrt[1/(1 - x^2 + x^4)]*EllipticF[I*ArcSinh[(-1)^(1/3)*x], (-1)^(2/3)]
The product of the algebraic factors involving x should be just 1, but
I was not able to get Mathematica to give that:
In[8]:=
FullSimplify[Sqrt[1 - (-1)^(1/3)*x^2]*Sqrt[1 + (-1)^(2/3)*x^2]*
Sqrt[1/(1 - x^2 + x^4)], x > 0]
Out[8]=
Sqrt[1 - (-1)^(1/3)*x^2]*Sqrt[(1 + (-1)^(2/3)*x^2)/(1 - x^2 + x^4)]
Maybe I didn't try hard enough. Can someone suggest how I should get
the above to simplify to 1, perhaps assuming x > 0? Anyway, I then
simplified Out[7] myself to
In[9]:= -(-1)^(1/6)*EllipticF[I*ArcSinh[(-1)^(1/3)*x], (-1)^(2/3)]
and then plotted the real part of that from x = 0 to 10. A jump
discontinuity is apparent, at roughly x = 2.45. If we ignore that
discontinuity and try to evaluate the integral from x = 0 to 10 by the
Fundamental Theorem, we get the same wrong answer that Mathematica
got. So now I think I know what Mathematica did wrong. But I'm a bit
surprised. Even when Mathematica returns antiderivatives having such
discontinuities, it normally takes such discontinuities correctly into
account when doing _definite_ integrals. But it didn't do so here, so
there is indeed a bug to be exterminated.
Part 2:
I thought I would help Mathematica by taking the discontinuity, caused
by a branch cut, into account by hand. But unfortunately it wasn't
obvious exactly where the cut was. So I went to the Wolfram Functions
site and looked under EllipticF. It was not helpful. The only relevant
statement was "Branch cut locations: complicated." I then gave up on
determining the exact location. But if I could have done that, it
would have been nicer than what follows. (Anyone know exactly where
the discontinuity is?)
Part 3:
I sought to transform the original integral into another one which
Mathematica might handle better. My first thought was to use the
substitution u = x^2, which led to
In[16]:= 1/2*Integrate[Sqrt[(1 + u)/(u*(1 + u^3))], u]
Out[16]=
(-1)^(1/6)*Sqrt[1 - (-1)^(1/3)/u]*Sqrt[1 + (-1)^(2/3)/u]*u^(3/2)*
Sqrt[1/(u - u^2 + u^3)]*EllipticF[I*ArcSinh[(-1)^(1/3)/Sqrt[u]],(-1)^(2/3)]
Simplifying by hand and putting the result back in terms of x gave
In[17]:= (-1)^(1/6)*EllipticF[I*ArcSinh[(-1)^(1/3)/x], (-1)^(2/3)].
Although different from In[9], this result still has a discontinuity.
But thankfully, it's in a different place, roughly x = 0.4. That
allows us to combine In[9] and In[17], using the former to integrate
from x0 to 1 and the latter to integrate from 1 to x1. Thus, assuming
x0 < 1 < x1, we can state that
Integrate[Sqrt[(1 + x^2)/(1 + x^6)], {x, x0, x1}]
is
In[21]:=
int = (-1)^(1/6)*(-2*EllipticF[I*ArcSinh[(-1)^(1/3)], (-1)^(2/3)] +
EllipticF[I*ArcSinh[(-1)^(1/3)*x0], (-1)^(2/3)] +
EllipticF[I*ArcSinh[(-1)^(1/3)/x1], (-1)^(2/3)]);
In[22]:= int /. {x0 -> 0, x1 -> 10}
Out[22]=
(-1)^(1/6)*(EllipticF[I*ArcSinh[(1/10)*(-1)^(1/3)], (-1)^(2/3)] -
2*EllipticF[I*ArcSinh[(-1)^(1/3)], (-1)^(2/3)])
which is the precise symbolic result desired.
In[23]:= N[%]
Out[23]= 2.0563492371150103 + 3.3306690738754696*^-16*I
David W. Cantrell
- Follow-Ups:
- Re: Re: integrate
- From: Andrzej Kozlowski <akoz@mimuw.edu.pl>
- Re: Re: integrate