A bug in nutationlib's equation_of_the_equinoxes_complimentary_terms()

[trufanov-nok]

I was trying to debug my TEME<->ITRF conversion by comparing to skyfield results and found out a point there my results are starting to differ. It's here: https://github.com/skyfielders/python-skyfield/blob/master/skyfield/nutationlib.py#L177

    fa[4] = (( 450160.398036 +
              (-482890.5431 +
              (      7.4722 +
              (      0.007702 +
              (     -0.00005939)
              * t) * t) * t) * t) * ASEC2RAD
              + (-5.0*t % 1.0) * tau)

In particular at -5.0*t % 1.0.
I was surprised that modulo operation with negative value behaves weird in python. For example:

>>> t =  0.16322491910627326
>>> -5.0*t % 1.0
0.1838754044686337

while I would expect

>>> -5.0*t 
-0.8161245955313663

And that's how it would work in C with fmod(-5.0*0.16322491910627326, 1.0).

So I tried to find an article or another source for this method. And in SOFA's C codebase (http://iausofa.org/2021_0512_C/sofa_c-20210512.tar.gz) at nut80.c::224 I found:

/* Longitude of the mean ascending node of the lunar orbit on the */
/* ecliptic, measured from the mean equinox of date. */
   om = iauAnpm(
        (450160.280 + (-482890.539 + (7.455 + 0.008 * t) * t) * t)
        * DAS2R + fmod(-5.0 * t, 1.0) * D2PI);

I assume that C code was a source and then it was ported to python. In this case Skyfield has a bug because C fmod(-5.0 * t, 1.0) is not equivalent to python's -5.0*t % 1.0.
I don't know python lang but perhaps it should be?:

>>> 5.0*t % 1.0 * -1
-0.8161245955313663

[trufanov-nok]

Or better, keeping in mind that the module has from numpy import ... fmod to use fmod(-5.0*t, 1.0)
Also the fa %= tau at https://github.com/skyfielders/python-skyfield/blob/master/skyfield/nutationlib.py#L189 will produce incorrect results. It should be replaced with fa = fmod(fa, tau)

[brandon-rhodes]

As we discovered in the PR, this doesn't seem to affect Skyfield's results measurably, so for the moment I'm included to leave that code alone. Please do let me know if you find a case in which Skyfield returns an incorrect angle, though!