From c0639c609a78e9208b83a912e32f8dd0e8638e21 Mon Sep 17 00:00:00 2001 From: Kat Richey Date: Thu, 27 Oct 2022 08:57:08 -0500 Subject: [PATCH] Wow, maths functions! --- math/log.c | 21 +++++++++++++++++++++ math/pow.c | 5 +++++ math/sqrt.c | 26 ++++++++++++++++++++++++++ 3 files changed, 52 insertions(+) create mode 100644 math/log.c create mode 100644 math/pow.c create mode 100644 math/sqrt.c diff --git a/math/log.c b/math/log.c new file mode 100644 index 0000000..a0a0c15 --- /dev/null +++ b/math/log.c @@ -0,0 +1,21 @@ +#include + +double _agm(double x, double y) { + double a = x; + double g = y; + double temp; + + for(int i = 0; i < 15; i++) { + temp = 0.5 * (a + g); + g = sqrt(a * g); + a = temp; + } + + return a; +} + +double log(double x) { + double s = x * exp2(24); + double nat_log = (M_PI / (2 * _agm(1, 4/s))) - (24 * M_LN2); + return nat_log; +} \ No newline at end of file diff --git a/math/pow.c b/math/pow.c new file mode 100644 index 0000000..9c88a26 --- /dev/null +++ b/math/pow.c @@ -0,0 +1,5 @@ +#include + +double pow(double x, double y) { + return exp(y * log(x)); +} \ No newline at end of file diff --git a/math/sqrt.c b/math/sqrt.c new file mode 100644 index 0000000..9d2b077 --- /dev/null +++ b/math/sqrt.c @@ -0,0 +1,26 @@ +#include + +typedef union { + double d; + struct { + unsigned long int mantissa : 52; + unsigned int exponent : 11; + unsigned int sign_bit : 1; + } parts; +} _dbl_cast; + +double sqrt(double x) { + _dbl_cast S = { .d = x }; + double estimate = (0.5 + 0.5 * S.parts.mantissa) * + exp2(S.parts.exponent - 1023); + + double a_n, b_n; + + for(int i = 0; i < 15; i++) { + a_n = (x - (estimate * estimate)) / (2 * estimate); + b_n = estimate + a_n; + estimate = b_n - ((a_n * a_n) / (2 * b_n)); + } + + return estimate; +}