hare

[hare] The Hare programming language
git clone https://git.torresjrjr.com/hare.git
Log | Files | Refs | README | LICENSE

commit 6ec62fc6657f2010540da02574697c257b14e102
parent 6389bca04359e19eb3f628cb8ade9ed355fbe25c
Author: Lassi Pulkkinen <lassi@pulk.fi>
Date:   Thu,  6 Apr 2023 20:02:19 +0300

math::modfracf*: Return integer part as float

i64 can't represent all integer values that f64 can (and vice versa).
Integer overflow in modfracf* caused floorf64, etc. to behave incorrectly
on large inputs.

Add tests for affected functions.

Signed-off-by: Lassi Pulkkinen <lassi@pulk.fi>

Diffstat:
Mmath/+test/data.ha | 22+++++++++++-----------
Mmath/+test/floats.ha | 30++++++++++++++++++------------
Mmath/+test/math.ha | 6++++++
Mmath/floats.ha | 17++++++++---------
Mmath/math.ha | 26++++++++++----------------
5 files changed, 53 insertions(+), 48 deletions(-)

diff --git a/math/+test/data.ha b/math/+test/data.ha @@ -308,17 +308,17 @@ const TEST_LOG1P: []f64 = [ 1.8088493239630770262045333e-02, -9.0865245631588989681559268e-02, ]; -const TEST_MODFRAC: [_](i64, f64) = [ - (4i64, 9.7901192488367350108546816e-01), - (7i64, 7.3887247457810456552351752e-01), - (-0i64, -2.7688005719200159404635997e-01), - (-5i64, -1.060361827107492160848778e-02), - (9i64, 6.3629370719841737980004837e-01), - (2i64, 9.2637723924396464525443662e-01), - (5i64, 2.2908343145930665230025625e-01), - (2i64, 7.2793991043601025126008608e-01), - (1i64, 8.2530809168085506044576505e-01), - (-8i64, -6.8592476857560136238589621e-01), +const TEST_MODFRAC: [_](f64, f64) = [ + (4f64, 9.7901192488367350108546816e-01), + (7f64, 7.3887247457810456552351752e-01), + (-0f64, -2.7688005719200159404635997e-01), + (-5f64, -1.060361827107492160848778e-02), + (9f64, 6.3629370719841737980004837e-01), + (2f64, 9.2637723924396464525443662e-01), + (5f64, 2.2908343145930665230025625e-01), + (2f64, 7.2793991043601025126008608e-01), + (1f64, 8.2530809168085506044576505e-01), + (-8f64, -6.8592476857560136238589621e-01), ]; const TEST_POW: [_]f64 = [ 9.5282232631648411840742957e+04, diff --git a/math/+test/floats.ha b/math/+test/floats.ha @@ -174,47 +174,53 @@ assert(isclose(res.1, TEST_MODFRAC[idx].1)); }; let res = modfracf64(1.75f64); - assert(res.0 == 1i64); + assert(res.0 == 1f64); assert(res.1 == 0.75f64); res = modfracf64(0.75f64); - assert(res.0 == 0i64); + assert(res.0 == 0f64); assert(res.1 == 0.75f64); res = modfracf64(-0.75f64); - assert(res.0 == -0i64); + assert(res.0 == -0f64); assert(res.1 == -0.75f64); res = modfracf64(0f64); - assert(res.0 == 0i64); + assert(res.0 == 0f64); assert(res.1 == 0f64); assert(signf(res.1) > 0); res = modfracf64(-0f64); - assert(res.0 == -0i64); + assert(res.0 == -0f64); assert(res.1 == -0f64); assert(signf(res.1) < 0); res = modfracf64(23.50f64); - assert(res.0 == 23i64); + assert(res.0 == 23f64); assert(res.1 == 0.50f64); + res = modfracf64(F64_MAX_NORMAL); + assert(res.0 == F64_MAX_NORMAL); + assert(res.1 == 0f64); // 32 let res = modfracf32(1.75f32); - assert(res.0 == 1i32); + assert(res.0 == 1f32); assert(res.1 == 0.75f32); res = modfracf32(0.75f32); - assert(res.0 == 0i32); + assert(res.0 == 0f32); assert(res.1 == 0.75f32); res = modfracf32(-0.75f32); - assert(res.0 == -0i32); + assert(res.0 == -0f32); assert(res.1 == -0.75f32); res = modfracf32(0.0f32); - assert(res.0 == 0i32); + assert(res.0 == 0f32); assert(res.1 == 0.0f32); assert(signf(res.1) > 0); res = modfracf32(-0.0f32); - assert(res.0 == -0i32); + assert(res.0 == -0f32); assert(res.1 == -0.0f32); assert(signf(res.1) < 0); res = modfracf32(23.50f32); - assert(res.0 == 23i32); + assert(res.0 == 23f32); assert(res.1 == 0.50f32); + res = modfracf32(F32_MAX_NORMAL); + assert(res.0 == F32_MAX_NORMAL); + assert(res.1 == 0f32); }; @test fn nextafter() void = { diff --git a/math/+test/math.ha b/math/+test/math.ha @@ -209,8 +209,10 @@ TEST_CEIL[idx])); }; assert(ceilf64(-INF) == -INF); + assert(ceilf64(-F64_MAX_NORMAL) == -F64_MAX_NORMAL); assert(ceilf64(-0f64) == -0f64); assert(ceilf64(0f64) == 0f64); + assert(ceilf64(F64_MAX_NORMAL) == F64_MAX_NORMAL); assert(ceilf64(INF) == INF); assert(isnan(ceilf64(NAN))); }; @@ -222,8 +224,10 @@ TEST_TRUNC[idx])); }; assert(truncf64(-INF) == -INF); + assert(truncf64(-F64_MAX_NORMAL) == -F64_MAX_NORMAL); assert(truncf64(-0f64) == -0f64); assert(truncf64(0f64) == 0f64); + assert(truncf64(F64_MAX_NORMAL) == F64_MAX_NORMAL); assert(truncf64(INF) == INF); assert(isnan(truncf64(NAN))); }; @@ -235,8 +239,10 @@ TEST_ROUND[idx])); }; assert(roundf64(-INF) == -INF); + assert(roundf64(-F64_MAX_NORMAL) == -F64_MAX_NORMAL); assert(roundf64(-0f64) == -0f64); assert(roundf64(0f64) == 0f64); + assert(roundf64(F64_MAX_NORMAL) == F64_MAX_NORMAL); assert(roundf64(INF) == INF); assert(isnan(roundf64(NAN))); }; diff --git a/math/floats.ha b/math/floats.ha @@ -478,16 +478,16 @@ export fn ldexpf32(mantissa: f32, exp: i64) f32 = { }; // Returns the integer and fractional parts of an f64. -export fn modfracf64(n: f64) (i64, f64) = { +export fn modfracf64(n: f64) (f64, f64) = { if (n < 1f64) { if (n < 0f64) { let positive_parts = modfracf64(-n); return (-positive_parts.0, -positive_parts.1); }; if (n == 0f64) { - return ((n: i64), n); + return (n, n); }; - return (0i64, n); + return (0f64, n); }; let bits = f64bits(n); const exp = (((bits >> F64_MANTISSA_BITS) & F64_EXPONENT_MASK): i64) - @@ -501,20 +501,20 @@ export fn modfracf64(n: f64) (i64, f64) = { }; const int_part = f64frombits(bits); const frac_part = n - int_part; - return ((int_part: i64), frac_part); + return (int_part, frac_part); }; // Returns the integer and fractional parts of an f32. -export fn modfracf32(n: f32) (i32, f32) = { +export fn modfracf32(n: f32) (f32, f32) = { if (n < 1.0f32) { if (n < 0.0f32) { let positive_parts = modfracf32(-n); return (-positive_parts.0, -positive_parts.1); }; if (n == 0.0f32) { - return ((n: i32), n); + return (n, n); }; - return (0i32, n); + return (0f32, n); }; let bits = f32bits(n); const exp = (((bits >> F32_MANTISSA_BITS) & F32_EXPONENT_MASK): i32) - @@ -528,10 +528,9 @@ export fn modfracf32(n: f32) (i32, f32) = { }; const int_part = f32frombits(bits); const frac_part = n - int_part; - return ((int_part: i32), frac_part); + return (int_part, frac_part); }; - // Returns the f32 that is closest to 'x' in direction of 'y'. Returns NaN // if either parameter is NaN. Returns 'x' if both parameters are same. export fn nextafterf32(x: f32, y: f32) f32 = { diff --git a/math/math.ha b/math/math.ha @@ -680,11 +680,9 @@ export fn sqrtf64(x: f64) f64 = { }; fn is_f64_odd_int(x: f64) bool = { - const x_int_frac = modfracf64(x); - const x_int = x_int_frac.0; - const x_frac = x_int_frac.1; + const (x_int, x_frac) = modfracf64(x); const has_no_frac = (x_frac == 0f64); - const is_odd = ((x_int & 1i64) == 1i64); + const is_odd = ((x_int: i64 & 1i64) == 1i64); return has_no_frac && is_odd; }; @@ -750,7 +748,7 @@ export fn powf64(x: f64, p: f64) f64 = { if (p_frac != 0f64) { if (p_frac > 0.5f64) { p_frac -= 1f64; - p_int += 1i64; + p_int += 1f64; }; res_mantissa = expf64(p_frac * logf64(x)); }; @@ -758,11 +756,11 @@ export fn powf64(x: f64, p: f64) f64 = { // Repeatedly square our number x, for each bit in our power p. // If the current bit is 1 in p, add the respective power of x to our // result. - for (let i = p_int; i != 0; i >>= 1) { + for (let i = p_int: i64; i != 0; i >>= 1) { // Check for over/underflow. if (x_exp <= -1i64 << (F64_EXPONENT_BITS: i64)) { return 0f64; - }; + }; if (x_exp >= 1i64 << (F64_EXPONENT_BITS: i64)) { return INF; }; @@ -794,16 +792,13 @@ export fn floorf64(x: f64) f64 = { return x; }; if (x < 0f64) { - const modfrac_res = modfracf64(-x); - let int_part = modfrac_res.0; - const frac_part = modfrac_res.1; + let (int_part, frac_part) = modfracf64(-x); if (frac_part != 0f64) { - int_part += 1; + int_part += 1f64; }; - return -(int_part: f64); + return -int_part; }; - const modfrac_res = modfracf64(x); - return (modfrac_res.0: f64); + return modfracf64(x).0; }; // Returns the least integer value greater than or equal to x. @@ -814,8 +809,7 @@ export fn truncf64(x: f64) f64 = { if (x == 0f64 || isnan(x) || isinf(x)) { return x; }; - const modfrac_res = modfracf64(x); - return (modfrac_res.0: f64); + return modfracf64(x).0; }; // Returns the nearest integer, rounding half away from zero.