Merge branch 'main' into insertion_gate

src/bin/bench_recursion.rs

@@ -21,8 +21,6 @@ fn main() {
 }
 
 fn bench_prove<F: Field + Extendable<D>, const D: usize>() {
-    let gmimc_gate = GMiMCGate::<F, D, GMIMC_ROUNDS>::with_automatic_constants();
-
     let config = CircuitConfig {
         num_wires: 134,
         num_routed_wires: 27,
@@ -32,7 +30,7 @@ fn bench_prove<F: Field + Extendable<D>, const D: usize>() {
         fri_config: FriConfig {
             proof_of_work_bits: 1,
             rate_bits: 3,
-            reduction_arity_bits: vec![1],
+            reduction_arity_bits: vec![2, 2, 2, 2, 2],
             num_query_rounds: 1,
         },
     };

src/circuit_builder.rs

@@ -36,7 +36,7 @@ pub struct CircuitBuilder<F: Extendable<D>, const D: usize> {
     /// The next available index for a public input.
     public_input_index: usize,
 
-    /// The next available index for a VirtualAdviceTarget.
+    /// The next available index for a `VirtualTarget`.
     virtual_target_index: usize,
 
     copy_constraints: Vec<(Target, Target)>,
@@ -63,6 +63,10 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         }
     }
 
+    pub fn num_gates(&self) -> usize {
+        self.gate_instances.len()
+    }
+
     pub fn add_public_input(&mut self) -> Target {
         let index = self.public_input_index;
         self.public_input_index += 1;
@@ -73,22 +77,18 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         (0..n).map(|_i| self.add_public_input()).collect()
     }
 
-    /// Adds a new "virtual" advice target. This is not an actual wire in the witness, but just a
+    /// Adds a new "virtual" target. This is not an actual wire in the witness, but just a target
-    /// target that help facilitate witness generation. In particular, a generator can assign a
+    /// that help facilitate witness generation. In particular, a generator can assign a values to a
-    /// values to a virtual target, which can then be copied to other (virtual or concrete) targets
+    /// virtual target, which can then be copied to other (virtual or concrete) targets. When we
-    /// via `generate_copy`. When we generate the final witness (a grid of wire values), these
+    /// generate the final witness (a grid of wire values), these virtual targets will go away.
-    /// virtual targets will go away.
+    pub fn add_virtual_target(&mut self) -> Target {
-    ///
-    /// Since virtual targets are not part of the actual permutation argument, they cannot be used
-    /// with `assert_equal`.
-    pub fn add_virtual_advice_target(&mut self) -> Target {
         let index = self.virtual_target_index;
         self.virtual_target_index += 1;
-        Target::VirtualAdviceTarget { index }
+        Target::VirtualTarget { index }
     }
 
-    pub fn add_virtual_advice_targets(&mut self, n: usize) -> Vec<Target> {
+    pub fn add_virtual_targets(&mut self, n: usize) -> Vec<Target> {
-        (0..n).map(|_i| self.add_virtual_advice_target()).collect()
+        (0..n).map(|_i| self.add_virtual_target()).collect()
     }
 
     pub fn add_gate_no_constants(&mut self, gate_type: GateRef<F, D>) -> usize {
@@ -97,6 +97,11 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
 
     /// Adds a gate to the circuit, and returns its index.
     pub fn add_gate(&mut self, gate_type: GateRef<F, D>, constants: Vec<F>) -> usize {
+        assert_eq!(
+            gate_type.0.num_constants(),
+            constants.len(),
+            "Number of constants doesn't match."
+        );
         // If we haven't seen a gate of this type before, check that it's compatible with our
         // circuit configuration, then register it.
         if !self.gates.contains(&gate_type) {
@@ -250,7 +255,7 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
     /// polynomials (which are opened at only one location) and for the Z polynomials (which are
     /// opened at two).
     fn blinding_counts(&self) -> (usize, usize) {
-        let num_gates = self.gates.len();
+        let num_gates = self.gate_instances.len();
         let mut degree_estimate = 1 << log2_ceil(num_gates);
 
         loop {
@@ -272,6 +277,10 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
 
     fn blind_and_pad(&mut self) {
         let (regular_poly_openings, z_openings) = self.blinding_counts();
+        info!(
+            "Adding {} blinding terms for witness polynomials, and {}*2 for Z polynomials",
+            regular_poly_openings, z_openings
+        );
 
         let num_routed_wires = self.config.num_routed_wires;
         let num_wires = self.config.num_wires;
@@ -319,12 +328,11 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         }
     }
 
-    fn constant_polys(&self, gates: &[PrefixedGate<F, D>]) -> Vec<PolynomialValues<F>> {
+    fn constant_polys(
-        let num_constants = gates
+        &self,
-            .iter()
+        gates: &[PrefixedGate<F, D>],
-            .map(|gate| gate.gate.0.num_constants() + gate.prefix.len())
+        num_constants: usize,
-            .max()
+    ) -> Vec<PolynomialValues<F>> {
-            .unwrap();
         let constants_per_gate = self
             .gate_instances
             .iter()
@@ -360,7 +368,11 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         }
 
         for index in 0..self.public_input_index {
-            target_partitions.add_partition(Target::PublicInput { index })
+            target_partitions.add_partition(Target::PublicInput { index });
+        }
+
+        for index in 0..self.virtual_target_index {
+            target_partitions.add_partition(Target::VirtualTarget { index });
         }
 
         for &(a, b) in &self.copy_constraints {
@@ -375,22 +387,21 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
     pub fn build(mut self) -> CircuitData<F, D> {
         let start = Instant::now();
         info!(
-            "degree before blinding & padding: {}",
+            "Degree before blinding & padding: {}",
             self.gate_instances.len()
         );
         self.blind_and_pad();
         let degree = self.gate_instances.len();
-        info!("degree after blinding & padding: {}", degree);
+        info!("Degree after blinding & padding: {}", degree);
 
         let gates = self.gates.iter().cloned().collect();
-        let gate_tree = Tree::from_gates(gates);
+        let (gate_tree, max_filtered_constraint_degree, num_constants) = Tree::from_gates(gates);
         let prefixed_gates = PrefixedGate::from_tree(gate_tree);
 
         let degree_bits = log2_strict(degree);
         let subgroup = F::two_adic_subgroup(degree_bits);
 
-        let constant_vecs = self.constant_polys(&prefixed_gates);
+        let constant_vecs = self.constant_polys(&prefixed_gates, num_constants);
-        let num_constants = constant_vecs.len();
 
         let k_is = get_unique_coset_shifts(degree, self.config.num_routed_wires);
         let sigma_vecs = self.sigma_vecs(&k_is, &subgroup);
@@ -438,7 +449,7 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
             config: self.config,
             degree_bits,
             gates: prefixed_gates,
-            max_filtered_constraint_degree_bits: 3, // TODO: compute this correctly once filters land.
+            max_filtered_constraint_degree,
             num_gate_constraints,
             num_constants,
             k_is,

src/circuit_data.rs

@@ -39,7 +39,7 @@ impl Default for CircuitConfig {
             fri_config: FriConfig {
                 proof_of_work_bits: 1,
                 rate_bits: 1,
-                reduction_arity_bits: vec![1],
+                reduction_arity_bits: vec![1, 1, 1, 1],
                 num_query_rounds: 1,
             },
         }
@@ -54,14 +54,14 @@ impl CircuitConfig {
     pub(crate) fn large_config() -> Self {
         Self {
             num_wires: 134,
-            num_routed_wires: 12,
+            num_routed_wires: 28,
             security_bits: 128,
             rate_bits: 3,
             num_challenges: 3,
             fri_config: FriConfig {
                 proof_of_work_bits: 1,
                 rate_bits: 3,
-                reduction_arity_bits: vec![1],
+                reduction_arity_bits: vec![1, 1, 1, 1],
                 num_query_rounds: 1,
             },
         }
@@ -146,7 +146,7 @@ pub struct CommonCircuitData<F: Extendable<D>, const D: usize> {
     pub(crate) gates: Vec<PrefixedGate<F, D>>,
 
     /// The maximum degree of a filter times a constraint by any gate.
-    pub(crate) max_filtered_constraint_degree_bits: usize,
+    pub(crate) max_filtered_constraint_degree: usize,
 
     /// The largest number of constraints imposed by any gate.
     pub(crate) num_gate_constraints: usize,
@@ -184,7 +184,7 @@ impl<F: Extendable<D>, const D: usize> CommonCircuitData<F, D> {
     }
 
     pub fn quotient_degree(&self) -> usize {
-        ((1 << self.max_filtered_constraint_degree_bits) - 1) * self.degree()
+        (self.max_filtered_constraint_degree - 1) * self.degree()
     }
 
     pub fn total_constraints(&self) -> usize {

src/field/crandall_field.rs

@@ -136,6 +136,8 @@ impl Debug for CrandallField {
 }
 
 impl Field for CrandallField {
+    type PrimeField = Self;
+
     const ZERO: Self = Self(0);
     const ONE: Self = Self(1);
     const TWO: Self = Self(2);
@@ -143,6 +145,7 @@ impl Field for CrandallField {
 
     const ORDER: u64 = 18446744071293632513;
     const TWO_ADICITY: usize = 28;
+    const CHARACTERISTIC: u64 = Self::ORDER;
 
     const MULTIPLICATIVE_GROUP_GENERATOR: Self = Self(5);
     const POWER_OF_TWO_GENERATOR: Self = Self(10281950781551402419);

src/field/extension_field/quadratic.rs

@@ -43,11 +43,14 @@ impl From<<Self as FieldExtension<2>>::BaseField> for QuadraticCrandallField {
 }
 
 impl Field for QuadraticCrandallField {
+    type PrimeField = CrandallField;
+
     const ZERO: Self = Self([CrandallField::ZERO; 2]);
     const ONE: Self = Self([CrandallField::ONE, CrandallField::ZERO]);
     const TWO: Self = Self([CrandallField::TWO, CrandallField::ZERO]);
     const NEG_ONE: Self = Self([CrandallField::NEG_ONE, CrandallField::ZERO]);
 
+    const CHARACTERISTIC: u64 = CrandallField::ORDER;
     // Does not fit in 64-bits.
     const ORDER: u64 = 0;
     const TWO_ADICITY: usize = 29;

src/field/extension_field/quartic.rs

@@ -50,6 +50,8 @@ impl From<<Self as FieldExtension<4>>::BaseField> for QuarticCrandallField {
 }
 
 impl Field for QuarticCrandallField {
+    type PrimeField = CrandallField;
+
     const ZERO: Self = Self([CrandallField::ZERO; 4]);
     const ONE: Self = Self([
         CrandallField::ONE,
@@ -70,6 +72,7 @@ impl Field for QuarticCrandallField {
         CrandallField::ZERO,
     ]);
 
+    const CHARACTERISTIC: u64 = CrandallField::ORDER;
     // Does not fit in 64-bits.
     const ORDER: u64 = 0;
     const TWO_ADICITY: usize = 30;

src/field/extension_field/target.rs

@@ -5,7 +5,6 @@ use crate::circuit_builder::CircuitBuilder;
 use crate::field::extension_field::algebra::ExtensionAlgebra;
 use crate::field::extension_field::{Extendable, FieldExtension, OEF};
 use crate::field::field::Field;
-use crate::gates::mul_extension::MulExtensionGate;
 use crate::target::Target;
 
 /// `Target`s representing an element of an extension field.
@@ -110,160 +109,6 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         self.constant_ext_algebra(ExtensionAlgebra::ZERO)
     }
 
-    pub fn add_extension(
-        &mut self,
-        mut a: ExtensionTarget<D>,
-        b: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        for i in 0..D {
-            a.0[i] = self.add(a.0[i], b.0[i]);
-        }
-        a
-    }
-
-    pub fn add_ext_algebra(
-        &mut self,
-        mut a: ExtensionAlgebraTarget<D>,
-        b: ExtensionAlgebraTarget<D>,
-    ) -> ExtensionAlgebraTarget<D> {
-        for i in 0..D {
-            a.0[i] = self.add_extension(a.0[i], b.0[i]);
-        }
-        a
-    }
-
-    pub fn add_many_extension(&mut self, terms: &[ExtensionTarget<D>]) -> ExtensionTarget<D> {
-        let mut sum = self.zero_extension();
-        for term in terms {
-            sum = self.add_extension(sum, *term);
-        }
-        sum
-    }
-
-    /// TODO: Change this to using an `arithmetic_extension` function once `MulExtensionGate` supports addend.
-    pub fn sub_extension(
-        &mut self,
-        mut a: ExtensionTarget<D>,
-        b: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        for i in 0..D {
-            a.0[i] = self.sub(a.0[i], b.0[i]);
-        }
-        a
-    }
-
-    pub fn sub_ext_algebra(
-        &mut self,
-        mut a: ExtensionAlgebraTarget<D>,
-        b: ExtensionAlgebraTarget<D>,
-    ) -> ExtensionAlgebraTarget<D> {
-        for i in 0..D {
-            a.0[i] = self.sub_extension(a.0[i], b.0[i]);
-        }
-        a
-    }
-
-    pub fn mul_extension_with_const(
-        &mut self,
-        const_0: F,
-        multiplicand_0: ExtensionTarget<D>,
-        multiplicand_1: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        let gate = self.add_gate(MulExtensionGate::new(), vec![const_0]);
-
-        let wire_multiplicand_0 =
-            ExtensionTarget::from_range(gate, MulExtensionGate::<D>::wires_multiplicand_0());
-        let wire_multiplicand_1 =
-            ExtensionTarget::from_range(gate, MulExtensionGate::<D>::wires_multiplicand_1());
-        let wire_output = ExtensionTarget::from_range(gate, MulExtensionGate::<D>::wires_output());
-
-        self.route_extension(multiplicand_0, wire_multiplicand_0);
-        self.route_extension(multiplicand_1, wire_multiplicand_1);
-        wire_output
-    }
-
-    pub fn mul_extension(
-        &mut self,
-        multiplicand_0: ExtensionTarget<D>,
-        multiplicand_1: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        self.mul_extension_with_const(F::ONE, multiplicand_0, multiplicand_1)
-    }
-
-    pub fn mul_ext_algebra(
-        &mut self,
-        a: ExtensionAlgebraTarget<D>,
-        b: ExtensionAlgebraTarget<D>,
-    ) -> ExtensionAlgebraTarget<D> {
-        let mut res = [self.zero_extension(); D];
-        let w = self.constant(F::Extension::W);
-        for i in 0..D {
-            for j in 0..D {
-                let ai_bi = self.mul_extension(a.0[i], b.0[j]);
-                res[(i + j) % D] = if i + j < D {
-                    self.add_extension(ai_bi, res[(i + j) % D])
-                } else {
-                    let w_ai_bi = self.scalar_mul_ext(w, ai_bi);
-                    self.add_extension(w_ai_bi, res[(i + j) % D])
-                }
-            }
-        }
-        ExtensionAlgebraTarget(res)
-    }
-
-    pub fn mul_many_extension(&mut self, terms: &[ExtensionTarget<D>]) -> ExtensionTarget<D> {
-        let mut product = self.one_extension();
-        for term in terms {
-            product = self.mul_extension(product, *term);
-        }
-        product
-    }
-
-    /// Like `mul_add`, but for `ExtensionTarget`s. Note that, unlike `mul_add`, this has no
-    /// performance benefit over separate muls and adds.
-    /// TODO: Change this to using an `arithmetic_extension` function once `MulExtensionGate` supports addend.
-    pub fn mul_add_extension(
-        &mut self,
-        a: ExtensionTarget<D>,
-        b: ExtensionTarget<D>,
-        c: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        let product = self.mul_extension(a, b);
-        self.add_extension(product, c)
-    }
-
-    /// Like `mul_sub`, but for `ExtensionTarget`s. Note that, unlike `mul_sub`, this has no
-    /// performance benefit over separate muls and subs.
-    /// TODO: Change this to using an `arithmetic_extension` function once `MulExtensionGate` supports addend.
-    pub fn scalar_mul_sub_extension(
-        &mut self,
-        a: Target,
-        b: ExtensionTarget<D>,
-        c: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        let product = self.scalar_mul_ext(a, b);
-        self.sub_extension(product, c)
-    }
-
-    /// Returns `a * b`, where `b` is in the extension field and `a` is in the base field.
-    pub fn scalar_mul_ext(&mut self, a: Target, b: ExtensionTarget<D>) -> ExtensionTarget<D> {
-        let a_ext = self.convert_to_ext(a);
-        self.mul_extension(a_ext, b)
-    }
-
-    /// Returns `a * b`, where `b` is in the extension of the extension field, and `a` is in the
-    /// extension field.
-    pub fn scalar_mul_ext_algebra(
-        &mut self,
-        a: ExtensionTarget<D>,
-        mut b: ExtensionAlgebraTarget<D>,
-    ) -> ExtensionAlgebraTarget<D> {
-        for i in 0..D {
-            b.0[i] = self.mul_extension(a, b.0[i]);
-        }
-        b
-    }
-
     pub fn convert_to_ext(&mut self, t: Target) -> ExtensionTarget<D> {
         let zero = self.zero();
         let mut arr = [zero; D];

src/field/fft.rs

@@ -1,142 +1,304 @@
+use std::option::Option;
+
 use crate::field::field::Field;
 use crate::polynomial::polynomial::{PolynomialCoeffs, PolynomialValues};
-use crate::util::{log2_ceil, log2_strict};
+use crate::util::{log2_strict, reverse_index_bits};
 
-/// Permutes `arr` such that each index is mapped to its reverse in binary.
+// TODO: Should really do some "dynamic" dispatch to handle the
-fn reverse_index_bits<T: Copy>(arr: Vec<T>) -> Vec<T> {
+// different FFT algos rather than C-style enum dispatch.
-    let n = arr.len();
+enum FftStrategy { Classic, Unrolled }
-    let n_power = log2_strict(n);
 
-    let mut result = Vec::with_capacity(n);
+const FFT_STRATEGY: FftStrategy = FftStrategy::Classic;
-    for i in 0..n {
+
-        result.push(arr[reverse_bits(i, n_power)]);
+type FftRootTable<F: Field> = Vec<Vec<F>>;
+
+fn fft_classic_root_table<F: Field>(n: usize) -> FftRootTable<F> {
+    let lg_n = log2_strict(n);
+    // bases[i] = g^2^i, for i = 0, ..., lg_n - 1
+    let mut bases = Vec::with_capacity(lg_n);
+    let mut base = F::primitive_root_of_unity(lg_n);
+    bases.push(base);
+    for _ in 1..lg_n {
+        base = base.square(); // base = g^2^_
+        bases.push(base);
     }
-    result
-}
 
-fn reverse_bits(n: usize, num_bits: usize) -> usize {
+    let mut root_table = Vec::with_capacity(lg_n);
-    let mut result = 0;
+    for lg_m in 1..=lg_n {
-    for i in 0..num_bits {
+        let half_m = 1 << (lg_m - 1);
-        let i_rev = num_bits - i - 1;
+        let base = bases[lg_n - lg_m];
-        result |= (n >> i & 1) << i_rev;
+        let root_row = base.powers().take(half_m.max(2)).collect();
+        root_table.push(root_row);
     }
-    result
+    root_table
 }
 
-pub(crate) struct FftPrecomputation<F: Field> {
+
-    /// For each layer index i, stores the cyclic subgroup corresponding to the evaluation domain of
+fn fft_unrolled_root_table<F: Field>(n: usize) -> FftRootTable<F> {
-    /// layer i. The indices within these subgroup vectors are bit-reversed.
+    // Precompute a table of the roots of unity used in the main
-    subgroups_rev: Vec<Vec<F>>,
+    // loops.
+
+    // Suppose n is the size of the outer vector and g is a primitive nth
+    // root of unity. Then the [lg(m) - 1][j] element of the table is
+    // g^{ n/2m * j } for j = 0..m-1
+
+    let lg_n = log2_strict(n);
+    // bases[i] = g^2^i, for i = 0, ..., lg_n - 2
+    let mut bases = Vec::with_capacity(lg_n);
+    let mut base = F::primitive_root_of_unity(lg_n);
+    bases.push(base);
+    // NB: If n = 1, then lg_n is zero, so we can't do 1..(lg_n-1) here
+    for _ in 2..lg_n {
+        base = base.square(); // base = g^2^(_-1)
+        bases.push(base);
+    }
+
+    let mut root_table = Vec::with_capacity(lg_n);
+    for lg_m in 1..lg_n {
+        let m = 1 << lg_m;
+        let base = bases[lg_n - lg_m - 1];
+        let root_row = base.powers().take(m.max(2)).collect();
+        root_table.push(root_row);
+    }
+    root_table
 }
 
-impl<F: Field> FftPrecomputation<F> {
+#[inline]
-    pub fn size(&self) -> usize {
+fn fft_dispatch<F: Field>(
-        self.subgroups_rev.last().unwrap().len()
+    input: Vec<F>,
+    zero_factor: Option<usize>,
+    root_table: Option<FftRootTable<F>>
+) -> Vec<F> {
+    let n = input.len();
+    match FFT_STRATEGY {
+        FftStrategy::Classic
+            => fft_classic(input,
+                           zero_factor.unwrap_or(0),
+                           root_table.unwrap_or_else(|| fft_classic_root_table(n))),
+        FftStrategy::Unrolled
+            => fft_unrolled(input,
+                            zero_factor.unwrap_or(0),
+                            root_table.unwrap_or_else(|| fft_unrolled_root_table(n)))
     }
 }
 
+#[inline]
 pub fn fft<F: Field>(poly: PolynomialCoeffs<F>) -> PolynomialValues<F> {
-    let precomputation = fft_precompute(poly.len());
+    fft_with_options(poly, None, None)
-    fft_with_precomputation_power_of_2(poly, &precomputation)
 }
 
-pub(crate) fn fft_precompute<F: Field>(degree: usize) -> FftPrecomputation<F> {
+#[inline]
-    let degree_log = log2_ceil(degree);
+pub fn fft_with_options<F: Field>(
-
+    poly: PolynomialCoeffs<F>,
-    let mut subgroups_rev = Vec::new();
+    zero_factor: Option<usize>,
-    let mut subgroup = F::two_adic_subgroup(degree_log);
+    root_table: Option<FftRootTable<F>>
-    for _i in 0..=degree_log {
+) -> PolynomialValues<F> {
-        let subsubgroup = subgroup.iter().step_by(2).copied().collect();
+    let PolynomialCoeffs { coeffs } = poly;
-        let subgroup_rev = reverse_index_bits(subgroup);
+    PolynomialValues { values: fft_dispatch(coeffs, zero_factor, root_table) }
-        subgroups_rev.push(subgroup_rev);
-        subgroup = subsubgroup;
-    }
-    subgroups_rev.reverse();
-
-    FftPrecomputation { subgroups_rev }
 }
 
-pub(crate) fn ifft_with_precomputation_power_of_2<F: Field>(
+#[inline]
+pub fn ifft<F: Field>(poly: PolynomialValues<F>) -> PolynomialCoeffs<F> {
+    ifft_with_options(poly, None, None)
+}
+
+pub fn ifft_with_options<F: Field>(
     poly: PolynomialValues<F>,
-    precomputation: &FftPrecomputation<F>,
+    zero_factor: Option<usize>,
+    root_table: Option<FftRootTable<F>>
 ) -> PolynomialCoeffs<F> {
     let n = poly.len();
-    let n_inv = F::from_canonical_usize(n).try_inverse().unwrap();
+    let lg_n = log2_strict(n);
+    let n_inv = F::inverse_2exp(lg_n);
 
     let PolynomialValues { values } = poly;
-    let PolynomialValues { values: mut result } =
+    let mut coeffs = fft_dispatch(values, zero_factor, root_table);
-        fft_with_precomputation_power_of_2(PolynomialCoeffs { coeffs: values }, precomputation);
 
     // We reverse all values except the first, and divide each by n.
-    result[0] *= n_inv;
+    coeffs[0] *= n_inv;
-    result[n / 2] *= n_inv;
+    coeffs[n / 2] *= n_inv;
     for i in 1..(n / 2) {
         let j = n - i;
-        let result_i = result[j] * n_inv;
+        let coeffs_i = coeffs[j] * n_inv;
-        let result_j = result[i] * n_inv;
+        let coeffs_j = coeffs[i] * n_inv;
-        result[i] = result_i;
+        coeffs[i] = coeffs_i;
-        result[j] = result_j;
+        coeffs[j] = coeffs_j;
     }
-    PolynomialCoeffs { coeffs: result }
+    PolynomialCoeffs { coeffs }
 }
 
-pub(crate) fn fft_with_precomputation_power_of_2<F: Field>(
+/// FFT implementation based on Section 32.3 of "Introduction to
-    poly: PolynomialCoeffs<F>,
+/// Algorithms" by Cormen et al.
-    precomputation: &FftPrecomputation<F>,
+///
-) -> PolynomialValues<F> {
+/// The parameter r signifies that the first 1/2^r of the entries of
-    debug_assert_eq!(
+/// input may be non-zero, but the last 1 - 1/2^r entries are
-        poly.len(),
+/// definitely zero.
-        precomputation.subgroups_rev.last().unwrap().len(),
+pub(crate) fn fft_classic<F: Field>(
-        "Number of coefficients does not match size of subgroup in precomputation"
+    input: Vec<F>,
-    );
+    r: usize,
+    root_table: FftRootTable<F>
+) -> Vec<F> {
+    let mut values = reverse_index_bits(input);
 
-    let half_degree = poly.len() >> 1;
+    let n = values.len();
-    let degree_log = poly.log_len();
+    let lg_n = log2_strict(n);
 
-    // In the base layer, we're just evaluating "degree 0 polynomials", i.e. the coefficients
+    if root_table.len() != lg_n {
-    // themselves.
+        panic!("Expected root table of length {}, but it was {}.", lg_n, root_table.len());
-    let PolynomialCoeffs { coeffs } = poly;
+    }
-    let mut evaluations = reverse_index_bits(coeffs);
 
-    for i in 1..=degree_log {
+    // After reverse_index_bits, the only non-zero elements of values
-        // In layer i, we're evaluating a series of polynomials, each at 2^i points. In practice
+    // are at indices i*2^r for i = 0..n/2^r.  The loop below copies
-        // we evaluate a pair of points together, so we have 2^(i - 1) pairs.
+    // the value at i*2^r to the positions [i*2^r + 1, i*2^r + 2, ...,
-        let points_per_poly = 1 << i;
+    // (i+1)*2^r - 1]; i.e. it replaces the 2^r - 1 zeros following
-        let pairs_per_poly = 1 << (i - 1);
+    // element i*2^r with the value at i*2^r.  This corresponds to the
-
+    // first r rounds of the FFT when there are 2^r zeros at the end
-        let mut new_evaluations = Vec::new();
+    // of the original input.
-        for pair_index in 0..half_degree {
+    if r > 0 { // if r == 0 then this loop is a noop.
-            let poly_index = pair_index / pairs_per_poly;
+        let mask = !((1 << r) - 1);
-            let pair_index_within_poly = pair_index % pairs_per_poly;
+        for i in 0..n {
-
+            values[i] = values[i & mask];
-            let child_index_0 = poly_index * points_per_poly + pair_index_within_poly;
-            let child_index_1 = child_index_0 + pairs_per_poly;
-
-            let even = evaluations[child_index_0];
-            let odd = evaluations[child_index_1];
-
-            let point_0 = precomputation.subgroups_rev[i][pair_index_within_poly * 2];
-            let product = point_0 * odd;
-            new_evaluations.push(even + product);
-            new_evaluations.push(even - product);
         }
-        evaluations = new_evaluations;
     }
 
-    // Reorder so that evaluations' indices correspond to (g_0, g_1, g_2, ...)
+    let mut m = 1 << (r + 1);
-    let values = reverse_index_bits(evaluations);
+    for lg_m in (r+1)..=lg_n {
-    PolynomialValues { values }
+        let half_m = m / 2;
+        for k in (0..n).step_by(m) {
+            for j in 0..half_m {
+                let omega = root_table[lg_m - 1][j];
+                let t = omega * values[k + half_m + j];
+                let u = values[k + j];
+                values[k + j] = u + t;
+                values[k + half_m + j] = u - t;
+            }
+        }
+        m *= 2;
+    }
+    values
 }
 
-pub(crate) fn ifft<F: Field>(poly: PolynomialValues<F>) -> PolynomialCoeffs<F> {
+/// FFT implementation inspired by Barretenberg's (but with extra unrolling):
-    let precomputation = fft_precompute(poly.len());
+/// https://github.com/AztecProtocol/barretenberg/blob/master/barretenberg/src/aztec/polynomials/polynomial_arithmetic.cpp#L58
-    ifft_with_precomputation_power_of_2(poly, &precomputation)
+/// https://github.com/AztecProtocol/barretenberg/blob/master/barretenberg/src/aztec/polynomials/evaluation_domain.cpp#L30
+///
+/// The parameter r signifies that the first 1/2^r of the entries of
+/// input may be non-zero, but the last 1 - 1/2^r entries are
+/// definitely zero.
+fn fft_unrolled<F: Field>(
+    input: Vec<F>,
+    r_orig: usize,
+    root_table: FftRootTable<F>
+) -> Vec<F> {
+    let n = input.len();
+    let lg_n = log2_strict(input.len());
+
+    let mut values = reverse_index_bits(input);
+
+    // FFT of a constant polynomial (including zero) is itself.
+    if n < 2 {
+        return values
+    }
+
+    // The 'm' corresponds to the specialisation from the 'm' in the
+    // main loop (m >= 4) below.
+
+    // (See comment in fft_classic near same code.)
+    let mut r = r_orig;
+    let mut m = 1 << r;
+    if r > 0 { // if r == 0 then this loop is a noop.
+        let mask = !((1 << r) - 1);
+        for i in 0..n {
+            values[i] = values[i & mask];
+        }
+    }
+
+    // m = 1
+    if m == 1 {
+        for k in (0..n).step_by(2) {
+            let t = values[k + 1];
+            values[k + 1] = values[k] - t;
+            values[k] += t;
+        }
+        r += 1;
+        m *= 2;
+    }
+
+    if n == 2 {
+        return values
+    }
+
+    if root_table.len() != (lg_n - 1) {
+        panic!("Expected root table of length {}, but it was {}.", lg_n, root_table.len());
+    }
+
+    // m = 2
+    if m <= 2 {
+        for k in (0..n).step_by(4) {
+            // NB: Grouping statements as is done in the main loop below
+            // does not seem to help here (worse by a few millis).
+            let omega_0 = root_table[0][0];
+            let tmp_0 = omega_0 * values[k + 2 + 0];
+            values[k + 2 + 0] = values[k + 0] - tmp_0;
+            values[k + 0] += tmp_0;
+
+            let omega_1 = root_table[0][1];
+            let tmp_1 = omega_1 * values[k + 2 + 1];
+            values[k + 2 + 1] = values[k + 1] - tmp_1;
+            values[k + 1] += tmp_1;
+        }
+        r += 1;
+        m *= 2;
+    }
+
+    // m >= 4
+    for lg_m in r..lg_n {
+        for k in (0..n).step_by(2*m) {
+            // Unrolled the commented loop by groups of 4 and
+            // rearranged the lines. Improves runtime by about
+            // 10%.
+            /*
+            for j in (0..m) {
+                let omega = root_table[lg_m - 1][j];
+                let tmp = omega * values[k + m + j];
+                values[k + m + j] = values[k + j] - tmp;
+                values[k + j] += tmp;
+            }
+            */
+            for j in (0..m).step_by(4) {
+                let off1 = k + j;
+                let off2 = k + m + j;
+
+                let omega_0 = root_table[lg_m - 1][j + 0];
+                let omega_1 = root_table[lg_m - 1][j + 1];
+                let omega_2 = root_table[lg_m - 1][j + 2];
+                let omega_3 = root_table[lg_m - 1][j + 3];
+
+                let tmp_0 = omega_0 * values[off2 + 0];
+                let tmp_1 = omega_1 * values[off2 + 1];
+                let tmp_2 = omega_2 * values[off2 + 2];
+                let tmp_3 = omega_3 * values[off2 + 3];
+
+                values[off2 + 0] = values[off1 + 0] - tmp_0;
+                values[off2 + 1] = values[off1 + 1] - tmp_1;
+                values[off2 + 2] = values[off1 + 2] - tmp_2;
+                values[off2 + 3] = values[off1 + 3] - tmp_3;
+                values[off1 + 0] += tmp_0;
+                values[off1 + 1] += tmp_1;
+                values[off1 + 2] += tmp_2;
+                values[off1 + 3] += tmp_3;
+            }
+        }
+        m *= 2;
+    }
+    values
 }
 
+
 #[cfg(test)]
 mod tests {
     use crate::field::crandall_field::CrandallField;
-    use crate::field::fft::{fft, ifft};
+    use crate::field::fft::{fft, ifft, fft_with_options};
     use crate::field::field::Field;
     use crate::polynomial::polynomial::{PolynomialCoeffs, PolynomialValues};
     use crate::util::{log2_ceil, log2_strict};
@@ -162,6 +324,12 @@ mod tests {
         for i in degree..degree_padded {
             assert_eq!(interpolated_coefficients.coeffs[i], F::ZERO);
         }
+
+        for r in 0..4 {
+            // expand ceofficients by factor 2^r by filling with zeros
+            let zero_tail = coefficients.clone().lde(r);
+            assert_eq!(fft(zero_tail.clone()), fft_with_options(zero_tail, Some(r), None));
+        }
     }
 
     fn evaluate_naive<F: Field>(coefficients: &PolynomialCoeffs<F>) -> PolynomialValues<F> {

src/field/field.rs

@@ -32,11 +32,14 @@ pub trait Field:
     + Send
     + Sync
 {
+    type PrimeField: Field;
+
     const ZERO: Self;
     const ONE: Self;
     const TWO: Self;
     const NEG_ONE: Self;
 
+    const CHARACTERISTIC: u64;
     const ORDER: u64;
     const TWO_ADICITY: usize;
 
@@ -101,6 +104,31 @@ pub trait Field:
         x_inv
     }
 
+    /// Compute the inverse of 2^exp in this field.
+    #[inline]
+    fn inverse_2exp(exp: usize) -> Self {
+        let p = Self::CHARACTERISTIC;
+
+        if exp <= Self::PrimeField::TWO_ADICITY {
+            // The inverse of 2^exp is p-(p-1)/2^exp when char(F) = p and exp is
+            // at most the TWO_ADICITY of the prime field.
+            //
+            // NB: PrimeFields fit in 64 bits => TWO_ADICITY < 64 =>
+            // exp < 64 => this shift amount is legal.
+            Self::from_canonical_u64(p - ((p - 1) >> exp))
+        } else {
+            // In the general case we compute 1/2 = (p+1)/2 and then exponentiate
+            // by exp to get 1/2^exp. Costs about log_2(exp) operations.
+            let half = Self::from_canonical_u64((p + 1) >> 1);
+            half.exp(exp as u64)
+
+            // TODO: Faster to combine several high powers of 1/2 using multiple
+            // applications of the trick above. E.g. if the 2-adicity is v, then
+            // compute 1/2^(v^2 + v + 13) with 1/2^((v + 1) * v + 13), etc.
+            // (using the v-adic expansion of m). Costs about log_v(exp) operations.
+        }
+    }
+
     fn primitive_root_of_unity(n_log: usize) -> Self {
         assert!(n_log <= Self::TWO_ADICITY);
         let mut base = Self::POWER_OF_TWO_GENERATOR;

src/field/field_testing.rs

@@ -315,6 +315,20 @@ macro_rules! test_arithmetic {
                 assert_eq!(x, F::ONE);
                 assert_eq!(F::ZERO - x, F::NEG_ONE);
             }
+
+            #[test]
+            fn inverse_2exp() {
+                // Just check consistency with try_inverse()
+                type F = $field;
+
+                let v = <F as Field>::PrimeField::TWO_ADICITY;
+
+                for e in [0, 1, 2, 3, 4, v - 2, v - 1, v, v + 1, v + 2, 123*v] {
+                    let x = F::TWO.exp(e as u64).inverse();
+                    let y = F::inverse_2exp(e);
+                    assert_eq!(x, y);
+                }
+            }
         }
     };
 }

src/gadgets/arithmetic.rs

@@ -1,15 +1,6 @@
-use std::convert::TryInto;
-
 use crate::circuit_builder::CircuitBuilder;
-use crate::field::extension_field::target::ExtensionTarget;
+use crate::field::extension_field::Extendable;
-use crate::field::extension_field::{Extendable, FieldExtension};
-use crate::field::field::Field;
-use crate::gates::arithmetic::ArithmeticGate;
-use crate::gates::mul_extension::MulExtensionGate;
-use crate::generator::SimpleGenerator;
 use crate::target::Target;
-use crate::wire::Wire;
-use crate::witness::PartialWitness;
 
 impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
     /// Computes `-x`.
@@ -43,30 +34,18 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         {
             return result;
         }
+        let multiplicand_0_ext = self.convert_to_ext(multiplicand_0);
+        let multiplicand_1_ext = self.convert_to_ext(multiplicand_1);
+        let addend_ext = self.convert_to_ext(addend);
 
-        let gate = self.add_gate(ArithmeticGate::new(), vec![const_0, const_1]);
+        self.arithmetic_extension(
-
+            const_0,
-        let wire_multiplicand_0 = Wire {
+            const_1,
-            gate,
+            multiplicand_0_ext,
-            input: ArithmeticGate::WIRE_MULTIPLICAND_0,
+            multiplicand_1_ext,
-        };
+            addend_ext,
-        let wire_multiplicand_1 = Wire {
+        )
-            gate,
+        .0[0]
-            input: ArithmeticGate::WIRE_MULTIPLICAND_1,
-        };
-        let wire_addend = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_ADDEND,
-        };
-        let wire_output = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_OUTPUT,
-        };
-
-        self.route(multiplicand_0, Target::Wire(wire_multiplicand_0));
-        self.route(multiplicand_1, Target::Wire(wire_multiplicand_1));
-        self.route(addend, Target::Wire(wire_addend));
-        Target::Wire(wire_output)
     }
 
     /// Checks for special cases where the value of
@@ -144,6 +123,7 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
         self.arithmetic(F::ONE, x, one, F::ONE, y)
     }
 
+    // TODO: Can be made `2*D` times more efficient by using all wires of an `ArithmeticExtensionGate`.
     pub fn add_many(&mut self, terms: &[Target]) -> Target {
         let mut sum = self.zero();
         for term in terms {
@@ -174,21 +154,31 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
     }
 
     // TODO: Optimize this, maybe with a new gate.
+    // TODO: Test
     /// Exponentiate `base` to the power of `exponent`, where `exponent < 2^num_bits`.
     pub fn exp(&mut self, base: Target, exponent: Target, num_bits: usize) -> Target {
         let mut current = base;
-        let one = self.one();
+        let one_ext = self.one_extension();
-        let mut product = one;
+        let mut product = self.one();
         let exponent_bits = self.split_le(exponent, num_bits);
 
         for bit in exponent_bits.into_iter() {
-            product = self.mul_many(&[bit, current, product]);
+            let current_ext = self.convert_to_ext(current);
+            let multiplicand = self.select(bit, current_ext, one_ext);
+            product = self.mul(product, multiplicand.0[0]);
             current = self.mul(current, current);
         }
 
         product
     }
 
+    /// Exponentiate `base` to the power of a known `exponent`.
+    // TODO: Test
+    pub fn exp_u64(&mut self, base: Target, exponent: u64) -> Target {
+        let base_ext = self.convert_to_ext(base);
+        self.exp_u64_extension(base_ext, exponent).0[0]
+    }
+
     /// Computes `q = x / y` by witnessing `q` and requiring that `q * y = x`. This can be unsafe in
     /// some cases, as it allows `0 / 0 = <anything>`.
     pub fn div_unsafe(&mut self, x: Target, y: Target) -> Target {
@@ -207,224 +197,8 @@ impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
             return self.constant(x_const / y_const);
         }
 
-        // Add an `ArithmeticGate` to compute `q * y`.
+        let x_ext = self.convert_to_ext(x);
-        let gate = self.add_gate(ArithmeticGate::new(), vec![F::ONE, F::ZERO]);
+        let y_ext = self.convert_to_ext(y);
-
+        self.div_unsafe_extension(x_ext, y_ext).0[0]
-        let wire_multiplicand_0 = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_MULTIPLICAND_0,
-        };
-        let wire_multiplicand_1 = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_MULTIPLICAND_1,
-        };
-        let wire_addend = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_ADDEND,
-        };
-        let wire_output = Wire {
-            gate,
-            input: ArithmeticGate::WIRE_OUTPUT,
-        };
-
-        let q = Target::Wire(wire_multiplicand_0);
-        self.add_generator(QuotientGenerator {
-            numerator: x,
-            denominator: y,
-            quotient: q,
-        });
-
-        self.route(y, Target::Wire(wire_multiplicand_1));
-
-        // This can be anything, since the whole second term has a weight of zero.
-        self.route(zero, Target::Wire(wire_addend));
-
-        let q_y = Target::Wire(wire_output);
-        self.assert_equal(q_y, x);
-
-        q
-    }
-
-    /// Computes `q = x / y` by witnessing `q` and requiring that `q * y = x`. This can be unsafe in
-    /// some cases, as it allows `0 / 0 = <anything>`.
-    pub fn div_unsafe_extension(
-        &mut self,
-        x: ExtensionTarget<D>,
-        y: ExtensionTarget<D>,
-    ) -> ExtensionTarget<D> {
-        // Add an `ArithmeticGate` to compute `q * y`.
-        let gate = self.add_gate(MulExtensionGate::new(), vec![F::ONE]);
-
-        let multiplicand_0 =
-            Target::wires_from_range(gate, MulExtensionGate::<D>::wires_multiplicand_0());
-        let multiplicand_0 = ExtensionTarget(multiplicand_0.try_into().unwrap());
-        let multiplicand_1 =
-            Target::wires_from_range(gate, MulExtensionGate::<D>::wires_multiplicand_1());
-        let multiplicand_1 = ExtensionTarget(multiplicand_1.try_into().unwrap());
-        let output = Target::wires_from_range(gate, MulExtensionGate::<D>::wires_output());
-        let output = ExtensionTarget(output.try_into().unwrap());
-
-        self.add_generator(QuotientGeneratorExtension {
-            numerator: x,
-            denominator: y,
-            quotient: multiplicand_0,
-        });
-
-        self.route_extension(y, multiplicand_1);
-
-        self.assert_equal_extension(output, x);
-
-        multiplicand_0
-    }
-}
-
-struct QuotientGenerator {
-    numerator: Target,
-    denominator: Target,
-    quotient: Target,
-}
-
-impl<F: Field> SimpleGenerator<F> for QuotientGenerator {
-    fn dependencies(&self) -> Vec<Target> {
-        vec![self.numerator, self.denominator]
-    }
-
-    fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
-        let num = witness.get_target(self.numerator);
-        let den = witness.get_target(self.denominator);
-        PartialWitness::singleton_target(self.quotient, num / den)
-    }
-}
-
-struct QuotientGeneratorExtension<const D: usize> {
-    numerator: ExtensionTarget<D>,
-    denominator: ExtensionTarget<D>,
-    quotient: ExtensionTarget<D>,
-}
-
-impl<F: Extendable<D>, const D: usize> SimpleGenerator<F> for QuotientGeneratorExtension<D> {
-    fn dependencies(&self) -> Vec<Target> {
-        let mut deps = self.numerator.to_target_array().to_vec();
-        deps.extend(&self.denominator.to_target_array());
-        deps
-    }
-
-    fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
-        let num = witness.get_extension_target(self.numerator);
-        let dem = witness.get_extension_target(self.denominator);
-        let quotient = num / dem;
-        let mut pw = PartialWitness::new();
-        for i in 0..D {
-            pw.set_target(
-                self.quotient.to_target_array()[i],
-                quotient.to_basefield_array()[i],
-            );
-        }
-        pw
-    }
-}
-
-/// An iterator over the powers of a certain base element `b`: `b^0, b^1, b^2, ...`.
-#[derive(Clone)]
-pub struct PowersTarget<const D: usize> {
-    base: ExtensionTarget<D>,
-    current: ExtensionTarget<D>,
-}
-
-impl<const D: usize> PowersTarget<D> {
-    pub fn next<F: Extendable<D>>(
-        &mut self,
-        builder: &mut CircuitBuilder<F, D>,
-    ) -> ExtensionTarget<D> {
-        let result = self.current;
-        self.current = builder.mul_extension(self.base, self.current);
-        result
-    }
-
-    pub fn repeated_frobenius<F: Extendable<D>>(
-        self,
-        k: usize,
-        builder: &mut CircuitBuilder<F, D>,
-    ) -> Self {
-        let Self { base, current } = self;
-        Self {
-            base: base.repeated_frobenius(k, builder),
-            current: current.repeated_frobenius(k, builder),
-        }
-    }
-}
-
-impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
-    pub fn powers(&mut self, base: ExtensionTarget<D>) -> PowersTarget<D> {
-        PowersTarget {
-            base,
-            current: self.one_extension(),
-        }
-    }
-}
-
-#[cfg(test)]
-mod tests {
-    use crate::circuit_builder::CircuitBuilder;
-    use crate::circuit_data::CircuitConfig;
-    use crate::field::crandall_field::CrandallField;
-    use crate::field::extension_field::quartic::QuarticCrandallField;
-    use crate::field::field::Field;
-    use crate::fri::FriConfig;
-    use crate::gates::arithmetic::ArithmeticGate;
-    use crate::target::Target;
-    use crate::witness::PartialWitness;
-
-    #[test]
-    fn test_div() {
-        type F = CrandallField;
-        type FF = QuarticCrandallField;
-        const D: usize = 4;
-
-        let config = CircuitConfig::large_config();
-
-        let mut builder = CircuitBuilder::<F, D>::new(config);
-
-        let x = F::rand();
-        let y = F::rand();
-        let mut pw = PartialWitness::new();
-        /// Computes x*x + 0*y = x^2.
-        let square_gate = builder.add_gate(ArithmeticGate::new(), vec![F::ONE, F::ZERO]);
-        pw.set_target(Target::wire(square_gate, 0), x);
-        pw.set_target(Target::wire(square_gate, 1), x);
-        let x2t = Target::wire(square_gate, ArithmeticGate::WIRE_OUTPUT);
-        let yt = Target::wire(square_gate, ArithmeticGate::WIRE_ADDEND);
-        pw.set_target(yt, y);
-        // Constant for x*x/y.
-        let zt = builder.constant(x * x / y);
-        // Computed division for x*x/y using the division gadget.
-        let comp_zt = builder.div_unsafe(x2t, yt);
-        builder.assert_equal(zt, comp_zt);
-
-        let data = builder.build();
-        let proof = data.prove(pw);
-    }
-
-    #[test]
-    fn test_div_extension() {
-        type F = CrandallField;
-        type FF = QuarticCrandallField;
-        const D: usize = 4;
-
-        let config = CircuitConfig::large_config();
-
-        let mut builder = CircuitBuilder::<F, D>::new(config);
-
-        let x = FF::rand();
-        let y = FF::rand();
-        let z = x / y;
-        let xt = builder.constant_extension(x);
-        let yt = builder.constant_extension(y);
-        let zt = builder.constant_extension(z);
-        let comp_zt = builder.div_unsafe_extension(xt, yt);
-        builder.assert_equal_extension(zt, comp_zt);
-
-        let data = builder.build();
-        let proof = data.prove(PartialWitness::new());
     }
 }

src/gadgets/arithmetic_extension.rs

@@ -0,0 +1,465 @@
+use std::convert::TryInto;
+use std::ops::Range;
+
+use itertools::Itertools;
+use num::Integer;
+
+use crate::circuit_builder::CircuitBuilder;
+use crate::field::extension_field::target::{ExtensionAlgebraTarget, ExtensionTarget};
+use crate::field::extension_field::{Extendable, OEF};
+use crate::field::field::Field;
+use crate::gates::arithmetic::ArithmeticExtensionGate;
+use crate::generator::SimpleGenerator;
+use crate::target::Target;
+use crate::util::bits_u64;
+use crate::witness::PartialWitness;
+
+impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
+    pub fn double_arithmetic_extension(
+        &mut self,
+        const_0: F,
+        const_1: F,
+        fixed_multiplicand: ExtensionTarget<D>,
+        multiplicand_0: ExtensionTarget<D>,
+        addend_0: ExtensionTarget<D>,
+        multiplicand_1: ExtensionTarget<D>,
+        addend_1: ExtensionTarget<D>,
+    ) -> (ExtensionTarget<D>, ExtensionTarget<D>) {
+        let gate = self.add_gate(ArithmeticExtensionGate::new(), vec![const_0, const_1]);
+
+        let wire_fixed_multiplicand = ExtensionTarget::from_range(
+            gate,
+            ArithmeticExtensionGate::<D>::wires_fixed_multiplicand(),
+        );
+        let wire_multiplicand_0 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_multiplicand_0());
+        let wire_addend_0 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_addend_0());
+        let wire_multiplicand_1 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_multiplicand_1());
+        let wire_addend_1 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_addend_1());
+        let wire_output_0 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_output_0());
+        let wire_output_1 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_output_1());
+
+        self.route_extension(fixed_multiplicand, wire_fixed_multiplicand);
+        self.route_extension(multiplicand_0, wire_multiplicand_0);
+        self.route_extension(addend_0, wire_addend_0);
+        self.route_extension(multiplicand_1, wire_multiplicand_1);
+        self.route_extension(addend_1, wire_addend_1);
+        (wire_output_0, wire_output_1)
+    }
+
+    pub fn arithmetic_extension(
+        &mut self,
+        const_0: F,
+        const_1: F,
+        multiplicand_0: ExtensionTarget<D>,
+        multiplicand_1: ExtensionTarget<D>,
+        addend: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let zero = self.zero_extension();
+        self.double_arithmetic_extension(
+            const_0,
+            const_1,
+            multiplicand_0,
+            multiplicand_1,
+            addend,
+            zero,
+            zero,
+        )
+        .0
+    }
+
+    pub fn add_extension(
+        &mut self,
+        a: ExtensionTarget<D>,
+        b: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let one = self.one_extension();
+        self.arithmetic_extension(F::ONE, F::ONE, one, a, b)
+    }
+
+    pub fn add_two_extension(
+        &mut self,
+        a0: ExtensionTarget<D>,
+        b0: ExtensionTarget<D>,
+        a1: ExtensionTarget<D>,
+        b1: ExtensionTarget<D>,
+    ) -> (ExtensionTarget<D>, ExtensionTarget<D>) {
+        let one = self.one_extension();
+        self.double_arithmetic_extension(F::ONE, F::ONE, one, a0, b0, a1, b1)
+    }
+
+    pub fn add_ext_algebra(
+        &mut self,
+        a: ExtensionAlgebraTarget<D>,
+        b: ExtensionAlgebraTarget<D>,
+    ) -> ExtensionAlgebraTarget<D> {
+        // We run two additions in parallel. So `[a0,a1,a2,a3] + [b0,b1,b2,b3]` is computed with two
+        // `add_two_extension`, first `[a0,a1]+[b0,b1]` then `[a2,a3]+[b2,b3]`.
+        let mut res = Vec::with_capacity(D);
+        // We need some extra logic if D is odd.
+        let d_even = D & (D ^ 1); // = 2 * (D/2)
+        for mut chunk in &(0..d_even).chunks(2) {
+            let i = chunk.next().unwrap();
+            let j = chunk.next().unwrap();
+            let (o0, o1) = self.add_two_extension(a.0[i], b.0[i], a.0[j], b.0[j]);
+            res.extend([o0, o1]);
+        }
+        if D.is_odd() {
+            res.push(self.add_extension(a.0[D - 1], b.0[D - 1]));
+        }
+        ExtensionAlgebraTarget(res.try_into().unwrap())
+    }
+
+    pub fn add_many_extension(&mut self, terms: &[ExtensionTarget<D>]) -> ExtensionTarget<D> {
+        let zero = self.zero_extension();
+        let mut terms = terms.to_vec();
+        if terms.len().is_odd() {
+            terms.push(zero);
+        }
+        // We maintain two accumulators, one for the sum of even elements, and one for odd elements.
+        let mut acc0 = zero;
+        let mut acc1 = zero;
+        for chunk in terms.chunks_exact(2) {
+            (acc0, acc1) = self.add_two_extension(acc0, chunk[0], acc1, chunk[1]);
+        }
+        // We sum both accumulators to get the final result.
+        self.add_extension(acc0, acc1)
+    }
+
+    pub fn sub_extension(
+        &mut self,
+        a: ExtensionTarget<D>,
+        b: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let one = self.one_extension();
+        self.arithmetic_extension(F::ONE, F::NEG_ONE, one, a, b)
+    }
+
+    pub fn sub_two_extension(
+        &mut self,
+        a0: ExtensionTarget<D>,
+        b0: ExtensionTarget<D>,
+        a1: ExtensionTarget<D>,
+        b1: ExtensionTarget<D>,
+    ) -> (ExtensionTarget<D>, ExtensionTarget<D>) {
+        let one = self.one_extension();
+        self.double_arithmetic_extension(F::ONE, F::NEG_ONE, one, a0, b0, a1, b1)
+    }
+
+    pub fn sub_ext_algebra(
+        &mut self,
+        a: ExtensionAlgebraTarget<D>,
+        b: ExtensionAlgebraTarget<D>,
+    ) -> ExtensionAlgebraTarget<D> {
+        // See `add_ext_algebra`.
+        let mut res = Vec::with_capacity(D);
+        let d_even = D & (D ^ 1); // = 2 * (D/2)
+        for mut chunk in &(0..d_even).chunks(2) {
+            let i = chunk.next().unwrap();
+            let j = chunk.next().unwrap();
+            let (o0, o1) = self.sub_two_extension(a.0[i], b.0[i], a.0[j], b.0[j]);
+            res.extend([o0, o1]);
+        }
+        if D.is_odd() {
+            res.push(self.sub_extension(a.0[D - 1], b.0[D - 1]));
+        }
+        ExtensionAlgebraTarget(res.try_into().unwrap())
+    }
+
+    pub fn mul_extension_with_const(
+        &mut self,
+        const_0: F,
+        multiplicand_0: ExtensionTarget<D>,
+        multiplicand_1: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let zero = self.zero_extension();
+        self.double_arithmetic_extension(
+            const_0,
+            F::ZERO,
+            multiplicand_0,
+            multiplicand_1,
+            zero,
+            zero,
+            zero,
+        )
+        .0
+    }
+
+    pub fn mul_extension(
+        &mut self,
+        multiplicand_0: ExtensionTarget<D>,
+        multiplicand_1: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        self.mul_extension_with_const(F::ONE, multiplicand_0, multiplicand_1)
+    }
+
+    /// Computes `x^2`.
+    pub fn square_extension(&mut self, x: ExtensionTarget<D>) -> ExtensionTarget<D> {
+        self.mul_extension(x, x)
+    }
+
+    pub fn mul_ext_algebra(
+        &mut self,
+        a: ExtensionAlgebraTarget<D>,
+        b: ExtensionAlgebraTarget<D>,
+    ) -> ExtensionAlgebraTarget<D> {
+        let mut res = [self.zero_extension(); D];
+        let w = self.constant(F::Extension::W);
+        for i in 0..D {
+            for j in 0..D {
+                res[(i + j) % D] = if i + j < D {
+                    self.mul_add_extension(a.0[i], b.0[j], res[(i + j) % D])
+                } else {
+                    let ai_bi = self.mul_extension(a.0[i], b.0[j]);
+                    self.scalar_mul_add_extension(w, ai_bi, res[(i + j) % D])
+                }
+            }
+        }
+        ExtensionAlgebraTarget(res)
+    }
+
+    pub fn mul_many_extension(&mut self, terms: &[ExtensionTarget<D>]) -> ExtensionTarget<D> {
+        let mut product = self.one_extension();
+        for term in terms {
+            product = self.mul_extension(product, *term);
+        }
+        product
+    }
+
+    /// Like `mul_add`, but for `ExtensionTarget`s.
+    pub fn mul_add_extension(
+        &mut self,
+        a: ExtensionTarget<D>,
+        b: ExtensionTarget<D>,
+        c: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        self.arithmetic_extension(F::ONE, F::ONE, a, b, c)
+    }
+
+    /// Like `mul_add`, but for `ExtensionTarget`s.
+    pub fn scalar_mul_add_extension(
+        &mut self,
+        a: Target,
+        b: ExtensionTarget<D>,
+        c: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let a_ext = self.convert_to_ext(a);
+        self.arithmetic_extension(F::ONE, F::ONE, a_ext, b, c)
+    }
+
+    /// Like `mul_sub`, but for `ExtensionTarget`s.
+    pub fn scalar_mul_sub_extension(
+        &mut self,
+        a: Target,
+        b: ExtensionTarget<D>,
+        c: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        let a_ext = self.convert_to_ext(a);
+        self.arithmetic_extension(F::ONE, F::NEG_ONE, a_ext, b, c)
+    }
+
+    /// Returns `a * b`, where `b` is in the extension field and `a` is in the base field.
+    pub fn scalar_mul_ext(&mut self, a: Target, b: ExtensionTarget<D>) -> ExtensionTarget<D> {
+        let a_ext = self.convert_to_ext(a);
+        self.mul_extension(a_ext, b)
+    }
+
+    /// Returns `a * b`, where `b` is in the extension of the extension field, and `a` is in the
+    /// extension field.
+    pub fn scalar_mul_ext_algebra(
+        &mut self,
+        a: ExtensionTarget<D>,
+        mut b: ExtensionAlgebraTarget<D>,
+    ) -> ExtensionAlgebraTarget<D> {
+        for i in 0..D {
+            b.0[i] = self.mul_extension(a, b.0[i]);
+        }
+        b
+    }
+
+    /// Exponentiate `base` to the power of a known `exponent`.
+    // TODO: Test
+    pub fn exp_u64_extension(
+        &mut self,
+        base: ExtensionTarget<D>,
+        exponent: u64,
+    ) -> ExtensionTarget<D> {
+        let mut current = base;
+        let mut product = self.one_extension();
+
+        for j in 0..bits_u64(exponent as u64) {
+            if (exponent >> j & 1) != 0 {
+                product = self.mul_extension(product, current);
+            }
+            current = self.square_extension(current);
+        }
+        product
+    }
+
+    /// Computes `q = x / y` by witnessing `q` and requiring that `q * y = x`. This can be unsafe in
+    /// some cases, as it allows `0 / 0 = <anything>`.
+    pub fn div_unsafe_extension(
+        &mut self,
+        x: ExtensionTarget<D>,
+        y: ExtensionTarget<D>,
+    ) -> ExtensionTarget<D> {
+        // Add an `ArithmeticExtensionGate` to compute `q * y`.
+        let gate = self.add_gate(ArithmeticExtensionGate::new(), vec![F::ONE, F::ZERO]);
+
+        let multiplicand_0 = ExtensionTarget::from_range(
+            gate,
+            ArithmeticExtensionGate::<D>::wires_fixed_multiplicand(),
+        );
+        let multiplicand_1 =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_multiplicand_0());
+        let output =
+            ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_output_0());
+
+        self.add_generator(QuotientGeneratorExtension {
+            numerator: x,
+            denominator: y,
+            quotient: multiplicand_0,
+        });
+        // We need to zero out the other wires for the `ArithmeticExtensionGenerator` to hit.
+        self.add_generator(ZeroOutGenerator {
+            gate_index: gate,
+            ranges: vec![
+                ArithmeticExtensionGate::<D>::wires_addend_0(),
+                ArithmeticExtensionGate::<D>::wires_multiplicand_1(),
+                ArithmeticExtensionGate::<D>::wires_addend_1(),
+            ],
+        });
+
+        self.route_extension(y, multiplicand_1);
+        self.assert_equal_extension(output, x);
+
+        multiplicand_0
+    }
+}
+
+/// Generator used to zero out wires at a given gate index and ranges.
+pub struct ZeroOutGenerator {
+    gate_index: usize,
+    ranges: Vec<Range<usize>>,
+}
+
+impl<F: Field> SimpleGenerator<F> for ZeroOutGenerator {
+    fn dependencies(&self) -> Vec<Target> {
+        Vec::new()
+    }
+
+    fn run_once(&self, _witness: &PartialWitness<F>) -> PartialWitness<F> {
+        let mut pw = PartialWitness::new();
+        for t in self
+            .ranges
+            .iter()
+            .flat_map(|r| Target::wires_from_range(self.gate_index, r.clone()))
+        {
+            pw.set_target(t, F::ZERO);
+        }
+
+        pw
+    }
+}
+
+struct QuotientGeneratorExtension<const D: usize> {
+    numerator: ExtensionTarget<D>,
+    denominator: ExtensionTarget<D>,
+    quotient: ExtensionTarget<D>,
+}
+
+impl<F: Extendable<D>, const D: usize> SimpleGenerator<F> for QuotientGeneratorExtension<D> {
+    fn dependencies(&self) -> Vec<Target> {
+        let mut deps = self.numerator.to_target_array().to_vec();
+        deps.extend(&self.denominator.to_target_array());
+        deps
+    }
+
+    fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
+        let num = witness.get_extension_target(self.numerator);
+        let dem = witness.get_extension_target(self.denominator);
+        let quotient = num / dem;
+        let mut pw = PartialWitness::new();
+        pw.set_extension_target(self.quotient, quotient);
+
+        pw
+    }
+}
+
+/// An iterator over the powers of a certain base element `b`: `b^0, b^1, b^2, ...`.
+#[derive(Clone)]
+pub struct PowersTarget<const D: usize> {
+    base: ExtensionTarget<D>,
+    current: ExtensionTarget<D>,
+}
+
+impl<const D: usize> PowersTarget<D> {
+    pub fn next<F: Extendable<D>>(
+        &mut self,
+        builder: &mut CircuitBuilder<F, D>,
+    ) -> ExtensionTarget<D> {
+        let result = self.current;
+        self.current = builder.mul_extension(self.base, self.current);
+        result
+    }
+
+    pub fn repeated_frobenius<F: Extendable<D>>(
+        self,
+        k: usize,
+        builder: &mut CircuitBuilder<F, D>,
+    ) -> Self {
+        let Self { base, current } = self;
+        Self {
+            base: base.repeated_frobenius(k, builder),
+            current: current.repeated_frobenius(k, builder),
+        }
+    }
+}
+
+impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
+    pub fn powers(&mut self, base: ExtensionTarget<D>) -> PowersTarget<D> {
+        PowersTarget {
+            base,
+            current: self.one_extension(),
+        }
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use crate::circuit_builder::CircuitBuilder;
+    use crate::circuit_data::CircuitConfig;
+    use crate::field::crandall_field::CrandallField;
+    use crate::field::extension_field::quartic::QuarticCrandallField;
+    use crate::field::field::Field;
+    use crate::fri::FriConfig;
+    use crate::witness::PartialWitness;
+
+    #[test]
+    fn test_div_extension() {
+        type F = CrandallField;
+        type FF = QuarticCrandallField;
+        const D: usize = 4;
+
+        let config = CircuitConfig::large_config();
+
+        let mut builder = CircuitBuilder::<F, D>::new(config);
+
+        let x = FF::rand();
+        let y = FF::rand();
+        let z = x / y;
+        let xt = builder.constant_extension(x);
+        let yt = builder.constant_extension(y);
+        let zt = builder.constant_extension(z);
+        let comp_zt = builder.div_unsafe_extension(xt, yt);
+        builder.assert_equal_extension(zt, comp_zt);
+
+        let data = builder.build();
+        let proof = data.prove(PartialWitness::new());
+    }
+}

src/gadgets/interpolation.rs

@@ -72,13 +72,10 @@ mod tests {
     fn test_interpolate() {
         type F = CrandallField;
         type FF = QuarticCrandallField;
-        let config = CircuitConfig {
+        let config = CircuitConfig::large_config();
-            num_routed_wires: 18,
-            ..CircuitConfig::large_config()
-        };
         let mut builder = CircuitBuilder::<F, 4>::new(config);
 
-        let len = 2;
+        let len = 4;
         let points = (0..len)
             .map(|_| (F::rand(), FF::rand()))
             .collect::<Vec<_>>();

src/gadgets/mod.rs

@@ -1,4 +1,5 @@
 pub mod arithmetic;
+pub mod arithmetic_extension;
 pub mod hash;
 pub mod insert;
 pub mod interpolation;

src/gadgets/split_join.rs

@@ -9,14 +9,14 @@ use crate::wire::Wire;
 use crate::witness::PartialWitness;
 
 impl<F: Extendable<D>, const D: usize> CircuitBuilder<F, D> {
-    /// Split the given integer into a list of virtual advice targets, where each one represents a
+    /// Split the given integer into a list of virtual targets, where each one represents a bit of
-    /// bit of the integer, with little-endian ordering.
+    /// the integer, with little-endian ordering.
     ///
     /// Note that this only handles witness generation; it does not enforce that the decomposition
     /// is correct. The output should be treated as a "purported" decomposition which must be
     /// enforced elsewhere.
     pub(crate) fn split_le_virtual(&mut self, integer: Target, num_bits: usize) -> Vec<Target> {
-        let bit_targets = self.add_virtual_advice_targets(num_bits);
+        let bit_targets = self.add_virtual_targets(num_bits);
         self.add_generator(SplitGenerator {
             integer,
             bits: bit_targets.clone(),

src/gates/arithmetic.rs

@@ -1,34 +1,47 @@
+use std::ops::Range;
+
 use crate::circuit_builder::CircuitBuilder;
 use crate::field::extension_field::target::ExtensionTarget;
 use crate::field::extension_field::Extendable;
-use crate::field::field::Field;
 use crate::gates::gate::{Gate, GateRef};
 use crate::generator::{SimpleGenerator, WitnessGenerator};
 use crate::target::Target;
 use crate::vars::{EvaluationTargets, EvaluationVars};
-use crate::wire::Wire;
 use crate::witness::PartialWitness;
 
-/// A gate which can be configured to perform various arithmetic. In particular, it computes
+/// A gate which can a linear combination `c0*x*y+c1*z` twice with the same `x`.
-///
-/// ```text
-/// output := const_0 * multiplicand_0 * multiplicand_1 + const_1 * addend
-/// ```
 #[derive(Debug)]
-pub struct ArithmeticGate;
+pub struct ArithmeticExtensionGate<const D: usize>;
 
-impl ArithmeticGate {
+impl<const D: usize> ArithmeticExtensionGate<D> {
-    pub fn new<F: Extendable<D>, const D: usize>() -> GateRef<F, D> {
+    pub fn new<F: Extendable<D>>() -> GateRef<F, D> {
-        GateRef::new(ArithmeticGate)
+        GateRef::new(ArithmeticExtensionGate)
     }
 
-    pub const WIRE_MULTIPLICAND_0: usize = 0;
+    pub fn wires_fixed_multiplicand() -> Range<usize> {
-    pub const WIRE_MULTIPLICAND_1: usize = 1;
+        0..D
-    pub const WIRE_ADDEND: usize = 2;
+    }
-    pub const WIRE_OUTPUT: usize = 3;
+    pub fn wires_multiplicand_0() -> Range<usize> {
+        D..2 * D
+    }
+    pub fn wires_addend_0() -> Range<usize> {
+        2 * D..3 * D
+    }
+    pub fn wires_multiplicand_1() -> Range<usize> {
+        3 * D..4 * D
+    }
+    pub fn wires_addend_1() -> Range<usize> {
+        4 * D..5 * D
+    }
+    pub fn wires_output_0() -> Range<usize> {
+        5 * D..6 * D
+    }
+    pub fn wires_output_1() -> Range<usize> {
+        6 * D..7 * D
+    }
 }
 
-impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticGate {
+impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticExtensionGate<D> {
     fn id(&self) -> String {
         format!("{:?}", self)
     }
@@ -36,12 +49,23 @@ impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticGate {
     fn eval_unfiltered(&self, vars: EvaluationVars<F, D>) -> Vec<F::Extension> {
         let const_0 = vars.local_constants[0];
         let const_1 = vars.local_constants[1];
-        let multiplicand_0 = vars.local_wires[Self::WIRE_MULTIPLICAND_0];
+
-        let multiplicand_1 = vars.local_wires[Self::WIRE_MULTIPLICAND_1];
+        let fixed_multiplicand = vars.get_local_ext_algebra(Self::wires_fixed_multiplicand());
-        let addend = vars.local_wires[Self::WIRE_ADDEND];
+        let multiplicand_0 = vars.get_local_ext_algebra(Self::wires_multiplicand_0());
-        let output = vars.local_wires[Self::WIRE_OUTPUT];
+        let addend_0 = vars.get_local_ext_algebra(Self::wires_addend_0());
-        let computed_output = const_0 * multiplicand_0 * multiplicand_1 + const_1 * addend;
+        let multiplicand_1 = vars.get_local_ext_algebra(Self::wires_multiplicand_1());
-        vec![computed_output - output]
+        let addend_1 = vars.get_local_ext_algebra(Self::wires_addend_1());
+        let output_0 = vars.get_local_ext_algebra(Self::wires_output_0());
+        let output_1 = vars.get_local_ext_algebra(Self::wires_output_1());
+
+        let computed_output_0 =
+            fixed_multiplicand * multiplicand_0 * const_0.into() + addend_0 * const_1.into();
+        let computed_output_1 =
+            fixed_multiplicand * multiplicand_1 * const_0.into() + addend_1 * const_1.into();
+
+        let mut constraints = (output_0 - computed_output_0).to_basefield_array().to_vec();
+        constraints.extend((output_1 - computed_output_1).to_basefield_array());
+        constraints
     }
 
     fn eval_unfiltered_recursively(
@@ -51,15 +75,30 @@ impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticGate {
     ) -> Vec<ExtensionTarget<D>> {
         let const_0 = vars.local_constants[0];
         let const_1 = vars.local_constants[1];
-        let multiplicand_0 = vars.local_wires[Self::WIRE_MULTIPLICAND_0];
-        let multiplicand_1 = vars.local_wires[Self::WIRE_MULTIPLICAND_1];
-        let addend = vars.local_wires[Self::WIRE_ADDEND];
-        let output = vars.local_wires[Self::WIRE_OUTPUT];
 
-        let product_term = builder.mul_many_extension(&[const_0, multiplicand_0, multiplicand_1]);
+        let fixed_multiplicand = vars.get_local_ext_algebra(Self::wires_fixed_multiplicand());
-        let addend_term = builder.mul_extension(const_1, addend);
+        let multiplicand_0 = vars.get_local_ext_algebra(Self::wires_multiplicand_0());
-        let computed_output = builder.add_many_extension(&[product_term, addend_term]);
+        let addend_0 = vars.get_local_ext_algebra(Self::wires_addend_0());
-        vec![builder.sub_extension(computed_output, output)]
+        let multiplicand_1 = vars.get_local_ext_algebra(Self::wires_multiplicand_1());
+        let addend_1 = vars.get_local_ext_algebra(Self::wires_addend_1());
+        let output_0 = vars.get_local_ext_algebra(Self::wires_output_0());
+        let output_1 = vars.get_local_ext_algebra(Self::wires_output_1());
+
+        let computed_output_0 = builder.mul_ext_algebra(fixed_multiplicand, multiplicand_0);
+        let computed_output_0 = builder.scalar_mul_ext_algebra(const_0, computed_output_0);
+        let scaled_addend_0 = builder.scalar_mul_ext_algebra(const_1, addend_0);
+        let computed_output_0 = builder.add_ext_algebra(computed_output_0, scaled_addend_0);
+
+        let computed_output_1 = builder.mul_ext_algebra(fixed_multiplicand, multiplicand_1);
+        let computed_output_1 = builder.scalar_mul_ext_algebra(const_0, computed_output_1);
+        let scaled_addend_1 = builder.scalar_mul_ext_algebra(const_1, addend_1);
+        let computed_output_1 = builder.add_ext_algebra(computed_output_1, scaled_addend_1);
+
+        let diff_0 = builder.sub_ext_algebra(output_0, computed_output_0);
+        let diff_1 = builder.sub_ext_algebra(output_1, computed_output_1);
+        let mut constraints = diff_0.to_ext_target_array().to_vec();
+        constraints.extend(diff_1.to_ext_target_array());
+        constraints
     }
 
     fn generators(
@@ -67,16 +106,21 @@ impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticGate {
         gate_index: usize,
         local_constants: &[F],
     ) -> Vec<Box<dyn WitnessGenerator<F>>> {
-        let gen = ArithmeticGenerator {
+        let gen0 = ArithmeticExtensionGenerator0 {
             gate_index,
             const_0: local_constants[0],
             const_1: local_constants[1],
         };
-        vec![Box::new(gen)]
+        let gen1 = ArithmeticExtensionGenerator1 {
+            gate_index,
+            const_0: local_constants[0],
+            const_1: local_constants[1],
+        };
+        vec![Box::new(gen0), Box::new(gen1)]
     }
 
     fn num_wires(&self) -> usize {
-        4
+        7 * D
     }
 
     fn num_constants(&self) -> usize {
@@ -88,70 +132,96 @@ impl<F: Extendable<D>, const D: usize> Gate<F, D> for ArithmeticGate {
     }
 
     fn num_constraints(&self) -> usize {
-        1
+        2 * D
     }
 }
 
-struct ArithmeticGenerator<F: Field> {
+struct ArithmeticExtensionGenerator0<F: Extendable<D>, const D: usize> {
     gate_index: usize,
     const_0: F,
     const_1: F,
 }
 
-impl<F: Field> SimpleGenerator<F> for ArithmeticGenerator<F> {
+struct ArithmeticExtensionGenerator1<F: Extendable<D>, const D: usize> {
+    gate_index: usize,
+    const_0: F,
+    const_1: F,
+}
+
+impl<F: Extendable<D>, const D: usize> SimpleGenerator<F> for ArithmeticExtensionGenerator0<F, D> {
     fn dependencies(&self) -> Vec<Target> {
-        vec![
+        ArithmeticExtensionGate::<D>::wires_fixed_multiplicand()
-            Target::Wire(Wire {
+            .chain(ArithmeticExtensionGate::<D>::wires_multiplicand_0())
-                gate: self.gate_index,
+            .chain(ArithmeticExtensionGate::<D>::wires_addend_0())
-                input: ArithmeticGate::WIRE_MULTIPLICAND_0,
+            .map(|i| Target::wire(self.gate_index, i))
-            }),
+            .collect()
-            Target::Wire(Wire {
-                gate: self.gate_index,
-                input: ArithmeticGate::WIRE_MULTIPLICAND_1,
-            }),
-            Target::Wire(Wire {
-                gate: self.gate_index,
-                input: ArithmeticGate::WIRE_ADDEND,
-            }),
-        ]
     }
 
     fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
-        let multiplicand_0_target = Wire {
+        let extract_extension = |range: Range<usize>| -> F::Extension {
-            gate: self.gate_index,
+            let t = ExtensionTarget::from_range(self.gate_index, range);
-            input: ArithmeticGate::WIRE_MULTIPLICAND_0,
+            witness.get_extension_target(t)
-        };
-        let multiplicand_1_target = Wire {
-            gate: self.gate_index,
-            input: ArithmeticGate::WIRE_MULTIPLICAND_1,
-        };
-        let addend_target = Wire {
-            gate: self.gate_index,
-            input: ArithmeticGate::WIRE_ADDEND,
-        };
-        let output_target = Wire {
-            gate: self.gate_index,
-            input: ArithmeticGate::WIRE_OUTPUT,
         };
 
-        let multiplicand_0 = witness.get_wire(multiplicand_0_target);
+        let fixed_multiplicand =
-        let multiplicand_1 = witness.get_wire(multiplicand_1_target);
+            extract_extension(ArithmeticExtensionGate::<D>::wires_fixed_multiplicand());
-        let addend = witness.get_wire(addend_target);
+        let multiplicand_0 =
+            extract_extension(ArithmeticExtensionGate::<D>::wires_multiplicand_0());
+        let addend_0 = extract_extension(ArithmeticExtensionGate::<D>::wires_addend_0());
 
-        let output = self.const_0 * multiplicand_0 * multiplicand_1 + self.const_1 * addend;
+        let output_target_0 = ExtensionTarget::from_range(
+            self.gate_index,
+            ArithmeticExtensionGate::<D>::wires_output_0(),
+        );
 
-        PartialWitness::singleton_wire(output_target, output)
+        let computed_output_0 = fixed_multiplicand * multiplicand_0 * self.const_0.into()
+            + addend_0 * self.const_1.into();
+
+        PartialWitness::singleton_extension_target(output_target_0, computed_output_0)
+    }
+}
+
+impl<F: Extendable<D>, const D: usize> SimpleGenerator<F> for ArithmeticExtensionGenerator1<F, D> {
+    fn dependencies(&self) -> Vec<Target> {
+        ArithmeticExtensionGate::<D>::wires_fixed_multiplicand()
+            .chain(ArithmeticExtensionGate::<D>::wires_multiplicand_1())
+            .chain(ArithmeticExtensionGate::<D>::wires_addend_1())
+            .map(|i| Target::wire(self.gate_index, i))
+            .collect()
+    }
+
+    fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
+        let extract_extension = |range: Range<usize>| -> F::Extension {
+            let t = ExtensionTarget::from_range(self.gate_index, range);
+            witness.get_extension_target(t)
+        };
+
+        let fixed_multiplicand =
+            extract_extension(ArithmeticExtensionGate::<D>::wires_fixed_multiplicand());
+        let multiplicand_1 =
+            extract_extension(ArithmeticExtensionGate::<D>::wires_multiplicand_1());
+        let addend_1 = extract_extension(ArithmeticExtensionGate::<D>::wires_addend_1());
+
+        let output_target_1 = ExtensionTarget::from_range(
+            self.gate_index,
+            ArithmeticExtensionGate::<D>::wires_output_1(),
+        );
+
+        let computed_output_1 = fixed_multiplicand * multiplicand_1 * self.const_0.into()
+            + addend_1 * self.const_1.into();
+
+        PartialWitness::singleton_extension_target(output_target_1, computed_output_1)
     }
 }
 
 #[cfg(test)]
 mod tests {
     use crate::field::crandall_field::CrandallField;
-    use crate::gates::arithmetic::ArithmeticGate;
+    use crate::gates::arithmetic::ArithmeticExtensionGate;
     use crate::gates::gate_testing::test_low_degree;
 
     #[test]
     fn low_degree() {
-        test_low_degree(ArithmeticGate::new::<CrandallField, 4>())
+        test_low_degree(ArithmeticExtensionGate::<4>::new::<CrandallField>())
     }
 }

src/gates/gate_tree.rs

@@ -57,8 +57,9 @@ impl<F: Extendable<D>, const D: usize> Tree<GateRef<F, D>> {
     /// For this construction, we use the greedy algorithm in `Self::find_tree`.
     /// This latter function greedily adds gates at the depth where
     /// `filtered_deg(gate)=D, constant_wires(gate)=C` to ensure no space is wasted.
-    /// We return the first tree found in this manner.
+    /// We return the first tree found in this manner, along with it's maximum filtered degree
-    pub fn from_gates(mut gates: Vec<GateRef<F, D>>) -> Self {
+    /// and the number of constant wires needed when using this tree.
+    pub fn from_gates(mut gates: Vec<GateRef<F, D>>) -> (Self, usize, usize) {
         let timer = std::time::Instant::now();
         gates.sort_unstable_by_key(|g| (-(g.0.degree() as isize), -(g.0.num_constants() as isize)));
 
@@ -67,14 +68,30 @@ impl<F: Extendable<D>, const D: usize> Tree<GateRef<F, D>> {
             // So we can restrict our search space by setting `max_degree` to a power of 2.
             let max_degree = 1 << max_degree_bits;
             for max_constants in 1..100 {
-                if let Some(mut tree) = Self::find_tree(&gates, max_degree, max_constants) {
+                if let Some(mut best_tree) = Self::find_tree(&gates, max_degree, max_constants) {
-                    tree.shorten();
+                    let mut best_num_constants = best_tree.num_constants();
+                    let mut best_degree = max_degree;
+                    // Iterate backwards from `max_degree` to try to find a tree with a lower degree
+                    // but the same number of constants.
+                    'optdegree: for degree in (0..max_degree).rev() {
+                        if let Some(mut tree) = Self::find_tree(&gates, degree, max_constants) {
+                            let num_constants = tree.num_constants();
+                            if num_constants > best_num_constants {
+                                break 'optdegree;
+                            } else {
+                                best_degree = degree;
+                                best_num_constants = num_constants;
+                                best_tree = tree;
+                            }
+                        }
+                    }
                     info!(
-                        "Found tree with max degree {} in {}s.",
+                        "Found tree with max degree {} and {} constants wires in {:.4}s.",
-                        max_degree,
+                        best_degree,
+                        best_num_constants,
                         timer.elapsed().as_secs_f32()
                     );
-                    return tree;
+                    return (best_tree, best_degree, best_num_constants);
                 }
             }
         }
@@ -89,6 +106,7 @@ impl<F: Extendable<D>, const D: usize> Tree<GateRef<F, D>> {
         for g in gates {
             tree.try_add_gate(g, max_degree, max_constants)?;
         }
+        tree.shorten();
         Some(tree)
     }
 
@@ -180,18 +198,35 @@ impl<F: Extendable<D>, const D: usize> Tree<GateRef<F, D>> {
             }
         }
     }
+
+    /// Returns the tree's maximum filtered constraint degree.
+    fn max_filtered_degree(&self) -> usize {
+        self.traversal()
+            .into_iter()
+            .map(|(g, p)| g.0.degree() + p.len())
+            .max()
+            .expect("Empty tree.")
+    }
+
+    /// Returns the number of constant wires needed to fit all prefixes and gate constants.
+    fn num_constants(&self) -> usize {
+        self.traversal()
+            .into_iter()
+            .map(|(g, p)| g.0.num_constants() + p.len())
+            .max()
+            .expect("Empty tree.")
+    }
 }
 
 #[cfg(test)]
 mod tests {
     use super::*;
     use crate::field::crandall_field::CrandallField;
-    use crate::gates::arithmetic::ArithmeticGate;
+    use crate::gates::arithmetic::ArithmeticExtensionGate;
     use crate::gates::base_sum::BaseSumGate;
     use crate::gates::constant::ConstantGate;
     use crate::gates::gmimc::GMiMCGate;
     use crate::gates::interpolation::InterpolationGate;
-    use crate::gates::mul_extension::MulExtensionGate;
     use crate::gates::noop::NoopGate;
     use crate::hash::GMIMC_ROUNDS;
 
@@ -204,15 +239,14 @@ mod tests {
         let gates = vec![
             NoopGate::get::<F, D>(),
             ConstantGate::get(),
-            ArithmeticGate::new(),
+            ArithmeticExtensionGate::new(),
             BaseSumGate::<4>::new(4),
             GMiMCGate::<F, D, GMIMC_ROUNDS>::with_automatic_constants(),
             InterpolationGate::new(4),
-            MulExtensionGate::new(),
         ];
         let len = gates.len();
 
-        let tree = Tree::from_gates(gates.clone());
+        let (tree, _, _) = Tree::from_gates(gates.clone());
         let mut gates_with_prefix = tree.traversal();
         for (g, p) in &gates_with_prefix {
             info!(

src/gates/mod.rs

@@ -1,11 +1,10 @@
-pub(crate) mod arithmetic;
+pub mod arithmetic;
 pub mod base_sum;
 pub mod constant;
 pub(crate) mod gate;
 pub mod gate_tree;
 pub mod gmimc;
 pub mod interpolation;
-pub mod mul_extension;
 pub(crate) mod noop;
 
 #[cfg(test)]

src/gates/mul_extension.rs

@@ -1,145 +0,0 @@
-use std::ops::Range;
-
-use crate::circuit_builder::CircuitBuilder;
-use crate::field::extension_field::target::ExtensionTarget;
-use crate::field::extension_field::{Extendable, FieldExtension};
-use crate::gates::gate::{Gate, GateRef};
-use crate::generator::{SimpleGenerator, WitnessGenerator};
-use crate::target::Target;
-use crate::vars::{EvaluationTargets, EvaluationVars};
-use crate::wire::Wire;
-use crate::witness::PartialWitness;
-
-/// A gate which can multiply two field extension elements.
-/// TODO: Add an addend if `NUM_ROUTED_WIRES` is large enough.
-#[derive(Debug)]
-pub struct MulExtensionGate<const D: usize>;
-
-impl<const D: usize> MulExtensionGate<D> {
-    pub fn new<F: Extendable<D>>() -> GateRef<F, D> {
-        GateRef::new(MulExtensionGate)
-    }
-
-    pub fn wires_multiplicand_0() -> Range<usize> {
-        0..D
-    }
-    pub fn wires_multiplicand_1() -> Range<usize> {
-        D..2 * D
-    }
-    pub fn wires_output() -> Range<usize> {
-        2 * D..3 * D
-    }
-}
-
-impl<F: Extendable<D>, const D: usize> Gate<F, D> for MulExtensionGate<D> {
-    fn id(&self) -> String {
-        format!("{:?}", self)
-    }
-
-    fn eval_unfiltered(&self, vars: EvaluationVars<F, D>) -> Vec<F::Extension> {
-        let const_0 = vars.local_constants[0];
-        let multiplicand_0 = vars.get_local_ext_algebra(Self::wires_multiplicand_0());
-        let multiplicand_1 = vars.get_local_ext_algebra(Self::wires_multiplicand_1());
-        let output = vars.get_local_ext_algebra(Self::wires_output());
-        let computed_output = multiplicand_0 * multiplicand_1 * const_0.into();
-        (output - computed_output).to_basefield_array().to_vec()
-    }
-
-    fn eval_unfiltered_recursively(
-        &self,
-        builder: &mut CircuitBuilder<F, D>,
-        vars: EvaluationTargets<D>,
-    ) -> Vec<ExtensionTarget<D>> {
-        let const_0 = vars.local_constants[0];
-        let multiplicand_0 = vars.get_local_ext_algebra(Self::wires_multiplicand_0());
-        let multiplicand_1 = vars.get_local_ext_algebra(Self::wires_multiplicand_1());
-        let output = vars.get_local_ext_algebra(Self::wires_output());
-        let computed_output = builder.mul_ext_algebra(multiplicand_0, multiplicand_1);
-        let computed_output = builder.scalar_mul_ext_algebra(const_0, computed_output);
-        let diff = builder.sub_ext_algebra(output, computed_output);
-        diff.to_ext_target_array().to_vec()
-    }
-
-    fn generators(
-        &self,
-        gate_index: usize,
-        local_constants: &[F],
-    ) -> Vec<Box<dyn WitnessGenerator<F>>> {
-        let gen = MulExtensionGenerator {
-            gate_index,
-            const_0: local_constants[0],
-        };
-        vec![Box::new(gen)]
-    }
-
-    fn num_wires(&self) -> usize {
-        12
-    }
-
-    fn num_constants(&self) -> usize {
-        1
-    }
-
-    fn degree(&self) -> usize {
-        3
-    }
-
-    fn num_constraints(&self) -> usize {
-        D
-    }
-}
-
-struct MulExtensionGenerator<F: Extendable<D>, const D: usize> {
-    gate_index: usize,
-    const_0: F,
-}
-
-impl<F: Extendable<D>, const D: usize> SimpleGenerator<F> for MulExtensionGenerator<F, D> {
-    fn dependencies(&self) -> Vec<Target> {
-        MulExtensionGate::<D>::wires_multiplicand_0()
-            .chain(MulExtensionGate::<D>::wires_multiplicand_1())
-            .map(|i| {
-                Target::Wire(Wire {
-                    gate: self.gate_index,
-                    input: i,
-                })
-            })
-            .collect()
-    }
-
-    fn run_once(&self, witness: &PartialWitness<F>) -> PartialWitness<F> {
-        let multiplicand_0_target = ExtensionTarget::from_range(
-            self.gate_index,
-            MulExtensionGate::<D>::wires_multiplicand_0(),
-        );
-        let multiplicand_0 = witness.get_extension_target(multiplicand_0_target);
-
-        let multiplicand_1_target = ExtensionTarget::from_range(
-            self.gate_index,
-            MulExtensionGate::<D>::wires_multiplicand_1(),
-        );
-        let multiplicand_1 = witness.get_extension_target(multiplicand_1_target);
-
-        let output_target =
-            ExtensionTarget::from_range(self.gate_index, MulExtensionGate::<D>::wires_output());
-
-        let computed_output =
-            F::Extension::from_basefield(self.const_0) * multiplicand_0 * multiplicand_1;
-
-        let mut pw = PartialWitness::new();
-        pw.set_extension_target(output_target, computed_output);
-        pw
-    }
-}
-
-#[cfg(test)]
-mod tests {
-    use crate::field::crandall_field::CrandallField;
-    use crate::gates::gate_testing::test_low_degree;
-    use crate::gates::mul_extension::MulExtensionGate;
-
-    #[test]
-    fn low_degree() {
-        test_low_degree(MulExtensionGate::<4>::new::<CrandallField>())
-    }
-}

src/lib.rs

@@ -1,3 +1,5 @@
+#![feature(destructuring_assignment)]
+
 pub mod circuit_builder;
 pub mod circuit_data;
 pub mod field;

src/permutation_argument.rs

@@ -57,7 +57,7 @@ impl TargetPartitions {
     }
 
     pub fn to_wire_partitions(&self) -> WirePartitions {
-        // Here we just drop all CircuitInputs, leaving all GateInputs.
+        // Here we keep just the Wire targets, filtering out everything else.
         let mut partitions = Vec::new();
         let mut indices = HashMap::new();
 

src/polynomial/commitment.rs

@@ -327,7 +327,7 @@ mod tests {
             },
             degree_bits: 0,
             gates: vec![],
-            max_filtered_constraint_degree_bits: 0,
+            max_filtered_constraint_degree: 0,
             num_gate_constraints: 0,
             num_constants: 4,
             k_is: vec![F::ONE; 6],

src/polynomial/polynomial.rs

@@ -5,7 +5,7 @@ use std::ops::{Add, AddAssign, Mul, MulAssign, Sub, SubAssign};
 use anyhow::{ensure, Result};
 
 use crate::field::extension_field::Extendable;
-use crate::field::fft::{fft, ifft};
+use crate::field::fft::{fft, fft_with_options, ifft};
 use crate::field::field::Field;
 use crate::util::log2_strict;
 
@@ -55,7 +55,7 @@ impl<F: Field> PolynomialValues<F> {
 
     pub fn lde(self, rate_bits: usize) -> Self {
         let coeffs = ifft(self).lde(rate_bits);
-        fft(coeffs)
+        fft_with_options(coeffs, Some(rate_bits), None)
     }
 
     pub fn degree(&self) -> usize {

src/prover.rs

@@ -12,7 +12,7 @@ use crate::polynomial::commitment::ListPolynomialCommitment;
 use crate::polynomial::polynomial::{PolynomialCoeffs, PolynomialValues};
 use crate::proof::Proof;
 use crate::timed;
-use crate::util::transpose;
+use crate::util::{log2_ceil, transpose};
 use crate::vars::EvaluationVarsBase;
 use crate::witness::{PartialWitness, Witness};
 
@@ -219,23 +219,22 @@ fn compute_quotient_polys<'a, F: Extendable<D>, const D: usize>(
     alphas: &[F],
 ) -> Vec<PolynomialCoeffs<F>> {
     let num_challenges = common_data.config.num_challenges;
+    let max_filtered_constraint_degree_bits = log2_ceil(common_data.max_filtered_constraint_degree);
     assert!(
-        common_data.max_filtered_constraint_degree_bits <= common_data.config.rate_bits,
+        max_filtered_constraint_degree_bits <= common_data.config.rate_bits,
         "Having constraints of degree higher than the rate is not supported yet. \
         If we need this in the future, we can precompute the larger LDE before computing the `ListPolynomialCommitment`s."
     );
 
     // We reuse the LDE computed in `ListPolynomialCommitment` and extract every `step` points to get
     // an LDE matching `max_filtered_constraint_degree`.
-    let step =
+    let step = 1 << (common_data.config.rate_bits - max_filtered_constraint_degree_bits);
-        1 << (common_data.config.rate_bits - common_data.max_filtered_constraint_degree_bits);
     // When opening the `Z`s polys at the "next" point in Plonk, need to look at the point `next_step`
     // steps away since we work on an LDE of degree `max_filtered_constraint_degree`.
-    let next_step = 1 << common_data.max_filtered_constraint_degree_bits;
+    let next_step = 1 << max_filtered_constraint_degree_bits;
 
-    let points = F::two_adic_subgroup(
+    let points =
-        common_data.degree_bits + common_data.max_filtered_constraint_degree_bits,
+        F::two_adic_subgroup(common_data.degree_bits + max_filtered_constraint_degree_bits);
-    );
     let lde_size = points.len();
 
     // Retrieve the LDE values at index `i`.
@@ -243,10 +242,8 @@ fn compute_quotient_polys<'a, F: Extendable<D>, const D: usize>(
         comm.get_lde_values(i * step)
     };
 
-    let z_h_on_coset = ZeroPolyOnCoset::new(
+    let z_h_on_coset =
-        common_data.degree_bits,
+        ZeroPolyOnCoset::new(common_data.degree_bits, max_filtered_constraint_degree_bits);
-        common_data.max_filtered_constraint_degree_bits,
-    );
 
     let quotient_values: Vec<Vec<F>> = points
         .into_par_iter()

src/recursive_verifier.rs

@@ -5,7 +5,7 @@ use crate::gates::gate::GateRef;
 use crate::proof::ProofTarget;
 
 const MIN_WIRES: usize = 120; // TODO: Double check.
-const MIN_ROUTED_WIRES: usize = 8; // TODO: Double check.
+const MIN_ROUTED_WIRES: usize = 28; // TODO: Double check.
 
 /// Recursively verifies an inner proof.
 pub fn add_recursive_verifier<F: Extendable<D>, const D: usize>(

src/target.rs

@@ -8,7 +8,10 @@ use crate::wire::Wire;
 pub enum Target {
     Wire(Wire),
     PublicInput { index: usize },
-    VirtualAdviceTarget { index: usize },
+    /// A target that doesn't have any inherent location in the witness (but it can be copied to
+    /// another target that does). This is useful for representing intermediate values in witness
+    /// generation.
+    VirtualTarget { index: usize },
 }
 
 impl Target {
@@ -20,7 +23,7 @@ impl Target {
         match self {
             Target::Wire(wire) => wire.is_routable(config),
             Target::PublicInput { .. } => true,
-            Target::VirtualAdviceTarget { .. } => false,
+            Target::VirtualTarget { .. } => true,
         }
     }
 

src/util/mod.rs

@@ -49,13 +49,13 @@ pub(crate) fn transpose<T: Clone>(matrix: &[Vec<T>]) -> Vec<Vec<T>> {
 }
 
 /// Permutes `arr` such that each index is mapped to its reverse in binary.
-pub(crate) fn reverse_index_bits<T: Clone>(arr: Vec<T>) -> Vec<T> {
+pub(crate) fn reverse_index_bits<T: Copy>(arr: Vec<T>) -> Vec<T> {
     let n = arr.len();
     let n_power = log2_strict(n);
 
     let mut result = Vec::with_capacity(n);
     for i in 0..n {
-        result.push(arr[reverse_bits(i, n_power)].clone());
+        result.push(arr[reverse_bits(i, n_power)]);
     }
     result
 }
@@ -73,12 +73,11 @@ pub(crate) fn reverse_index_bits_in_place<T>(arr: &mut Vec<T>) {
 }
 
 pub(crate) fn reverse_bits(n: usize, num_bits: usize) -> usize {
-    let mut result = 0;
+    // NB: The only reason we need overflowing_shr() here as opposed
-    for i in 0..num_bits {
+    // to plain '>>' is to accommodate the case n == num_bits == 0,
-        let i_rev = num_bits - i - 1;
+    // which would become `0 >> 64`. Rust thinks that any shift of 64
-        result |= (n >> i & 1) << i_rev;
+    // bits causes overflow, even when the argument is zero.
-    }
+    n.reverse_bits().overflowing_shr(usize::BITS - num_bits as u32).0
-    result
 }
 
 #[cfg(test)]

src/util/scaling.rs

@@ -1,7 +1,12 @@
 use std::borrow::Borrow;
 
-use crate::field::extension_field::Frobenius;
+use num::Integer;
+
+use crate::circuit_builder::CircuitBuilder;
+use crate::field::extension_field::target::ExtensionTarget;
+use crate::field::extension_field::{Extendable, Frobenius};
 use crate::field::field::Field;
+use crate::gates::arithmetic::ArithmeticExtensionGate;
 use crate::polynomial::polynomial::PolynomialCoeffs;
 
 /// When verifying the composition polynomial in FRI we have to compute sums of the form
@@ -73,3 +78,142 @@ impl<F: Field> ReducingFactor<F> {
         }
     }
 }
+
+#[derive(Debug, Copy, Clone)]
+pub struct ReducingFactorTarget<const D: usize> {
+    base: ExtensionTarget<D>,
+    count: u64,
+}
+
+impl<const D: usize> ReducingFactorTarget<D> {
+    pub fn new(base: ExtensionTarget<D>) -> Self {
+        Self { base, count: 0 }
+    }
+
+    /// Reduces a length `n` vector of `ExtensionTarget`s using `n/2` `ArithmeticExtensionGate`s.
+    /// It does this by batching two steps of Horner's method in each gate.
+    /// Here's an example with `n=4, alpha=2, D=1`:
+    /// 1st gate: 2  0 4  4 3  4 11 <- 2*0+4=4, 2*4+3=11
+    /// 2nd gate: 2 11 2 24 1 24 49 <- 2*11+2=24, 2*24+1=49
+    /// which verifies that `2.reduce([1,2,3,4]) = 49`.
+    pub fn reduce<F>(
+        &mut self,
+        terms: &[ExtensionTarget<D>], // Could probably work with a `DoubleEndedIterator` too.
+        builder: &mut CircuitBuilder<F, D>,
+    ) -> ExtensionTarget<D>
+    where
+        F: Extendable<D>,
+    {
+        let zero = builder.zero_extension();
+        let l = terms.len();
+        self.count += l as u64;
+
+        let mut terms_vec = terms.to_vec();
+        // If needed, we pad the original vector so that it has even length.
+        if terms_vec.len().is_odd() {
+            terms_vec.push(zero);
+        }
+        terms_vec.reverse();
+
+        let mut acc = zero;
+        for pair in terms_vec.chunks(2) {
+            // We will route the output of the first arithmetic operation to the multiplicand of the
+            // second, i.e. we compute the following:
+            //     out_0 = alpha acc + pair[0]
+            //     acc' = out_1 = alpha out_0 + pair[1]
+            let gate = builder.num_gates();
+            let out_0 =
+                ExtensionTarget::from_range(gate, ArithmeticExtensionGate::<D>::wires_output_0());
+            acc = builder
+                .double_arithmetic_extension(
+                    F::ONE,
+                    F::ONE,
+                    self.base,
+                    acc,
+                    pair[0],
+                    out_0,
+                    pair[1],
+                )
+                .1;
+        }
+        acc
+    }
+
+    pub fn shift<F>(
+        &mut self,
+        x: ExtensionTarget<D>,
+        builder: &mut CircuitBuilder<F, D>,
+    ) -> ExtensionTarget<D>
+    where
+        F: Extendable<D>,
+    {
+        let exp = builder.exp_u64_extension(self.base, self.count);
+        let tmp = builder.mul_extension(exp, x);
+        self.count = 0;
+        tmp
+    }
+
+    pub fn reset(&mut self) {
+        self.count = 0;
+    }
+
+    pub fn repeated_frobenius<F>(&self, count: usize, builder: &mut CircuitBuilder<F, D>) -> Self
+    where
+        F: Extendable<D>,
+    {
+        Self {
+            base: self.base.repeated_frobenius(count, builder),
+            count: self.count,
+        }
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::circuit_data::CircuitConfig;
+    use crate::field::crandall_field::CrandallField;
+    use crate::field::extension_field::quartic::QuarticCrandallField;
+    use crate::witness::PartialWitness;
+
+    fn test_reduce_gadget(n: usize) {
+        type F = CrandallField;
+        type FF = QuarticCrandallField;
+        const D: usize = 4;
+
+        let config = CircuitConfig::large_config();
+
+        let mut builder = CircuitBuilder::<F, D>::new(config);
+
+        let alpha = FF::rand();
+        let alpha = FF::ONE;
+        let vs = (0..n)
+            .map(|i| FF::from_canonical_usize(i))
+            .collect::<Vec<_>>();
+
+        let manual_reduce = ReducingFactor::new(alpha).reduce(vs.iter());
+        let manual_reduce = builder.constant_extension(manual_reduce);
+
+        let mut alpha_t = ReducingFactorTarget::new(builder.constant_extension(alpha));
+        let vs_t = vs
+            .iter()
+            .map(|&v| builder.constant_extension(v))
+            .collect::<Vec<_>>();
+        let circuit_reduce = alpha_t.reduce(&vs_t, &mut builder);
+
+        builder.assert_equal_extension(manual_reduce, circuit_reduce);
+
+        let data = builder.build();
+        let proof = data.prove(PartialWitness::new());
+    }
+
+    #[test]
+    fn test_reduce_gadget_even() {
+        test_reduce_gadget(10);
+    }
+
+    #[test]
+    fn test_reduce_gadget_odd() {
+        test_reduce_gadget(11);
+    }
+}

src/util/timing.rs

@@ -7,7 +7,7 @@ macro_rules! timed {
 
         let timer = Instant::now();
         let res = $a;
-        info!("{:.3}s {}", timer.elapsed().as_secs_f32(), $msg);
+        info!("{:.4}s {}", timer.elapsed().as_secs_f32(), $msg);
         res
     }};
 }

src/witness.rs

@@ -30,7 +30,7 @@ impl<F: Field> Witness<F> {
         F: Extendable<D>,
     {
         for &(a, b) in copy_constraints {
-            // TODO: Take care of public inputs once they land.
+            // TODO: Take care of public inputs once they land, and virtual targets.
             if let (
                 Target::Wire(Wire {
                     gate: a_gate,
@@ -78,6 +78,18 @@ impl<F: Field> PartialWitness<F> {
         witness
     }
 
+    pub fn singleton_extension_target<const D: usize>(
+        et: ExtensionTarget<D>,
+        value: F::Extension,
+    ) -> Self
+    where
+        F: Extendable<D>,
+    {
+        let mut witness = PartialWitness::new();
+        witness.set_extension_target(et, value);
+        witness
+    }
+
     pub fn is_empty(&self) -> bool {
         self.target_values.is_empty()
     }