Re: Re: Re: Integrate vs NIntegrate
- To: mathgroup at smc.vnet.net
- Subject: [mg89211] Re: [mg89182] Re: [mg89147] Re: [mg89129] Integrate vs NIntegrate
- From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
- Date: Thu, 29 May 2008 07:05:33 -0400 (EDT)
- References: <9489155.1211802216995.JavaMail.root@m08> <op.ubr1yjh52c6ksp@bobbys-imac> <2b8d8f40805261936u2ec81454r1b6ba16be0897662@mail.gmail.com> <200805271114.HAA01380@smc.vnet.net> <29744585.1211975587467.JavaMail.root@m08> <op.ubvkuekr2c6ksp@bobbys-imac>
Perhaps you should pay a little more attention to other posts in the
same thread (not just mine, however flatterning I find this). In case
you have not noticed, I strongly recommend David Cantrell's posts,
particularly on such kind of topics and, of course, Daniel Lichtblau.
Although I really do not see why you need my help in finding them,
here are the relevant quotes.
Here is David:
> Res1 will not be the desired value. Reason: f is not an
> antiderivative of h
> on the interval [0, 3] due to the fact that f has a discontinuity at
> x = 2.
>
> ...........
>
> We can, nonetheless, use f to get a reasonably accurate value for the
> integral:
>
> In[20]:= N[f[b] - f[2 + 10^-7] + f[2 - 10^-7] - f[a]]
> Out[20]= 1.65267 + 0. I
and here is Daniel:
> The antiderivative, f[x], has a jump discontinuity due to a branch cut
> that is either at or very near to x=2 (probably at). Have a look at
>
> Plot[{Re[f[x]], Im[f[x]]}, {x, 0, 3}]
>
> Thus evaluation at the endpoints will fail to recover the definite
> integral. Integrate[h[x], {x, a, b}] will give a more plausible
> result, in
> accord with the NIntegrate value.
By the way, I would have been very happy to have staked some real
money on this but I feel it would have been, kind of, unfair.
Andrzej Kozlowski
On 29 May 2008, at 02:43, DrMajorBob wrote:
>> Not a proof, of course, but I would bet a fair amount on it
>> [Integrate] being right.
>
> I believe someone recently said we should ignore predictions from
> those who haven't staked actual money on it.
>
> But anyway, is Integrate correct in BOTH these calculations?
>
> q = -1/2; a = 0; b = 3;
> h[x_] = (1 + x^3)^q;
> f[x_] = Integrate[h[x], x];
> f[b] - f[a] // N
>
> -2.55387 - 2.42865 I
>
> f[x_] = Integrate[h[t], {t, 0, x}];
> f[b] - f[a] // N
>
> 1.65267
>
> Along with the OP, I'm betting the first result from Integrate is
> wrong, in view of
>
> NIntegrate[h[x], {x, a, b}]
>
> 1.65267
>
> If NIntegrate is wrong instead, that hardly improves matters.
>
> Bobby
>
> On Wed, 28 May 2008 03:46:27 -0500, Andrzej Kozlowski <akoz at mimuw.edu.pl
> > wrote:
>
>> I would say that Integrate is correct. Evidence:
>>
>> In[1]:= q = -2^(-1);
>>
>> In[2]:= f = Integrate[(1 + x^3)^q, x];
>>
>> In[3]:= g = Chop[N[D[Normal[f + O[x]^12], x]]]
>>
>> Out[3]= -0.3125000000000013*x^9 + 0.37500000000000105*x^6 -
>> 0.5000000000001258*x^3 +
>> 0.9999999999999998
>>
>> and, on the other hand:
>>
>> In[4]:= N[Normal[(1 + x^3)^q + O[x]^10]]
>> Out[4]= -0.3125*x^9 + 0.375*x^6 - 0.5*x^3 + 1.
>>
>>
>> Not a proof, of course, but I would bet a fair amount on it being
>> right.
>>
>> Andrzej Kozlowski
>>
>>
>> On 27 May 2008, at 20:14, DrMajorBob wrote:
>>
>>> Because they're completely different. One is an antiderivative in
>>> Complex
>>> variables; the other is a definite integral from 0 to x (so you know
>>> that
>>> it is zero at zero, among other things). Antiderivatives tend to
>>> have
>>> branch cuts, poles, etc., and Mathematica's chosen branch is complex
>>> in
>>> part of your region.
>>>
>>> Or... it's entirely possible Mathematica's answer for
>>> Integrate[h[x], x]
>>> is simply WRONG.
>>>
>>> Let's check whether it has the correct derivative:
>>>
>>> q = -1/2; a = 0; b = 3;
>>> h[x_] = (1 + x^3)^q;
>>> f[x_] = Integrate[h[x], x];
>>> FullSimplify[f'[x] == h[x], a < x < b]
>>>
>>> Sqrt[-(1 + x) (-1 + (-1)^(1/3) x)] == (
>>> 2 (-1)^(1/12) Sqrt[-(-1)^(1/6) ((-1)^(2/3) + x) (1 + x^3)])/(
>>> Sqrt[1 - I/Sqrt[3] + (-1 - I/Sqrt[3]) x] Sqrt[
>>> 3 + I Sqrt[3] + I (3 I + Sqrt[3]) x])
>>>
>>> As you see, FullSimplify doesn't think Integrate got a correct
>>> result. Or
>>> at least, it couldn't verify it.
>>>
>>> Here's an attempt to provide FullSimplify some adult supervision:
>>>
>>> simple= FullSimplify[ #, a<x<b]&;
>>> test= f'[x]== h[x]//simple
>>>
>>> Sqrt[(-(1 + x))*(-1 + (-1)^(1/3)*x)] ==
>>> (2*(-1)^(1/12)*Sqrt[(-(-1)^(1/6))*((-1)^(2/3) + x)*(1 + x^3)])/
>>> (Sqrt[1 -
>>> I/Sqrt[3] + (-1 - I/Sqrt[3])*x]*Sqrt[3 + I*Sqrt[3] + I*(3*I +
>>> Sqrt[3])*x])
>>>
>>> test2 = simple[test /. {(-1)^(1/3) -> -1, (-1)^(2/3) ->
>>> ComplexExpand[Exp[I*(Pi/3)]], (-1)^(1/6) -> ComplexExpand[Exp[I*(Pi/
>>> 6)]],
>>> (-1)^(1/12) -> ComplexExpand[Exp[I*(Pi/12)]]}]
>>>
>>> 1 + x == (((1/2 + I/2)*(3 - I*Sqrt[3])*Sqrt[(-(2*I + (I +
>>> Sqrt[3])*x))*(1
>>> + x^3)])/Sqrt[3 - I*Sqrt[3] + (-3 - I*Sqrt[3])*x])*Sqrt[3 +
>>> I*Sqrt[3]
>>> + I*(3*I + Sqrt[3])*x]
>>>
>>> test3 = test2 /. (z_)^(r_) :> ComplexExpand[z]^r
>>>
>>> 1 + x == (((1/2 + I/2)*(3 - I*Sqrt[3])*Sqrt[(-Sqrt[3])*x -
>>> Sqrt[3]*x^4
>>> + I*(-2 - x - 2*x^3 - x^4)])/Sqrt[3 - 3*x + I*(-Sqrt[3] -
>>> Sqrt[3]*x)])*Sqrt[3 - 3*x + I*(Sqrt[3] + Sqrt[3]*x)]
>>>
>>> simple[(#1^2 & ) /@ test3]
>>>
>>> 2*I + (3*I + Sqrt[3])*x == 0
>>>
>>> I don't think Integrate has the right answer.
>>>
>>> Bobby
>>>
>>> On Mon, 26 May 2008 21:36:31 -0500, Armen Kocharyan
>>> <armen.kocharyan at gmail.com> wrote:
>>>
>>>> Dear Bob,
>>>>
>>>> Sorry for the asterisks, it wasn't intentional. I just cut and
>>>> paste
>>>> from my
>>>> notebook.
>>>>
>>>> I understand. But why would Mathematica give different results for
>>>> the
>>>> following?
>>>>
>>>> Integrate[h[t], {t, 0, x}] and Integrate[h[x], x]
>>>>
>>>>
>>>> Regards,
>>>> Armen
>>>>
>>>> 2008/5/27 DrMajorBob <drmajorbob at att.net>:
>>>>
>>>>> The fundamental theorem of calculus frequently does not apply
>>>>> when, as
>>>>> in
>>>>> your example, Integrate is free to return a Complex-valued
>>>>> function.
>>>>>
>>>>> (A better question is why you put in all those extraneous stars to
>>>>> inhibit
>>>>> copy/pasting the code.)
>>>>>
>>>>> Despite that, here's a slightly different approach where the FTC
>>>>> does
>>>>> apply:
>>>>>
>>>>> q = -1/2; a = 0; b = 3;
>>>>> h[x_] = (1 + x^3)^q;
>>>>> f[x_] = Integrate[h[t], {t, 0, x}]
>>>>>
>>>>> If[x >= -1, (x^3)^(1/3) Hypergeometric2F1[1/3, 1/2, 4/3, -x^3],
>>>>> Integrate[1/Sqrt[1 + t^3], {t, 0, x}, Assumptions -> x < -1]]
>>>>>
>>>>> Res1 = N[f[b] - f[a]]
>>>>> Res2 = NIntegrate[h[x], {x, a, b}]
>>>>>
>>>>> 1.65267
>>>>>
>>>>> 1.65267
>>>>>
>>>>> Contrast that definition of f with yours:
>>>>>
>>>>> Integrate[h[x], x]
>>>>>
>>>>> (1/(3^(1/4) Sqrt[1 + x^3]))2 (-1)^(
>>>>> 1/6) Sqrt[-(-1)^(1/6) ((-1)^(2/3) + x)] Sqrt[
>>>>> 1 + (-1)^(1/3) x + (-1)^(2/3) x^2]
>>>>> EllipticF[ArcSin[Sqrt[-(-1)^(5/6) (1 + x)]/3^(1/4)], (-1)^(1/3)]
>>>>>
>>>>> If you plot both versions from a to b (use PlotRange->All to be
>>>>> sure),
>>>>> you'll find that your version is complex on part of the range
>>>>> (between
>>>>> 2 and
>>>>> 3), but mine is real.
>>>>>
>>>>> Bobby
>>>>>
>>>>>
>>>>> On Mon, 26 May 2008 05:23:35 -0500, Armen Kocharyan <
>>>>> armen.kocharyan at gmail.com> wrote:
>>>>>
>>>>> *q=-1/2; a=0; b=3;*
>>>>>> *h[x_]=(1+x^3)^q;*
>>>>>> *f[x_]=Integrate[h[x],x];*
>>>>>> *Res1=N[f[b]-f[a]]*
>>>>>> *Res2=NIntegrate[h[x],{x,a,b}]*
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> DrMajorBob at longhorns.com
>>>>>
>>>
>>>
>>>
>>> --
>>>
>>> DrMajorBob at longhorns.com
>>>
>>
>>
>>
>
>
>
> --
> DrMajorBob at longhorns.com
- References:
- Re: Integrate vs NIntegrate
- From: DrMajorBob <drmajorbob@att.net>
- Re: Integrate vs NIntegrate