From 45d430b3343b82c8cb60fbd4a51a04913ea4d361 Mon Sep 17 00:00:00 2001 From: EliteAsian Date: Thu, 10 Nov 2022 19:37:03 -0500 Subject: [PATCH] `pow` and `exp` --- Math/Core.scope | 93 ++++++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 89 insertions(+), 4 deletions(-) diff --git a/Math/Core.scope b/Math/Core.scope index fbede9a..67fa4a9 100644 --- a/Math/Core.scope +++ b/Math/Core.scope @@ -61,6 +61,25 @@ func dec powInt(dec base, int exponent) { } } +/% +@approx: This function could have some inaccuracies at higher values. + Use @"powInt" if the @"exponent" is an integer for a more accurate result. +@todo: Add support for negative bases. + +Computes @"base" ^ @"exponent". + +This function uses @"exp" and @"ln" to compute the power. +%/ +func dec pow(dec base, dec exponent) { + // log10(base ^ exponent) = log10(base) * exponent + + if (base < 0.0) { + ret nan; + } + + ret exp(ln(base) * exponent); +} + /% Computes the log (base 10) of @"x". @@ -90,8 +109,8 @@ func dec log10(dec x) { ret infinity; } - // The Taylor series below only works between 0 and 1 - // Transform x to fit in the range with a summand (b) + // The Taylor series below only works between 0 and 1. + // Transform x to fit in the range with a summand (b). // log10(x : 0..inf) = log10(x : 0..=1) + b // Get "b" @@ -122,7 +141,7 @@ func dec log10(dec x) { // Finally, return // log10(x) = ln(x) * log10(e) + b - // log10(e) = 0.434294481903251 + // log10(e) = 0.434294481903251 ret seriesSum * 0.434294481903251 + b; } @@ -149,7 +168,73 @@ func dec log(dec x, dec base) { // TODO: prevent x <= 0.0 ret nan; } - + // log_base(x) = log10(x) / log10(base) ret log10(x) / log10(base); +} + +/% +Computes e ^ @"x". + +This function uses the Taylor series to compute the result. +%/ +func dec exp(dec x) { + // Temp until global constants + dec E = 2.718281828459045; + + // Special cases + + if (x == -infinity) { + ret 0.0; + } + + if (x == 0.0) { + ret 1.0; + } + + if (x == infinity) { + ret infinity; + } + + // Deal with negatives later (see below) + + bool neg = x < 0.0; + if (neg) { + x = -x; + } + + // The Taylor series below is only accurate between -2 and 2 + // Transform x to fit in the range. + // e ^ 8.2 = e^3.2 * e^6.0 + // Note that the exponent in the second power is an integer + // which we can compute already. Let b be the second power. + // This does not work if x is negative. We can deal with that + // by doing `1.0 / result` at the end. + + // Get "b" + + dec b = 1.0; + while (x > 2.0) { + x -= 1.0; + b *= E; + } + + // Taylor series (for seriesSum) + + // exp(x) = (x^0 / 0!) + (x^1 / 1!) + (x^2 / 2!) + ... + dec seriesSum = 1.0; + int factorial = 1; + for (int k : 1..23) { + factorial *= k; + seriesSum += powInt(x, k) / factorial -> dec; + } + + // Finally, return + + dec result = seriesSum * b; + if (neg) { + ret 1.0 / result; + } else { + ret result; + } } \ No newline at end of file