maths/fast: actually account for symmetries in sin

This commit is contained in:
Danny Robson 2021-05-13 17:01:30 +10:00
parent 5cedd22d7a
commit e1f864cd88

View File

@ -21,9 +21,11 @@ namespace cruft::maths::fast {
return std::fma (b, t, a) + std::fma (d, t, c) * t2; return std::fma (b, t, a) + std::fma (d, t, c) * t2;
} }
// approximates sin over the range [0; pi/4] using a polynomial. // approximates sin over the range [-2pi, 2pi] using a polynomial
// approximation over [0; pi/4]
// //
// the relative error should be on the order of 1e-6. // the relative error should be on the order of 1e-6.
// the absolute error should be maximally around 1e-5
// //
// sollya: // sollya:
// display=hexadecimal; // display=hexadecimal;
@ -38,15 +40,21 @@ namespace cruft::maths::fast {
constexpr float constexpr float
sin (float x) noexcept sin (float x) noexcept
{ {
CHECK_LT (x, 2 * pi<float>); CHECK_LE (x, 2 * pi<float>);
CHECK_GT (x, -2 * pi<float>); CHECK_GE (x, -2 * pi<float>);
// flip negative angles due to having an even function // Push values into the range [0, 2pi]
float neg = x < 0 ? x *= -1, -1 : 1; if (x < 0)
x += 2 * pi<float>;
// reduce the angle due to the period // Values in the range [pi, 2pi] are the negation of [0, pi]
int periods = static_cast<int> (x / 2 / pi<float>); float neg = x > pi<float> ? -1 : 1;
x -= periods * 2 * pi<float>; if (x > pi<float>)
x -= pi<float>;
// Values in [0,pi] are symmetric about pi/2
if (x > pi<float> / 2)
x = pi<float> - x;
// 1,3,5,7: 1, -0x1.555546p-3f, 0x1.11077ap-7f, -0x1.9954e8p-13f; // 1,3,5,7: 1, -0x1.555546p-3f, 0x1.11077ap-7f, -0x1.9954e8p-13f;
// 1,3,5,7,9: 1, -0x1.555556p-3f, 0x1.11115cp-7f, -0x1.a0403ap-13f, 0x1.75dc26p-19f; // 1,3,5,7,9: 1, -0x1.555556p-3f, 0x1.11115cp-7f, -0x1.a0403ap-13f, 0x1.75dc26p-19f;