diff --git a/library.c b/library.c index 8f768dd..3d5a30e 100644 --- a/library.c +++ b/library.c @@ -1,9 +1,14 @@ #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++) { @@ -112,7 +117,7 @@ inline float cot(const float x) { return tan(PI * 0.5f - x); } -float arctan(const float x) { +float atan(const float x) { const float a = 0.6204500718160764f; const float b = 0.6204500718160764f; const float c = 0.9031580043522688f; @@ -120,14 +125,35 @@ float arctan(const float x) { return a * x / (b + sqrt(c + x * x)); } -inline float arccot(const float x) { - return PI * 0.5f - arctan(x); +inline float acot(const float x) { + return PI * 0.5f - atan(x); } -float atan2(const float x, const float y) { - return arctan(y/x); +inline float acsc(const float x) { + return asin(1/x); } -float arcsin(const float 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; } \ No newline at end of file diff --git a/library.h b/library.h index ba4cd14..15fb1ec 100644 --- a/library.h +++ b/library.h @@ -3,6 +3,7 @@ #define TAYLOR_SERIES_DEGREE 34 #define SQRT_NEWTON_ITERATIONS 26 +#define INV_SQRT_NEWTON_STEPS 2 const float PI = 3.1415926535897932384626433f; @@ -14,6 +15,14 @@ const float PI = 3.1415926535897932384626433f; */ float sgn(float x); +/** + * Computes the signum of the parameter x. + * Returns 1 if x > 0, 0 if x equals 0 and -1 if x < 0 + * @param x + * @return + */ +float signum(float x); + /** * Compute the factorial of parameter x. * This function only properly works on natural numbers. @@ -155,8 +164,50 @@ float cot(float x); * @param x * @return */ -float arctan(float x); +float atan(float x); +/** + * Approximates the inverse of the sine of x + * @param x + * @return + */ +float asin(float x); +/** + * Approximates the inverse of the cosine of x + * @param x + * @return + */ +float acos(float x); + +/** + * Approximates the inverse of the cotangent of x + * @param x + * @return + */ +float acot(float x); + +/** + * Approximates the inverse of the cosecant of x + * @param x + * @return + */ +float acsc(float x); + +/** + * Approximates the inverse of the secant of x + * @param x + * @return + */ +float asec(float x); + +/** + * Approximate the inverse square root. + * This implementation is based on the algorithm introduced in the Quake game with + * various improvements. + * @param x + * @return 1/sqrt(x) + */ +float inv_sqrt(float x); #endif //COMMON_MATH_LIBRARY_H