4

Consider this code:

import numpy
numpy.seterr(under='warn')
x1 = 1 + 1j / (1 << 533)
x2 = 1 - 1j / (1 << 533)
y1 = x1 * 1.1
y2 = x2 * 1.1
z1 = x1 / 1.1
z2 = x2 / 1.1
print(numpy.divide(1, x1))  #              1-3.55641399918e-161j  # OK
print(numpy.divide(1, x2))  #              1+3.55641399918e-161j  # OK
print(numpy.divide(1, y1))  # 0.909090909091-3.23310363561e-161j  # underflow
print(numpy.divide(1, y2))  # 0.909090909091+3.23310363561e-161j  # underflow
print(numpy.divide(1, z1))  #            1.1-3.91205539909e-161j  # underflow
print(numpy.divide(1, z2))  #            1.1+3.91205539909e-161j  # underflow

The underflow doesn't seem to make sense no matter how I look at it. Like Wikipedia says,

Underflow is a condition in a computer program where the result of a calculation is a number of smaller absolute value than the computer can actually store in memory on its CPU.

But clearly the computer is capable of storing numbers in the general proximity of the values in question, so the definition doesn't seem at all consistent with the behavior I'm seeing here.

Could someone explain why exactly some of these are giving underflows but others are not?
Is this correct behavior or is it a bug?

user541686
  • 205,094
  • 128
  • 528
  • 886
  • An intermediate calculation in the C code that implements complex division computes `y1.imag*(y1.imag/y1.real)`. With `y1.imag = 3.912055399093737e-161`, that calculation results in a [denormal number](https://en.wikipedia.org/wiki/Denormal_number), which should trigger the underflow FPE. This is a comment, not an answer, because I don't know why the other calculations don't trigger the underflow. – Warren Weckesser Aug 16 '17 at 04:27
  • I think this is the same question as yours: Let `a = x1.imag`. Why does `np.multiply(a, 1.6*a)` trigger the underflow warning, but `np.multiply(a, 1.5*a)` and `np.multiply(a, 2.0*a)` do not? – Warren Weckesser Aug 16 '17 at 04:36
  • @WarrenWeckesser: Thanks for the comment! Yeah that might be equivalent too. – user541686 Aug 16 '17 at 04:54
  • 2
    This might be a numpy bug. Let's find out: https://github.com/numpy/numpy/issues/9568 – Warren Weckesser Aug 16 '17 at 05:09
  • [When does underflow occur?](https://stackoverflow.com/q/42277132/2410359) may help. Underflow is sometimes not just a range issue, but a precision one too. – chux - Reinstate Monica Aug 17 '17 at 15:25

1 Answers1

0

To do this math, you need to first "scale", so-to-speak, the 3.23310363561e-161. That is what triggers the underflow

  0.909090909091
- 3.23310363561e-161 =
--------------------

  0.909090909091 x 10^0
- 3.23310363561 x 10^-161 =
-------------------------

  0.909090909091 x 10^0
- 0.000...                 0323310363561 =
    ^^^^^^^^^^^^^^^^^^^^^^^^
    160 0s
  --------------------------------------
Andrew
  • 1
  • 4
  • 19