wip

math/src/circle/cfft.rs

@@ -6,6 +6,8 @@ pub fn inplace_cfft(
     input: &mut [FieldElement<Mersenne31Field>],
     twiddles: Vec<Vec<FieldElement<Mersenne31Field>>>,
 ) {
+    use super::twiddles::TwiddlesConfig;
+
     let mut group_count = 1;
     let mut group_size = input.len();
     let mut round = 0;
@@ -51,109 +53,3 @@ pub fn reverse_cfft_index(index: usize, log_2_size: u32) -> usize {
     }
     new_index.reverse_bits() >> (usize::BITS - log_2_size)
 }
-
-#[cfg(test)]
-mod tests {
-    use super::*;
-    use crate::circle::{cosets::Coset, twiddles::get_twiddles};
-    type FpE = FieldElement<Mersenne31Field>;
-
-    fn evaluate_poly(coef: &[FpE; 8], x: FpE, y: FpE) -> FpE {
-        coef[0]
-            + coef[1] * y
-            + coef[2] * x
-            + coef[3] * x * y
-            + coef[4] * (x.square().double() - FpE::one())
-            + coef[5] * (x.square().double() - FpE::one()) * y
-            + coef[6] * ((x.square() * x).double() - x)
-            + coef[7] * ((x.square() * x).double() - x) * y
-    }
-
-    fn evaluate_poly_16(coef: &[FpE; 16], x: FpE, y: FpE) -> FpE {
-        let mut a = x;
-        let mut v = Vec::new();
-        v.push(FpE::one());
-        v.push(x);
-        for _ in 2..4 {
-            a = a.square().double() - FpE::one();
-            v.push(a);
-        }
-
-        coef[0] * v[0]
-            + coef[1] * y * v[0]
-            + coef[2] * v[1]
-            + coef[3] * y * v[1]
-            + coef[4] * v[2]
-            + coef[5] * y * v[2]
-            + coef[6] * v[1] * v[2]
-            + coef[7] * y * v[1] * v[2]
-            + coef[8] * v[3]
-            + coef[9] * y * v[3]
-            + coef[10] * v[1] * v[3]
-            + coef[11] * y * v[1] * v[3]
-            + coef[12] * v[2] * v[3]
-            + coef[13] * y * v[2] * v[3]
-            + coef[14] * v[1] * v[2] * v[3]
-            + coef[15] * y * v[1] * v[2] * v[3]
-    }
-
-    #[test]
-    fn cfft_test() {
-        let coset = Coset::new_standard(3);
-        let points = Coset::get_coset_points(&coset);
-        let twiddles = get_twiddles(coset);
-        let mut input = [
-            FpE::from(1),
-            FpE::from(2),
-            FpE::from(3),
-            FpE::from(4),
-            FpE::from(5),
-            FpE::from(6),
-            FpE::from(7),
-            FpE::from(8),
-        ];
-        let mut expected_result: Vec<FpE> = Vec::new();
-        for point in points {
-            let point_eval = evaluate_poly(&input, point.x, point.y);
-            expected_result.push(point_eval);
-        }
-        inplace_cfft(&mut input, twiddles);
-        inplace_order_cfft_values(&mut input);
-        let result: &[FpE] = &input;
-        assert_eq!(result, expected_result);
-    }
-
-    #[test]
-    fn cfft_test_16() {
-        let coset = Coset::new_standard(4);
-        let points = Coset::get_coset_points(&coset);
-        let twiddles = get_twiddles(coset);
-        let mut input = [
-            FpE::from(1),
-            FpE::from(2),
-            FpE::from(3),
-            FpE::from(4),
-            FpE::from(5),
-            FpE::from(6),
-            FpE::from(7),
-            FpE::from(8),
-            FpE::from(9),
-            FpE::from(10),
-            FpE::from(11),
-            FpE::from(12),
-            FpE::from(13),
-            FpE::from(14),
-            FpE::from(15),
-            FpE::from(16),
-        ];
-        let mut expected_result: Vec<FpE> = Vec::new();
-        for point in points {
-            let point_eval = evaluate_poly_16(&input, point.x, point.y);
-            expected_result.push(point_eval);
-        }
-        inplace_cfft(&mut input, twiddles);
-        inplace_order_cfft_values(&mut input);
-        let result: &[FpE] = &input;
-        assert_eq!(result, expected_result);
-    }
-}

math/src/circle/mod.rs

@@ -3,3 +3,4 @@ pub mod cosets;
 pub mod errors;
 pub mod point;
 pub mod twiddles;
+pub mod polynomial;

math/src/circle/point.rs

