flux-ffi 0.1.0

use std::os::raw::c_int;

/// Eisenstein integer norm: |a + bω|² = a² - ab + b²
#[no_mangle]
pub extern "C" fn flux_eisenstein_norm(a: c_int, b: c_int) -> i64 {
    let a = a as i64;
    let b = b as i64;
    a * a - a * b + b * b
}

/// Laman graph edges: 2n - 3 for n >= 2
#[no_mangle]
pub extern "C" fn flux_laman_edges(vertices: c_int) -> c_int {
    if vertices < 2 { return 0; }
    2 * vertices - 3
}

/// Check if (vertices, edges) satisfies Laman rigidity condition
#[no_mangle]
pub extern "C" fn flux_is_rigid(vertices: c_int, edges: c_int) -> bool {
    if vertices < 2 {
        return edges == 0;
    }
    edges == 2 * vertices - 3
}

/// Holonomy check: sum of transforms should be identity (zero for our encoding)
#[no_mangle]
pub extern "C" fn flux_holonomy_check(transforms: *const i64, len: usize) -> bool {
    if transforms.is_null() || len == 0 {
        return true;
    }
    let slice = unsafe { std::slice::from_raw_parts(transforms, len) };
    // Holonomy trivial if sum is zero
    slice.iter().sum::<i64>() == 0
}

/// Manhattan (L1) distance between two integer vectors
#[no_mangle]
pub extern "C" fn flux_manhattan_distance(a: *const c_int, b: *const c_int, len: usize) -> i64 {
    if a.is_null() || b.is_null() || len == 0 {
        return 0;
    }
    let sa = unsafe { std::slice::from_raw_parts(a, len) };
    let sb = unsafe { std::slice::from_raw_parts(b, len) };
    sa.iter().zip(sb.iter())
        .map(|(&x, &y)| (x as i64 - y as i64).abs())
        .sum()
}

/// Pythagorean 48-cell encoding: maps (x,y) to a lattice index
#[no_mangle]
pub extern "C" fn flux_pythagorean48_encode(x: f64, y: f64) -> c_int {
    // Quantize to integer lattice, compute sector and radius
    let ix = x.round() as i32;
    let iy = y.round() as i32;
    let radius = (ix * ix + iy * iy) as i32;
    if radius == 0 {
        return 0;
    }
    // 48-fold symmetry: atan2 maps to 0..47 sectors
    let angle = (iy as f64).atan2(ix as f64);
    let sector = ((angle + std::f64::consts::PI) / (2.0 * std::f64::consts::PI / 48.0)).round() as i32;
    radius * 48 + sector
}

/// Count how many constraints are violated: values[i] > bounds[i]
#[no_mangle]
pub extern "C" fn flux_constraint_check(values: *const f64, bounds: *const f64, len: usize) -> c_int {
    if values.is_null() || bounds.is_null() || len == 0 {
        return 0;
    }
    let sv = unsafe { std::slice::from_raw_parts(values, len) };
    let sb = unsafe { std::slice::from_raw_parts(bounds, len) };
    sv.iter().zip(sb.iter())
        .filter(|(&v, &b)| v > b)
        .count() as c_int
}

/// Simple linear interpolation along control points at parameter t ∈ [0,1]
#[no_mangle]
pub extern "C" fn flux_spline_interpolate(control_points: *const f64, t: f64, n_points: usize) -> f64 {
    if control_points.is_null() || n_points == 0 {
        return 0.0;
    }
    if n_points == 1 {
        return unsafe { *control_points };
    }
    let slice = unsafe { std::slice::from_raw_parts(control_points, n_points) };
    let t_clamped = t.clamp(0.0, 1.0);
    let idx = t_clamped * (n_points - 1) as f64;
    let lo = idx.floor() as usize;
    let hi = (lo + 1).min(n_points - 1);
    let frac = idx - lo as f64;
    slice[lo] * (1.0 - frac) + slice[hi] * frac
}

