ComputationalCryoEM / ASPIRE-Python

Algorithms for Single Particle Reconstruction
http://spr.math.princeton.edu
GNU General Public License v3.0
46 stars 21 forks source link

Possible flaky `Volume` unit test. #1012

Closed garrettwrong closed 11 months ago

garrettwrong commented 1 year ago

This test seems to fail reliably in a clean test environment on caf at v0.12.0 (examples below). I think its recent code, and I haven't seen it fail in our CI yet. Tried a few quick things (downgrades), but no magic bullet yet... so for now, just a report. This should be confirmed first then resolved/closed.

(v0.12.0_release_test) ➜  ASPIRE-Python.git_v0.12.0 git:(v0.12.0) ✗ pwd
/scratch/gbwright/releases/ASPIRE-Python.git_v0.12.0
(v0.12.0_release_test) ➜  ASPIRE-Python.git_v0.12.0 git:(v0.12.0) ✗ pytest --lf
============================= test session starts ==============================
platform linux -- Python 3.8.17, pytest-7.4.1, pluggy-1.3.0
Test order randomisation NOT enabled. Enable with --random-order or --random-order-bucket=<bucket_type>
rootdir: /scratch/gbwright/releases/ASPIRE-Python.git_v0.12.0
configfile: tox.ini
testpaths: tests
plugins: anyio-4.0.0, xdist-3.3.1, random-order-1.1.0, cov-4.1.0
collected 1 item                                                               
run-last-failure: rerun previous 1 failure (skipped 64 files)

tests/test_volume.py F                                                                                                                             [100%]

======================================================================== FAILURES ========================================================================
__________________________________________ test_rotate_broadcast_unicast[dtype=<class 'numpy.float32'>-res=43] ___________________________________________

vols_1 = 3 float32 volumes arranged as a (3,) stack each of size 43x43x43., dtype = <class 'numpy.float32'>

    def test_rotate_broadcast_unicast(vols_1, dtype):
        # Build `Rotation` objects. A singleton for broadcasting and a stack for unicasting.
        # The stack consists of copies of the singleton.
        angles = np.array([pi, pi / 2, 0], dtype=dtype)
        angles = np.tile(angles, (3, 1))
        rot_mat = Rotation.from_euler(angles, dtype=dtype).matrices
        rot = Rotation(rot_mat[0])
        rots = Rotation(rot_mat)

        # Broadcast the singleton `Rotation` across the `Volume` stack.
        vols_broadcast = vols_1.rotate(rot)

        # Unicast the `Rotation` stack across the `Volume` stack.
        vols_unicast = vols_1.rotate(rots)

        for i in range(N):
>           assert np.allclose(vols_broadcast[i], vols_unicast[i])
E           assert False
E            +  where False = <function allclose at 0x7fadb11fe940>(1 float32 volumes arranged as a (1,) stack each of size 43x43x43., 1 float32 volumes arranged as a (1,) stack each of size 43x43x43.)
E            +    where <function allclose at 0x7fadb11fe940> = np.allclose

tests/test_volume.py:368: AssertionError
------------------------------------------------------------------- Captured log call --------------------------------------------------------------------
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend cufinufft
DEBUG    aspire.nufft:__init__.py:83 NFFT backend cufinufft not usable:
    No module named 'pycuda'
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend finufft
DEBUG    aspire.nufft:__init__.py:85 NFFT backend finufft usable.
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend pynfft
DEBUG    aspire.nufft:__init__.py:83 NFFT backend pynfft not usable:
    No module named 'pynfft'
