//! `wai.det.fixed64` — Q32.32 signed fixed-point (extension §6.2).
//!
//! The floor numeric profile. Replay-equivalence carries no tolerance
//! parameter, so the floor arithmetic must agree bit-for-bit on every
//! machine that exists; integer arithmetic does, unconditionally —
//! the same reasoning that put PNG (not AVIF) in the core floor, and
//! the simulation analog of `wai.neural.int_synth`.
//!
//! Normative semantics, fixed forever once registered:
//! - representation: `i64`, 32 integer bits, 32 fractional bits
//! - add/sub: saturating
//! - mul/div: computed in 128-bit, rounded **toward negative infinity**
//!   (floor), then saturated to the i64 range
//! - div by zero: saturates (`MAX` for non-negative dividend, else `MIN`)

/// Q32.32 fixed-point value.
#[derive(Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Debug, Default, Hash)]
pub struct Fx(i64);

impl Fx {
    pub const FRAC_BITS: u32 = 32;
    pub const ZERO: Fx = Fx(0);
    pub const ONE: Fx = Fx(1i64 << 32);
    pub const MAX: Fx = Fx(i64::MAX);
    pub const MIN: Fx = Fx(i64::MIN);

    /// Raw bit pattern — the wire form (i64 little-endian everywhere).
    #[inline]
    pub const fn from_bits(bits: i64) -> Fx {
        Fx(bits)
    }

    #[inline]
    pub const fn to_bits(self) -> i64 {
        self.0
    }

    #[inline]
    pub const fn from_int(v: i32) -> Fx {
        Fx((v as i64) << 32)
    }

    /// Test/authoring convenience; never used in the simulation path
    /// (floats stay out of state by construction).
    pub fn approx_f64(self) -> f64 {
        self.0 as f64 / (1u64 << 32) as f64
    }

    /// Deterministic square root, rounded down (floor), per the
    /// profile's rounding rule. `sqrt(v)` in Q32.32 is
    /// `isqrt(raw · 2^32)`: the 128-bit integer square root of the
    /// raw bits shifted up by the fractional width. Negative input
    /// (which a distance computation cannot produce, but a malformed
    /// payload could) is defined as 0 — every input must map to one
    /// answer on every machine.
    pub fn sqrt(self) -> Fx {
        if self.0 <= 0 {
            return Fx::ZERO;
        }
        let n = (self.0 as u128) << Self::FRAC_BITS;
        Fx(isqrt_u128(n) as i64)
    }
}

/// Floor integer square root via Newton's method — converges from
/// above, exits the first time the iterate stops decreasing, which
/// for integer Newton is exactly floor(√n).
fn isqrt_u128(n: u128) -> u128 {
    if n == 0 {
        return 0;
    }
    // Initial guess: a power of two >= sqrt(n).
    let mut x = 1u128 << ((128 - n.leading_zeros()).div_ceil(2));
    loop {
        let y = (x + n / x) >> 1;
        if y >= x {
            return x;
        }
        x = y;
    }
}

#[inline]
fn saturate(wide: i128) -> i64 {
    if wide > i64::MAX as i128 {
        i64::MAX
    } else if wide < i64::MIN as i128 {
        i64::MIN
    } else {
        wide as i64
    }
}

impl std::ops::Add for Fx {
    type Output = Fx;
    #[inline]
    fn add(self, rhs: Fx) -> Fx {
        Fx(self.0.saturating_add(rhs.0))
    }
}
impl std::ops::Sub for Fx {
    type Output = Fx;
    #[inline]
    fn sub(self, rhs: Fx) -> Fx {
        Fx(self.0.saturating_sub(rhs.0))
    }
}
impl std::ops::Mul for Fx {
    type Output = Fx;
    #[inline]
    fn mul(self, rhs: Fx) -> Fx {
        // i128 product, arithmetic shift right 32 = floor for negatives.
        Fx(saturate((self.0 as i128 * rhs.0 as i128) >> Fx::FRAC_BITS))
    }
}
impl std::ops::Div for Fx {
    type Output = Fx;
    #[inline]
    fn div(self, rhs: Fx) -> Fx {
        if rhs.0 == 0 {
            return if self.0 >= 0 { Fx::MAX } else { Fx::MIN };
        }
        // Floor division (toward negative infinity), per the profile.
        Fx(saturate(((self.0 as i128) << Fx::FRAC_BITS).div_euclid(rhs.0 as i128)))
    }
}
impl std::ops::Neg for Fx {
    type Output = Fx;
    #[inline]
    fn neg(self) -> Fx {
        Fx(self.0.checked_neg().unwrap_or(i64::MAX))
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn identities() {
        let three = Fx::from_int(3);
        assert_eq!(three * Fx::ONE, three);
        assert_eq!(three + Fx::ZERO, three);
        assert_eq!(Fx::from_int(6) / Fx::from_int(2), three);
    }

    #[test]
    fn floor_rounding_is_toward_negative_infinity() {
        // -1 / 2 in Q32.32 floors to -0.5 exactly; -1/3 floors below it.
        let minus_one = Fx::from_int(-1);
        let half = Fx::from_int(1) / Fx::from_int(2);
        assert_eq!(minus_one / Fx::from_int(2), Fx::ZERO - half);
        // (-1) * 0.333.. : floor, so result <= exact value.
        let third = Fx::from_int(1) / Fx::from_int(3);
        let prod = minus_one * third;
        assert!(prod.to_bits() <= -(third.to_bits()) , "floor must not round toward zero");
    }

    #[test]
    fn sqrt_floor_semantics() {
        assert_eq!(Fx::from_int(4).sqrt(), Fx::from_int(2));
        assert_eq!(Fx::from_int(9).sqrt(), Fx::from_int(3));
        assert_eq!(Fx::ZERO.sqrt(), Fx::ZERO);
        assert_eq!(Fx::from_int(-4).sqrt(), Fx::ZERO, "negative input is defined as 0");
        // floor: sqrt(2) in Q32.32 must satisfy s*s <= 2 < (s+1ulp)^2
        let s = Fx::from_int(2).sqrt();
        assert!(s * s <= Fx::from_int(2));
        let s1 = Fx::from_bits(s.to_bits() + 1);
        assert!(s1 * s1 >= Fx::from_int(2) || s1 * s1 >= s * s, "monotone above floor");
        assert!((s.approx_f64() - std::f64::consts::SQRT_2).abs() < 1e-9);
    }

    #[test]
    fn saturation() {
        assert_eq!(Fx::MAX + Fx::ONE, Fx::MAX);
        assert_eq!(Fx::MIN - Fx::ONE, Fx::MIN);
        assert_eq!(Fx::MAX * Fx::from_int(2), Fx::MAX);
        assert_eq!(Fx::ONE / Fx::ZERO, Fx::MAX);
        assert_eq!(Fx::from_int(-1) / Fx::ZERO, Fx::MIN);
    }
}
