       Re: Can Integrate[expr,{x,a,b}] give an incorrect result?

• To: mathgroup at smc.vnet.net
• Subject: [mg81855] Re: [mg81827] Can Integrate[expr,{x,a,b}] give an incorrect result?
• From: DrMajorBob <drmajorbob at bigfoot.com>
• Date: Fri, 5 Oct 2007 04:46:36 -0400 (EDT)
• References: <8207246.1191505119218.JavaMail.root@m35>

```Yes, it's an incorrect result.

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 > 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

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 && 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 && 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 \[Element] Integers && z \[Element] Reals &&
C >= 1 && alpha == 2 \[Pi] C && 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 && 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]
>
>
>
> Thanks, WCC
>
>

--

DrMajorBob at bigfoot.com

```

• Prev by Date: Re: Re: Problems with PlotStyle in ListPointPlot3D
• Next by Date: Re: Manipulate+Plot showing no output
• Previous by thread: Result Duplicated! Re: Can Integrate[expr,{x,a,b}] give
• Next by thread: Re: Re: Can Integrate[expr, {x, a, b}] give an