DEBUG    aspire.nufft:__init__.py:92 Selected NFFT backend = finufft.
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
================================================================ short test summary info =================================================================
FAILED tests/test_volume.py::test_rotate_broadcast_unicast[dtype=<class 'numpy.float32'>-res=43] - assert False
=================================================================== 1 failed in 51.24s ===================================================================
(v0.12.0_release_test) ➜  ASPIRE-Python.git_v0.12.0 git:(v0.12.0) ✗ e tests/test_volume.py
(v0.12.0_release_test) ➜  ASPIRE-Python.git_v0.12.0 git:(v0.12.0) ✗ pytest --lf           
============================================================================ test session starts =============================================================================
platform linux -- Python 3.8.17, pytest-7.4.1, pluggy-1.3.0
Test order randomisation NOT enabled. Enable with --random-order or --random-order-bucket=<bucket_type>
rootdir: /scratch/gbwright/releases/ASPIRE-Python.git_v0.12.0
configfile: tox.ini
testpaths: tests
plugins: anyio-4.0.0, xdist-3.3.1, random-order-1.1.0, cov-4.1.0
collected 1 item                                                                                                                                                             
run-last-failure: rerun previous 1 failure (skipped 64 files)

tests/test_volume.py F                                                                                                                                                 [100%]

================================================================================== FAILURES ==================================================================================
____________________________________________________ test_rotate_broadcast_unicast[dtype=<class 'numpy.float32'>-res=43] _____________________________________________________

vols_1 = 3 float32 volumes arranged as a (3,) stack each of size 43x43x43., dtype = <class 'numpy.float32'>

    def test_rotate_broadcast_unicast(vols_1, dtype):
        # Build `Rotation` objects. A singleton for broadcasting and a stack for unicasting.
        # The stack consists of copies of the singleton.
        angles = np.array([pi, pi / 2, 0], dtype=dtype)
        angles = np.tile(angles, (3, 1))
        rot_mat = Rotation.from_euler(angles, dtype=dtype).matrices
        rot = Rotation(rot_mat[0])
        rots = Rotation(rot_mat)

        # Broadcast the singleton `Rotation` across the `Volume` stack.
        vols_broadcast = vols_1.rotate(rot)

        # Unicast the `Rotation` stack across the `Volume` stack.
        vols_unicast = vols_1.rotate(rots)

        for i in range(N):
>           np.testing.assert_allclose(vols_broadcast[i], vols_unicast[i])

tests/test_volume.py:368: 
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 

args = (<function assert_allclose.<locals>.compare at 0x7fcc5cf505e0>, array([[[[7.95078828e+04, 7.76586641e+04, 7.58096641e+....58108594e+04, 7.39615703e+04, ...,
          3.69899365e+03, 1.85006433e+03, 7.26836920e-01]]]],
      dtype=float32))
kwds = {'equal_nan': True, 'err_msg': '', 'header': 'Not equal to tolerance rtol=1e-07, atol=0', 'verbose': True}

    @wraps(func)
    def inner(*args, **kwds):
        with self._recreate_cm():
