Wrong Limit[...] for inexact expressions

Question

I noticed that Limit can return nonsense when using inexact parameters. In the following code, a is my (exact/inexact) parameter, and I want to know the limit for b -> 1 (or b -> 1., this does not matter).

expr := Limit[(-b*Cos[a*b] Sin[a/2] + Sin[a*b] Cos[a/2])/(b^2 - 1), b -> 1];
a = 2. Pi;
expr    (* Returns -Infinity *)
a = 2 Pi;
expr    (* Returns -Pi *)

I observe this behavior in Mathematica 8, 9, and 10.

I then turned to the documentation of Limit, which states:

Limit may return an incorrect answer for an inexact input:
Limit[Log[1 - (Log[Exp[z]/z - 1] + Log[z])/z]/z, z -> 100.] 
(* -Infinity *)
The result is correct when an exact input is used:
Limit[Log[1 - (Log[Exp[z]/z - 1] + Log[z])/z]/z, z -> 100]
(* 1/100 Log[1 - 1/100 Log[-100 + E^100]] *)

I can see why Limit fails here: The subexpression 1 - 1/100*Log[-100 + E^100] evaluates to approximately 4*10^-44, which is prohibitively small (in single precision, but not in double - this confuses me a bit).

In my example, however, I can't see tiny or huge numbers that could cause this sort of problem. Can you help me see my error?

Edit: I found this problem when I was calculating the Fourier transform of an RF-pulse to find the maximum power density of its spectrum.

Trace[expr] should help to clarify why the approximate case gives what it gives. — Daniel Lichtblau, Jan 27 '16 at 18:05
Have a look at a = 2. [Pi] and a = 2 Pi without the semicolon 2 Pi ist not the same as 2*Pi, see also Trace as stated by @Daniel Lichtblau — , Jan 27 '16 at 18:46
@Louis: I am thoroughly confused. What do the square brackets in a = 2. [Pi] mean? I see that the output is different, but I am unfamiliar with that syntax (and could not find it here or here). Further, I don't see different output when not using the semicolon, or why 2 Pi is not equal to 2*Pi. — Martin J.H., Jan 27 '16 at 18:59
I think @Louis is saying that "2Pi" and "2.Pi" are different, but you already knew that. Try setting a to 2*Pi+e to see an odd result: apparently, there's a sign change in the limit near 2*Pi. — , Jan 27 '16 at 19:04

score 4 · Accepted Answer · answered Jan 27 '16 at 19:29

Mathematica knows how to simplify when a is exactly 2*Pi:


(-b*Cos[a*b] Sin[a/2] + Sin[a*b] Cos[a/2])/(b^2 - 1) /. a -> 2*Pi // InputForm 
-(Sin[2*b*Pi]/(-1 + b^2))

It then applies the numerical limit for b=1.

For the approximate number 2.*Pi, Mathematica can't make this simplification, and it turns out the limit is +Infinity for a<2*Pi and -Infinity for a>2*Pi, so even a small variation from 2*Pi completely changes the problem:

 
expr := Limit[(-b*Cos[a*b] Sin[a/2] + Sin[a*b] Cos[a/2])/(b^2 - 1), b -> 1]; 
a = 2*Pi-10^-10 
expr (* yields Infinity *) 
a = 2*Pi+10^-10 
expr (* yields -Infinity *) 
a = 2*Pi 
expr (* yields -Pi as expected *)

And if we let e represent an arbitrary distance from 2*Pi:


a = 2*Pi + e 
expr // InputForm 
(* yields DirectedInfinity[(I*Sign[-1 + E^(I*e)])/E^((I/2)*Re[e])] *)

So why does 2.*Pi return -Infinity? Because it's actually a little different from 2*Pi:

 
a = 2.*Pi; 
SetPrecision[a, Infinity] (* result: 884279719003555/140737488355328 *) 
expr (* -Infinity *)

Barry, maybe I am misunderstanding, but the expression has a removable discontinuity at $b=1$. The limits approaching from either side of $1$ are the same, and equal to $\pi$, not infinity. You can convince yourself of this by plotting this expression in the neighbourhood of $b=1$ after substituting the $a=2\pi$ value in it (plot). — MarcoB, Jan 27 '16 at 20:37
I meant setting a equal to 2*Pi-e and 2*Pi+e and then taking the limit. If Mathematica is accurate above, a small change in a makes the limit infinite. — , Jan 27 '16 at 20:47

Wrong Limit[...] for inexact expressions

1 Answers1