vsrg 0.3.0

use std::f64::consts::{LN_2, PI};

use super::{Easable, Easing, split_and_calc};

pub fn evaluate_simple<T, const N: usize>(ease: Easing, from: T, to: T, t: T) -> T
where
    T: Easable<N>,
{
    split_and_calc::<T, N>(from, to, t, |from, to, t| {
        if t < 0.0 {
            from * t
        } else if t > 1.0 {
            let i1 = integrate_simple(ease, 1.0);
            from * (t - i1) + to * (i1 + t - 1.0)
        } else {
            let i = integrate_simple(ease, t);
            from * (t - i) + to * i
        }
    })
}

#[cfg(any(test, feature = "lut_integration"))]
use crate::math::AdvancedEasing;

#[cfg(any(test, feature = "lut_integration"))]
pub fn evaluate_advanced<T, const N: usize>(ease: AdvancedEasing, from: T, to: T, t: T) -> T
where
    T: Easable<N>,
{
    if let Some(simple_ease) = ease.to_monotonic() {
        evaluate_simple(simple_ease, from, to, t)
    } else {
        split_and_calc::<T, N>(from, to, t, |from, to, t| {
            if t < 0.0 {
                from * t
            } else if t > 1.0 {
                let i1 = integrate_advanced(ease, 1.0);
                from * (t - i1) + to * (i1 + t - 1.0)
            } else {
                let i = integrate_advanced(ease, t);
                from * (t - i) + to * i
            }
        })
    }
}

pub(crate) fn integrate_simple(ease: Easing, t: f64) -> f64 {
    match ease {
        Easing::Linear => t * t / 2.0,

        Easing::InConst => 0.0,
        Easing::OutConst => t,
        Easing::InOutConst => (t - 0.5).max(0.0),

        Easing::InSine => t - 2.0 * f64::sin(PI * t / 2.0) / PI,
        Easing::OutSine => -2.0 * f64::cos(PI * t / 2.0) / PI + 2.0 / PI,
        Easing::InOutSine => 0.5 * t - 0.5 / PI * f64::sin(PI * t),

        Easing::InQuad => t * t * t / 3.0,
        Easing::OutQuad => t * t - t * t * t / 3.0,
        Easing::InOutQuad => {
            if t < 0.5 {
                2.0 * t * t * t / 3.0
            } else {
                -2.0 * t * t * t / 3.0 + 2.0 * t * t - t + 1.0 / 6.0
            }
        }
        Easing::InCubic => f64::powi(t, 4) / 4.0,
        Easing::OutCubic => f64::powi(t, 4) / 4.0 - t * t * t + 1.5 * t * t,
        Easing::InOutCubic => {
            if t < 0.5 {
                f64::powi(t, 4)
            } else {
                f64::powi(t, 4) - 4.0 * t * t * t + 6.0 * t * t - 3.0 * t + 0.5
            }
        }

        Easing::InQuart => f64::powi(t, 5) / 5.0,
        Easing::OutQuart => t - 0.2 * f64::powi(t - 1.0, 5) - 0.2,
        Easing::InOutQuart => {
            if t < 0.5 {
                1.6 * t * t * t * t * t
            } else {
                t - 1.6 * f64::powi(t - 1.0, 5) - 0.5
            }
        }

        Easing::InQuint => f64::powi(t, 6) / 6.0,
        Easing::OutQuint => f64::powi(t - 1.0, 6) / 6.0 + t - 1.0 / 6.0,
        Easing::InOutQuint => {
            if t < 0.5 {
                8.0 * f64::powi(t, 6) / 3.0
            } else {
                t + 8.0 * f64::powi(t - 1.0, 6) / 3.0 - 0.5
            }
        }
        Easing::InExpo => (f64::powf(2.0, 10.0 * t - 11.0) - 1.0 / 2048.0) / (5.0 * LN_2),
        Easing::OutExpo => t + (f64::powf(2.0, -10.0 * t) - 1.0) / (10.0 * LN_2),
        Easing::InOutExpo => {
            if t < 0.5 {
                (f64::powf(2.0, 20.0 * t - 13.0) - 1.0 / 8192.0) / (5.0 * LN_2)
            } else {
                t + f64::powf(2.0, -20.0 * t + 7.0) / (5.0 * LN_2)
                    - 0.5
                    - 1.0 / (8192.0 * 5.0 * LN_2)
            }
        }
    }
}

