1
Fork 0
mirror of https://github.com/RGBCube/serenity synced 2025-07-24 21:57:35 +00:00

LibM: Implement various functions.

Path from Anonymous.
This commit is contained in:
Andreas Kling 2019-10-17 09:32:24 +02:00
parent e51c1aa3d4
commit c13be2c7ea
3 changed files with 147 additions and 53 deletions

View file

@ -14,6 +14,24 @@ TEST_CASE(trig)
EXPECT_CLOSE(cos(-1), 0.594715); EXPECT_CLOSE(cos(-1), 0.594715);
EXPECT_CLOSE(tan(-1), -1.459446); EXPECT_CLOSE(tan(-1), -1.459446);
EXPECT(isnan(sqrt(-1))); EXPECT(isnan(sqrt(-1)));
EXPECT(isnan(asin(1.1)));
EXPECT(isnan(asin(-1.1)));
EXPECT_CLOSE(asin(0), 0.0);
EXPECT_CLOSE(asin(0.01), 0.01);
EXPECT_CLOSE(asin(0.1), 0.100167);
EXPECT_CLOSE(asin(0.3), 0.304693);
EXPECT_CLOSE(asin(0.499), 0.522444);
EXPECT_CLOSE(asin(0.5), 0.523599);
EXPECT_CLOSE(asin(0.501), 0.524754);
EXPECT_CLOSE(asin(0.9), 1.119770);
EXPECT_CLOSE(asin(0.99), 1.429246);
EXPECT_CLOSE(asin(1.0), 1.570750);
EXPECT_CLOSE(atan(0), 0.0)
EXPECT_CLOSE(atan(0.5), 0.463648)
EXPECT_CLOSE(atan(-0.5), -0.463648)
EXPECT_CLOSE(atan(5.5), 1.390943)
EXPECT_CLOSE(atan(-5.5), -1.390943)
EXPECT_CLOSE(atan(555.5), 1.568996)
} }
TEST_CASE(other) TEST_CASE(other)
@ -51,4 +69,17 @@ TEST_CASE(exponents)
EXPECT_EQ(exp(1000), std::numeric_limits<double>::infinity()); EXPECT_EQ(exp(1000), std::numeric_limits<double>::infinity());
} }
TEST_CASE(logarithms)
{
EXPECT(isnan(log(-1)));
EXPECT(log(0) < -1000000);
EXPECT_CLOSE(log(0.5), -0.693233)
EXPECT_CLOSE(log(1.1), 0.095310)
EXPECT_CLOSE(log(5), 1.609480)
EXPECT_CLOSE(log(5.5), 1.704842)
EXPECT_CLOSE(log(500), 6.214104)
EXPECT_CLOSE(log2(5), 2.321989)
EXPECT_CLOSE(log10(5), 0.698988)
}
TEST_MAIN(Math) TEST_MAIN(Math)

View file

