Make day + fraction separation more robust

[olebole]

This should fix astropy/astropy#5425.
The code behaved diferently depending whether one uses the flags -O2 -fcaller-save (which gives correct results) or -O2 -fno-callers-save (which gives wrong results).

[mhvk]

Great! This definitely needs to go upstream to SOFA too... Not sure what the process is; cc @scottransom, @eteq.

For ease of reference, let me copy your example calls from astropy/astropy#5425 (comment)

Checking the conversion (on the -O2 lib):

>>> import astrop._erfa as erfa
>>> erfa.taiutc(2456925.0, 0.4876736111111111)
(array(2456925.0), array(1.4872685185185186))

which is different without optimization:

>>> erfa.taiutc(2456925.0, 0.4876736111111111)
(array(2456925.0), array(0.4872685185185185))

[olebole]

@eteq, @scottransom, could I ask for a comment here and in #41? I could add a bit on the README and/or contact the SOFA team. What would you propose?

[eteq]

@olebole The usual process is:

1. Report it to SOFA
2. Wait for a response from SOFA
3. If they say "yes it's a bug that we'll fix", we merge the PR immediately w/ the README note indicating the change, as it is sometimes months before SOFA gets to it and issues a new release.
4. If they don't, we have a further discussion to decide if ERFA has call to diverge on this issue
5. Do an erfa release and then integrate that into Astropy

I've been waiting to hear from @scottransom, as he has often conveyed it on, but he might be busy. If we haven't heard from him, by, say, Tuesday (it's been a holiday in the US since ~Wed), I'll try some of the other channels. But definitely want to get this fixed and into astropy by the next release (~5 weeks), so if we get stuck on waiting for SOFA we might just merge.

[scottransom]

Hi All... sorry for the delay. Yes, I've been doing US holiday and avoiding email. I'll send this to the SOFA board members and do a bit of testing myself. I suspect that it will go into the next SOFA release (assuming that it has no strange side-effects).

[eteq]

Ok, great, thanks @scottransom! Let us know as soon as you know here if it looks like it will go in - if so we can go ahead and merge this after updating the README

[scottransom]

Quick update on this. This is a confirmed bug in SOFA and we are working on getting a new release out with a bug fix. The fix will be slightly different than the one proposed here as we will likely use the "dnint" macro in sofam.h which replicates the behavior of the original Fortran function ANINT which was used in SLALIB.

If anyone has specific examples as to how the bug could show up and under what circumstances (and how likely that was), we would like to put that on the SOFA website when we release a new version. Thanks!

[olebole]

@scottransom Here is an example for erfa (in Python). Translated to C:

#include <stdlib.h>
#include <stdio.h>
#include <erfa.h>
int main(void) {
  double u1, u2;
  erfTaiutc(2456925.0, 0.4876736111111111, &u1, &u2);
  printf("%f %f\n", u1, u2);
  exit(0);
}

This problem appears under i386 with the 387 FPU, with -O2 when compiling the library. Then the output is

2456925.000000 1.487269

while it should be

2456925.000000 0.487269

[eteq]

@scottransom - out of curiosity, do you have a sense of when the next SOFA release will be? Basically I'm wondering if we should merge this vs waiting for the new SOFA (if it'll be soon).

Also, a side note (just for my own curiosity): @olebole, do you actually have a 386 w/387, or do you have an emulator of some kind?

[olebole]

@eteq: I just run it in a 32-bit Debian chroot on a very recent CPU. What we call "386" is actually a 686; everything below Pentium Pro doesn't work anymore.

[scottransom]

@eteq I suspect that a new release will be made within the next week or so. I just reported #41 as well, but that one seems like a no brainer, so I'm hoping it will get approved and added very quickly.

[eteq]

OK, in that case I'd say lets wait on this until we see that (same with #41). Assuming they end up on a SOFA that is then released soon, we can merge these (just to keep @olebole's contribution in the changelog), but then immediately pull in the new SOFA to make sure we're in-line with that.

[scottransom]

SOFA has already acknowledged that #41 is a good fix and we have a version for testing already. So I think that is a good plan. I suspect a new release quite soon.

[mhvk]

Supposedly this is included in 12c, which we should rebase on (#44)

[eteq]

Merging this so as to keep a record, but it well be superceded shortly by a PR to upgrade to the latest SOFA