#[cfg(any(test, feature = "lut_integration"))]
pub(super) fn integrate_advanced(ease: AdvancedEasing, t: f64) -> f64 {
    match ease {
        AdvancedEasing::InCirc => lut::sample_lut(lut::IN_CIRC, t),
        AdvancedEasing::OutCirc => lut::sample_lut(lut::OUT_CIRC, t),
        AdvancedEasing::InOutCirc => lut::sample_lut(lut::INOUT_CIRC, t),
        AdvancedEasing::InBack => lut::sample_lut(lut::IN_BACK, t),
        AdvancedEasing::OutBack => lut::sample_lut(lut::OUT_BACK, t),
        AdvancedEasing::InOutBack => lut::sample_lut(lut::INOUT_BACK, t),
        AdvancedEasing::InElastic => lut::sample_lut(lut::IN_ELASTIC, t),
        AdvancedEasing::OutElastic => lut::sample_lut(lut::OUT_ELASTIC, t),
        AdvancedEasing::InOutElastic => lut::sample_lut(lut::INOUT_ELASTIC, t),
        AdvancedEasing::InBounce => lut::sample_lut(lut::IN_BOUNCE, t),
        AdvancedEasing::OutBounce => lut::sample_lut(lut::OUT_BOUNCE, t),
        AdvancedEasing::InOutBounce => lut::sample_lut(lut::INOUT_BOUNCE, t),
        x => integrate_simple(x.to_monotonic().expect("All non-simple cases handled"), t),
    }
}

#[cfg(any(test, feature = "lut_integration"))]
mod lut {
    const INTEG_LUT_SIZE: usize = 10000;
    const INTEG_LUT_BYTES: usize = INTEG_LUT_SIZE * 4;

    const fn lut(bytes: &[u8]) -> &[u8; INTEG_LUT_BYTES] {
        bytes.as_array().unwrap()
    }
    pub static IN_CIRC: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_in_circ.bin"));
    pub static OUT_CIRC: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_out_circ.bin"));
    pub static INOUT_CIRC: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_inout_circ.bin"));
    pub static IN_BACK: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_in_back.bin"));
    pub static OUT_BACK: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_out_back.bin"));
    pub static INOUT_BACK: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_inout_back.bin"));
    pub static IN_ELASTIC: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_in_elastic.bin"));
    pub static OUT_ELASTIC: &[u8; INTEG_LUT_BYTES] =
        lut(include_bytes!("lut/integ_out_elastic.bin"));
    pub static INOUT_ELASTIC: &[u8; INTEG_LUT_BYTES] =
        lut(include_bytes!("lut/integ_inout_elastic.bin"));
    pub static IN_BOUNCE: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_in_bounce.bin"));
    pub static OUT_BOUNCE: &[u8; INTEG_LUT_BYTES] = lut(include_bytes!("lut/integ_out_bounce.bin"));
    pub static INOUT_BOUNCE: &[u8; INTEG_LUT_BYTES] =
        lut(include_bytes!("lut/integ_inout_bounce.bin"));

    pub fn sample_lut(lut: &[u8; INTEG_LUT_BYTES], t: f64) -> f64 {
        if !t.is_finite() {
            return f64::NAN;
        }

        let t = t.clamp(0.0, 1.0) * (INTEG_LUT_SIZE - 1) as f64;
        let lut_index_0 = t.floor() as usize;
        let lut_index_1 = lut_index_0 + 1;

        let v0 = f32::from_le_bytes(
            *lut[lut_index_0 * 4..lut_index_1 * 4]
                .as_array()
                .expect("4 bytes taken"),
        );

        if lut_index_0 >= INTEG_LUT_SIZE - 1 {
            return v0 as f64;
        }

        let v1 = f32::from_le_bytes(
            *lut[lut_index_1 * 4..(lut_index_1 + 1) * 4]
                .as_array()
                .expect("4 bytes taken"),
        );

        v0 as f64 + (v1 as f64 - v0 as f64) * t.fract()
    }
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::math::easing::eval::evaluate_easing;

    const FROM: f64 = 0.0;
    const TO: f64 = 100.0;
    const H: f64 = 1e-4;

