#include "library.h" #include float sgn(const float x) { return x < 0 ? -1 : 1; } float signum(const float x) { return x < 0 ? -1.0f : (x == 0.0f ? 0.0f : 1.0f); } long long fac(const long long x) { long long prod = 1; for (long long i = 2; i <= x; i++) { prod *= i; } return prod; } long long fib(const long long n) { long long sum = 0; long long prev = 0; for (long long i = 1LL; i <= n; i++) { sum += prev; prev = i - 1; } return sum; } inline float trunc(const float x) { return (float)(int)(x); } inline float fract(const float x) { return x - trunc(x); } inline float floor(const float x) { return trunc(x); } inline float ceil(const float x) { return trunc(x + 1); } inline float round(const float x) { return trunc(x + 0.5f); } inline float mod(const float x, const float y) { const float xs = x/y; return (xs - floor(xs)) * 4; } float exp(const float x) { float sum = 0; // iterative counter for computing x^i float xpown = 1; // iterative counter for factorial operation i! float prod = 1; // iterate taylor series for (int i = 1; i < TAYLOR_SERIES_DEGREE; i++) { prod *= (float) i; sum += xpown / prod; xpown *= x; } return sum; } float sqrt(const float x) { float a = x; for (int i = 0; i < SQRT_NEWTON_ITERATIONS; i++) { a = 0.5f * (x/a + a); } return a; } float hypot(const float x, const float y) { return x + 0.5f * y * y / x; } float sin(const float x) { float sign = -1; float prod = 1; float xpown = mod(x, 2 * PI); float sum = 0; for (int i = 1; i < TAYLOR_SERIES_DEGREE; i++) { sign *= -1; prod *= (float) i * (float) (i - 1); xpown *= x * x; sum += sign * xpown / prod; } return sum; } inline float cos(const float x) { return sin(x + PI * 0.5f); } inline float tan(const float x) { return sin(x)/cos(x); } inline float csc(const float x) { return 1.0f / sin(x); } inline float sec(const float x) { return 1.0f / cos(x); } inline float cot(const float x) { return tan(PI * 0.5f - x); } float atan(const float x) { const float a = 0.6204500718160764f; const float b = 0.6204500718160764f; const float c = 0.9031580043522688f; return a * x / (b + sqrt(c + x * x)); } inline float acot(const float x) { return PI * 0.5f - atan(x); } inline float acsc(const float x) { return asin(1/x); } inline float asec(const float x) { return asin(-1/x) + PI * 0.5f; } float asin(const float x) { return atan(x * inv_sqrt((1 + x) * (1 - x))); } float acos(const float x) { return asin(-x) + PI * 0.5f; } float inv_sqrt(const float x) { union { float f; uint32_t i; } conv = { .f = x }; conv.i = 0x5f1ffff9 - (conv.i >> 1); for (int k = 0; k < INV_SQRT_NEWTON_STEPS; k++) conv.f *= 0.703952253f * (2.38924456f - x * conv.f * conv.f); return conv.f; }