tests/test_Galactic.py::TestGalactic::test_proper_motion_identity failing

[dlakaplan]

This test evaluates a position at the POSEPOCH and compares against the baseline position. The threshold is set to 1e-11 arcsec. Normally this passes but I've been having intermittent failures. I don't know why it's intermittent (it should be deterministic) and only on some architectures. I also don't know if this is due to #1748.

Any ideas @abhisrkckl ?

[dlakaplan]

Errors look like:

=================================== FAILURES ===================================
___________________ TestGalactic.test_proper_motion_identity ___________________

self = <test_Galactic.TestGalactic object at 0x7f0ba784c380>

    def test_proper_motion_identity(self):
        # sanity check that evaluation at POSEPOCH returns something very close to 0
        J0613_icrs = self.modelJ0613.coords_as_ICRS()
        J0613_icrs_alt = self.modelJ0613.coords_as_ICRS(
            epoch=self.modelJ0613.POSEPOCH.quantity.mjd
        )
        sep = J0613_icrs_alt.separation(J0613_icrs)
>       assert sep < 1e-11 * u.arcsec
E       assert <Angle 1.33203881e-14 deg> < (1e-11 * Unit("arcsec"))
E        +  where Unit("arcsec") = u.arcsec

[dlakaplan]

And 1.3e-14 deg = 5e-11 arcsec

[abhisrkckl]

I don't think this is due to #1748. This test works on my system on the following setup:

>>> pint.print_info()
# Created: 2024-05-07T11:10:34.685880
# PINT_version: 1.0+59.g86f243c5
# User: Abhimanyu Susobhanan (abhimanyu)
# Host: abhimanyu-HP-Envy-x360-2-in-1-Laptop-15-fh0xxx
# OS: Linux-6.5.0-28-generic-x86_64-with-glibc2.35
# Python: 3.12.0 | packaged by Anaconda, Inc. | (main, Oct  2 2023, 17:29:18) [GCC 11.2.0]
# endian: little
# numpy_version: 1.26.0
# numpy_longdouble_precision: float128
# scipy_version: 1.11.3
# astropy_version: 5.3.4
# pyerfa_version: 2.0.0
# jplephem_version: 2.20
# matplotlib_version: 3.8.0
# loguru_version: 0.7.2
# Python_prefix: /home/abhimanyu/miniconda3/envs/pint-devel
# PINT_file: /home/abhimanyu/Work/PINT/src/pint/__init__.py
# Environment: conda
# conda_prefix: /home/abhimanyu/miniconda3/envs/pint-devel

But it fails on the following setup:

>>> pint.print_info()
# Created: 2024-05-07T11:13:44.844470
# PINT_version: 1.0+59.g86f243c5
# User: Abhimanyu Susobhanan (abhimanyu)
# Host: abhimanyu-HP-Envy-x360-2-in-1-Laptop-15-fh0xxx
# OS: Linux-6.5.0-28-generic-x86_64-with-glibc2.35
# Python: 3.12.3 | packaged by Anaconda, Inc. | (main, May  6 2024, 19:46:43) [GCC 11.2.0]
# endian: little
# numpy_version: 1.26.4
# numpy_longdouble_precision: float128
# scipy_version: 1.13.0
# astropy_version: 6.1.0
# pyerfa_version: 2.0.1.4
# jplephem_version: 2.22
# matplotlib_version: 3.8.4
# loguru_version: 0.7.2
# Python_prefix: /home/abhimanyu/miniconda3/envs/pint-new
# PINT_file: /home/abhimanyu/Work/PINT/src/pint/__init__.py
# Environment: conda
# conda_prefix: /home/abhimanyu/miniconda3/envs/pint-new

So probably a version issue.

[abhisrkckl]

astropy v6.1.0 was tagged two days ago. So it may be an astropy-related issue.

[abhisrkckl]

The test goes through when I restrict astropy<6.1.0.

[abhisrkckl]

The issue seems to go away if I remove the add_dummy_distance and remove_dummy_distance calls in the get_psr_coords function.

But that causes other tests to fail.

[abhisrkckl]

@dlakaplan What do you think?

[dlakaplan]

Thanks for investigating. Maybe the reason it was showing up only in some tests and not others was just a coincidence about when they started wrt when the astropy bump happened.

I think the simplest fix for now would be to increase the tolerance in the failing test. 5e-11 arcsec is still very small and it would pass. Or we could restrict the astropy version until we understand the origin of the difference. If those dummy functions are the issue there might be some change to the API.

[dlakaplan]

Oddly, this test compares two copies of what should be the same coordinate against each other. But I'm finding that the distance returned is not symmetric:
A.separation(B) = 1.3320388147503583e-14
but
B.separation(A) = 3.975693351829395e-16
So one looks good, while the other isn't.

I think we could fix this test by reversing the order of the calculation, but I'll try to understand the origin.

[dlakaplan]

OK, narrowing things down. I think that somehow this issue relates to units & precision. It compares the RA as hourangles, not degrees. The difference there is ~machine precision (8e-16). But then it gets multiplied by 15 to get degrees, and we see a difference of ~1e-14. I don't know exactly why this computation is different here, but I think that it isn't meant to be. Whether this is an astropy bug or what it will only crop up in rare situations. So I think we can fix just by reversing the operation order.