Merge pull request #991 from powdr-labs/ec-add

Arithmetic Machine: Add equations 1-4 (Witgen not working yet)

.github/workflows/pr-tests.yml

@@ -102,7 +102,7 @@ jobs:
     strategy:
       matrix:
         test:
-        - "'test(=test_keccak) | test(=test_vec_median) | test(=instruction_tests::addi)'"
+        - "'test(=test_keccak) | test(=test_vec_median) | test(=instruction_tests::addi) | test(=arith_test)'"
         - "'test(=test_many_chunks)'"
     needs: build
     runs-on: ubuntu-latest

pipeline/tests/pil.rs

@@ -1,6 +1,5 @@
 use powdr_number::GoldilocksField;
 use powdr_pipeline::test_util::{gen_estark_proof, gen_halo2_proof, verify_test_file};
-use std::path::Path;
 use test_log::test;
 
 pub fn verify_pil(file_name: &str, inputs: Vec<GoldilocksField>) {
@@ -175,16 +174,6 @@ fn conditional_fixed_constraints() {
     verify_pil("pil/conditional_fixed_constraints.pil", Default::default());
 }
 
-#[test]
-fn arith_improved() {
-    let f = "pil/arith_improved.pil";
-    powdr_pil_analyzer::analyze_file::<GoldilocksField>(
-        &Path::new(&format!("{}/../test_data/{f}", env!("CARGO_MANIFEST_DIR")))
-            .canonicalize()
-            .unwrap(),
-    );
-}
-
 #[test]
 fn referencing_arrays() {
     let f = "pil/referencing_array.pil";

pipeline/tests/powdr_std.rs

@@ -34,6 +34,7 @@ fn split_gl_test() {
 }
 
 #[test]
+#[ignore = "Too slow"]
 fn arith_test() {
     let f = "std/arith_test.asm";
     verify_test_file::<GoldilocksField>(f, Default::default(), vec![]);

std/arith.asm

@@ -17,18 +17,18 @@ machine Arith(CLK32_31, operation_id){
     
     // Performs elliptic curve addition of points (x1, y2) and (x2, y2).
     // Operation ID is 10 = 0b1010, i.e., we activate equations 1, 3, and 4.
-    // TODO: Implement these equations
+    // TODO: Witgen doesn't work for that yet.
     operation ec_add<10> x1_0, x1_1, x1_2, x1_3, x1_4, x1_5, x1_6, x1_7, y1_0, y1_1, y1_2, y1_3, y1_4, y1_5, y1_6, y1_7, x2_0, x2_1, x2_2, x2_3, x2_4, x2_5, x2_6, x2_7, y2_0, y2_1, y2_2, y2_3, y2_4, y2_5, y2_6, y2_7 -> x3_0, x3_1, x3_2, x3_3, x3_4, x3_5, x3_6, x3_7, y3_0, y3_1, y3_2, y3_3, y3_4, y3_5, y3_6, y3_7;
     
     // Performs elliptic curve doubling of point (x1, y2).
     // Operation ID is 12 = 0b1100, i.e., we activate equations 2, 3, and 4.
-    // TODO: Implement these equations
+    // TODO: Witgen doesn't work for that yet.
     operation ec_double<12> x1_0, x1_1, x1_2, x1_3, x1_4, x1_5, x1_6, x1_7, y1_0, y1_1, y1_2, y1_3, y1_4, y1_5, y1_6, y1_7 -> x3_0, x3_1, x3_2, x3_3, x3_4, x3_5, x3_6, x3_7, y3_0, y3_1, y3_2, y3_3, y3_4, y3_5, y3_6, y3_7;
     
     let BYTE = |i| i & 0xff;
     let BYTE2 = |i| i & 0xffff;
 
-    pol commit x1[16], y1[16], x2[16], y2[16], x3[16], y3[16];
+    pol commit x1[16], y1[16], x2[16], y2[16], x3[16], y3[16], s[16], q0[16], q1[16], q2[16];
 
     // Intermediate polynomials, 32-Bit each
     pol x1_0 = x1[1] * 2**16 + x1[0];
@@ -88,29 +88,60 @@ machine Arith(CLK32_31, operation_id){
     let CLK32: col[32] = array::new(32, |i| |row| if row % 32 == i { 1 } else { 0 });
     let CLK32_31 = CLK32[31];
 
+    // TODO: Add the equivalent of these constraints for soundness: https://github.com/0xPolygonHermez/zkevm-proverjs/blob/main/pil/arith.pil#L43-L243
+
     /****
     *
-    * LATCH POLS: x1,y1,x2,y2,y3
+    * LATCH POLS: x1,y1,x2,y2,x3,y3,s,q0,q1,q2
     *
     *****/
 
-    array::map(x1, |e| unchanged_until(e, CLK32[31]));
-    array::map(y1, |e| unchanged_until(e, CLK32[31]));
-    array::map(x2, |e| unchanged_until(e, CLK32[31]));
-    array::map(y2, |e| unchanged_until(e, CLK32[31]));
-    array::map(x3, |e| unchanged_until(e, CLK32[31]));
-    array::map(y3, |e| unchanged_until(e, CLK32[31]));
+    let fixed_inside_32_block = [|e| unchanged_until(e, CLK32[31])][0];
+
+    array::map(x1, fixed_inside_32_block);
+    array::map(y1, fixed_inside_32_block);
+    array::map(x2, fixed_inside_32_block);
+    array::map(y2, fixed_inside_32_block);
+    array::map(x3, fixed_inside_32_block);
+    array::map(y3, fixed_inside_32_block);
+    array::map(s, fixed_inside_32_block);
+    array::map(q0, fixed_inside_32_block);
+    array::map(q1, fixed_inside_32_block);
+    array::map(q2, fixed_inside_32_block);
 
     /****
     *
-    * RANGE CHECK x1,y1,x2,y2,y3
+    * RANGE CHECK x1,y1,x2,y2,x3,y3,s,q0,q1,q2
     *
     *****/
 
     sum(16, |i| x1[i] * CLK32[i]) + sum(16, |i| y1[i] * CLK32[16 + i]) in BYTE2;
     sum(16, |i| x2[i] * CLK32[i]) + sum(16, |i| y2[i] * CLK32[16 + i]) in BYTE2;
     sum(16, |i| x3[i] * CLK32[i]) + sum(16, |i| y3[i] * CLK32[16 + i]) in BYTE2;
+    // Note that for q0-q2, we only range-constrain the first 15 limbs here
+    sum(16, |i| s[i] * CLK32[i]) + sum(15, |i| q0[i] * CLK32[16 + i]) in BYTE2;
+    sum(15, |i| q1[i] * CLK32[i]) + sum(15, |i| q2[i] * CLK32[16 + i]) in BYTE2;
+
+    // The most significant limbs of q0-q2 are constrained to be 32 bits
+    // In Polygon's version they are 19 bits, but that requires increasing the minimum degree
+    // to fit the lookup.
+    // Instead, we decompose the most significant limb into two 16-Bit limbs.
+    // Having a larger range-constraint is fine, because we're only multiplying it with 16-bit
+    // limbs of the prime, so the result is within 48 bits, still far from overflowing the
+    // Goldilocks field.
+    pol witness q0_15_high, q0_15_low, q1_15_high, q1_15_low, q2_15_high, q2_15_low;
+    q0_15_high * CLK32[0] + q0_15_low * CLK32[1] + q1_15_high * CLK32[2] + q1_15_low * CLK32[3] + q2_15_high * CLK32[4] + q2_15_low * CLK32[5] in BYTE2;
+
+    fixed_inside_32_block(q0_15_high);
+    fixed_inside_32_block(q0_15_low);
+    fixed_inside_32_block(q1_15_high);
+    fixed_inside_32_block(q1_15_low);
+    fixed_inside_32_block(q2_15_high);
+    fixed_inside_32_block(q2_15_low);
 
+    q0[15] = 2**16 * q0_15_high + q0_15_low;
+    q1[15] = 2**16 * q1_15_high + q1_15_low;
+    q2[15] = 2**16 * q2_15_high + q2_15_low;
 
     /*******
     *
@@ -135,7 +166,12 @@ machine Arith(CLK32_31, operation_id){
     let y1f = array_as_fun(y1);
     let x2f = array_as_fun(x2);
     let y2f = array_as_fun(y2);
+    let x3f = array_as_fun(x3);
     let y3f = array_as_fun(y3);
+    let sf = array_as_fun(s);
+    let q0f = array_as_fun(q0);
+    let q1f = array_as_fun(q1);
+    let q2f = array_as_fun(q2);
 
     // Defined for arguments from 0 to 31 (inclusive)
     let eq0: int -> expr = |nr|
@@ -144,17 +180,76 @@ machine Arith(CLK32_31, operation_id){
         - shift_right(y2f, 16)(nr)
         - y3f(nr);
 
+    /*******
+    *
+    * EQ1: s * x2 - s * x1 - y2 + y1 + (q0 * p)
+    *
+    *******/
+
+    // 0xffffffffffffffffffffffffffffffffffffffffffffffffffff fffe ffff fc2f
+    let p = array_as_fun([
+        0xfc2f, 0xffff, 0xfffe, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff,
+        0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff
+    ]);
+
+    // The "- 4 * shift_right(p, 16)" effectively subtracts 4 * (p << 16 * 16) = 2 ** 258 * p
+    // As a result, the term computes `(x - 2 ** 258) * p`.
+    let product_with_p = (|| |x| |nr| product(p, x)(nr) - 4 * shift_right(p, 16)(nr))();
+
+    let eq1 = (|| |nr| product(sf, x2f)(nr) - product(sf, x1f)(nr) - y2f(nr) + y1f(nr) + product_with_p(q0f)(nr))();
+
+    /*******
+    *
+    * EQ2:  2 * s * y1 - 3 * x1 * x1 + (q0 * p)
+    *
+    *******/
+
+    let eq2 = (|| |nr| 2 * product(sf, y1f)(nr) - 3 * product(x1f, x1f)(nr) + product_with_p(q0f)(nr))();
+
+    /*******
+    *
+    * EQ3:  s * s - x1 - x2 - x3 + (q1 * p)
+    *
+    *******/
+
+    // If we're doing the ec_double operation (selEq[2] == 1), x2 is so far unconstrained and should be set to x1
+    array::new(16, |i| selEq[2] * (x1[i] - x2[i]) = 0);
+
+    let eq3 = (|| |nr| product(sf, sf)(nr) - x1f(nr) - x2f(nr) - x3f(nr) + product_with_p(q1f)(nr))();
+
+
+    /*******
+    *
+    * EQ4:  s * x1 - s * x3 - y1 - y3 + (q2 * p)
+    *
+    *******/
+
+    let eq4 = (|| |nr| product(sf, x1f)(nr) - product(sf, x3f)(nr) - y1f(nr) - y3f(nr) + product_with_p(q2f)(nr))();
+
+
+    /*******
+    *
+    * Equation Selectors
+    *
+    *******/
+
     // Binary selectors for the equations that are activated. Determined from the operation ID via bit-decomposition.
     // Note that there are only 4 selectors because equation 4 is activated iff. equation 3 is activated, so we can
     // re-use the same selector.
     pol commit selEq[4];
     // Note that this is not necessary, because the operation ID is already constant within the block
-    // array::map(selEq, |e| unchanged_until(e, CLK32[31]));
+    // array::map(selEq, fixed_inside_32_block);
     array::map(selEq, |c| force_bool(c));
     sum(4, |i| 2 ** i * selEq[i]) = operation_id;
+
+    /*******
+    *
+    * Carry
+    *
+    *******/
     
     // Note that Polygon uses a single 22-Bit column. However, this approach allows for a lower degree (2**16)
-    // while still preventing overflows: The 32-bit carry gets added to 32 16-Bit values, which can't overflow
+    // while still preventing overflows: The 32-bit carry gets added to 32 48-Bit values, which can't overflow
     // the Goldilocks field.
     pol witness carry_low[3], carry_high[3];
     { carry_low[0] } in { BYTE2 };
@@ -172,10 +267,21 @@ machine Arith(CLK32_31, operation_id){
     
     array::map(carry, |c| c * CLK32[0] = 0);
 
-    // TODO: Somehow witgen fails if a carry is unconstrained
-    selEq[0] * carry[1] = 0;
-    selEq[0] * carry[2] = 0;
-
+    /*******
+    *
+    * Putting everything together
+    *
+    *******/
+    
     col eq0_sum = sum(32, |i| eq0(i) * CLK32[i]);
+    col eq1_sum = sum(32, |i| eq1(i) * CLK32[i]);
+    col eq2_sum = sum(32, |i| eq2(i) * CLK32[i]);
+    col eq3_sum = sum(32, |i| eq3(i) * CLK32[i]);
+    col eq4_sum = sum(32, |i| eq4(i) * CLK32[i]);
+
     selEq[0] * (eq0_sum + carry[0]) = selEq[0] * carry[0]' * 2**16;
+    selEq[1] * (eq1_sum + carry[0]) = selEq[1] * carry[0]' * 2**16;
+    selEq[2] * (eq2_sum + carry[0]) = selEq[2] * carry[0]' * 2**16;
+    selEq[3] * (eq3_sum + carry[1]) = selEq[3] * carry[1]' * 2**16;
+    selEq[3] * (eq4_sum + carry[2]) = selEq[3] * carry[2]' * 2**16;
 }