>           return func(*args, **kwds)
E           AssertionError: 
E           Not equal to tolerance rtol=1e-07, atol=0
E           
E           Mismatched elements: 56976 / 79507 (71.7%)
E           Max absolute difference: 0.04296875
E           Max relative difference: 0.0139712
E            x: array([[[[7.950788e+04, 7.765866e+04, 7.580966e+04, ..., 5.547084e+03,
E                     3.697868e+03, 1.848958e+03],
E                    [7.946497e+04, 7.761577e+04, 7.576673e+04, ..., 5.504084e+03,...
E            y: array([[[[7.950788e+04, 7.765866e+04, 7.580966e+04, ..., 5.547079e+03,
E                     3.697866e+03, 1.848948e+03],
E                    [7.946497e+04, 7.761577e+04, 7.576671e+04, ..., 5.504078e+03,...

../../miniconda3/envs/v0.12.0_release_test/lib/python3.8/contextlib.py:75: AssertionError
----------------------------------------------------------------------------- Captured log call ------------------------------------------------------------------------------
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend cufinufft
DEBUG    aspire.nufft:__init__.py:83 NFFT backend cufinufft not usable:
    No module named 'pycuda'
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend finufft
DEBUG    aspire.nufft:__init__.py:85 NFFT backend finufft usable.
DEBUG    aspire.nufft:__init__.py:48 Trying NFFT backend pynfft
DEBUG    aspire.nufft:__init__.py:83 NFFT backend pynfft not usable:
    No module named 'pynfft'
DEBUG    aspire.nufft:__init__.py:92 Selected NFFT backend = finufft.
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
DEBUG    aspire.nufft:__init__.py:209 nufft passed real_type for signal, converting
DEBUG    aspire.nufft.finufft:finufft.py:44 FinufftPlan adjusted eps=1.1920928955078125e-07 from requested 1e-08.
========================================================================== short test summary info ===========================================================================
FAILED tests/test_volume.py::test_rotate_broadcast_unicast[dtype=<class 'numpy.float32'>-res=43] - AssertionError: 
======================================================================= 1 failed in 114.67s (0:01:54) ========================================================================
(v0.12.0_release_test) ➜  ASPIRE-Python.git_v0.12.0 git:(v0.12.0) ✗ python --version 
Python 3.8.17
aiosignal==1.3.1
alabaster==0.7.13
anyio==4.0.0
argon2-cffi==23.1.0
argon2-cffi-bindings==21.2.0
arrow==1.2.3
aspire==0.12.0
asttokens==2.4.0
async-lru==2.0.4
attrs==23.1.0
Babel==2.12.1
backcall==0.2.0
beautifulsoup4==4.12.2
black==23.7.0
bleach==6.0.0
build==1.0.0
bump2version==1.0.1
bumpversion==0.6.0
cachetools==5.3.1
certifi==2023.7.22
cffi==1.15.1
chardet==5.2.0
charset-normalizer==3.2.0
check-manifest==0.49
click==8.1.7
colorama==0.4.6
comm==0.1.4
confuse==2.0.1
contourpy==1.1.0
coverage==7.3.0
cryptography==41.0.3
cycler==0.11.0
debugpy==1.6.7.post1
decorator==5.1.1
defusedxml==0.7.1
distlib==0.3.7
docutils==0.17.1
exceptiongroup==1.1.3
execnet==2.0.2
executing==1.2.0
fastjsonschema==2.18.0
filelock==3.12.3
finufft==2.1.0
flake8==6.1.0
fonttools==4.42.1
fqdn==1.5.1
frozenlist==1.4.0
gemmi==0.6.2
grpcio==1.57.0
idna==3.4
imageio==2.31.3
imagesize==1.4.1
importlib-metadata==6.8.0
importlib-resources==6.0.1
iniconfig==2.0.0
ipykernel==6.25.2
ipython==8.12.2
ipython-genutils==0.2.0
ipywidgets==8.1.0
isoduration==20.11.0
isort==5.12.0
jaraco.classes==3.3.0
jedi==0.19.0
jeepney==0.8.0
Jinja2==3.1.2
joblib==1.3.2
json5==0.9.14
jsonpointer==2.4
jsonschema==4.19.0
jsonschema-specifications==2023.7.1
jupyter==1.0.0
jupyter-console==6.6.3
jupyter-events==0.7.0
jupyter-lsp==2.2.0
jupyter_client==8.3.1
jupyter_core==5.3.1
jupyter_server==2.7.3
jupyter_server_terminals==0.4.4
jupyterlab==4.0.5
jupyterlab-pygments==0.2.2
jupyterlab-widgets==3.0.8
jupyterlab_server==2.24.0
keyring==24.2.0
kiwisolver==1.4.5
latexcodec==2.0.1
lazy_loader==0.3
markdown-it-py==3.0.0
MarkupSafe==2.1.3
matplotlib==3.7.2
matplotlib-inline==0.1.6
mccabe==0.7.0
mdurl==0.1.2
mistune==3.0.1
more-itertools==10.1.0
mrcfile==1.4.3
msgpack==1.0.5
mypy-extensions==1.0.0
nbclient==0.8.0
nbconvert==7.8.0
nbformat==5.9.2
nest-asyncio==1.5.7
networkx==3.1
notebook==7.0.3
notebook_shim==0.2.3
numpy==1.24.4
overrides==7.4.0
packaging==23.1
pandocfilters==1.5.0
parameterized==0.9.0
parso==0.8.3
pathspec==0.11.2
pexpect==4.8.0
pickleshare==0.7.5
Pillow==10.0.0
pkginfo==1.9.6
pkgutil_resolve_name==1.3.10
platformdirs==3.10.0
pluggy==1.3.0
pooch==1.7.0.post5
prometheus-client==0.17.1
prompt-toolkit==3.0.39
protobuf==4.24.2
psutil==5.9.5
ptyprocess==0.7.0
pure-eval==0.2.2
pybtex==0.24.0
pybtex-docutils==1.0.3
pycodestyle==2.11.0
pycparser==2.21
pydocstyle==6.3.0
pyFFTW==0.13.1
pyflakes==3.1.0
Pygments==2.16.1
pymanopt==2.1.1
pyparsing==3.0.9
pyproject-api==1.6.1
pyproject_hooks==1.0.0
pytest==7.4.1
pytest-cov==4.1.0
pytest-random-order==1.1.0
pytest-xdist==3.3.1
python-dateutil==2.8.2
python-json-logger==2.0.7
pytz==2023.3.post1
PyWavelets==1.4.1
PyYAML==6.0.1
pyzmq==25.1.1
qtconsole==5.4.4
QtPy==2.4.0
ray==2.6.3
readme-renderer==41.0
referencing==0.30.2
requests==2.31.0
requests-toolbelt==1.0.0
rfc3339-validator==0.1.4
rfc3986==2.0.0
rfc3986-validator==0.1.1
rich==13.5.2
rpds-py==0.10.2
scikit-image==0.21.0
scikit-learn==1.3.0
scipy==1.10.1
SecretStorage==3.3.3
Send2Trash==1.8.2
six==1.16.0
snakeviz==2.2.0
sniffio==1.3.0
snowballstemmer==2.2.0
soupsieve==2.5
Sphinx==5.3.0
sphinx-gallery==0.14.0
sphinx-rtd-theme==1.3.0
sphinxcontrib-applehelp==1.0.4
sphinxcontrib-bibtex==2.6.1
sphinxcontrib-devhelp==1.0.2
sphinxcontrib-htmlhelp==2.0.1
sphinxcontrib-jquery==4.1
sphinxcontrib-jsmath==1.0.1
sphinxcontrib-mermaid==0.9.2
sphinxcontrib-qthelp==1.0.3
sphinxcontrib-serializinghtml==1.1.5
stack-data==0.6.2
terminado==0.17.1
threadpoolctl==3.2.0
tifffile==2023.7.10
tinycss2==1.2.1
tomli==2.0.1
tornado==6.3.3
tox==4.11.1
tqdm==4.66.1
traitlets==5.9.0
twine==4.0.2
typing_extensions==4.7.1
uri-template==1.3.0
urllib3==2.0.4
virtualenv==20.24.4
wcwidth==0.2.6
webcolors==1.13
webencodings==0.5.1
websocket-client==1.6.2
widgetsnbextension==4.0.8
zipp==3.16.2
garrettwrong commented 1 year ago

I did make the following change locally, but just to get an idea how far off things were.

diff --git a/tests/test_volume.py b/tests/test_volume.py
index aa3eaf3..4120dca 100644
--- a/tests/test_volume.py
+++ b/tests/test_volume.py
@@ -365,7 +365,7 @@ def test_rotate_broadcast_unicast(vols_1, dtype):
     vols_unicast = vols_1.rotate(rots)

     for i in range(N):
-        assert np.allclose(vols_broadcast[i], vols_unicast[i])
+        np.testing.assert_allclose(vols_broadcast[i], vols_unicast[i])

 def to_vec(vols_1, vec):
(END)
garrettwrong commented 12 months ago

Confirmed this is still occurring, this time on decaf with fresh develop 118d1a7c0f98b206b5f262d55ad2158c1f9bdb39 python 3.8.17 env.

For reference, freeze is below.

Appears to happen with float32, but the absolute amounts are a little concerning (not somethings I'd just atol over).


aiosignal==1.3.1
alabaster==0.7.13
anyio==4.0.0
argon2-cffi==23.1.0
argon2-cffi-bindings==21.2.0
arrow==1.2.3
-e git+ssh://git@github.com/ComputationalCryoEM/ASPIRE-Python.git@118d1a7c0f98b206b5f262d55ad2158c1f9bdb39#egg=aspire
asttokens==2.4.0
async-lru==2.0.4
attrs==23.1.0
Babel==2.12.1
backcall==0.2.0
beautifulsoup4==4.12.2
black==23.9.1
bleach==6.0.0
build==1.0.3
bump2version==1.0.1
bumpversion==0.6.0
cachetools==5.3.1
certifi==2023.7.22
cffi==1.15.1
chardet==5.2.0
charset-normalizer==3.2.0
check-manifest==0.49
click==8.1.7
colorama==0.4.6
comm==0.1.4
confuse==2.0.1
contourpy==1.1.1
coverage==7.3.1
cryptography==41.0.4
cycler==0.11.0
debugpy==1.8.0
decorator==5.1.1
defusedxml==0.7.1
distlib==0.3.7
docutils==0.17.1
exceptiongroup==1.1.3
execnet==2.0.2
executing==1.2.0
fastjsonschema==2.18.0
filelock==3.12.4
finufft==2.1.0
flake8==6.1.0
fonttools==4.42.1
fqdn==1.5.1
frozenlist==1.4.0
gemmi==0.6.3
grpcio==1.58.0
idna==3.4
imageio==2.31.4
imagesize==1.4.1
importlib-metadata==6.8.0
importlib-resources==6.1.0
iniconfig==2.0.0
ipykernel==6.25.2
ipython==8.12.2
ipython-genutils==0.2.0
ipywidgets==8.1.1
isoduration==20.11.0
isort==5.12.0
jaraco.classes==3.3.0
jedi==0.19.0
jeepney==0.8.0
Jinja2==3.1.2
joblib==1.3.2
json5==0.9.14
jsonpointer==2.4
jsonschema==4.19.1
jsonschema-specifications==2023.7.1
jupyter==1.0.0
jupyter-console==6.6.3
jupyter-events==0.7.0
jupyter-lsp==2.2.0
jupyter_client==8.3.1
jupyter_core==5.3.2
jupyter_server==2.7.3
jupyter_server_terminals==0.4.4
jupyterlab==4.0.6
jupyterlab-pygments==0.2.2
jupyterlab-widgets==3.0.9
jupyterlab_server==2.25.0
keyring==24.2.0
kiwisolver==1.4.5
latexcodec==2.0.1
lazy_loader==0.3
markdown-it-py==3.0.0
MarkupSafe==2.1.3
matplotlib==3.7.3
matplotlib-inline==0.1.6
mccabe==0.7.0
mdurl==0.1.2
mistune==3.0.1
more-itertools==10.1.0
mrcfile==1.4.3
msgpack==1.0.7
mypy-extensions==1.0.0
nbclient==0.8.0
nbconvert==7.8.0
nbformat==5.9.2
nest-asyncio==1.5.8
networkx==3.1
nh3==0.2.14
notebook==7.0.4
notebook_shim==0.2.3
numpy==1.24.4
overrides==7.4.0
packaging==23.1
pandocfilters==1.5.0
parso==0.8.3
pathspec==0.11.2
pexpect==4.8.0
pickleshare==0.7.5
Pillow==10.0.1
pkginfo==1.9.6
pkgutil_resolve_name==1.3.10
platformdirs==3.10.0
pluggy==1.3.0
pooch==1.7.0
prometheus-client==0.17.1
prompt-toolkit==3.0.39
protobuf==4.24.3
psutil==5.9.5
ptyprocess==0.7.0
pure-eval==0.2.2
pybtex==0.24.0
pybtex-docutils==1.0.3
pycodestyle==2.11.0
pycparser==2.21
pydocstyle==6.3.0
pyFFTW==0.13.1
pyflakes==3.1.0
Pygments==2.16.1
pymanopt==2.1.1
pyparsing==3.1.1
pyproject-api==1.6.1
pyproject_hooks==1.0.0
pytest==7.4.2
pytest-cov==4.1.0
pytest-random-order==1.1.0
pytest-xdist==3.3.1
python-dateutil==2.8.2
python-json-logger==2.0.7
pytz==2023.3.post1
PyWavelets==1.4.1
PyYAML==6.0.1
pyzmq==25.1.1
qtconsole==5.4.4
QtPy==2.4.0
ray==2.7.0
readme-renderer==42.0
referencing==0.30.2
requests==2.31.0
requests-toolbelt==1.0.0
rfc3339-validator==0.1.4
rfc3986==2.0.0
rfc3986-validator==0.1.1
rich==13.5.3
rpds-py==0.10.3
scikit-image==0.21.0
scikit-learn==1.3.1
scipy==1.10.1
SecretStorage==3.3.3
Send2Trash==1.8.2
six==1.16.0
snakeviz==2.2.0
sniffio==1.3.0
snowballstemmer==2.2.0
soupsieve==2.5
Sphinx==5.3.0
sphinx-gallery==0.14.0
sphinx-rtd-theme==1.3.0
sphinxcontrib-applehelp==1.0.4
sphinxcontrib-bibtex==2.6.1
sphinxcontrib-devhelp==1.0.2
sphinxcontrib-htmlhelp==2.0.1
sphinxcontrib-jquery==4.1
sphinxcontrib-jsmath==1.0.1
sphinxcontrib-mermaid==0.9.2
sphinxcontrib-qthelp==1.0.3
sphinxcontrib-serializinghtml==1.1.5
stack-data==0.6.2
terminado==0.17.1
threadpoolctl==3.2.0
tifffile==2023.7.10
tinycss2==1.2.1
tomli==2.0.1
tornado==6.3.3
tox==4.11.3
tqdm==4.66.1
traitlets==5.10.1
twine==4.0.2
typing_extensions==4.8.0
uri-template==1.3.0
urllib3==2.0.5
virtualenv==20.24.5
wcwidth==0.2.6
webcolors==1.13
webencodings==0.5.1
websocket-client==1.6.3
widgetsnbextension==4.0.9
zipp==3.17.0
garrettwrong commented 11 months ago

Josh has confirmed he can reproduce and is narrowing down the issue. Seems like a poor volume was chosen for reconstruction, but should understand if/why any results from finufft are different under broadcasting before we commit to changing the test.

j-c-c commented 11 months ago

Below is a minimal script that isolates a call to finufft using a stack of volumes vs singleton volumes. I have the volume from our unit test in there but have confirmed this still fails for singles using a random volume.

import numpy as np                                                                                                       
import finufft                                                                                                           

singles = np.float32, np.complex64                                                                                       
doubles =  np.float64, np.complex128                                                                                     

def run_finufft(dtypes):                                                                                                 
    # Generate a stack of volumes.                                                                                       
    dtype, complex_dtype = dtypes                                                                                        
    res = 23

    N = 3  # number of volumes                                                                                           
    vols = np.arange(1, 1 + N * res**3, dtype=dtype).reshape(N, res, res, res)                                           
    vols = np.array(vols,  copy=False, dtype=complex_dtype, order="C")                                                   

    # Generate set of fourier grid points.                                                                               
    pts_1d = np.linspace(-np.pi, np.pi, res, dtype=dtype)                                                                
    x, y, z = np.meshgrid(pts_1d, pts_1d, pts_1d)                                                                        
    fourier_pts = np.vstack([x.flatten(), y.flatten(), z.flatten()])                                                     

    # Make finufft plan for stack of volumes.                                                                            
    eps = 1e-8                                                                                                           
    if np.finfo(dtype).eps > eps:                                                                                        
        eps = np.finfo(dtype).eps                                                                                        
    transform_plan = finufft.Plan(nufft_type=2, n_modes_or_dim=(res,)*3, eps=eps, n_trans=3, dtype=dtype)                
    transform_plan.setpts(*fourier_pts)                                                                                  

    # Execute transform on volume stack.                                                                                 
    f_vols = transform_plan.execute(vols)                                                                                
    f_vols = f_vols.reshape(-1, res, res, res)                                                                           

    # Make finufft plan for single volumes.                                                                              
    transform_plan = finufft.Plan(nufft_type=2, n_modes_or_dim=(res,)*3, eps=eps, n_trans=1, dtype=dtype)                
    transform_plan.setpts(*fourier_pts)                                                                                  

    # Execute transform on single volumes.                                                                               
    f_vol_0 = transform_plan.execute(vols[0]).reshape(res, res, res)                                                     
    f_vol_1 = transform_plan.execute(vols[1]).reshape(res, res, res)                                                     
    f_vol_2 = transform_plan.execute(vols[2],).reshape(res, res, res)                                                    

    # Test that we get the same transforms with stack and singletons.                                                    
    print(dtypes)                                                                                                        
    print(np.allclose(f_vols[0], f_vol_0))                                                                               
    print(np.allclose(f_vols[1], f_vol_1))                                                                               
    print(np.allclose(f_vols[2], f_vol_2))                                                                               

run_finufft(doubles)                                                                                                     
run_finufft(singles)

finufft_check.txt

j-c-c commented 11 months ago

Results of running this script:

(flaky_vol_test) [jc5485@decaf ASPIRE-python.flaky_vol_test]$ python finnuft_check.py
(<class 'numpy.float64'>, <class 'numpy.complex128'>)
True
True
True
(<class 'numpy.float32'>, <class 'numpy.complex64'>)
False
False
False
garrettwrong commented 11 months ago

Cool! Let's see how far apart via max(abs()) and rmse.

It is possible that highly optimized code is not deterministic, which should be fine, but that should be pretty small if it is the case.

j-c-c commented 11 months ago

Some follow up here. I've tested finufft's latest upstream code with this script and got the same results. Here they are with max absolute difference and rmse values added:

(<class 'numpy.float64'>, <class 'numpy.complex128'>)
allclose: True,  max_abs_diff: 5.587935447692871e-09,  rmse: 1.389204659730386e-10
allclose: True,  max_abs_diff: 7.475015174206003e-09,  rmse: 3.1540061897167237e-10
allclose: True,  max_abs_diff: 5.962675727405252e-08,  rmse: 7.363349092859463e-10
(<class 'numpy.float32'>, <class 'numpy.complex64'>)
allclose: False,  max_abs_diff: 2.0,  rmse: 0.055754758417606354
allclose: False,  max_abs_diff: 16.03447914123535,  rmse: 0.21659733355045319
allclose: False,  max_abs_diff: 8.000244140625,  rmse: 0.24948608875274658
j-c-c commented 11 months ago

Adding a note that in order to compile the latest upstream finnuft code I had to alter the g++ compiler flag --std from c++14 to c++1y since our version of g++ (4.8.5, finnuft has >=5.0 dependency) does not support c++14.

garrettwrong commented 11 months ago

If finufft has a gcc 5 or greater dependency, we should be using a compatible compiler for this effort.

I doubt it will make a difference, but It could, and the whole point is to end up witha viable piece of evidence for a bug report. If we ignore a major build requirement the efficacy of the report is dubious. It really depends on the reason(s) for that requirement, which I am not aware of at this time. Best to follow the instructions.

Totally hypothetical, but consider the following. What if the issue is only occurs on older (somehow incompatible) compilers and this issue is caused by a similar bad build being packaged for publication.

j-c-c commented 11 months ago

I agree. That's why I made the note. I might be mistaken but I thought I needed sudo to upgrade the compiler.

garrettwrong commented 11 months ago

Nope, you can 100% build a compiler if you needed to (that's a right of passage sort of thing), but you don't need to. Several are available on the machine already.

j-c-c commented 11 months ago

Ah ok. I'll rebuild with the proper compiler and see if the bug still occurs

j-c-c commented 11 months ago

Closing this issue. The flaky test was resolved in #1024 and a bug report was made to finufft https://github.com/flatironinstitute/finufft/issues/363