@@ -117,6 +117,7 @@ impl<F: IsField + HasCircleParams<F>> CirclePoint<F> {
 
     /// Computes the inverse of the point.
     /// We are using -(x, y) = (x, -y), i.e. the inverse of the group opertion is conjugation.
+    /// because the norm of every point in the circle is one.
     pub fn conjugate(self) -> Self {
         Self {
             x: self.x,

math/src/circle/polynomial.rs

@@ -0,0 +1,161 @@
+use crate::field::{element::FieldElement, fields::mersenne31::field::Mersenne31Field};
+
+use super::{cfft::{inplace_cfft, inplace_order_cfft_values}, cosets::Coset, twiddles::{get_twiddles, TwiddlesConfig}};
+
+/// Given the 2^n coefficients of a two-variables polynomial in the basis {1, y, x, xy, 2xˆ2 -1, 2xˆ2y-y, 2xˆ3-x, 2xˆ3y-xy,...}
+/// returns the evaluation of the polynomianl on the points of the standard coset of size 2^n.
+/// Note that coeff has to be a vector with length a power of two 2^n.
+pub fn evaluate_cfft(mut coeff: Vec<FieldElement<Mersenne31Field>>) -> Vec<FieldElement<Mersenne31Field>>{
+    let domain_log_2_size: u32 = coeff.len().trailing_zeros();
+    let coset = Coset::new_standard(domain_log_2_size);
+    let config = TwiddlesConfig::Evaluation;
+    let twiddles = get_twiddles(coset, config);
+
+    inplace_cfft(&mut coeff, twiddles);
+    inplace_order_cfft_values(&mut coeff);
+    coeff
+}
+
+/// Interpolates the 2^n evaluations of a two-variables polynomial on the points of the standard coset of size 2^n.
+/// As a result we obtain the coefficients of the polynomial in the basis: {1, y, x, xy, 2xˆ2 -1, 2xˆ2y-y, 2xˆ3-x, 2xˆ3y-xy,...}
+/// Note that eval has to be a vector of length a power of two 2^n.
+pub fn interpolate_cfft(mut eval: Vec<FieldElement<Mersenne31Field>>) -> Vec<FieldElement<Mersenne31Field>>{
+    let domain_log_2_size: u32 = eval.len().trailing_zeros();
+    let coset = Coset::new_standard(domain_log_2_size);
+    let config = TwiddlesConfig::Interpolation;
+    let twiddles = get_twiddles(coset, config);
+
+    inplace_cfft(&mut eval, twiddles);
+    inplace_order_cfft_values(&mut eval);
+    eval
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::circle::cosets::Coset;
+    type FpE = FieldElement<Mersenne31Field>;
+
+    fn evaluate_poly(coef: &[FpE; 8], x: FpE, y: FpE) -> FpE {
+        coef[0]
+            + coef[1] * y
+            + coef[2] * x
+            + coef[3] * x * y
+            + coef[4] * (x.square().double() - FpE::one())
+            + coef[5] * (x.square().double() - FpE::one()) * y
+            + coef[6] * ((x.square() * x).double() - x)
+            + coef[7] * ((x.square() * x).double() - x) * y
+    }
+
+    fn evaluate_poly_16(coef: &[FpE; 16], x: FpE, y: FpE) -> FpE {
+        let mut a = x;
+        let mut v = Vec::new();
+        v.push(FpE::one());
+        v.push(x);
+        for _ in 2..4 {
+            a = a.square().double() - FpE::one();
+            v.push(a);
+        }
+
+        coef[0] * v[0]
+            + coef[1] * y * v[0]
+            + coef[2] * v[1]
+            + coef[3] * y * v[1]
+            + coef[4] * v[2]
+            + coef[5] * y * v[2]
+            + coef[6] * v[1] * v[2]
+            + coef[7] * y * v[1] * v[2]
+            + coef[8] * v[3]
+            + coef[9] * y * v[3]
+            + coef[10] * v[1] * v[3]
+            + coef[11] * y * v[1] * v[3]
+            + coef[12] * v[2] * v[3]
+            + coef[13] * y * v[2] * v[3]
+            + coef[14] * v[1] * v[2] * v[3]
+            + coef[15] * y * v[1] * v[2] * v[3]
+    }
+
+    #[test]
+    fn cfft_evaluation_8_points() {
+        // We create the coset points and evaluate them without the fft.
+        let coset = Coset::new_standard(3);
+        let points = Coset::get_coset_points(&coset);
+        let mut input = [
+            FpE::from(1),
+            FpE::from(2),
+            FpE::from(3),
+            FpE::from(4),
+            FpE::from(5),
+            FpE::from(6),
+            FpE::from(7),
+            FpE::from(8),
+        ];
+        let mut expected_result: Vec<FpE> = Vec::new();
+        for point in points {
+            let point_eval = evaluate_poly(&input, point.x, point.y);
+            expected_result.push(point_eval);
+        }
+
+        let result = evaluate_cfft(input.to_vec());
+        let slice_result: &[FpE] = &result;
+        assert_eq!(slice_result, expected_result);
+    }
+
+    #[test]
+    fn cfft_evaluation_16_points() {
+        let coset = Coset::new_standard(4);
+        let points = Coset::get_coset_points(&coset);
+        let mut input = [
+            FpE::from(1),
+            FpE::from(2),
+            FpE::from(3),
+            FpE::from(4),
+            FpE::from(5),
+            FpE::from(6),
+            FpE::from(7),
+            FpE::from(8),
+            FpE::from(9),
+            FpE::from(10),
+            FpE::from(11),
+            FpE::from(12),
+            FpE::from(13),
+            FpE::from(14),
+            FpE::from(15),
+            FpE::from(16),
+        ];
+        let mut expected_result: Vec<FpE> = Vec::new();
+        for point in points {
+            let point_eval = evaluate_poly_16(&input, point.x, point.y);
+            expected_result.push(point_eval);
+        }
+
+        let result = evaluate_cfft(input.to_vec());
+        let slice_result: &[FpE] = &result;
+        assert_eq!(slice_result, expected_result);
+    }
+
+    #[test]
+    fn evaluate_and_interpolate_8() {
+        // 
+        let coeff = vec![
+            FpE::from(1),
+            FpE::from(2),
+            FpE::from(3),
+            FpE::from(4),
+            FpE::from(5),
+            FpE::from(6),
+            FpE::from(7),
+            FpE::from(8),
+        ];
+
+
+        let evals = evaluate_cfft(coeff.clone());
+        let factor = FpE::from(8).inv().unwrap();
+        let mut new_coeff = interpolate_cfft(evals);
+        new_coeff = new_coeff.iter()
+            .map(|coeff| factor * coeff)
+            .collect();
+        assert_eq!(new_coeff, coeff);
+    }
+}
+

math/src/circle/twiddles.rs

@@ -6,17 +6,20 @@ use crate::{
 };
 use alloc::vec::Vec;
 
+#[derive(PartialEq)]
+pub enum TwiddlesConfig {
+    Evaluation,
+    Interpolation,
+}
 #[cfg(feature = "alloc")]
-pub fn get_twiddles(domain: Coset) -> Vec<Vec<FieldElement<Mersenne31Field>>> {
+pub fn get_twiddles(domain: Coset, config: TwiddlesConfig) -> Vec<Vec<FieldElement<Mersenne31Field>>> {
     let mut half_domain_points = Coset::get_coset_points(&Coset::half_coset(domain.clone()));
     in_place_bit_reverse_permute::<CirclePoint<Mersenne31Field>>(&mut half_domain_points[..]);
 
-    let mut twiddles: Vec<Vec<FieldElement<Mersenne31Field>>> =
-        vec![half_domain_points.iter().map(|p| p.y).collect()];
+    let mut twiddles: Vec<Vec<FieldElement<Mersenne31Field>>> =  vec![half_domain_points.iter().map(|p| p.y).collect()];
 
     if domain.log_2_size >= 2 {
         twiddles.push(half_domain_points.iter().step_by(2).map(|p| p.x).collect());
-
         for _ in 0..(domain.log_2_size - 2) {
             let prev = twiddles.last().unwrap();
             let cur = prev
@@ -28,6 +31,12 @@ pub fn get_twiddles(domain: Coset) -> Vec<Vec<FieldElement<Mersenne31Field>>> {
         }
     }
     twiddles.reverse();
+
+    if config == TwiddlesConfig::Interpolation {
+        twiddles.iter_mut().for_each(|x| {
+            FieldElement::<Mersenne31Field>::inplace_batch_inverse(x).unwrap();
+        });
+    }
     twiddles
 }
 
@@ -38,7 +47,8 @@ mod tests {
     #[test]
     fn twiddles_vectors_lenght() {
         let domain = Coset::new_standard(3);
-        let twiddles = get_twiddles(domain);
+        let config = TwiddlesConfig::Evaluation;
+        let twiddles = get_twiddles(domain, config);
         for i in 0..twiddles.len() - 1 {
             assert_eq!(twiddles[i].len(), 2 * twiddles[i + 1].len())
         }