    const INOUT_EASINGS_SYMMETRIC: &[AdvancedEasing] = &[
        AdvancedEasing::InOutSine,
        AdvancedEasing::InOutQuad,
        AdvancedEasing::InOutCubic,
        AdvancedEasing::InOutQuart,
        AdvancedEasing::InOutQuint,
        AdvancedEasing::InOutExpo,
        AdvancedEasing::InOutCirc,
        AdvancedEasing::InOutBack,
        AdvancedEasing::InOutElastic,
        AdvancedEasing::InOutBounce,
    ];

    const SAMPLE_TS: &[f64] = &[0.1, 0.25, 0.5, 0.75, 0.9];

    fn approx_eq(a: f64, b: f64, tol: f64) -> bool {
        (a - b).abs() <= tol
    }

    #[test]
    fn t0_returns_zero() {
        for &ease in AdvancedEasing::ALL {
            let result = evaluate_advanced::<f64, 1>(ease, FROM, TO, 0.0);
            assert!(
                approx_eq(result, 0.0, 1e-3),
                "{ease:?}: expected cumulate(0)=0, got {result}"
            );
        }
    }

    #[test]
    fn derivative_of_cumulate_matches_eval() {
        // d/dt cumulate(from,to,t) = from + (to-from)*e(t) = eval(from,to,t)
        let mut failures: Vec<(AdvancedEasing, f64, f64, f64)> = Vec::new();
        for &ease in AdvancedEasing::ALL {
            if matches!(
                ease,
                AdvancedEasing::InConst | AdvancedEasing::OutConst | AdvancedEasing::InOutConst
            ) {
                continue;
            }

            for &t in SAMPLE_TS {
                let vp = evaluate_advanced::<f64, 1>(ease, FROM, TO, t + H);
                let vm = evaluate_advanced::<f64, 1>(ease, FROM, TO, t - H);
                let numerical = (vp - vm) / (2.0 * H);
                let eval = evaluate_easing::<f64, 1>(ease, FROM, TO, t);
                // Simpson error + finite-diff error compound; allow modest tolerance
                if !approx_eq(numerical, eval, 1.0) {
                    failures.push((ease, t, numerical, eval));
                }
            }
        }
        assert!(
            failures.is_empty(),
            "d/dt cumulate did not match eval: {}",
            failures
                .iter()
                .map(|(e, t, n, v)| format!("{e:?}@t={t}: d/dt={n}, eval={v}"))
                .collect::<Vec<_>>()
                .join("; ")
        );
    }

    #[test]
    fn linear_cumulate_at_1_is_midpoint() {
        // ∫₀¹ [from + (to-from)*s] ds = from/2 + to/2 = midpoint
        let result = evaluate_advanced::<f64, 1>(AdvancedEasing::Linear, FROM, TO, 1.0);
        let midpoint = (FROM + TO) / 2.0;
        assert!(
            approx_eq(result, midpoint, 1e-2),
            "Linear cumulate(1) should be {midpoint}, got {result}"
        );
    }

    #[test]
    fn inout_cumulate_at_1_is_midpoint() {
        // For symmetric InOut easings e(s)+e(1-s)=1, so ∫₀¹ e(s)ds = 0.5.
        // Thus cumulate(1) = from + (to-from)*0.5 = midpoint.
        let midpoint = (FROM + TO) / 2.0;
        let mut failures: Vec<(AdvancedEasing, f64)> = Vec::new();
        for &ease in INOUT_EASINGS_SYMMETRIC {
            let result = evaluate_advanced::<f64, 1>(ease, FROM, TO, 1.0);
            if !approx_eq(result, midpoint, 0.1) {
                failures.push((ease, result));
            }
        }
        assert!(
            failures.is_empty(),
            "InOut cumulate(1) should be midpoint ({}), but got: {}",
            midpoint,
            failures
                .iter()
                .map(|(e, r)| format!("{e:?}={r}"))
                .collect::<Vec<_>>()
                .join(", ")
        );
    }