/// Deadband filter: return value if |value - baseline| > threshold, else baseline
#[no_mangle]
pub extern "C" fn flux_deadband_filter(value: f64, baseline: f64, threshold: f64) -> f64 {
    if (value - baseline).abs() > threshold {
        value
    } else {
        baseline
    }
}

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

    // --- Eisenstein norm ---
    #[test]
    fn test_eisenstein_norm_zero() {
        assert_eq!(flux_eisenstein_norm(0, 0), 0);
    }

    #[test]
    fn test_eisenstein_norm_unit() {
        assert_eq!(flux_eisenstein_norm(1, 0), 1);
        assert_eq!(flux_eisenstein_norm(0, 1), 1);
        assert_eq!(flux_eisenstein_norm(1, 1), 1);
    }

    #[test]
    fn test_eisenstein_norm_general() {
        // a=2, b=3: 4 - 6 + 9 = 7
        assert_eq!(flux_eisenstein_norm(2, 3), 7);
        // a=-1, b=2: 1 + 2 + 4 = 7
        assert_eq!(flux_eisenstein_norm(-1, 2), 7);
    }

    // --- Laman edges ---
    #[test]
    fn test_laman_edges_small() {
        assert_eq!(flux_laman_edges(0), 0);
        assert_eq!(flux_laman_edges(1), 0);
        assert_eq!(flux_laman_edges(2), 1);
        assert_eq!(flux_laman_edges(3), 3);
        assert_eq!(flux_laman_edges(4), 5);
    }

    // --- Is rigid ---
    #[test]
    fn test_is_rigid_laman() {
        assert!(flux_is_rigid(4, 5));   // 2*4-3 = 5
        assert!(flux_is_rigid(2, 1));    // 2*2-3 = 1
        assert!(!flux_is_rigid(4, 4));   // under-constrained
        assert!(!flux_is_rigid(4, 6));   // over-constrained
    }

    #[test]
    fn test_is_rigid_trivial() {
        assert!(flux_is_rigid(0, 0));
        assert!(flux_is_rigid(1, 0));
        assert!(!flux_is_rigid(1, 1));
    }

    // --- Holonomy check ---
    #[test]
    fn test_holonomy_trivial() {
        assert!(flux_holonomy_check(std::ptr::null(), 0));
    }

    #[test]
    fn test_holonomy_zero_sum() {
        let transforms: [i64; 4] = [1, 2, -3, 0];
        assert!(flux_holonomy_check(transforms.as_ptr(), 4));
    }

    #[test]
    fn test_holonomy_nonzero() {
        let transforms: [i64; 3] = [1, 2, 3];
        assert!(!flux_holonomy_check(transforms.as_ptr(), 3));
    }

    // --- Manhattan distance ---
    #[test]
    fn test_manhattan_identity() {
        let a: [i32; 3] = [1, 2, 3];
        let b: [i32; 3] = [1, 2, 3];
        assert_eq!(flux_manhattan_distance(a.as_ptr(), b.as_ptr(), 3), 0);
    }

    #[test]
    fn test_manhattan_basic() {
        let a: [i32; 3] = [0, 0, 0];
        let b: [i32; 3] = [1, 2, 3];
        assert_eq!(flux_manhattan_distance(a.as_ptr(), b.as_ptr(), 3), 6);
    }

    // --- Pythagorean 48 ---
    #[test]
    fn test_pythagorean48_origin() {
        assert_eq!(flux_pythagorean48_encode(0.0, 0.0), 0);
    }

    #[test]
    fn test_pythagorean48_unit_x() {
        let idx = flux_pythagorean48_encode(1.0, 0.0);
        // radius=1, sector at angle 0 => PI offset => sector ~24
        assert!(idx > 0);
    }

    // --- Constraint check ---
    #[test]
    fn test_constraint_check_none_violated() {
        let values: [f64; 3] = [1.0, 2.0, 3.0];
        let bounds: [f64; 3] = [5.0, 5.0, 5.0];
        assert_eq!(flux_constraint_check(values.as_ptr(), bounds.as_ptr(), 3), 0);
    }

    #[test]
    fn test_constraint_check_some_violated() {
        let values: [f64; 4] = [1.0, 6.0, 3.0, 10.0];
        let bounds: [f64; 4] = [5.0, 5.0, 5.0, 5.0];
        assert_eq!(flux_constraint_check(values.as_ptr(), bounds.as_ptr(), 4), 2);
    }

    // --- Spline interpolate ---
    #[test]
    fn test_spline_single_point() {
        let pts: [f64; 1] = [42.0];
        assert!((flux_spline_interpolate(pts.as_ptr(), 0.5, 1) - 42.0).abs() < 1e-10);
    }

    #[test]
    fn test_spline_two_points() {
        let pts: [f64; 2] = [0.0, 10.0];
        assert!((flux_spline_interpolate(pts.as_ptr(), 0.0, 2) - 0.0).abs() < 1e-10);
        assert!((flux_spline_interpolate(pts.as_ptr(), 1.0, 2) - 10.0).abs() < 1e-10);
        assert!((flux_spline_interpolate(pts.as_ptr(), 0.5, 2) - 5.0).abs() < 1e-10);
    }

    #[test]
    fn test_spline_three_points() {
        let pts: [f64; 3] = [0.0, 10.0, 20.0];
        assert!((flux_spline_interpolate(pts.as_ptr(), 0.25, 3) - 5.0).abs() < 1e-10);
    }

    // --- Deadband filter ---
    #[test]
    fn test_deadband_within_threshold() {
        assert!((flux_deadband_filter(10.0, 10.0, 0.5) - 10.0).abs() < 1e-10);
        assert!((flux_deadband_filter(10.3, 10.0, 0.5) - 10.0).abs() < 1e-10);
    }

    #[test]
    fn test_deadband_outside_threshold() {
        assert!((flux_deadband_filter(11.0, 10.0, 0.5) - 11.0).abs() < 1e-10);
        assert!((flux_deadband_filter(8.0, 10.0, 1.0) - 8.0).abs() < 1e-10);
    }

    #[test]
    fn test_deadband_negative() {
        assert!((flux_deadband_filter(-1.0, 0.0, 0.5) - (-1.0)).abs() < 1e-10);
    }
}