1.1.1: TestPhysicsSphericalRepresentation failure with numpy 1.11.0 beta 2

[olebole]

There are two failures in TestPhysicsSphericalRepresentation with the new numy beta release:

____________ TestPhysicsSphericalRepresentation.test_init_quantity _____________
self = <astropy.coordinates.tests.test_representation.TestPhysicsSphericalRepresentation object at 0x7f507be01630>

    def test_init_quantity(self):

        s3 = PhysicsSphericalRepresentation(phi=8 * u.hourangle, theta=5 * u.deg, r=10 * u.kpc)
>       assert s3.phi == 8. * u.hourangle
E       assert <Angle 7.999999999999999 hourangle> == (8.0 * Unit("hourangle"))
E        +  where <Angle 7.999999999999999 hourangle> = <PhysicsSphericalRepresentation (phi, theta, r) in (hourangle, deg, kpc)\n    (8.0, 5.0, 10.0)>.phi
E        +  and   Unit("hourangle") = u.hourangle

/usr/lib/python3/dist-packages/astropy/coordinates/tests/test_representation.py:331: AssertionError
____________ TestPhysicsSphericalRepresentation.test_init_phitheta _____________

self = <astropy.coordinates.tests.test_representation.TestPhysicsSphericalRepresentation object at 0x7f50817cb6a0>

    def test_init_phitheta(self):

        s2 = PhysicsSphericalRepresentation(Angle(8, u.hour),
                                            Angle(5, u.deg),
                                            Distance(10, u.kpc))

>       assert s2.phi == 8. * u.hourangle
E       assert <Angle 7.999999999999999 hourangle> == (8.0 * Unit("hourangle"))
E        +  where <Angle 7.999999999999999 hourangle> = <PhysicsSphericalRepresentation (phi, theta, r) in (hourangle, deg, kpc)\n    (8.0, 5.0, 10.0)>.phi
E        +  and   Unit("hourangle") = u.hourangle

/usr/lib/python3/dist-packages/astropy/coordinates/tests/test_representation.py:345: AssertionError

This code worked well with numpy-1.10. wcslib version for both (just in case that it matters) is 5.13.
Full test log here

[olebole]

There is another failure (from the same log that may have the same cause:

________________________________ test_longitude ________________________________

    def test_longitude():
[...]
        # Test setting a mutable value and having it wrap
        lon[1] = -10 * u.deg
>       assert np.all(lon == Angle(['10d', '350d']))
E       assert <function all at 0x7f508e131d08>(<Longitude [  10., 350.] deg> == <Angle [  10., 350.] deg>
E        +  where <function all at 0x7f508e131d08> = np.all
E         Use -v to get the full diff)

/usr/lib/python3/dist-packages/astropy/coordinates/tests/test_angles.py:643: AssertionError

[mhvk]

Hmm, similar failures in recent travis tests against master. Looks like rounding issues, but unclear why those did not turn up before.

[njsmith]

Filed numpy/numpy#7161 to track -- please do report these things upstream! We want to know! :-)

Nothing comes to mind as an obvious cause, but it seems eminently bisectable if someone wants to try. (Basically: write a little script that exits successfully if the test passes and exits with an error code if the test fails, then clone numpy and do git bisect start; git bisect bad master; git bisect good origin/maintenance/1.10.x; git bisect run ./runtests.py my-script.py.)

[njsmith]

The two TestPhysicsSphericalRepresentation failures are due to a change to make divmod rounding more consistent: numpy/numpy@34c2369
My guess is that for those the change is OK, and these are just over-picky tests?

The test_longitude test also bisects to that same commit. This test I'm a little more nervous about, because of the line in the test output:

<Longitude [  10., 350.] deg> == <Angle [  10., 350.] deg>

looks like somehow there's a type error? This feels like the sort of bug that @mhvk always comes hassling us about as soon as we release a new version, like maybe we aren't preserving types properly somewhere :-)

[mhvk]

@njsmith - thanks very much for investigating! I'm still overwhelmed with graduate admission, but will try to have a look. The type change is fine in this case -- the intent is just to test the values, not the class -- though the test should probably be changed just so there is no confusion.

[charris]

The type error is probably because scalar math remainder is now done using c operators instead of calling the ufunc fmod. This could be a problem with other scalar functions that don't use ufuncs...

[njsmith]

I just installed numpy 1.9.3 + astropy 1.1.1, and that configuration shows the same funny types as are reported above -- so unless someone tells me otherwise, I'm going to assume that the only change in numpy 1.11 is that the division rounding is slightly different (and more consistent, see numpy/numpy@34c2369 for details), and that the only issue here is that astropy has some tests that are a bit too strict. (test_longitude uses == to compare floats...)

If that's wrong and there's something in numpy causing a genuine problem, please let us know and comment here or re-open numpy/numpy#7161; otherwise we'll assume that this is not a blocker for the numpy 1.11 release.

[mhvk]

@njsmith - Agreed with your analysis, and funny that it an issue I raised myself that has returned to bite us (though, as you said, in a way that just exposes that our tests were to strict; am still a bit surprised about that, as I thought I wrote them on purpose in integer degrees, so rounding errors should not be a problem; but anyway, that would be on our side too).

[njsmith]

I haven't looked at the actual change in detail -- so it's possible that
the old rounding was better :-) not intended that way obviously, but
software is tricky. If you look into it and decide that you really do think
that the exact tests ought to work then please let us know and we can talk
about it. I'm tentatively not worrying about following up on that as a
potential release critical bug though because AFAIK the worst case is
1ulp-ish errors.
On Feb 9, 2016 5:04 PM, "Marten van Kerkwijk" notifications@github.com
wrote:

@njsmith https://github.com/njsmith - Agreed with your analysis, and
funny that it an issue I raised myself that has returned to bite us
(though, as you said, in a way that just exposes that our tests were to
strict; am still a bit surprised about that, as I thought I wrote them on
purpose in integer degrees, so rounding errors should not be a problem; but
anyway, that would be on our side too).

—
Reply to this email directly or view it on GitHub
#4548 (comment).

[adrn]

The reason they are not identical is that the input value angles are always wrapped with .wrap_at(360*u.degree), which actually uses numpy.mod() even if the input angle was already within the correct bounds. Two options here: (1) I can change these tests to use assert_allclose(), or (2) we can avoid using mod() if the angle is already properly bounded. Thoughts? @mhvk @eteq

[njsmith]

Check that it still happens with numpy 1.11rc1 or so before making any changes -- the handling of mod for reworked a bit in response to issues like this.

[mhvk]

While we may want to consider using assert_allclose here independently of this issue, is it still a problem with the more recent 1.11? My guess is that numpy/numpy#7258 has fixed it. ... ... ... Indeed, running a quick test, I don't see a problem with maintenance/1.11.x.

[mhvk]

In fact, now that the tests worked, I think we should consider this fixed. I do like the idea that one can use equality on integer numbers, and in this case finding this bug has led to what I think is a very nice improvement in numpy, so all would seem well!

[mhvk]

p.s. thanks for reporting, @olebole