Which is just a support function for cos (which doesn't do argument reduction). See [1] for the actual function, and [2] & [3] for the argument reduction code.
You really only need to compute one eighth of a circle segment (so, zero to pi/4), everything else you can do by symmetry by switching signs and flipping x, y. (And if you supported a larger range, you’d waste precision).
The function supports doubledouble precision (which is a trick, different from long double or quadruple precision, in which you express a number as an unevaluated sum of two doubles, the head (large) and the tail (much smaller)).
So, the “real” cosine function first takes its argument, and then reduces it modulo pi/2 to the desired range -pi/4, +pi/4 (with minimal loss of precision - you utilise the sign bit, so that’s better than using the range 0, +pi/2) in double double precision (rem_pio2 returns an int denoting which quadrant the input was in, and two doubles), and those two doubles are then passed to this function (or its cousin __sin, with the sign flipped as appropriate, depending on the quadrant).
Many people really have put in quite some care and thought to come up with this neat and precise solution.
You really only need to compute one eighth of a circle segment (so, zero to pi/4), everything else you can do by symmetry by switching signs and flipping x, y. (And if you supported a larger range, you’d waste precision).
In high school I was writing a paint program in BASIC on the Apple ][. Drawing circles was mind-numbingly slow until I read an article somewhere that pointed this out. It was like magic.
[1] https://github.com/ifduyue/musl/blob/master/src/math/cos.c
[2] https://github.com/ifduyue/musl/blob/master/src/math/__rem_p...
[3] https://github.com/ifduyue/musl/blob/master/src/math/__rem_p...