libss 0.1.0

use crate::Field;

// TODO: Find a faster algorithm than naive
/// Perform lagrange interpolation on the slice of points given as an argument
pub fn interpolate<Fq: Field>(points: &[(Fq, Fq)]) -> Fq {
    let mut result = Fq::zero();
    for i in 0..points.len() {
        let mut term = points[i].1;
        for j in 0..points.len() {
            if i != j {
                let num = -points[j].0;
                let denom = points[i].0 - points[j].0;
                term *= num / denom;
            }
        }
        result += term;
    }

    result
}

// http://terpconnect.umd.edu/~petersd/666/BarycentricLagrange1.pdf
/// Perform barycentric lagrange interpolation on the slice of points given as an argument
pub fn barycentric_interpolate<Fq: Field>(points: &[(Fq, Fq)]) -> Fq {
    let mut result = Fq::zero();
    let mut l_x = Fq::one();
    let mut w = Vec::with_capacity(points.len());
    for j in 0..points.len() {
        let x_j = points[j].0;
        w.push(Fq::one());
        l_x *= -x_j;
        for (k, &(x_k, _)) in points.iter().enumerate() {
            if j != k {
                let denom = x_j - x_k;
                w[j] *= Fq::one() / denom;
            }
        }
    }

    for (j, &(x_j, y_j)) in points.iter().enumerate() {
        result += (w[j] * y_j) / -x_j;
    }

    result *= l_x;

    result
}

#[cfg(test)]
mod test {
    use crate::gf256::GF256;
    use crate::interpolate::{barycentric_interpolate, interpolate};
    use crate::polynomial::Polynomial;

    quickcheck! {
        fn test_interpolate_with_intercept(intercept: GF256) -> bool {
            let size = 32;
            let poly: Polynomial<GF256> = Polynomial::from_intercept(intercept, size);
            let points = (1..size+1)
                .map(|i| (GF256(i as u8), poly.evaluate_at(GF256(i as u8))))
                .collect::<Vec<(GF256,GF256)>>();
            let interpolated_point = interpolate(&points);
            intercept == interpolated_point
        }

        fn test_barycentric_interpolate_with_intercept(intercept: GF256) -> bool {
            let size = 32;
            let poly: Polynomial<GF256> = Polynomial::from_intercept(intercept, size);
            let points = (1..size+1)
                .map(|i| (GF256(i as u8), poly.evaluate_at(GF256(i as u8))))
                .collect::<Vec<(GF256,GF256)>>();
            let interpolated_point = barycentric_interpolate(&points);
            intercept == interpolated_point
        }
    }
}