runmat-runtime 0.6.0

Core runtime for RunMat with builtins, BLAS/LAPACK integration, and execution APIs
Documentation
const MIN_ALWAYS_INTEGRAL: f64 = 4_503_599_627_370_496.0;

fn half_turn_mod4(value: f64) -> Option<i32> {
    if !value.is_finite() {
        return None;
    }
    let doubled = value * 2.0;
    if !doubled.is_finite() {
        return (value.abs() >= MIN_ALWAYS_INTEGRAL).then_some(0);
    }
    let rounded = doubled.round();
    if doubled != rounded {
        return None;
    }
    match rounded.rem_euclid(4.0) {
        0.0 => Some(0),
        1.0 => Some(1),
        2.0 => Some(2),
        3.0 => Some(3),
        _ => None,
    }
}

pub(crate) fn sinpi_real(value: f64) -> f64 {
    if !value.is_finite() {
        return f64::NAN;
    }
    match half_turn_mod4(value) {
        Some(0 | 2) => 0.0,
        Some(1) => 1.0,
        Some(3) => -1.0,
        _ => (std::f64::consts::PI * value).sin(),
    }
}

pub(crate) fn cospi_real(value: f64) -> f64 {
    if !value.is_finite() {
        return f64::NAN;
    }
    match half_turn_mod4(value) {
        Some(0) => 1.0,
        Some(1 | 3) => 0.0,
        Some(2) => -1.0,
        _ => (std::f64::consts::PI * value).cos(),
    }
}

fn mul_preserving_exact_zero(factor: f64, magnitude: f64) -> f64 {
    if factor == 0.0 && !magnitude.is_nan() {
        factor
    } else {
        factor * magnitude
    }
}

pub(crate) fn sinpi_complex(re: f64, im: f64) -> (f64, f64) {
    let scaled_im = std::f64::consts::PI * im;
    let sin_re = sinpi_real(re);
    let cos_re = cospi_real(re);
    (
        mul_preserving_exact_zero(sin_re, scaled_im.cosh()),
        mul_preserving_exact_zero(cos_re, scaled_im.sinh()),
    )
}

pub(crate) fn cospi_complex(re: f64, im: f64) -> (f64, f64) {
    let scaled_im = std::f64::consts::PI * im;
    let cos_re = cospi_real(re);
    let neg_sin_re = -sinpi_real(re);
    (
        mul_preserving_exact_zero(cos_re, scaled_im.cosh()),
        mul_preserving_exact_zero(neg_sin_re, scaled_im.sinh()),
    )
}

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

    #[cfg_attr(target_arch = "wasm32", wasm_bindgen_test::wasm_bindgen_test)]
    #[cfg_attr(not(target_arch = "wasm32"), test)]
    fn sinpi_exact_half_turns() {
        assert_eq!(sinpi_real(0.0), 0.0);
        assert_eq!(sinpi_real(0.5), 1.0);
        assert_eq!(sinpi_real(1.0), 0.0);
        assert_eq!(sinpi_real(1.5), -1.0);
        assert_eq!(sinpi_real(-0.5), -1.0);
        assert_eq!(sinpi_real(-1.0), 0.0);
        assert_eq!(sinpi_real(9_007_199_254_740_991.0), 0.0);
        assert_eq!(sinpi_real(9_007_199_254_740_992.0), 0.0);
        assert_eq!(sinpi_real(1.0e300), 0.0);
    }

    #[cfg_attr(target_arch = "wasm32", wasm_bindgen_test::wasm_bindgen_test)]
    #[cfg_attr(not(target_arch = "wasm32"), test)]
    fn cospi_exact_half_turns() {
        assert_eq!(cospi_real(0.0), 1.0);
        assert_eq!(cospi_real(0.5), 0.0);
        assert_eq!(cospi_real(1.0), -1.0);
        assert_eq!(cospi_real(1.5), 0.0);
        assert_eq!(cospi_real(-0.5), 0.0);
        assert_eq!(cospi_real(-1.0), -1.0);
        assert_eq!(cospi_real(9_007_199_254_740_991.0), -1.0);
        assert_eq!(cospi_real(9_007_199_254_740_992.0), 1.0);
        assert_eq!(cospi_real(1.0e300), 1.0);
    }

    #[cfg_attr(target_arch = "wasm32", wasm_bindgen_test::wasm_bindgen_test)]
    #[cfg_attr(not(target_arch = "wasm32"), test)]
    fn complex_exact_zero_factors_survive_overflowing_imaginary_scale() {
        let (re, im) = sinpi_complex(0.0, f64::INFINITY);
        assert_eq!(re, 0.0);
        assert!(im.is_infinite() && im.is_sign_positive());

        let (re, im) = sinpi_complex(0.5, f64::INFINITY);
        assert!(re.is_infinite() && re.is_sign_positive());
        assert_eq!(im, 0.0);

        let (re, im) = cospi_complex(0.5, f64::INFINITY);
        assert_eq!(re, 0.0);
        assert!(im.is_infinite() && im.is_sign_negative());

        let (re, im) = cospi_complex(1.0, 1.0e300);
        assert!(re.is_infinite() && re.is_sign_negative());
        assert_eq!(im, -0.0);
    }

    #[cfg_attr(target_arch = "wasm32", wasm_bindgen_test::wasm_bindgen_test)]
    #[cfg_attr(not(target_arch = "wasm32"), test)]
    fn real_nonfinite_inputs_return_nan() {
        assert!(sinpi_real(f64::NAN).is_nan());
        assert!(sinpi_real(f64::INFINITY).is_nan());
        assert!(sinpi_real(f64::NEG_INFINITY).is_nan());
        assert!(cospi_real(f64::NAN).is_nan());
        assert!(cospi_real(f64::INFINITY).is_nan());
        assert!(cospi_real(f64::NEG_INFINITY).is_nan());
    }
}