Merge pull request #3 from mileschet/feature/lisp

Feature/lisp

lisp/core.rs

@@ -8,23 +8,17 @@ use crate::printer::pr_seq;
 use crate::reader::read_str;
 use crate::types::MalErr::ErrMalVal;
 use crate::types::MalVal::{
-    Atom, Bool, Func, Hash, Int, List, MalFunc, Nil, Str,
-    Sym, Vector
+    Atom, Bool, Func, Hash, Int, List, MalFunc, Nil, Str, Sym, Vector, ZKScalar,
 };
 use crate::types::{MalArgs, MalRet, MalVal, _assoc, _dissoc, atom, error, func, hash_map};
-use MalVal::ZKScalar;
-use bellman::{gadgets::Assignment, groth16, Circuit, ConstraintSystem, SynthesisError};
 
-use bls12_381::Scalar;
+use bls12_381;
 use ff::{Field, PrimeField};
 
 use sapvi::bls_extensions::BlsStringConversion;
 
-
-
 use std::ops::{AddAssign, MulAssign, SubAssign};
 
-
 macro_rules! fn_t_int_int {
     ($ret:ident, $fn:expr) => {{
         |a: MalArgs| match (a[0].clone(), a[1].clone()) {
@@ -171,11 +165,11 @@ fn unpack_bits(a: MalArgs) -> MalRet {
     let mut result = vec![];
     match a[0].clone() {
         Str(ref s) => {
-            let value = Scalar::from_string(s);
+            let value = bls12_381::Scalar::from_string(s);
             for (_, bit) in value.to_le_bits().into_iter().cloned().enumerate() {
                 match bit {
-                    true => result.push(Scalar::one()),
-                    false => result.push(Scalar::zero()),
+                    true => result.push(bls12_381::Scalar::one()),
+                    false => result.push(bls12_381::Scalar::zero()),
                 }
             }
             Ok(list!(result
@@ -261,7 +255,10 @@ fn conj(a: MalArgs) -> MalRet {
 fn sub_scalar(a: MalArgs) -> MalRet {
     match (a[0].clone(), a[1].clone()) {
         (Str(a0), Str(a1)) => {
-            let (mut s0, s1) = (Scalar::from_string(&a0), Scalar::from_string(&a1));
+            let (mut s0, s1) = (
+                bls12_381::Scalar::from_string(&a0),
+                bls12_381::Scalar::from_string(&a1),
+            );
             s0.sub_assign(s1);
             Ok(Str(std::string::ToString::to_string(&s0)[2..].to_string()))
         }
@@ -283,7 +280,10 @@ fn mul_scalar(a: MalArgs) -> MalRet {
 fn div_scalar(a: MalArgs) -> MalRet {
     match (a[0].clone(), a[1].clone()) {
         (Str(a0), Str(a1)) => {
-            let (s0, s1) = (Scalar::from_string(&a0), Scalar::from_string(&a1));
+            let (s0, s1) = (
+                bls12_381::Scalar::from_string(&a0),
+                bls12_381::Scalar::from_string(&a1),
+            );
             let ret = s1.invert().map(|other| *&s0 * other);
             Ok(Str(
                 std::string::ToString::to_string(&ret.unwrap())[2..].to_string()
@@ -298,66 +298,70 @@ fn range(a: MalArgs) -> MalRet {
     match (a[0].clone(), a[1].clone()) {
         (Int(a0), Int(a1)) => {
             for n in a0..a1 {
-               result.push(n);
-            };
-            Ok(list!(result
-                .iter()
-                .map(|_a| Nil) 
-                .collect::<Vec<MalVal>>()))
-        },
-        _ => error("expected int int")
+                result.push(n);
+            }
+            Ok(list!(result.iter().map(|_a| Nil).collect::<Vec<MalVal>>()))
+        }
+        _ => error("expected int int"),
     }
 }
 
-fn alloc_input(a: MalArgs) -> MalRet {
-    // TODO implement
-    Ok(ZKScalar(bls12_381::Scalar::zero()))
-}
-
 fn scalar_zero(a: MalArgs) -> MalRet {
-    Ok(ZKScalar(bls12_381::Scalar::zero()))
+    Ok(vector![vec![
+        ZKScalar(bls12_381::Scalar::zero()),
+        a[0].clone()
+    ]])
 }
 
 fn scalar_one(a: MalArgs) -> MalRet {
-    Ok(ZKScalar(bls12_381::Scalar::one()))
+    match a.len() {
+        0 => Ok(vector![vec![ZKScalar(bls12_381::Scalar::one())]]),
+        _ => Ok(vector![vec![
+            ZKScalar(bls12_381::Scalar::one()),
+            a[0].clone()
+        ]]),
+    }
+}
+
+fn cs_one(a: MalArgs) -> MalRet {
+    Ok(vector![vec![Sym("cs::one".to_string())]])
 }
 
 fn negate_from(a: MalArgs) -> MalRet {
-    match a[0].apply(vec![])? {
-        ZKScalar(a0) => {
-            Ok(ZKScalar(a0.neg()))
+    match a[0].clone() {
+        ZKScalar(a0) => Ok(ZKScalar(a0.neg())),
+        _ => match a[0].apply(vec![])? {
+            List(v, _) | Vector(v, _) => match v[0] {
+                ZKScalar(val) => Ok(vector![vec![ZKScalar(val.neg())]]),
+                _ => error("not scalar"),
+            },
+            _ => return error("non zkscalar passed to negate"),
         },
-        Nil => error("nil not supported"),
-        _ => error("negate error, expected (zkscalar)"),
     }
 }
+
 fn scalar_from(a: MalArgs) -> MalRet {
-    println!("{:?}", a);
     match a[0].clone() {
         Str(a0) => {
-            let s0 = Scalar::from_string(&a0.to_string());
+            let s0 = bls12_381::Scalar::from_string(&a0.to_string());
             Ok(ZKScalar(s0))
-        },
-        _ => error("expected (string)"),
+        }
+        Int(a0) => {
+            println!("{:?}", a0);
+            let s0 = bls12_381::Scalar::from(a0 as u64);
+            Ok(ZKScalar(s0))
+        }
+        _ => error("expected (string or int)"),
     }
 }
-fn alloc(a: MalArgs) -> MalRet {
-    println!("{:?}", a);
-    Ok(Nil)
-}
-fn cs_one(a: MalArgs) -> MalRet {
-    println!("{:?}", a);
-    Ok(Nil)
-}
-fn bellman_one(a: MalArgs) -> MalRet {
-    println!("{:?}", a);
-    Ok(Nil)
-}
 
 fn add_scalar(a: MalArgs) -> MalRet {
     match (a[0].clone(), a[1].clone()) {
         (Str(a0), Str(a1)) => {
-            let (mut s0, s1) = (Scalar::from_string(&a0), Scalar::from_string(&a1));
+            let (mut s0, s1) = (
+                bls12_381::Scalar::from_string(&a0),
+                bls12_381::Scalar::from_string(&a1),
+            );
             s0.add_assign(s1);
             Ok(Str(std::string::ToString::to_string(&s0)[2..].to_string()))
         }
@@ -476,14 +480,10 @@ pub fn ns() -> Vec<(&'static str, MalVal)> {
         ("swap!", func(|a| a[0].swap_bang(&a[1..].to_vec()))),
         ("unpack-bits", func(unpack_bits)),
         ("range", func(range)),
-        ("alloc", func(alloc)),
-        ("alloc-input", func(alloc_input)),
         ("scalar::one", func(scalar_one)),
         ("neg", func(negate_from)),
         ("scalar::zero", func(scalar_zero)),
-        // TODO add .neg maybe neg, add and sub
         ("scalar", func(scalar_from)),
         ("cs::one", func(cs_one)),
-        ("bellman::one", func(bellman_one)),
     ]
 }

lisp/lisp.rs

@@ -1,26 +1,27 @@
 #![allow(non_snake_case)]
 
-use crate::MalVal::Zk;
+use crate::MalVal::Enforce;
 use crate::groth16::VerifyingKey;
 use crate::types::LispCircuit;
-use sapvi::{ZKVMCircuit, ZKVirtualMachine};
+use crate::MalVal::Zk;
+use bellman::groth16::PreparedVerifyingKey;
 use sapvi::bls_extensions::BlsStringConversion;
+use sapvi::{ZKVMCircuit, ZKVirtualMachine};
 
 use simplelog::*;
 
-use bellman::{
-    gadgets::{
-        Assignment,
-    },
-    groth16, Circuit, ConstraintSystem, SynthesisError,
-};
+use bellman::{gadgets::Assignment, groth16, Circuit, ConstraintSystem, SynthesisError};
 use bls12_381::Bls12;
 use bls12_381::Scalar;
 use ff::{Field, PrimeField};
 use rand::rngs::OsRng;
-use std::ops::{AddAssign, MulAssign, SubAssign};
+use types::EnforceAllocation;
+use std::{borrow::BorrowMut, rc::Rc};
 use std::time::Instant;
-use std::rc::Rc;
+use std::{
+    cell::RefCell,
+    ops::{AddAssign, MulAssign, SubAssign},
+};
 
 //use std::collections::HashMap;
 use fnv::FnvHashMap;
@@ -38,11 +39,8 @@ extern crate regex;
 #[macro_use]
 mod types;
 use crate::types::MalErr::{ErrMalVal, ErrString};
-use crate::types::MalVal::{
-    Bool, Func, Hash, List, MalFunc, Nil, Str,
-    Sym, Vector,
-};
-use crate::types::{error, format_error, MalArgs, MalErr, MalRet, MalVal};
+use crate::types::MalVal::{Bool, Func, Hash, List, MalFunc, Nil, Str, Sym, Vector};
+use crate::types::{error, format_error, Allocation, MalArgs, MalErr, MalRet, MalVal};
 mod env;
 mod printer;
 mod reader;
@@ -50,7 +48,7 @@ use crate::env::{env_bind, env_find, env_get, env_new, env_set, env_sets, Env};
 #[macro_use]
 mod core;
 
-pub const ZK_CIRCUIT_ENV_KEY : &str = "ZKC";
+pub const ZK_CIRCUIT_ENV_KEY: &str = "ZKC";
 
 // read
 fn read(str: &str) -> MalRet {
@@ -287,30 +285,7 @@ fn eval(mut ast: MalVal, mut env: Env) -> MalRet {
                             _ => Ok(Nil),
                         }
                     }
-                    Sym(ref a0sym) if a0sym == "setup" => {
-                        let a1 = l[1].clone();
-                        let circuit = setup(&ast)?;
-                        println!("{:?}", a1); 
-                        env_sets(&env, ZK_CIRCUIT_ENV_KEY, circuit);
-                        eval(a1.clone(), env.clone())
-                    }
-                    Sym(ref a0sym) if a0sym == "prove" => {
-                        let a1 = l[1].clone();
-                        println!("{:?}", a1);
-                        prove(a1.clone(), env.clone())
-                    }
-                    //Sym(ref a0sym) if a0sym == "verify" => {
-                    Sym(ref a0sym) if a0sym == "enforce" => {
-                        let (a1, a2) = (l[0].clone(), l[1].clone());
-                        let value = eval_ast(&a2, &env)?;
-                        match value {
-                            List(ref el, _) => {
-                                println!("{:?}", el.to_vec());
-                            }
-                            _ => println!("invalid format"),
-                        }
-                        Ok(Nil)
-                    }
+
                     Sym(ref a0sym) if a0sym == "fn*" => {
                         let (a1, a2) = (l[1].clone(), l[2].clone());
                         Ok(MalFunc {
@@ -329,6 +304,103 @@ fn eval(mut ast: MalVal, mut env: Env) -> MalRet {
                         }
                         continue 'tco;
                     }
+                    Sym(ref a0sym) if a0sym == "setup" => {
+                        let a1 = l[1].clone();
+                        let pvk = setup(a1.clone(), env.clone())?;
+                        ast = eval(a1.clone(), env.clone())?;
+                        continue 'tco;
+                    }
+                    Sym(ref a0sym) if a0sym == "prove" => {
+                        let a1 = l[0].clone();
+                        println!("prove {:?}", a1);
+                        prove(a1.clone(), env.clone())
+                    }
+                    Sym(ref a0sym) if a0sym == "alloc-input" => {
+                        let a1 = l[1].clone();
+                        let value = eval(l[2].clone(), env.clone())?;
+                        let result = eval(value.clone(), env.clone())?;
+                        //                        let symbol = MalVal::Sym(a1.pr_str(false));
+                        //                       env_set(&env, Sym(a1.pr_str(false)), result.clone());
+                        if let Hash(allocs, _) = get_allocations(&env, "AllocationsInput")? {
+                            let mut new_hm: FnvHashMap<String, MalVal> = FnvHashMap::default();
+                            for (k, v) in allocs.iter() {
+                                new_hm.insert(k.to_string(), eval(v.clone(), env.clone())?);
+                            }
+                            new_hm.insert(a1.pr_str(false), result);
+                            env_set(
+                                &env,
+                                Sym("AllocationsInput".to_string()),
+                                Hash(Rc::new(new_hm), Rc::new(Nil)),
+                            );
+                        };
+                        Ok(Nil)
+                    }
+                    Sym(ref a0sym) if a0sym == "alloc" => {
+                        let a1 = l[1].clone();
+                        let value = eval(l[2].clone(), env.clone())?;
+                        let result = eval(value.clone(), env.clone())?;
+                        if let Hash(allocs, _) = get_allocations(&env, "Allocations")? {
+                            let mut new_hm: FnvHashMap<String, MalVal> = FnvHashMap::default();
+                            for (k, v) in allocs.iter() {
+                                new_hm.insert(k.to_string(), eval(v.clone(), env.clone())?);
+                            }
+                            new_hm.insert(a1.pr_str(false), result);
+                            env_set(
+                                &env,
+                                Sym("Allocations".to_string()),
+                                Hash(Rc::new(new_hm), Rc::new(Nil)),
+                            );
+                        };
+                        Ok(Nil)
+                    }
+                    //Sym(ref a0sym) if a0sym == "verify" => {
+                    Sym(ref a0sym) if a0sym == "enforce" => {
+                        // here i'm considering that we always have tuple with only two elements
+                        // also it's important to keep in mind for the sake of brevity of this v0
+                        // we will not allow calculation or any lisp evaluations inside the enforce 
+                        // it means that every symbol will be on allocations and we will do the 
+                        // find/replace on the bellman circuit, it's nasty v0
+                        let left = match l[1].clone() {
+                            List(v, _) | Vector(v, _) => {                                
+                                (v[0].pr_str(false), v[1].pr_str(false))
+                            }
+                            _ => {("".to_string(), "".to_string())}
+                        };
+                        let right = match l[1].clone() {
+                            List(v, _) | Vector(v, _) => {                                
+                                (v[0].pr_str(false), v[1].pr_str(false))
+                            }
+                            _ => {("".to_string(), "".to_string())}
+                        };
+                        let output = match l[1].clone() {
+                            List(v, _) | Vector(v, _) => {                                
+                                (v[0].pr_str(false), v[1].pr_str(false))
+                            }
+                            _ => {("".to_string(), "".to_string())}
+                        };
+                        let enforce = EnforceAllocation{left : left, right : right, output : output};
+                        let mut new_vec: Vec<EnforceAllocation> = vec![enforce];
+                        match get_enforce_allocs(&env)? {
+                            Enforce(v) => { 
+                                println!("---> {:?}", v);
+                                for value in v.iter() {
+                                    new_vec.push(value.to_owned());
+                                }
+                                },
+                                _ => {}
+                        };
+                        env_set(
+                            &env,
+                            Sym("AllocationsEnforce".to_string()),
+                            vector![vec![Enforce(Rc::new(new_vec))]],
+                        );
+
+                        // println!("allocations {:?}", get_allocations(&env, "Allocations"));
+                        // println!("allocations input {:?}", get_allocations(&env, "AllocationsInput"));
+                        println!("allocations enforce {:?}", get_enforce_allocs(&env));
+
+                        Ok(vector![vec![]])
+                    }
                     _ => match eval_ast(&ast, &env)? {
                         List(ref el, _) => {
                             let ref f = el[0].clone();
@@ -348,7 +420,8 @@ fn eval(mut ast: MalVal, mut env: Env) -> MalRet {
                                     continue 'tco;
                                 }
                                 _ => {
-                                    Ok(Nil)
+                                    Ok(vector![el.to_vec()])
+
                                     //error("call non-function")
                                 }
                             }
@@ -366,51 +439,79 @@ fn eval(mut ast: MalVal, mut env: Env) -> MalRet {
     ret
 }
 
-pub fn setup(ast: &MalVal) -> MalRet {
-    // TODO get params from ast 
+pub fn get_enforce_allocs(env: &Env) -> MalRet {
+    let found = match env_find(env, "AllocationsEnforce") {
+        Some(e) => match env_get(&e, &Sym("AllocationsEnforce".to_string())) {
+            Ok(f) => Ok(f),
+            _ => Ok(vector![vec![]])
+        },
+        _ => Ok(vector![vec![]])
+    };
+    found
+}
+pub fn get_allocations(env: &Env, key: &str) -> MalRet {
+    let mut alloc_hm: FnvHashMap<String, MalVal> = FnvHashMap::default();
+    match env_find(env, key) {
+        Some(e) => match env_get(&e, &Sym(key.to_string())) {
+            Ok(f) => Ok(f),
+            _ => Ok(Hash(Rc::new(alloc_hm), Rc::new(Nil))),
+        },
+        _ => Ok(Hash(Rc::new(alloc_hm), Rc::new(Nil))),
+    }
+}
+
+pub fn setup(ast: MalVal, mut env: Env) -> Result<PreparedVerifyingKey<Bls12>, MalErr> {
     let start = Instant::now();
     // Create parameters for our circuit. In a production deployment these would
     // be generated securely using a multiparty computation.
-    let c = LispCircuit { params: Rc::new(vector!(vec![])) };
-    // TODO move to another fn    
-    let random_parameters = groth16::generate_random_parameters::<Bls12, _, _>(c.clone(), &mut OsRng).unwrap();
+
+    // get all allocs from env
+
+    let mut c = LispCircuit {
+        params: vec![],
+        allocs: vec![],
+        alloc_inputs: vec![],
+        constraints: vec![],
+        env: env.clone(),
+    };
+    // TODO move to another fn
+    let random_parameters =
+        groth16::generate_random_parameters::<Bls12, _, _>(c, &mut OsRng).unwrap();
     let pvk = groth16::prepare_verifying_key(&random_parameters.vk);
     println!("Setup: [{:?}]", start.elapsed());
 
-    Ok(MalVal::Zk(Rc::new(c)))
+    Ok(pvk)
 }
 
 pub fn prove(mut ast: MalVal, mut env: Env) -> MalRet {
-    let c = match env_find(&env, ZK_CIRCUIT_ENV_KEY) {
-        Some(e) => match env_get(&e, &Sym(ZK_CIRCUIT_ENV_KEY.to_string()))? {
-            Zk(c) => {
-                MalVal::Zk(c)
-            }
-            _ => { MalVal::Nil }
-        }
-        None => { println!("circuit not found."); MalVal::Nil }
-    };
-
-    println!("{:?}", c);
-    
-    // Pick a preimage and compute its hash.
+    // TODO remove it
     let quantity = bls12_381::Scalar::from(3);
-    
+
     // Create an instance of our circuit (with the preimage as a witness).
-    let c = LispCircuit {
-        params: Rc::new(vector![vec![
-            ZKScalar(quantity),
-            ZKScalar(quantity * quantity),
-            ZKScalar(quantity * quantity * quantity),
-        ]]),
+    let params = {
+        let c = LispCircuit {
+            params: vec![],
+            allocs: vec![],
+            alloc_inputs: vec![],
+            constraints: vec![],
+            env: env.clone(),
+        };
+        groth16::generate_random_parameters::<Bls12, _, _>(c, &mut OsRng).unwrap()
     };
 
+    let circuit = LispCircuit {
+        params: vec![],
+        allocs: vec![],
+        alloc_inputs: vec![],
+        constraints: vec![],
+        env: env.clone(),
+    };
     let start = Instant::now();
     // Create a Groth16 proof with our parameters.
-    //let proof = groth16::create_random_proof(c, &params, &mut OsRng).unwrap();
+    let proof = groth16::create_random_proof(circuit, &params, &mut OsRng).unwrap();
     println!("Prove: [{:?}]", start.elapsed());
     Ok(MalVal::Nil)
-} 
+}
 
 pub fn verify(ast: &MalVal) -> MalRet {
     let public_input = vec![bls12_381::Scalar::from(27)];
@@ -476,7 +577,7 @@ fn repl_load(file: String) -> Result<(), ()> {
         "(def! load-file (fn* (f) (eval (read-string (str \"(do \" (slurp f) \"\nnil)\")))))",
         &repl_env,
     );
-    let _ = rep("(defmacro! cond (fn* (& xs) (if (> (count xs) 0) (list 'if (first xs) (if (> (count xs) 1) (nth xs 1) (throw \"odd number of forms to cond\")) (cons 'cond (rest (rest xs)))))))", &repl_env);
+    //let _ = rep("(defmacro! cond (fn* (& xs) (if (> (count xs) 0) (list 'if (first xs) (if (> (count xs) 1) (nth xs 1) (throw \"odd number of forms to cond\")) (cons 'cond (rest (rest xs)))))))", &repl_env);
     match rep(&format!("(load-file \"{}\")", file), &repl_env) {
         Ok(_) => std::process::exit(0),
         Err(e) => {

lisp/new-cs.lisp

@@ -1,41 +1,35 @@
-;; defzk!
-;; enforce LABEL 
-;; alloc
-;; alloc-input 
-;; scalar::one
-;; scalar::zero
-;; scalar
-;; cs::one
-;; bellman::zero
-;; setup
-;; prove
-;; verify
 (println "new-cs.lisp")
-(def! MyCircuit (fn* [aux]  
-(let* [x  (alloc "num" (first aux))
-     x2 (alloc "product num" (first (rest aux))) 
-     x3 (alloc "product num" (last aux))
-     input (alloc-input "input variable" (last aux))]
-   ;; Lc0: [(Scalar::one(), CS::one()), (Scalar::one().neg(), C)]
-;; Lc1: [(Scalar::one(), y)]
-;; Lc2: [(Scalar::one(), U), (Scalar::one().neg(), A), (Scalar::one().neg(), B)]
-(println
-    (enforce ((scalar::one x) ((neg scalar::one) x2) ((neg scalar::one) x3)))
-)
-)))
-(def! a (scalar "0000000000000000000000000000000000000000000000000000000000000003"))
-(setup (MyCircuit (a (* a a) (* (* a a) a))))
-(prove MyCircuit)
-;; (verify (prove MyCircuit) (scalar 27))
-;; (U - A - B) / (1 - C)
-;; [(1 - C)] * [y] = [U - A - B]
-;; Lc0: [(Scalar::one(), CS::one()), (Scalar::one().neg(), C)]
-;; Lc1: [(Scalar::one(), y)]
-;; Lc2: [(Scalar::one(), U), (Scalar::one().neg(), A), (Scalar::one().neg(), B)]
-;; assert (x1 + y1) * (x2 + y2) == U
-;; Lc0: [(Scalar::one(), x1), (Scalar::one(), y1)]
-;; Lc1: [(Scalar::one(), x2), (Scalar::one(), y2)]
-;; Lc2: [(Scalar::one(), U)]
-;; (enforce ((scalar::one x1) (scalar::one y1)) ((scalar::one x2) (scalar::one y2) ((scalar::one U)))
 
-     
+( (let* [aux (scalar 3)
+      x (alloc "x" aux)
+      x2 (alloc "x2" (* aux aux))
+      x3 (alloc "x3" (* aux (* aux aux)))
+      input (alloc-input "input" aux)
+      ]
+ (setup 
+   ;; (enforce left right output)
+  (
+  (enforce  
+    (scalar::one x)
+    ;;(scalar::one::neg x)
+    (scalar::one x)
+    (scalar::one x2)
+  )
+
+  (enforce 
+    (scalar::one x2)
+    (scalar::one x)
+    (scalar::one x3)
+  )
+
+  (enforce 
+    (scalar::one input)
+    (scalar::one cs::one)
+    (scalar::one x3)  
+  )
+  )
+  )
+ )
+(prove)
+)
+;; (println 'verify  (MyCircuit (scalar 27)))

lisp/types.rs

@@ -1,26 +1,42 @@
-use bellman::SynthesisError;
-use bellman::ConstraintSystem;
+use bellman::groth16::PreparedVerifyingKey;
 use bellman::Circuit;
+use bellman::ConstraintSystem;
+use bellman::SynthesisError;
+use bls12_381::Bls12;
 use std::cell::RefCell;
 use std::rc::Rc;
 //use std::collections::HashMap;
 use fnv::FnvHashMap;
 use itertools::Itertools;
 
+use crate::env;
 use crate::env::{env_bind, Env};
 use crate::types::MalErr::{ErrMalVal, ErrString};
-use crate::types::MalVal::{
-    Atom, Bool, Func, Hash, Int, List, MalFunc, Nil, Str, Sym, Vector,
-};
-
+use crate::types::MalVal::{Atom, Bool, Func, Hash, Int, List, MalFunc, Nil, Str, Sym, Vector};
 use bls12_381::Scalar;
-use sapvi::{
-    BlsStringConversion, ConstraintInstruction,
-};
+use sapvi::{BlsStringConversion, ConstraintInstruction};
+
+#[derive(Debug, Clone)]
+pub struct Allocation {
+    pub symbol: String,
+    pub value: Scalar,
+}
+
+#[derive(Debug, Clone)]
+pub struct EnforceAllocation {
+    pub left: (String, String),
+    pub right:  (String, String),
+    pub output: (String, String)
+}
 
-#[derive(Clone, Debug)]
+#[derive(Debug, Clone)]
 pub struct LispCircuit {
-    pub params: Rc<MalVal>,
+    // TODO refactor to vec
+    pub params: Vec<Option<Scalar>>,
+    pub allocs: Vec<Option<Allocation>>,
+    pub alloc_inputs: Vec<Option<Allocation>>,
+    pub constraints: Vec<Option<Scalar>>,
+    pub env: Env,
 }
 
 impl Circuit<bls12_381::Scalar> for LispCircuit {
@@ -28,6 +44,12 @@ impl Circuit<bls12_381::Scalar> for LispCircuit {
         self,
         cs: &mut CS,
     ) -> Result<(), SynthesisError> {
+        for alloc_value in &self.allocs {
+            //            let var = cs.alloc(|| "private alloc", ||)?;
+            // TODO use env
+            println!("{:?}", alloc_value);
+        }
+
         Ok(())
     }
 }
@@ -52,10 +74,9 @@ pub enum MalVal {
         meta: Rc<MalVal>,
     },
     Atom(Rc<RefCell<MalVal>>),
-    Zk(Rc<LispCircuit>),
-    Enforce(Rc<Vec<MalVal>>),
-    // TODO maybe change to bls scalar
-    ZKScalar(bls12_381::Scalar)
+    Zk(Rc<LispCircuit>), // TODO remote it
+    Enforce(Rc<Vec<EnforceAllocation>>),
+    ZKScalar(bls12_381::Scalar),
 }
 
 #[derive(Debug)]

scripts/finite_fields/polynomial.py

@@ -4,7 +4,7 @@ except ImportError:
     from itertools import izip_longest as zip_longest
 import fractions
 
-from numbertype import *
+from .numbertype import *
 
 # strip all copies of elt from the end of the list
 def strip(L, elt):

scripts/zk/3.3-encrypted-polynomial.py

@@ -0,0 +1,103 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+import random
+import numpy as np
+
+# Section 3.3.4 from "Why and How zk-SNARK Works"
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#x = rand_scalar()
+#y = ec.y_for_x(x)
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(Q, 0), Fq(Q, 1), True, bls12381)
+assert g1 + null == g1
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+s = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+d = 10
+encrypted_powers = [
+    g1 * (s**i) for i in range(d)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+target = (s - 1) * (s - 2)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+# E(p(s)) = p(s)G
+#         = c_d s^d G + ... + c_1 s^1 G + c_0 s^0 G
+#         = s^3 G - 3 s^2 G + 2 s G
+# E(h(s)) = sG
+# t(s) = s^2 - 3s + 2
+# E(h(s)) t(s) = s^3 G - 3 s^2 G + 2 s G
+
+# Lets test these manually:
+
+e_s = encrypted_powers
+e_p_s = e_s[3] - 3 * e_s[2] + 2 * e_s[1]
+e_h_s = e_s[1]
+t_s = s**2 - 3*s + 2
+assert t_s == target
+assert e_p_s == e_h_s * t_s
+
+#############################
+
+# x^3 - 3x^2 + 2x
+main_poly = np.poly1d([1, -3, 2, 0])
+# (x - 1)(x - 2)
+target_poly = np.poly1d([1, -1]) * np.poly1d([1, -2])
+
+# Calculates polynomial h(x) = p(x) / t(x)
+cofactor, remainder = main_poly / target_poly
+assert remainder == np.poly1d([0])
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers):
+    coeffs = list(poly.coef)[::-1]
+    result = null
+    for power, coeff in zip(encrypted_powers, coeffs):
+        #print(coeff, power)
+        coeff = int(coeff)
+        # I have to do this for some strange reason
+        # Because if coeff is negative and I do += power * coeff
+        # then it gives me a different result than what I expect
+        if coeff < 0:
+            result -= power * (-coeff)
+        else:
+            result += power * coeff
+    return result
+
+encrypted_poly = evaluate(main_poly, encrypted_powers)
+assert encrypted_poly == e_p_s
+encrypted_cofactor = evaluate(cofactor, encrypted_powers)
+
+# resulting g^p and g^h are provided to the verifier
+
+#################################
+# Verifier
+#################################
+
+# Last check that p = t(s) h
+
+assert encrypted_poly == encrypted_cofactor * target
+

scripts/zk/3.4-restricted-polynomial.py

@@ -0,0 +1,119 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+import random
+import numpy as np
+
+# Section 3.4 from "Why and How zk-SNARK Works"
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#x = rand_scalar()
+#y = ec.y_for_x(x)
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(Q, 0), Fq(Q, 1), True, bls12381)
+assert g1 + null == g1
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+s = rand_scalar()
+
+# calculate the shift
+a = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+d = 10
+encrypted_powers = [
+    g1 * (s**i) for i in range(d)
+]
+encrypted_shifted_powers = [
+    g1 * (a * s**i) for i in range(d)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+target = (s - 1) * (s - 2)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+# E(p(s)) = p(s)G
+#         = c_d s^d G + ... + c_1 s^1 G + c_0 s^0 G
+#         = s^3 G - 3 s^2 G + 2 s G
+# E(h(s)) = sG
+# t(s) = s^2 - 3s + 2
+# E(h(s)) t(s) = s^3 G - 3 s^2 G + 2 s G
+
+# Lets test these manually:
+
+e_s = encrypted_powers
+e_p_s = e_s[3] - 3 * e_s[2] + 2 * e_s[1]
+e_h_s = e_s[1]
+t_s = s**2 - 3*s + 2
+assert t_s == target
+assert e_p_s == e_h_s * t_s
+
+e_as = encrypted_shifted_powers
+e_p_as = e_as[3] - 3 * e_as[2] + 2 * e_as[1]
+assert e_p_s * a == e_p_as
+
+#############################
+
+# x^3 - 3x^2 + 2x
+main_poly = np.poly1d([1, -3, 2, 0])
+# (x - 1)(x - 2)
+target_poly = np.poly1d([1, -1]) * np.poly1d([1, -2])
+
+# Calculates polynomial h(x) = p(x) / t(x)
+cofactor, remainder = main_poly / target_poly
+assert remainder == np.poly1d([0])
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers):
+    coeffs = list(poly.coef)[::-1]
+    result = null
+    for power, coeff in zip(encrypted_powers, coeffs):
+        #print(coeff, power)
+        coeff = int(coeff)
+        # I have to do this for some strange reason
+        # Because if coeff is negative and I do += power * coeff
+        # then it gives me a different result than what I expect
+        if coeff < 0:
+            result -= power * (-coeff)
+        else:
+            result += power * coeff
+    return result
+
+encrypted_poly = evaluate(main_poly, encrypted_powers)
+assert encrypted_poly == e_p_s
+encrypted_cofactor = evaluate(cofactor, encrypted_powers)
+
+# Alpha shifted powers
+encrypted_shift_poly = evaluate(main_poly, encrypted_shifted_powers)
+
+# resulting g^p and g^h are provided to the verifier
+
+#################################
+# Verifier
+#################################
+
+# Last check that p = t(s) h
+
+assert encrypted_poly == encrypted_cofactor * target
+
+# Verify (g^p)^a == g^p'
+
+assert encrypted_poly * a == encrypted_shift_poly

scripts/zk/3.5-zero-knowledge.py

@@ -0,0 +1,129 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+import random
+import numpy as np
+
+# Section 3.5 from "Why and How zk-SNARK Works"
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#x = rand_scalar()
+#y = ec.y_for_x(x)
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(Q, 0), Fq(Q, 1), True, bls12381)
+assert g1 + null == g1
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+s = rand_scalar()
+
+# calculate the shift
+a = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+d = 10
+encrypted_powers = [
+    g1 * (s**i) for i in range(d)
+]
+encrypted_shifted_powers = [
+    g1 * (a * s**i) for i in range(d)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+target = (s - 1) * (s - 2)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+# delta shift
+delta = rand_scalar()
+
+# E(p(s)) = p(s)G
+#         = c_d s^d G + ... + c_1 s^1 G + c_0 s^0 G
+#         = s^3 G - 3 s^2 G + 2 s G
+# E(h(s)) = sG
+# t(s) = s^2 - 3s + 2
+# E(h(s)) t(s) = s^3 G - 3 s^2 G + 2 s G
+
+# Lets test these manually:
+
+e_s = encrypted_powers
+e_p_s = e_s[3] - 3 * e_s[2] + 2 * e_s[1]
+e_h_s = e_s[1]
+t_s = s**2 - 3*s + 2
+# exponentiate with delta
+e_p_s *= delta
+e_h_s *= delta
+assert t_s == target
+assert e_p_s == e_h_s * t_s
+
+e_as = encrypted_shifted_powers
+e_p_as = e_as[3] - 3 * e_as[2] + 2 * e_as[1]
+# exponentiate with delta
+e_p_as *= delta
+assert e_p_s * a == e_p_as
+
+#############################
+
+# x^3 - 3x^2 + 2x
+main_poly = np.poly1d([1, -3, 2, 0])
+# (x - 1)(x - 2)
+target_poly = np.poly1d([1, -1]) * np.poly1d([1, -2])
+
+# Calculates polynomial h(x) = p(x) / t(x)
+cofactor, remainder = main_poly / target_poly
+assert remainder == np.poly1d([0])
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers):
+    coeffs = list(poly.coef)[::-1]
+    result = null
+    for power, coeff in zip(encrypted_powers, coeffs):
+        #print(coeff, power)
+        coeff = int(coeff)
+        # I have to do this for some strange reason
+        # Because if coeff is negative and I do += power * coeff
+        # then it gives me a different result than what I expect
+        if coeff < 0:
+            result -= power * (-coeff)
+        else:
+            result += power * coeff
+    # Add delta to the result
+    # Free extra obfuscation to the polynomial
+    return result * delta
+
+encrypted_poly = evaluate(main_poly, encrypted_powers)
+assert encrypted_poly == e_p_s
+encrypted_cofactor = evaluate(cofactor, encrypted_powers)
+
+# Alpha shifted powers
+encrypted_shift_poly = evaluate(main_poly, encrypted_shifted_powers)
+
+# resulting g^p and g^h are provided to the verifier
+
+#################################
+# Verifier
+#################################
+
+# Last check that p = t(s) h
+
+assert encrypted_poly == encrypted_cofactor * target
+
+# Verify (g^p)^a == g^p'
+
+assert encrypted_poly * a == encrypted_shift_poly

scripts/zk/3.6-trusted-setup.py

@@ -0,0 +1,152 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+import random
+import numpy as np
+
+# Section 3.6 from "Why and How zk-SNARK Works"
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#x = rand_scalar()
+#y = ec.y_for_x(x)
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(Q, 0), Fq(Q, 1), True, bls12381)
+assert g1 + null == g1
+null2 = ec.AffinePoint(Fq2.zero(Q), Fq2.zero(Q), True, bls12381)
+assert null2 + g2 == g2
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+s = rand_scalar()
+
+# calculate the shift
+a = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+d = 10
+encrypted_powers = [
+    g1 * (s**i) for i in range(d)
+]
+encrypted_powers_g2 = [
+    g2 * (s**i) for i in range(d)
+]
+encrypted_shifted_powers = [
+    g1 * (a * s**i) for i in range(d)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+target = (s - 1) * (s - 2)
+# CRS = common reference string = trusted setup parameters
+target_crs = g1 * target
+alpha_crs = g2 * a
+
+# Proving key = (encrypted_powers, encrypted_shifted_powers)
+# Verify key = (target_crs, alpha_crs)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+# delta shift
+delta = rand_scalar()
+
+# E(p(s)) = p(s)G
+#         = c_d s^d G + ... + c_1 s^1 G + c_0 s^0 G
+#         = s^3 G - 3 s^2 G + 2 s G
+# E(h(s)) = sG
+# t(s) = s^2 - 3s + 2
+# E(h(s)) t(s) = s^3 G - 3 s^2 G + 2 s G
+
+# Lets test these manually:
+
+e_s = encrypted_powers
+e_p_s = e_s[3] - 3 * e_s[2] + 2 * e_s[1]
+e_h_s = e_s[1]
+t_s = s**2 - 3*s + 2
+# exponentiate with delta
+e_p_s *= delta
+e_h_s *= delta
+assert t_s == target
+assert e_p_s == e_h_s * t_s
+
+e_as = encrypted_shifted_powers
+e_p_as = e_as[3] - 3 * e_as[2] + 2 * e_as[1]
+# exponentiate with delta
+e_p_as *= delta
+assert e_p_s * a == e_p_as
+
+#############################
+
+# x^3 - 3x^2 + 2x
+main_poly = np.poly1d([1, -3, 2, 0])
+# (x - 1)(x - 2)
+target_poly = np.poly1d([1, -1]) * np.poly1d([1, -2])
+
+# Calculates polynomial h(x) = p(x) / t(x)
+cofactor, remainder = main_poly / target_poly
+assert remainder == np.poly1d([0])
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers, identity):
+    coeffs = list(poly.coef)[::-1]
+    result = identity
+    for power, coeff in zip(encrypted_powers, coeffs):
+        #print(coeff, power)
+        coeff = int(coeff)
+        # I have to do this for some strange reason
+        # Because if coeff is negative and I do += power * coeff
+        # then it gives me a different result than what I expect
+        if coeff < 0:
+            result -= power * (-coeff)
+        else:
+            result += power * coeff
+    # Add delta to the result
+    # Free extra obfuscation to the polynomial
+    return result * delta
+
+encrypted_poly = evaluate(main_poly, encrypted_powers, null)
+assert encrypted_poly == e_p_s
+encrypted_cofactor = evaluate(cofactor, encrypted_powers_g2, null2)
+
+# Alpha shifted powers
+encrypted_shift_poly = evaluate(main_poly, encrypted_shifted_powers, null)
+
+# resulting g^p and g^h are provided to the verifier
+
+# proof = (encrypted_poly, encrypted_cofactor, encrypted_shift_poly)
+
+#################################
+# Verifier
+#################################
+
+# Last check that p = t(s) h
+
+# Check polynomial cofactors:
+#assert encrypted_poly == encrypted_cofactor * target
+# e(g^p, g) == e(g^t, g^h)
+res1 = pairing.ate_pairing(encrypted_poly, g2)
+res2 = pairing.ate_pairing(target_crs, encrypted_cofactor)
+assert res1 == res2
+
+# Verify (g^p)^a == g^p'
+# Check polynomial restriction:
+
+res1 = pairing.ate_pairing(encrypted_shift_poly, g2)
+res2 = pairing.ate_pairing(encrypted_poly, alpha_crs)
+assert res1 == res2
+#assert encrypted_poly * a == encrypted_shift_poly
+

scripts/zk/4.4-proof-of-operation.py

@@ -0,0 +1,146 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+import random
+import numpy as np
+
+# Section 3.6 from "Why and How zk-SNARK Works"
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#x = rand_scalar()
+#y = ec.y_for_x(x)
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(Q, 0), Fq(Q, 1), True, bls12381)
+assert g1 + null == g1
+null2 = ec.AffinePoint(Fq2.zero(Q), Fq2.zero(Q), True, bls12381)
+assert null2 + g2 == g2
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+s = rand_scalar()
+
+# calculate the shift
+a = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+d = 10
+encrypted_powers = [
+    g1 * (s**i) for i in range(d)
+]
+encrypted_powers_g2 = [
+    g2 * (s**i) for i in range(d)
+]
+encrypted_shifted_powers = [
+    g1 * (a * s**i) for i in range(d)
+]
+encrypted_shifted_powers_g2 = [
+    g2 * (a * s**i) for i in range(d)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+target = (s - 1)
+# CRS = common reference string = trusted setup parameters
+target_crs = g1 * target
+alpha_crs = g2 * a
+alpha_crs_g1 = g1 * a
+
+# Proving key = (encrypted_powers, encrypted_shifted_powers)
+# Verify key = (target_crs, alpha_crs)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+left_poly = np.poly1d([3])
+right_poly = np.poly1d([2])
+out_poly = np.poly1d([6])
+
+# x^3 - 3x^2 + 2x
+main_poly = left_poly * right_poly - out_poly
+# (x - 1)
+target_poly = np.poly1d([1, -1])
+
+# Calculates polynomial h(x) = p(x) / t(x)
+cofactor, remainder = main_poly / target_poly
+assert remainder == np.poly1d([0])
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers, identity):
+    coeffs = list(poly.coef)[::-1]
+    result = identity
+    for power, coeff in zip(encrypted_powers, coeffs):
+        #print(coeff, power)
+        coeff = int(coeff)
+        # I have to do this for some strange reason
+        # Because if coeff is negative and I do += power * coeff
+        # then it gives me a different result than what I expect
+        if coeff < 0:
+            result -= power * (-coeff)
+        else:
+            result += power * coeff
+    return result
+
+assert left_poly * right_poly == out_poly
+
+encrypted_left_poly = evaluate(left_poly, encrypted_powers, null)
+encrypted_right_poly = evaluate(right_poly, encrypted_powers_g2, null2)
+encrypted_out_poly = evaluate(out_poly, encrypted_powers, null)
+
+#assert encrypted_poly == e_p_s
+encrypted_cofactor = evaluate(cofactor, encrypted_powers_g2, null2)
+
+# Alpha shifted powers
+encrypted_shift_left_poly = evaluate(left_poly, encrypted_shifted_powers, null)
+encrypted_shift_right_poly = evaluate(right_poly, encrypted_shifted_powers_g2, null2)
+encrypted_shift_out_poly = evaluate(out_poly, encrypted_shifted_powers, null)
+
+# resulting g^p and g^h are provided to the verifier
+
+# proof = (encrypted_poly, encrypted_cofactor, encrypted_shift_poly)
+
+#################################
+# Verifier
+#################################
+
+# Last check that p = t(s) h
+
+assert pairing.ate_pairing(2 * g1, g2) == pairing.ate_pairing(g1, g2) * pairing.ate_pairing(g1, g2)
+
+# Verify (g^p)^a == g^p'
+# Check polynomial restriction:
+
+def check_polynomial_restriction(encrypted_shift_poly, encrypted_poly):
+    res1 = pairing.ate_pairing(encrypted_shift_poly, g2)
+    res2 = pairing.ate_pairing(encrypted_poly, alpha_crs)
+    assert res1 == res2
+
+def check_polynomial_restriction_swapped(encrypted_shift_poly, encrypted_poly):
+    res1 = pairing.ate_pairing(g1, encrypted_shift_poly)
+    res2 = pairing.ate_pairing(alpha_crs_g1, encrypted_poly)
+    assert res1 == res2
+
+check_polynomial_restriction(encrypted_shift_left_poly, encrypted_left_poly)
+check_polynomial_restriction_swapped(encrypted_shift_right_poly, encrypted_right_poly)
+check_polynomial_restriction(encrypted_shift_out_poly, encrypted_out_poly)
+
+# Valid operation check
+# e(g^l, g^r) == e(g^t, g^h) * e(g^o, g)
+res1 = pairing.ate_pairing(encrypted_left_poly, encrypted_right_poly)
+res2 = pairing.ate_pairing(target_crs, encrypted_cofactor) * \
+       pairing.ate_pairing(encrypted_out_poly, g2)
+assert res1 == res2
+

scripts/zk/4.5.1-polynomial-interpolation.py

@@ -0,0 +1,30 @@
+import numpy as np
+
+def lagrange(points):
+    result = np.poly1d([0])
+    for i, (x_i, y_i) in enumerate(points):
+        poly = np.poly1d([y_i])
+        for j, (x_j, y_j) in enumerate(points):
+            if i == j:
+                continue
+            poly *= np.poly1d([1, -x_j]) / (x_i - x_j)
+        #print(poly)
+        #print(poly(1), poly(2), poly(3))
+        result += poly
+    return result
+
+left = lagrange([
+    (1, 2), (2, 2), (3, 6)
+])
+print(left)
+
+right = lagrange([
+    (1, 1), (2, 3), (3, 2)
+])
+print(right)
+
+out = lagrange([
+    (1, 2), (2, 6), (3, 12)
+])
+print(out)
+

scripts/zk/4.5.2-multi-operation-polynomials.py

@@ -0,0 +1,167 @@
+from bls_py import bls12381
+from bls_py import pairing
+from bls_py import ec
+from bls_py.fields import Fq, Fq2, Fq6, Fq12, bls12381_q as Q
+from finite_fields.modp import IntegersModP
+from finite_fields.polynomial import polynomialsOver
+import random
+
+n = bls12381.n
+
+g1 = ec.generator_Fq(bls12381)
+g2 = ec.generator_Fq2(bls12381)
+
+null = ec.AffinePoint(Fq(n, 0), Fq(n, 1), True, bls12381)
+assert null + g1 == g1
+null2 = ec.AffinePoint(Fq2.zero(n), Fq2.zero(n), True, bls12381)
+assert null2 + g2 == g2
+
+mod_field = IntegersModP(n)
+poly = polynomialsOver(mod_field).factory
+
+def lagrange(points):
+    result = poly([0])
+    for i, (x_i, y_i) in enumerate(points):
+        p = poly([y_i])
+        for j, (x_j, y_j) in enumerate(points):
+            if i == j:
+                continue
+            p *= poly([-x_j, 1]) / (x_i - x_j)
+        #print(poly)
+        #print(poly(1), poly(2), poly(3))
+        result += p
+    return result
+
+def poly_call(poly, x):
+    result = mod_field(0)
+    for degree, coeff in enumerate(poly):
+        result += coeff * (x**degree)
+    return result.n
+
+left_points = [
+    (1, 2), (2, 2), (3, 6)
+]
+left_poly = lagrange(left_points)
+#l = poly([2]) * poly([1, -1])
+print("Left:")
+print(left_poly)
+for x, y in left_points:
+    assert poly_call(left_poly, x) == y
+
+right_points = [
+    (1, 1), (2, 3), (3, 2)
+]
+right_poly = lagrange(right_points)
+print("Right:")
+print(right_poly)
+for x, y in right_points:
+    assert poly_call(right_poly, x) == y
+
+out_points = [
+    (1, 2), (2, 6), (3, 12)
+]
+out_poly = lagrange(out_points)
+print("Out:")
+print(out_poly)
+for x, y in out_points:
+    assert poly_call(out_poly, x) == y
+
+target_poly = poly([-1, 1]) * poly([-2, 1]) * poly([-3, 1])
+assert poly_call(target_poly, 1) == 0
+assert poly_call(target_poly, 2) == 0
+assert poly_call(target_poly, 3) == 0
+
+main_poly = left_poly * right_poly - out_poly
+cofactor_poly = main_poly / target_poly
+
+assert left_poly * right_poly - out_poly == target_poly * cofactor_poly
+
+def rand_scalar():
+    return random.randrange(1, bls12381.q)
+
+#################################
+# Verifier (trusted setup)
+#################################
+
+# samples a random value (a secret)
+toxic_scalar = rand_scalar()
+# calculate the shift
+alpha_shift = rand_scalar()
+
+# calculates encryptions of s for all powers i in 0 to d
+# E(s^i) = g^s^i
+degree = 10
+enc_s1 = [
+    g1 * (toxic_scalar**i) for i in range(degree)
+]
+enc_s2 = [
+    g2 * (toxic_scalar**i) for i in range(degree)
+]
+enc_s1_shift = [
+    g1 * (alpha_shift * toxic_scalar**i) for i in range(degree)
+]
+enc_s2_shift = [
+    g2 * (alpha_shift * toxic_scalar**i) for i in range(degree)
+]
+
+# evaluates unencrypted target polynomial with s: t(s)
+toxic_target = (toxic_scalar - 1) * (toxic_scalar - 2) * (toxic_scalar - 3)
+# CRS = common reference string = trusted setup parameters
+target_crs = g1 * toxic_target
+alpha_crs = g2 * alpha_shift
+alpha_crs_g1 = g1 * alpha_shift
+
+# Proving key = (encrypted_powers, encrypted_shifted_powers)
+# Verify key = (target_crs, alpha_crs)
+
+# encrypted values of s provided to the prover
+# Actual values of s are toxic waste and discarded
+
+#################################
+# Prover
+#################################
+
+# Using encrypted powers and coefficients, evaluates
+# E(p(s)) and E(h(s))
+def evaluate(poly, encrypted_powers, identity):
+    result = identity
+    for power, coeff in zip(encrypted_powers, poly):
+        result += power * coeff.n
+    return result
+
+enc_left = evaluate(left_poly, enc_s1, null)
+enc_right = evaluate(right_poly, enc_s2, null2)
+enc_out = evaluate(out_poly, enc_s1, null)
+
+enc_cofactor = evaluate(cofactor_poly, enc_s2, null2)
+
+# Alpha shifted powers
+enc_left_shift = evaluate(left_poly, enc_s1_shift, null)
+enc_right_shift = evaluate(right_poly, enc_s2_shift, null2)
+enc_out_shift = evaluate(out_poly, enc_s1_shift, null)
+
+#################################
+# Verifier
+#################################
+
+def restrict_polynomial_g1(encrypted_shift_poly, encrypted_poly):
+    res1 = pairing.ate_pairing(encrypted_shift_poly, g2)
+    res2 = pairing.ate_pairing(encrypted_poly, alpha_crs)
+    assert res1 == res2
+
+def restrict_polynomial_g2(encrypted_shift_poly, encrypted_poly):
+    res1 = pairing.ate_pairing(g1, encrypted_shift_poly)
+    res2 = pairing.ate_pairing(alpha_crs_g1, encrypted_poly)
+    assert res1 == res2
+
+restrict_polynomial_g1(enc_left_shift, enc_left)
+restrict_polynomial_g2(enc_right_shift, enc_right)
+restrict_polynomial_g1(enc_out_shift, enc_out)
+
+# Valid operation check
+# e(g^l, g^r) == e(g^t, g^h) * e(g^o, g)
+res1 = pairing.ate_pairing(enc_left, enc_right)
+res2 = pairing.ate_pairing(target_crs, enc_cofactor) * \
+       pairing.ate_pairing(enc_out, g2)
+assert res1 == res2
+