Re: Re: Can Integrate[expr,{x,a,b}] give an incorrect result?
- To: mathgroup at smc.vnet.net
- Subject: [mg81883] Re: [mg81855] Re: [mg81827] Can Integrate[expr,{x,a,b}] give an incorrect result?
- From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
- Date: Sat, 6 Oct 2007 04:29:33 -0400 (EDT)
- References: <8207246.1191505119218.JavaMail.root@m35> <200710050846.EAA00304@smc.vnet.net>
> So, from that, D and/or Simplify must be wrong.
Must they? Even without any investigation I would say that it is much
more likely that this integral should not be evaluated by means of
the Leibniz rule (substituting the limits into the anti-derivative
and subtracting). In fact, it is easy to see what happens in your
special case:
integrand = (r*(r - rho*Cos[alpha - t]))/
(3*(r^2 - 2*rho*Cos[alpha - t]*r + rho^2 +
(z - zeta)^2)^3);
int = integrand /. {r -> 1, rho -> 11/10,
zeta -> 0, alpha -> 0, z -> 12/10}
(1 - (11*Cos[t])/10)/(3*(73/20 - (11*Cos[t])/5)^3)
First, let us compute the integral using Integrate:
Integrate[int, {t, 0, 2*Pi}]
(1772800*Pi)/(11512449*Sqrt[377])
N[%]
0.0249156
Now, let's do it by numerical integration:
NIntegrate[int, {t, 0, 2*Pi}]
0.024915592618816908
so far, so good. Now let's compute the indefinte integral:
indef = Integrate[int, t];
Using the Liebniz rule gives clearly the wrong answer:
indef = Integrate[int, t];
but why should it give the right answer? The function indef is
clearly discontinuous in the interval 0 to 2Pi
Plot[indef, {t, 0, 2 Pi}]
so the Leibniz rule does not apply. There is no reason at all to
suspect the very reliable function Simplify and almost as reliabel
Interate.
Andrzej Kozlowski
Andrzej Kozlowski
On 5 Oct 2007, at 17:46, DrMajorBob wrote:
> Yes, it's an incorrect result.
>
> To paraphrase your post:
>
> integrand = (r*(r -
> rho*Cos[alpha - t]))/(3*(r^2 + rho^2 + (z - zeta)^2 -
> 2*r*rho*Cos[alpha - t])^3);
> assumptions = {r > 0, rho > 0, Element[zeta, Reals],
> Element[z, Reals], alpha >In[15]:= 0};
>
> Integrate[integrand, {t, 0, 2 Pi}, Assumptions -> assumptions]
>
> 0
>
> Plot[integrand /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 0,
> z -> 1.2}, {t, 0, 2 Pi}, PlotRange -> All]
>
> "Visual integration" is unconvincing, but you're right, the
> integral isn't
> zero:
>
> NIntegrate[
> integrand /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 0,
> z -> 1.2}, {t, 0, 2 Pi}]
>
> 0.0249156
>
> Your example violates "assumptions" (alpha -> 0 vs. alpha > 0), but
> that's
> not the problem:
>
> NIntegrate[
> integrand /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 1,
> z -> 1.2}, {t, 0, 2 Pi}]
>
> 0.0249156
>
> So, as often happens, you've found a case where Mathematica got a
> general
> solution that isn't correct for some values of the parameters. The
> quadratic formula is another example:
>
> Solve[a x^2 + b x + c == 0, x]
>
> {{x -> (-b - Sqrt[b^2 - 4 a c])/(
> 2 a)}, {x -> (-b + Sqrt[b^2 - 4 a c])/(2 a)}}
>
> Not true when a == 0, of course.
>
> But let's look further:
>
> indefinite = Integrate[integrand, t, Assumptions -> assumptions];
> definite = Subtract @@ (indefinite /. {{t -> 2 Pi}, {t -> 0}});
> Simplify[definite, assumptions]
>
> 0
>
> Is that a Simplify error, or an Integrate error?
>
> definite /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 1, z -> 1.2}
>
> -2.05998*10^-18 + 0. \[ImaginaryI]
>
> Ah. Simplify may be doing its job just fine (hard to be certain),
> but it
> appears the indefinite integral is incorrect. But look:
>
> D[indefinite, t] == integrand // Simplify
>
> True
>
> and
>
> D[indefinite, t] == integrand // Simplify[#, assumptions] &
>
> True
>
> and also
>
> definite // Together // Numerator // Simplify
>
> 0
>
> So, from that, D and/or Simplify must be wrong.
>
> If this is correct:
>
> numdefinite // Together // Numerator // TrigFactor
>
> -4 r^2 (r^4 + r^2 rho^2 - 2 rho^4 + 2 r^2 z^2 - rho^2 z^2 + z^4 -
> 4 r^2 z zeta + 2 rho^2 z zeta - 4 z^3 zeta + 2 r^2 zeta^2 -
> rho^2 zeta^2 + 6 z^2 zeta^2 - 4 z zeta^3 +
> zeta^4) (ArcTanh[((r^2 + 2 r rho + rho^2 + (z - zeta)^2) Tan
> [alpha/
> 2])/Sqrt[-r^4 +
> 2 r^2 (rho^2 - (z - zeta)^2) - (rho^2 + (z - zeta)^2)^2]] -
> ArcTanh[((r^2 + 2 r rho + rho^2 + (z - zeta)^2) Tan[
> 1/2 (alpha - 2 \[Pi])])/
> Sqrt[-r^4 +
> 2 r^2 (rho^2 - (z - zeta)^2) - (rho^2 + (z - zeta)^2)^2]])
> (r^2 +
> rho^2 + z^2 - 2 z zeta + zeta^2 - 2 r rho Cos[alpha])^2
>
> The next-to-last factor (the only one that depends on t) is
>
> num[[-2]]
>
> ArcTanh[((r^2 + 2 r rho + rho^2 + (z - zeta)^2) Tan[alpha/2])/
> Sqrt[-r^4 +
> 2 r^2 (rho^2 - (z - zeta)^2) - (rho^2 + (z - zeta)^2)^2]] -
> ArcTanh[((r^2 + 2 r rho + rho^2 + (z - zeta)^2) Tan[
> 1/2 (alpha - 2 \[Pi])])/
> Sqrt[-r^4 + 2 r^2 (rho^2 - (z - zeta)^2) - (rho^2 + (z - zeta)^2)
> ^2]
> ]
>
> so it all comes down (seemingly) to
>
> Tan[alpha/2] == Tan[1/2 (alpha - 2 \[Pi])] // Simplify
>
> True
>
> or
>
> Tan[any] == Tan[any - Pi]
>
> True
>
> BUT... what's likely happening is that Simplify doesn't properly
> take into
> account the branch-cut discontinuities of ArcTanh. In that case,
> there's a
> significant set of parameters for which the integral IS zero, such as:
>
> Reduce[Flatten@{num == 0, assumptions}]
>
> ors = (zeta \[Element]
> Reals && ((z < zeta &&
> rho > Sqrt[z^2 - 2 z zeta + zeta^2]/Sqrt[2] && alpha > 0 &&
> r == Root[-2 rho^4 - rho^2 z^2 + z^4 + 2 rho^2 z zeta -
> 4 z^3 zeta - rho^2 zeta^2 + 6 z^2 zeta^2 - 4 z zeta^3 +
> zeta^4 + (rho^2 + 2 z^2 - 4 z zeta +
> 2 zeta^2) #1^2 + #1^4 &, 2]) || (z == zeta && rho >
> 0 &&
> alpha > 0 &&
> r == Root[-2 rho^4 - rho^2 z^2 + z^4 + 2 rho^2 z zeta -
> 4 z^3 zeta - rho^2 zeta^2 + 6 z^2 zeta^2 - 4 z zeta^3 +
> zeta^4 + (rho^2 + 2 z^2 - 4 z zeta +
> 2 zeta^2) #1^2 + #1^4 &, 2]) || (z > zeta &&
> rho > Sqrt[z^2 - 2 z zeta + zeta^2]/Sqrt[2] && alpha > 0 &&
> r == Root[-2 rho^4 - rho^2 z^2 + z^4 + 2 rho^2 z zeta -
> 4 z^3 zeta - rho^2 zeta^2 + 6 z^2 zeta^2 - 4 z zeta^3 +
> zeta^4 + (rho^2 + 2 z^2 - 4 z zeta +
> 2 zeta^2) #1^2 + #1^4 &, 2]))) || ((z |
> zeta) \[Element] Reals && alpha > 0 && r > 0 &&
> rho > 0) || (C[1] \[Element] Integers && z \[Element] Reals &&
> C[1] >= 1 && alpha == 2 \[Pi] C[1] && r > 0 && rho == r &&
> zeta == z)
>
> But look:
>
> % /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 1, z -> 1.2}
>
> True
>
> Hold on, though... there are a lot of Ors in there. The least
> restrictive
> seems to be
>
> ors[[1, 2, 1]]
>
> z < zeta && rho > Sqrt[z^2 - 2 z zeta + zeta^2]/Sqrt[2] && alpha >
> 0 &&
> r == Root[-2 rho^4 - rho^2 z^2 + z^4 + 2 rho^2 z zeta - 4 z^3
> zeta -
> rho^2 zeta^2 + 6 z^2 zeta^2 - 4 z zeta^3 +
> zeta^4 + (rho^2 + 2 z^2 - 4 z zeta + 2 zeta^2) #1^2 + #1^4 &, 2]
>
> But that's a pretty strict condition (on r, especially), and
>
> % /. {r -> 1, rho -> 1.1, zeta -> 0, alpha -> 1, z -> 1.2}
>
> False
>
> SO. Integrate found a set of conditions for which the definite
> integral
>
> was zero; it just didn't tell you what they were! This is supposed
> to be
> the fix for that, I think:
>
> Integrate[integrand, {t, 0, 2 Pi}, Assumptions -> assumptions,
> GenerateConditions -> True]
>
> 0
>
> BUT, as you see, that didn't help one bit!
>
> Bobby
>
> On Thu, 04 Oct 2007 03:27:21 -0500, W. Craig Carter <ccarter at mit.edu>
> wrote:
>
>>
>> I believe I am getting an incorrect result from a definite
>> integration:
>>
>> InputForm[integrand] is
>> (R*(R - rho*Cos[alpha - t]))/(3*(R^2 + rho^2 + (z - zeta)^2 -
>> 2*R*rho*Cos[alpha - t])^3)
>>
>> InputForm[assumptions] is
>> {R > 0, L > 0, rho > 0, Element[zeta, Reals], Element[z, Reals],
>> alpha>
>> 0}
>>
>> Integrate[integrand,{t,0,2Pi},Assumptions->assumptions]
>> returns 0
>>
>> But compare this to:
>> (visually integrate...)
>> Plot[integrand/.{R -> 1, rho -> 1.1, zeta -> 0, alpha -> 0, z ->
>> 1.2},{t,0,2 Pi},PlotRange->All]
>>
>> (numerically integrate...)
>> Plot[NIntegrate[integrand/.{{R -> 1, rho -> 1.1, zeta -> 0, alpha -
>> > 0,
>> z -> 1.2},{t,0,tau}],{tau,0,2Pi}, PlotRange->All]
>>
>>
>> Something isn't adding up??
>>
>> Thanks, WCC
>>
>>
>
>
>
> --
>
> DrMajorBob at bigfoot.com
>
- References:
- Re: Can Integrate[expr,{x,a,b}] give an incorrect result?
- From: DrMajorBob <drmajorbob@bigfoot.com>
- Re: Can Integrate[expr,{x,a,b}] give an incorrect result?