28

fixed in 10.1 (windows)


I'm running Mathematica 10.0.0 and encountered a disturbing error in the symbolic integration of a rather simple function

Integrate[(1 - x)*(1 + 2*x)^6/Sqrt[1 - x^2], {x, -1, 1}]/Pi

The correct value for this integral is 15 (and NIntegrate gives that correctly) but Mathematica evaluates it symbolically as 1/π+29/2. I tried Wolfram Alpha, and it also gives the wrong answer. Any idea what is going on?

incorrect Wolfram Alpha symbolic evaluation of this integral

the correct answer is 15π=47.1239


splitting the integrand into two parts does give the correct answer 15,

Integrate[(1 + 2*x)^6/Sqrt[1 - x^2], {x, -1, 1}]/Pi - 
 Integrate[x*(1 + 2*x)^6/Sqrt[1 - x^2], {x, -1, 1}]/Pi

somehow Mathematica has difficulty with square root singularities in the integrand?

Nasser
  • 143,286
  • 11
  • 154
  • 359
Carlo Beenakker
  • 935
  • 7
  • 15
  • 4
    Are you that Carlo Beenakker ? Anyway this seems better int = Integrate[(1 - x)*(1 + 2*x)^6/Sqrt[1 - x^2], x] and then Limit[int, x -> 1] - Limit[int, x -> -1]. – b.gates.you.know.what Sep 05 '14 at 12:11
  • 2
    thanks, actually splitting the integrand in two parts also gives the correct answer, but I do find it worrisome that I need to take these precautions for what is a simple square-root singularity. – Carlo Beenakker Sep 05 '14 at 12:18
  • Hmm...it gets this equivalent integral right: Integrate[(1 + 2*x)^6*Sqrt[(1 - x)/(1 + x)], {x, -1, 1}] – Michael E2 Sep 05 '14 at 12:34
  • 1
    This works, too: Integrate[(1 - x)*(1 + 2*x)^6/Sqrt[1 - x^2], {x, -1, 1}, GenerateConditions -> False] – Michael E2 Sep 05 '14 at 12:44
  • This also works: displace the endpoints of the integration, evaluate the integration, then move the endpoints back to their original locations: Integrate[1/\[Pi] ((1 - x) (1 + 2 x)^6)/Sqrt[1 - x^2], {x, -1 + a, 1 + a}] /. a -> 0. – Stephen Luttrell Sep 05 '14 at 12:57
  • 2
    Will investigate.. – Daniel Lichtblau Sep 05 '14 at 15:11
  • Also this bug I reported three years ago, never fixed: $Assumptions = w ∈ Reals; FEx[w_] = (I E^(-25 + 5 I w - w^2/4) Sinh[5 w])/Sqrt[2]; NIntegrate[(FEx[w])\[Conjugate]*FEx[w], {w, -∞, ∞}] NIntegrate[(FEx[w])\[Conjugate]* FEx[w], {w, -∞, ∞}, Method -> {Automatic, "SymbolicProcessing" -> 0}] – xslittlegrass Sep 05 '14 at 15:44
  • 1
    @xslittlegrass I am not able to figure out what might be the issue you are reporting. If it was this NIntegrate[FEx[w]*Conjugate[FEx[w]], {w, -Infinity, Infinity}] (or same, but with Integrate) then I do not see an incorrect result. In any case, sending a report via comment, without proper formatting of the input and without full explanation of the issue, is not all that sensible in terms of getting it addressed. – Daniel Lichtblau Sep 05 '14 at 17:51
  • @DanielLichtblau sorry for the confusion. With the same assumption $Assumptions = w ∈ Reals;, on my computer I get different result for these two cases: (1) NIntegrate[(FEx[w])\[Conjugate]*FEx[w], {w, -\[Infinity], \[Infinity]}] (2) NIntegrate[(FEx[w])\[Conjugate]* FEx[w], {w, -\[Infinity], \[Infinity]}, Method -> {Automatic, "SymbolicProcessing" -> 0}]. I reported to the technical support three years ago, and it was confirmed. I thought this is related to the bug here. I will send a proper report with explanation to the support. – xslittlegrass Sep 05 '14 at 21:50
  • @xslittlegrass As formatted (above), these will be unlikely to behave at all: the "conjugate" appears as a symbol multiplying FEx[w]. When I rewrite as Conjugate[FEx[w]] it seems to go better, and both cases agreed. – Daniel Lichtblau Sep 05 '14 at 22:05
  • @DanielLichtblau I have a screenshot here, also for your reference, the ticket number is "TS 33505" when I reported last time. – xslittlegrass Sep 05 '14 at 22:35
  • @xslittlegrass Thanks, the new post is helpful. I'll have a look. Definitely a bug of some sort. – Daniel Lichtblau Sep 05 '14 at 22:44
  • @DanielLichtblau Thanks for confirming that :) – xslittlegrass Sep 05 '14 at 22:48
  • For reference: 5.2 returns the correct answer. It seems tweaks to the internal algorithms inadvertently introduced this bug in later versions up until the fix. – J. M.'s missing motivation May 16 '16 at 17:13

