Function tan_eval

system/lib/llvm-libc/src/math/generic/tan.cpp:38–72

namespace {

LIBC_INLINE double tan_eval(const DoubleDouble &u, DoubleDouble &result) {
  // Evaluate tan(y) = tan(x - k * (pi/128))
  // We use the degree-9 Taylor approximation:
  //   tan(y) ~ P(y) = y + y^3/3 + 2*y^5/15 + 17*y^7/315 + 62*y^9/2835
  // Then the error is bounded by:
  //   |tan(y) - P(y)| < 2^-6 * |y|^11 < 2^-6 * 2^-66 = 2^-72.
  // For y ~ u_hi + u_lo, fully expanding the polynomial and drop any terms
  // < ulp(u_hi^3) gives us:
  //   P(y) = y + y^3/3 + 2*y^5/15 + 17*y^7/315 + 62*y^9/2835 = ...
  // ~ u_hi + u_hi^3 * (1/3 + u_hi^2 * (2/15 + u_hi^2 * (17/315 +
  //                                                     + u_hi^2 * 62/2835))) +
  //        + u_lo (1 + u_hi^2 * (1 + u_hi^2 * 2/3))
  double u_hi_sq = u.hi * u.hi; // Error < ulp(u_hi^2) < 2^(-6 - 52) = 2^-58.
  // p1 ~ 17/315 + u_hi^2 62 / 2835.
  double p1 =
      fputil::multiply_add(u_hi_sq, 0x1.664f4882c10fap-6, 0x1.ba1ba1ba1ba1cp-5);
  // p2 ~ 1/3 + u_hi^2 2 / 15.
  double p2 =
      fputil::multiply_add(u_hi_sq, 0x1.1111111111111p-3, 0x1.5555555555555p-2);
  // q1 ~ 1 + u_hi^2 * 2/3.
  double q1 = fputil::multiply_add(u_hi_sq, 0x1.5555555555555p-1, 1.0);
  double u_hi_3 = u_hi_sq * u.hi;
  double u_hi_4 = u_hi_sq * u_hi_sq;
  // p3 ~ 1/3 + u_hi^2 * (2/15 + u_hi^2 * (17/315 + u_hi^2 * 62/2835))
  double p3 = fputil::multiply_add(u_hi_4, p1, p2);
  // q2 ~ 1 + u_hi^2 * (1 + u_hi^2 * 2/3)
  double q2 = fputil::multiply_add(u_hi_sq, q1, 1.0);
  double tan_lo = fputil::multiply_add(u_hi_3, p3, u.lo * q2);
  // Overall, |tan(y) - (u_hi + tan_lo)| < ulp(u_hi^3) <= 2^-71.
  // And the relative errors is:
  // |(tan(y) - (u_hi + tan_lo)) / tan(y) | <= 2*ulp(u_hi^2) < 2^-64
  result = fputil::exact_add(u.hi, tan_lo);
  return fputil::multiply_add(fputil::FPBits<double>(u_hi_3).abs().get_val(),
                              0x1.0p-51, 0x1.0p-102);
}

#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
// Accurate evaluation of tan for small u.