@ -3,13 +3,33 @@
#include <stdint.h> #include <stdint.h>
#include <stdlib.h> #include <stdlib.h>
template<size_t> constexpr double e_to_power(); template<size_t>
template<> constexpr double e_to_power<0>() { return 1; } constexpr double e_to_power();
template<size_t exponent> constexpr double e_to_power() { return M_E * e_to_power<exponent - 1>(); } template<>
constexpr double e_to_power<0>() { return 1; }
template<size_t exponent>
constexpr double e_to_power() { return M_E * e_to_power<exponent - 1>(); }
template<size_t> constexpr size_t factorial(); template<size_t>
template<> constexpr size_t factorial<0>() { return 1; } constexpr size_t factorial();
template<size_t value> constexpr size_t factorial() { return value * factorial<value - 1>(); } template<>
constexpr size_t factorial<0>() { return 1; }
template<size_t value>
constexpr size_t factorial() { return value * factorial<value - 1>(); }
template<size_t>
constexpr size_t product_even();
template<>
constexpr size_t product_even<2>() { return 2; }
template<size_t value>
constexpr size_t product_even() { return value * product_even<value - 2>(); }
template<size_t>
constexpr size_t product_odd();
template<>
constexpr size_t product_odd<1>() { return 1; }
template<size_t value>
constexpr size_t product_odd() { return value * product_odd<value - 2>(); }
extern "C" { extern "C" {
double trunc(double x) double trunc(double x)
@ -67,7 +87,7 @@ double tanh(double x)
return (exponentiated - 1) / (exponentiated + 1); return (exponentiated - 1) / (exponentiated + 1);
} }
double plusX = exp(x); double plusX = exp(x);
double minusX = exp(-x); double minusX = 1 / plusX;
return (plusX - minusX) / (plusX + minusX); return (plusX - minusX) / (plusX + minusX);
} }
@ -79,29 +99,38 @@ double tan(double angle)
double sqrt(double x) double sqrt(double x)
{ {
double res; double res;
__asm__("fsqrt" : "=t"(res) : "0"(x)); __asm__("fsqrt"
: "=t"(res)
: "0"(x));
return res; return res;
} }
double sinh(double x) double sinh(double x)
{ {
if (x > 0) { double exponentiated = exp(x);
double exponentiated = exp(x); if (x > 0)
return (exponentiated * exponentiated - 1) / 2 / exponentiated; return (exponentiated * exponentiated - 1) / 2 / exponentiated;
} return (exponentiated - 1 / exponentiated) / 2;
return (exp(x) - exp(-x)) / 2;
} }
double log10(double) double log10(double x)
{ {
ASSERT_NOT_REACHED(); return log(x) / M_LN10;
return 0;
} }
double log(double) double log(double x)
{ {
ASSERT_NOT_REACHED(); if (x < 0)
return 0; return __builtin_nan("");
if (x == 0)
return -__builtin_huge_val();
double y = 1 + 2 * (x - 1) / (x + 1);
double exponentiated = exp(y);
y = y + 2 * (x - exponentiated) / (x + exponentiated);
exponentiated = exp(y);
y = y + 2 * (x - exponentiated) / (x + exponentiated);
exponentiated = exp(y);
return y + 2 * (x - exponentiated) / (x + exponentiated);
} }
double fmod(double index, double period) double fmod(double index, double period)
@ -114,14 +143,21 @@ double exp(double exponent)
double result = 1; double result = 1;
if (exponent >= 1) { if (exponent >= 1) {
size_t integer_part = (size_t)exponent; size_t integer_part = (size_t)exponent;
if (integer_part & 1) result *= e_to_power<1>(); if (integer_part & 1)
if (integer_part & 2) result *= e_to_power<2>(); result *= e_to_power<1>();
if (integer_part & 2)
result *= e_to_power<2>();
if (integer_part > 3) { if (integer_part > 3) {
if (integer_part & 4) result *= e_to_power<4>(); if (integer_part & 4)
if (integer_part & 8) result *= e_to_power<8>(); result *= e_to_power<4>();
if (integer_part & 16) result *= e_to_power<16>(); if (integer_part & 8)
if (integer_part & 32) result *= e_to_power<32>(); result *= e_to_power<8>();
if (integer_part >= 64) return __builtin_huge_val(); if (integer_part & 16)
result *= e_to_power<16>();
if (integer_part & 32)
result *= e_to_power<32>();
if (integer_part >= 64)
return __builtin_huge_val();
} }
exponent -= integer_part; exponent -= integer_part;
} else if (exponent < 0) } else if (exponent < 0)
@ -140,35 +176,64 @@ double exp(double exponent)
double cosh(double x) double cosh(double x)
{ {
if (x < 0) { double exponentiated = exp(-x);
double exponentiated = exp(-x); if (x < 0)
return (1 + exponentiated * exponentiated) / 2 / exponentiated; return (1 + exponentiated * exponentiated) / 2 / exponentiated;
return (1 / exponentiated + exponentiated) / 2;
}
double atan2(double y, double x)
{
if (x > 0)
return atan(y / x);
if (x == 0) {
if (y > 0)
return M_PI_2;
if (y < 0)
return -M_PI_2;
return 0;
} }
return (exp(x) + exp(-x)) / 2; if (y >= 0)
return atan(y / x) + M_PI;
return atan(y / x) - M_PI;
} }
double atan2(double, double) double atan(double x)
{ {
ASSERT_NOT_REACHED(); if (x < 0)
return 0; return -atan(-x);
if (x > 1)
return M_PI_2 - atan(1 / x);
double squared = x * x;
return x / (1 + 1 * 1 * squared / (3 + 2 * 2 * squared / (5 + 3 * 3 * squared / (7 + 4 * 4 * squared / (9 + 5 * 5 * squared / (11 + 6 * 6 * squared / (13 + 7 * 7 * squared)))))));
} }
double atan(double) double asin(double x)
{ {
ASSERT_NOT_REACHED(); if (x > 1 || x < -1)
return 0; return __builtin_nan("");
if (x > 0.5 || x < -0.5)
return 2 * atan(x / (1 + sqrt(1 - x * x)));
double squared = x * x;
double value = x;
double i = x * squared;
value += i * product_odd<1>() / product_even<2>() / 3;
i *= squared;
value += i * product_odd<3>() / product_even<4>() / 5;
i *= squared;
value += i * product_odd<5>() / product_even<6>() / 7;
i *= squared;
value += i * product_odd<7>() / product_even<8>() / 9;
i *= squared;
value += i * product_odd<9>() / product_even<10>() / 11;
i *= squared;
value += i * product_odd<11>() / product_even<12>() / 13;
return value;
} }
double asin(double) double acos(double x)
{ {
ASSERT_NOT_REACHED(); return M_PI_2 - asin(x);
return 0;
}
double acos(double)
{
ASSERT_NOT_REACHED();
return 0;
} }
double fabs(double value) double fabs(double value)
@ -176,22 +241,19 @@ double fabs(double value)
return value < 0 ? -value : value; return value < 0 ? -value : value;
} }
double log2(double) double log2(double x)
{ {
ASSERT_NOT_REACHED(); return log(x) / M_LN2;
return 0;
} }
float log2f(float) float log2f(float x)
{ {
ASSERT_NOT_REACHED(); return log2(x);
return 0;
} }
long double log2l(long double) long double log2l(long double x)
{ {
ASSERT_NOT_REACHED(); return log2(x);
return 0;
} }
double frexp(double, int*) double frexp(double, int*)
@ -211,5 +273,4 @@ long double frexpl(long double, int*)
ASSERT_NOT_REACHED(); ASSERT_NOT_REACHED();
return 0; return 0;
} }
} }

View file

@ -9,6 +9,8 @@ __BEGIN_DECLS
#define M_PI 3.141592653589793 #define M_PI 3.141592653589793
#define M_PI_2 (M_PI / 2) #define M_PI_2 (M_PI / 2)
#define M_TAU (M_PI * 2) #define M_TAU (M_PI * 2)
#define M_LN2 0.69314718055995
#define M_LN10 2.30258509299405
double acos(double); double acos(double);
float acosf(float); float acosf(float);