Issue 1:
Describe the issue:

This may not be a bug, but the behavior is mysterious and unclear to me. Part of the problem is that in trying to recreate an MRE I cannot.

The setup is I have two arrays I'm multiplying, one with very small numbers. I use pdb to drop into interactive mode at the exception. Poking around a bit, the error is happening at the multiplication step, and I can recreate it by taking the smallest non-zero value and multiplying it by the second array, e.g.

Pdb) x[x>0].min() * W
*** FloatingPointError: underflow encountered in multiply

## or

(Pdb) x[x>0].min() * 1e-10
*** FloatingPointError: underflow encountered in double_scalars

Ok, the behavior I would like is for this to just go to zero — I don't care about the loss of precision. But then I try to create an MRE to experiment with solutions. What is the minimum?

(Pdb) y = x[x>0]
(Pdb) y.min()
9.117510882823816e-308

Then I try in a new Python session,

>>> import numpy as np
>>> np.array([9.117510882823816e-308], dtype='float64')
array([9.11751088e-308])
>>> f = np.array([9.117510882823816e-308], dtype='float64')
>>> f * 1e-10
array([9.11751e-318])
# No error

I then think that this is because of the loss of precision in printing the float. I wasn't sure how to handle this, but I tried np.set_printoptions(precision=1000) but this seems to only impact how arrays are printed, not the np.float64 values as they're pulled out. So I next thought I can work with the raw data.

(Pdb) y.min().tobytes()
b'\xe4\x94\x10\xf0\xf6c0\x00'

Then in the second session,

$ python
Python 3.8.12 | packaged by conda-forge | (default, Jan 30 2022, 23:42:07)
[GCC 9.4.0] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import numpy as np
>>> np.frombuffer(b'\xe4\x94\x10\xf0\xf6c0\x00', dtype='float64')
array([9.11751088e-308])
>>> f = np.frombuffer(b'\xe4\x94\x10\xf0\xf6c0\x00', dtype='float64')
>>> f * 1e-10
array([9.11751e-318])
# No error again!

I find it odd that I can't create an MRE even with the buffer. What could be going on here?

I also tried using format strings,

(Pdb) f"{y:.10000000}"
'9.117510882823815665548977178863615628591324749590460596300089284407848522066006531932055383210337983249882939129197635693204683436837619342899485403616406753151023440802800570432302810469374755878262353843400205063415657530940618334968743282787718401379326462405248209351474798409859277445681925875157863523372638365924454046245513719836560774569949082153016343041264284375069892098945977461004851844859972461215974549935010425299501133344015649450602626874461701983236523291202172361897246010705343303184104772030271936525192961955598148050240545429201711268321687328347974353590798632145138753708909956222270239479037724670024066080647006736702136140242161352259516590652676898809691225616065391385488011421236917205938345143412249171888106502592563629150390625e-308'

Then, in the interactive python session I pasted this in

>>> ff = np.array([9.117510882823815665548977178863615628591324749590460596300089284407848522066006531932055383210337983249882939129197635693204683436837619342899485403616406753151023440802800570432302810469374755878262353843400205063415657530940618334968743282787718401379326462405248209351474798409859277445681925875157863523372638365924454046245513719836560774569949082153016343041264284375069892098945977461004851844859972461215974549935010425299501133344015649450602626874461701983236523291202172361897246010705343303184104772030271936525192961955598148050240545429201711268321687328347974353590798632145138753708909956222270239479037724670024066080647006736702136140242161352259516590652676898809691225616065391385488011421236917205938345143412249171888106502592563629150390625e-308], dtype='float64')
>>> ff
array([9.11751088e-308])
>>> ff * 1e-10
array([9.11751e-318])
# still no error!

So I'm at a loss as to how to create an MRE of this underflow (and how to get the behavior I see and want in the interactive session, that is, it just goes to zero). Apologies if I'm missing something.

Reproduce the code example:

# create the MRE is part of the bug, see above

Error message:

No response

Runtime information:

Pdb session:

(Pdb) import sys, numpy; print(numpy.__version__); print(sys.version)
1.24.1
3.8.12 | packaged by conda-forge | (default, Jan 30 2022, 23:42:07)
[GCC 9.4.0]

Second session:

>>> import sys, numpy; print(numpy.__version__); print(sys.version)
1.24.1
3.8.12 | packaged by conda-forge | (default, Jan 30 2022, 23:42:07)
[GCC 9.4.0]

Context for the issue:

This is important because I cannot create a working MRE due to the bug. The behavior is different even though I used the raw bytes to try to recreate the example. It would be helpful too if there was documentation on how to create a proper MRE given precision/underflow issues.

Issue 2:
While it is threadsafe, the with statement gives the false impression that it's always safe. It is not:

import asyncio
import numpy as np

def divide_by_zero():
np.float64(1) / 0

async def foo():
for j in range(3):
with np.errstate(divide='raise'):
for i in range(2):
try:
divide_by_zero()
except FloatingPointError:
pass
else:
print("foo: Failed to raise!")
await asyncio.sleep(0.15)

async def bar():
for j in range(3):
with np.errstate(divide='ignore'):
for i in range(2):
try:
divide_by_zero()
except FloatingPointError:
print("bar: raised anyway!")
await asyncio.sleep(0.11)

loop = asyncio.get_event_loop()
loop.run_until_complete(asyncio.gather(foo(), bar()))
loop.close()

gives:

foo: Failed to raise!
bar: raised anyway!
foo: Failed to raise!

because the error-state is thread-local, but asyncio tasks share a thread.

I don't know if there is any way this can be fixed - it seems like this might be a python bug, and context managers might need an __yield_enter__ and __yield_leave__ mechanism, to be notified when control is leaving the suite via a yield/await