    #[test]
    fn cumulate_matches_closed_form() {
        // Rigorous checks against exact antiderivatives for smooth curves.
        // cumulate(from,to,t) = from*t + (to-from)*I(t), with from=0, to=100.
        // Linear: I(t)=t²/2; InQuad: I(t)=t³/3; InCubic: I(t)=t⁴/4;
        // OutQuad: e=2t-t², I(t)=t²-t³/3; InSine: I(t)=t-(2/π)sin(πt/2).
        let cases: &[(AdvancedEasing, f64, f64, &str)] = &[
            (AdvancedEasing::Linear, 0.5, 50.0 * 0.5 * 0.5, "50 t^2"),
            (
                AdvancedEasing::InQuad,
                0.5,
                100.0 * 0.125 / 3.0,
                "100 t^3/3",
            ),
            (AdvancedEasing::InCubic, 0.5, 25.0 * 0.0625, "25 t^4"),
            (
                AdvancedEasing::OutQuad,
                0.5,
                100.0 * (0.25 - 0.125 / 3.0),
                "100(t^2-t^3/3)",
            ),
            (
                AdvancedEasing::InSine,
                0.5,
                100.0 * (0.5 - 2.0 / std::f64::consts::PI * (std::f64::consts::PI / 4.0).sin()),
                "100(t-(2/pi)sin(pi t/2))",
            ),
        ];
        let mut failures: Vec<(AdvancedEasing, f64, f64, &str)> = Vec::new();
        for &(ease, t, expected, label) in cases {
            let result = evaluate_advanced::<f64, 1>(ease, FROM, TO, t);
            if !approx_eq(result, expected, 1e-2) {
                failures.push((ease, result, expected, label));
            }
        }
        assert!(
            failures.is_empty(),
            "cumulate did not match closed form: {}",
            failures
                .iter()
                .map(|(e, r, exp, lbl)| format!("{e:?}: got {r}, expected {exp} ({lbl})"))
                .collect::<Vec<_>>()
                .join("; ")
        );
    }

    // Temporary LUT generator.
    // Produces 12 .bin files in src/easing/lut/ using the old Simpson's rule.
    // use crate::math::easing::eval::easing_value;
    //
    // fn simpson_integral(ease: Easing, t: f64, n: usize) -> f64 {
    //     if t <= 0.0 {
    //         return 0.0;
    //     }
    //     let n = if n % 2 == 1 { n + 1 } else { n };
    //     let h = t / n as f64;
    //     let mut sum = easing_value(ease, 0.0) + easing_value(ease, t);
    //     for k in 1..n {
    //         let s = k as f64 * h;
    //         let w = if k % 2 == 0 { 2.0 } else { 4.0 };
    //         sum += w * easing_value(ease, s);
    //     }
    //     sum * h / 3.0
    // }

    // #[test]
    // #[ignore]
    // fn generate_luts() {
    //     use std::fs::File;
    //     use std::io::{BufWriter, Write};

    //     let luts: &[(Easing, &str)] = &[
    //         (Easing::InCirc, "integ_in_circ.bin"),
    //         (Easing::OutCirc, "integ_out_circ.bin"),
    //         (Easing::InOutCirc, "integ_inout_circ.bin"),
    //         (Easing::InBack, "integ_in_back.bin"),
    //         (Easing::OutBack, "integ_out_back.bin"),
    //         (Easing::InOutBack, "integ_inout_back.bin"),
    //         (Easing::InElastic, "integ_in_elastic.bin"),
    //         (Easing::OutElastic, "integ_out_elastic.bin"),
    //         (Easing::InOutElastic, "integ_inout_elastic.bin"),
    //         (Easing::InBounce, "integ_in_bounce.bin"),
    //         (Easing::OutBounce, "integ_out_bounce.bin"),
    //         (Easing::InOutBounce, "integ_inout_bounce.bin"),
    //     ];

    //     const SIMPSON_N: usize = 20000;

    //     for &(ease, filename) in luts {
    //         let mut data = Vec::with_capacity(INTEG_LUT_SIZE * 4);

    //         for k in 0..INTEG_LUT_SIZE {
    //             let t = k as f64 / (INTEG_LUT_SIZE - 1) as f64;
    //             let integral = simpson_integral(ease, t, SIMPSON_N);
    //             data.extend_from_slice(&(integral as f32).to_le_bytes());
    //         }

    //         let path = format!("src/easing/lut/{filename}");
    //         let file = File::create(&path).expect("create lut file");
    //         let mut writer = BufWriter::new(file);
    //         writer.write_all(&data).expect("write lut data");
    //         writer.flush().expect("flush lut file");
    //         eprintln!("Wrote {path} ({INTEG_LUT_SIZE} entries)");
    //     }
    // }
}