lib/math: math::floats fold-2a frexpf64 decompose

lib/math/floats.ww

@@ -1,15 +1,11 @@
-// floats — f64 classification, sign, bit-reinterpret core, and the
-// subnormal-normalize step of the f64 decompose half. Ported from
+// floats — f64 classification, sign, bit-reinterpret core, and the f64
+// decompose half (subnormal-normalize + frexp). Ported from
 // ref/hare/math/floats.ha (fold-1: classify/sign/bits; fold-2a:
-// issubnormalf64/normalizef64). frexpf64 (floats.ha:278) is held back: its
-// Hare-exact zero guard `n == 0f64` miscompiles — a no-decimal `0f64`
-// literal in an f64 comparison is materialized into a GPR and never moved
-// to XMM, so the compare reads a stale operand (both stages identically,
-// so the byte-id gates are blind to it). `0.0` would dodge it, but that is
-// a workaround (CLAUDE.md rule 7); frexpf64 lands once the cgen bug is
-// fixed. f32 variants, NAN/INF + magnitude consts, and the
-// ldexp/modfrac/nextafter family (need f64 DIVIDE + the INF const) are
-// deferred to a later fold.
+// issubnormalf64/normalizef64/frexpf64). frexpf64's zero guard `n == 0f64`
+// rides the #103 fix (no-decimal f64 literal now materialized into XMM)
+// and its (f64, i64) tuple return rides the #105 fix (tuple f64-word read).
+// f32 variants stay deferred (#104 blocks the f32 call-arg narrowing), as
+// do the ldexp/modfrac/nextafter family (need f64 DIVIDE + the INF const).
 
 package math;
 
@@ -56,6 +52,18 @@ export def F64_EXPONENT_MASK: u64 = (1 << F64_EXPONENT_BITS) - 1;
 // ref/hare/math/floats.ha:75
 def F64_SIGN_MASK: u64 = 1u64 << 63;
 
+// Mask that clears an f64's exponent field, keeping sign + mantissa.
+// ref/hare/math/floats.ha:77. Hare hardcodes the 0x800FFFFFFFFFFFFF binary
+// literal because its lexer can't const-fold the expression; ww's #88
+// def-const-fold can, so the readable form is kept. floats.ha:79's NOTE
+// expression has an `0u64 &` upstream typo (it would yield 0); the value it
+// documents is exactly ~(F64_EXPONENT_MASK << F64_MANTISSA_BITS).
+def F64_EXP_REMOVAL_MASK: u64 = ~(F64_EXPONENT_MASK << F64_MANTISSA_BITS);
+
+// The f64 bit pattern whose exponent field evaluates to zero (0.5 scale).
+// ref/hare/math/floats.ha:84
+def F64_EXP_ZERO: u64 = (F64_EXPONENT_BIAS - 1) << F64_MANTISSA_BITS;
+
 // Returns true if the given floating-point number is NaN.
 // ref/hare/math/floats.ha:144 (Hare's expression body inlined into a
 // block: ww has no expression-bodied fn form, only brace blocks).
@@ -131,3 +139,22 @@ export fn normalizef64(n: f64) (f64, i64) = {
 	};
 	return (n, 0);
 };
+
+// Breaks a f64 down into its mantissa and exponent. The mantissa will be
+// between 0.5 and 1.
+// ref/hare/math/floats.ha:278
+export fn frexpf64(n: f64) (f64, i64) = {
+	if (isnan(n) || isinf(n) || n == 0f64) {
+		return (n, 0);
+	};
+	const normalized = normalizef64(n);
+	const normal_float = normalized.0;
+	const normalization_exp = normalized.1;
+	const bits = f64bits(normal_float);
+	const raw_exp: u64 = (bits >> F64_MANTISSA_BITS) & F64_EXPONENT_MASK;
+	const exp: i64 = normalization_exp +
+		(raw_exp: i64) - (F64_EXPONENT_BIAS: i64) + 1;
+	const mantissa: f64 =
+		f64frombits((bits & F64_EXP_REMOVAL_MASK) | F64_EXP_ZERO);
+	return (mantissa, exp);
+};