Is there special significance to 16331239353195370.0?

Question

Using import numpy as np I've noticed that

np.tan(np.pi/2)

gives the number in the title and not np.inf

16331239353195370.0

I'm curious about this number. Is it related to some system machine precision parameter? Could I have calculated it from something? (I'm thinking along the lines of something similar to sys.float_info)

EDIT: The same result is indeed reproducible in other environments such as Java, octace, matlab... The suggested dupe does not explain why, though.

http://stackoverflow.com/questions/12713790/tan-in-java-returning-a-strange-value — Padraic Cunningham, Jul 10 '16 at 19:04
I don't like that answer - it's entirely hand-wavy, not really explaining the cause. "Well, tan(pi/2) in radians is essentially infinite, isn't it?" doesn't explain anything about why - as the OP asked here - the answer _isn't_ in fact `np.inf`. But it's straightforward to not only explain why it isn't, but to also explain why the answer is exactly what was seen - and so I did ;-) — Tim Peters, Jul 11 '16 at 02:36

Tim Peters · Accepted Answer · 2016-07-14T18:18:56.990

pi isn't exactly representable as Python float (same as the platform C's double type). The closest representable approximation is used.

Here's the exact approximation in use on my box (probably the same as on your box):

>>> import math
>>> (math.pi / 2).as_integer_ratio()
(884279719003555, 562949953421312)

To find the tangent of that ratio, I'm going to switch to wxMaxima now:

(%i1) fpprec: 32;
(%o1) 32
(%i2) tan(bfloat(884279719003555) / 562949953421312);
(%o2) 1.6331239353195369755967737041529b16

So essentially identical to what you got. The binary approximation to pi/2 used is a little bit less than the mathematical ("infinite precision") value of pi/2. So you get a very large tangent instead of infinity. The computed tan() is appropriate for the actual input!

For exactly the same kinds of reasons, e.g.,

>>> math.sin(math.pi)
1.2246467991473532e-16

doesn't return 0. The approximation math.pi is a little bit less than pi, and the displayed result is correct given that truth.

OTHER WAYS OF SEEING math.pi

There are several ways to see the exact approximation in use:

>>> import math
>>> math.pi.as_integer_ratio()
(884279719003555, 281474976710656)

math.pi is exactly equal to the mathematical ("infinite precision") value of that ratio.

Or as an exact float in hex notation:

>>> math.pi.hex()
'0x1.921fb54442d18p+1'

Or in a way most easily understood by just about everyone:

>>> import decimal
>>> decimal.Decimal(math.pi)
Decimal('3.141592653589793115997963468544185161590576171875')

While it may not be immediately obvious, every finite binary float is exactly representable as a finite decimal float (the reverse is not true; e.g. the decimal 0.1 is not exactly representable as a finite binary float), and the Decimal(some_float) constructor produces the exact equivalent.

Here's the true value of pi followed by the exact decimal value of math.pi, and a caret on the third line points to the first digit where they differ:

true    3.14159265358979323846264338327950288419716939937510...
math.pi 3.141592653589793115997963468544185161590576171875
                         ^

math.pi is the same across "almost all" boxes now, because almost all boxes now use the same binary floating-point format (IEEE 754 double precision). You can use any of the ways above to confirm that on your box, or to find the precise approximation in use if your box is an exception.

@Tim Peters - This is perfectly clear. For completeness I'm guessing that this representation of `np.pi` is the closest rational representation to within system's epsilon? — Aguy, Jul 11 '16 at 03:19
Assuming `np.pi` has the same value as Python's `math.pi` (I didn't check, but you can ;-) ), it's the closest value to the mathematical pi representable in the platform's native `C double` floating-point format. Which means IEEE 754 double-precision on almost all boxes now, and so the closest binary float with 53 bits of (mantissa) precision. So the set of rationals is constrained to the form `+/- I * 2**J` where integer `I` is 0 or `2**52 <= I < 2**53`, and the range of integer `J` is way broad enough to cover all rationals of this form anywhere near `pi`. — Tim Peters, Jul 11 '16 at 03:32
And this is why I'd _love_ if "binary" trigonometric functions were more commonly implemented. Since pi can never be represented in a rational, it would be handy with a set of functions operating on angles from 0 to 1. — pipe, Jul 11 '16 at 09:27
@Έρικ Κωνσταντόπουλος `math.pi`, `np.pi`, and `scipy.pi` are all the same; they're duplicated just for naming convenience; http://stackoverflow.com/questions/12645547/should-i-use-scipy-pi-numpy-pi-or-math-pi — Tim Peters, Jul 11 '16 at 16:10

Is there special significance to 16331239353195370.0?

1 Answers1

OTHER WAYS OF SEEING math.pi