stillwater-sc / universal / 29844909844

85%

master: 84%

Build Type

push

github

Committed by web-flow

Commit Message

feat(elreal): Payne-Hanek range reduction for sin/cos/tan at large |x| (#1050) (#1178)

detail::sincos reduced x = n*(pi/2) + t using a host-double estimate
nd = round(to_double_approx(x)/(pi/2)). For |x| >= ~2^53 that double loses all
the fractional bits of x/(pi/2), so nd (and the octant n mod 4) were wrong and
the cancellation x - nd*(pi/2) left a garbage reduced argument -- sin/cos/tan
were unusable past ~1e15.

For |x| >= 2^50 sincos now uses a Payne-Hanek reduction that works entirely in
ZBCL space:
- reddepth blocks of pi/(2) chosen to cover |x|'s binade (from x.head().exponent(),
  the idiom log() already uses) plus the output precision;
- q = div_online(|x|, pi/2) -- the quotient at full precision;
- N = zbcl_round_to_int(q) -- the nearest integer, accumulated into a wide
  integer<256> by scaling each block's exact m*2^p contribution (frexp gives the
  signed k-bit significand), so N is exact even past 2^53 (n ~ 6.4e19 for 1e20);
- t = (q - N)*(pi/2), a well-conditioned multiply of the small fraction q - N
  (pooled + priestRenorm to keep 0-overlap, #1044), never a subtraction of
  near-equal large quantities.
The octant recombination and the Maclaurin series are unchanged; the reduction
works on |x| and folds the sign back (sin odd, cos even). |x| < 2^50 keeps the
fast host-double path (exact there).

Validated on the double host: sin/cos(1e3, 1e6, 1e12, 1e20) agree with mpmath
320-digit references to 180-194 digits, sin^2+cos^2 == 1 holds to ~193 digits at
1e20, tan == sin/cos, and the reduced argument stays 0-overlap. Large-|x| on a
narrow host (float, which cannot represent integers past 2^24) is out of scope,
noted alongside the existing asin/acos double-host caveat.

- trigonometry.hpp: zbcl_round_to_int + int_to_zbcl helpers, sincos rewrite,
  online_multiply/online_divide includes, header notes.
- reference_constants.hpp + generate_reference_constants.py: sin/cos(1e3..1e20)
  mpmath references (arguments are exact doubles).
- elastic... (continued)

Coverage Stats

61 of 62 new or added lines in 1 file covered. (98.39%)

33 existing lines in 4 files now uncovered.

49893 of 58859 relevant lines covered (84.77%)

7019651.64 hits per line