3 Answers3

17

Having experienced similar problematic issues with Mathematica I instantly thought that expanding the fraction in the integrand i.e. applying Appart could resolve the problem, and indeed it does:

Integrate[ Apart[(1 - x)(1 + 2x)^6/Sqrt[1 - x^2]], {x, -1, 1}]/Pi
15

These arguments apply to this case as well Bug in mathematica analytic integration? i.e.

... underlying complex functions behind the symbolic result are not defined in the whole complex plane...

Whenever an indetermined integral involves ArcSin, ArcCos, ArcTan etc. (here we have ArcSin[x]) one should carefully compare numerical and indetermined integrals remembering that Mathematica functions may have some arbitrary branch cuts in the complex plane. For a related issue see also How to calculate contour integrals with Mathematica? where various user defined branch cuts of analytic functions are compared.

Artes
  • 57,212
  • 12
  • 157
  • 245
  • oddly it gets it correct if we split up the integral as Integrate[ .. ,{x,-1,1/2}]+Integrate[..,{x,1/2,1}] (depending on where you make the cut you get either the right or the wrong result.. ) – george2079 Sep 05 '14 at 18:35
4

Main

The bug cannot only be attributed to the Sqrt in the integrand. It is trickier.

In fact, define for t=0,1,2,...

f[t_] := Integrate[(1 - x)*(1 + 2*x)^t/Sqrt[1 - x^2], {x, -1, 1}]/Pi

Then

{#, f[#]} & /@ {0, 1, 2, 3, 4, 5, 6, 7, 8}

(* Out[33]= {{0, 1}, {1, 0}, {2, 1}, {3, 1}, {4, 3}, {5, 6}, {6, (
  1 + (29 \[Pi])/2)/\[Pi]}, {7, (1 + (71 \[Pi])/2)/\[Pi]}, {8, (
  1 + (181 \[Pi])/2)/\[Pi]}} *)

Which is correct for t=0..5 but becomes wrong for t=6, 7, ... but in all cases the branch cut problem of the Sqrt in the integrand should appear. So what is so special about t=6, 7, ...?

I don't have an answer but have continued the study: writing in the integrand (for |x|<1)

Simplify[(1 - x)/Sqrt[1 - x^2 ] == Sqrt[(1 - x)/(1 + x)], -1 <= x <= 1]

(* Out[40]= True *)

the integral becomes

f1[t_] := Integrate[Sqrt[(1 - x)/(1 + x)] (1 + 2*x)^t, {x, -1, 1}]/Pi

At integer values we have

{#, f1[#]} & /@ {0, 1, 2, 3, 4, 5, 6, 7, 8}

(* Out[37]= {{0, 1}, {1, 0}, {2, 1}, {3, 1}, {4, 3}, {5, 6}, {6, 15}, {7, 36}, {8, 91}} *)

and no problem is encountered.

The integral can even be solved for real t (do the indefinite integral, then insert the limits x=1, x=-1), with the analytical result

f2[t_] := -3^t (-2 Hypergeometric2F1[1/2, -t, 1, 4/3] + Hypergeometric2F1[1/2, -t, 2, 4/3])

Except for non negative integers t this function is complex as can be seen by plotting it:

Plot[{Re[#], Im[#]} &[f2[t]], {t, -7, 7}]

(* graph now shown here *)

Best regards, Wolfgang

EDIT 07.09.14 00:55

Even stranger:

Integrate[(1 + 2 x)^6/Sqrt[1 - x^2] (1 - x), {x, -1, 1}]/Pi

(* Out[88]= (1 + (29 \[Pi])/2)/\[Pi] *)

Decomposing trivially (1-x) = 1 - x gives two integer parts

Integrate[(1 + 2 x)^6/Sqrt[1 - x^2] (1), {x, -1, 1}]/Pi

(* Out[91]= 141 *)

Integrate[(1 + 2 x)^6/Sqrt[1 - x^2] (-x), {x, -1, 1}]/Pi

(* Out[92]= -126 *)

and added = 15, as it should be

% + %%

(* Out[93]= 15 *)

Regards, Wolfgang

Dr. Wolfgang Hintze
  • 13,039
  • 17
  • 47
4

fixed in 10.1 (windows):

Mathematica graphics


Mathematica graphics

code:

Clear[x]
Integrate[(1 - x)*(1 + 2*x)^6/Sqrt[1 - x^2], {x, -1, 1}]/Pi
Nasser
  • 143,286
  • 11
  • 154
  • 359