some refactor, useful additions, tests; as a preparation for dynamic proofs

README.md

@@ -30,4 +30,5 @@ at your choice.
 - [ ] add Groth16 notes
 - [ ] document the `snarkjs` circuit-specific setup `H` points convention
 - [x] precalculate stuff for "partial" proofs
+- [ ] implement Dynark style "dynamic proofs" too
 - [ ] make it work for different curves

groth16/bn128/arrays.nim

@@ -32,6 +32,22 @@ proc scaleFrSeqInPlace*(s: Fr[BN254_Snarks], arr: var seq[Fr[BN254_Snarks]] ) =
   for i in 0..<N:
     arr[i] *= s
 
+proc pointwiseProdFr*(xs: seq[Fr[BN254_Snarks]], ys: seq[Fr[BN254_Snarks]]): seq[Fr[BN254_Snarks]] =
+  let N = xs.len
+  assert( N == ys.len )
+  var zs : seq[Fr[BN254_Snarks]] = newSeq[Fr[BN254_Snarks]]( N )
+  for i in 0..<N:
+    zs[i] = xs[i] * ys[i]
+  return zs
+
+func dotProdFr*(xs, ys: seq[Fr[BN254_Snarks]]): Fr[BN254_Snarks] =
+  let n = xs.len
+  assert( n == ys.len, "dotProdFr: incompatible vector lengths" )
+  var s : Fr[BN254_Snarks] = zeroFr
+  for i in 0..<n:
+    s += xs[i] * ys[i]
+  return s
+
 #-------------------------------------------------------------------------------
 # G1 arrays
 
@@ -58,4 +74,12 @@ proc scaleG1SeqInPlace*(s: Fr[BN254_Snarks], arr: var seq[G1] ) =
   for i in 0..<N:
     arr[i] = s ** arr[i]
 
+proc pointwiseScaleG1*(xs: seq[Fr[BN254_Snarks]], gs: seq[G1]): seq[G1] =
+  let N = xs.len
+  assert( N == gs.len )
+  var hs : seq[G1] = newSeq[G1]( N )
+  for i in 0..<N:
+    hs[i] = xs[i] ** gs[i]
+  return hs
+
 #-------------------------------------------------------------------------------

groth16/bn128/fields.nim

@@ -7,7 +7,6 @@
 #
 
 import sugar
-
 import std/bitops
 # import std/sequtils
 
@@ -34,9 +33,13 @@ const primeP*  : B = fromHex( B, "0x30644e72e131a029b85045b68181585d97816a916871
 
 const zeroFp*   = fromHex( Fp[BN254_Snarks], "0x00" )
 const zeroFr*   = fromHex( Fr[BN254_Snarks], "0x00" )
+
 const oneFp*    = fromHex( Fp[BN254_Snarks], "0x01" )
 const oneFr*    = fromHex( Fr[BN254_Snarks], "0x01" )
 
+const twoFp*    = fromHex( Fp[BN254_Snarks], "0x02" )
+const twoFr*    = fromHex( Fr[BN254_Snarks], "0x02" )
+
 const zeroFp2* = mkFp2( zeroFp, zeroFp )
 const oneFp2*  = mkFp2( oneFp , zeroFp )
 

groth16/dynamic/types.nim

@@ -0,0 +1,53 @@
+
+{.push raises:[].}
+
+#import constantine/named/properties_fields
+
+import groth16/bn128
+import groth16/math/domain
+#import groth16/math/ntt
+#import groth16/math/poly
+
+import groth16/partial/types
+export types
+
+#-------------------------------------------------------------------------------
+
+type
+
+  # note: already exported from partial/types, and Nim is stupid as usual
+  # F* = Fr[BN254_Snarks]
+
+  PowersOfTau* = object
+    elements* : seq[G1]
+
+  LagrangeOfTau* = object
+    domain*   : Domain
+    elements* : seq[G1]
+
+#-------------------------------------------------------------------------------
+
+type 
+
+  # things we can compute at circuit setup time
+  DynaSetupV1* = object
+    pointsDeltaL* : seq[G1]      # the points `delta^-1 * L_i(tau) * Z(tau) * g1`
+    weightVec*    : seq[F]       # the weights `W_k = 1/N/(omega^-k - 1)`
+    wConvDeltaL*  : seq[G1]      # the convolution of `W` and `deltaL`
+
+  # things we can compute from the partial witness
+  DynaPreprocess* = object
+    projA0*       : seq[G1]      # the points `delta^-1 * A0(tau) * L_i(tau) * g1` (well, "modulo Z(tau)")
+    projB0*       : seq[G1]      # the same for B0
+
+  DynaPreProof* = object
+    partialProof*   : PartialProof
+    dyanPreprocess* : DynaPreprocess
+
+  PartialAB* = object
+    valuesAz*    : seq[F]
+    valuesBz*    : seq[F]
+    complImageA* : seq[bool]     # image of the complement of the partial witness under A
+    complImageB* : seq[bool]     # same for B
+
+#-------------------------------------------------------------------------------

groth16/fake_setup.nim

@@ -13,6 +13,7 @@ import constantine/math/arithmetic
 import constantine/named/properties_fields
 
 import groth16/bn128
+#import groth16/bn128/arrays
 import groth16/math/domain
 import groth16/math/poly
 import groth16/zkey_types
@@ -189,16 +190,6 @@ func r1csToSparseMatrices*(r1cs: R1CS): SparseMatrices =
 
 #-------------------------------------------------------------------------------
 
-func dotProdFr(xs, ys: seq[Fr[BN254_Snarks]]): Fr[BN254_Snarks] =
-  let n = xs.len
-  assert( n == ys.len, "dotProdFr: incompatible vector lengths" )
-  var s : Fr[BN254_Snarks] = zeroFr
-  for i in 0..<n:
-    s += xs[i] * ys[i]
-  return s
-
-#-------------------------------------------------------------------------------
-
 func fakeCircuitSetup*(r1cs: R1CS, toxic: ToxicWaste, flavour=Snarkjs): ZKey = 
  
   let neqs       = r1cs.constraints.len

groth16/math/convert.nim

@@ -0,0 +1,151 @@
+
+#
+# convert between powers and Lagrange bases representations
+#
+
+{.push raises:[].}
+
+import std/sugar
+#import std/sequtils
+
+import constantine/math/arithmetic
+#import constantine/named/properties_fields
+
+import groth16/bn128
+import groth16/bn128/arrays
+
+import groth16/math/domain
+import groth16/math/ntt
+import groth16/math/group_fft
+import groth16/math/poly
+
+#import groth16/zkey_types
+import groth16/dynamic/types
+
+#-------------------------------------------------------------------------------
+# Lagrange basis version of powers-of-tau
+
+# convert the power basis to Lagrange basis
+func powersToLagrange*(inp: PowersOfTau): LagrangeOfTau = 
+  let D  = createDomain(inp.elements.len)
+  let hs = inverseGroupFFT(inp.elements , D)
+  return LagrangeOfTau(domain: D, elements:hs)
+
+#-------------------------------------------------------------------------------
+# Convert the "H" points between Jordi's and the Groth16 paper version
+
+#[
+
+the following should hold, because they have both detree 2N-1, 
+and they also take the same values on the 2N sized subgroup:
+
+  L_{2i+1}(x)  =  Z(x) * L_i( eta^-1 * x )
+
+from this fact we should be able to convert?
+
+]#
+
+func convertPointsFromJordi*( D: Domain, pointsJordi: seq[G1] ): seq[G1] =
+  let N   : int    = D.domainSize
+  let D2  : Domain = createDomain( 2*N )
+  let eta : F      = D2.domainGen
+
+  # do an FFT
+  var arr = forwardGroupFFT( pointsJordi , D )
+
+  # multiply by powers of eta + a constant -2
+  var s : F = negFr(twoFr) 
+  for i in 0..<N:
+    arr[i] = s ** arr[i] 
+    s *= eta
+
+  return arr
+
+func convertPointsToJordi*( D: Domain, pointsJens: seq[G1] ): seq[G1] =
+  let N   : int    = D.domainSize
+  let D2  : Domain = createDomain( 2*N )
+  let eta    : F   = D2.domainGen
+  let etaInv : F   = invFr(eta)
+
+  # multiply by powers of eta + a constant -1/2
+  var arr : seq[G1] = newSeq[G1]( N )
+  var s : F = negFr(oneHalfFr) 
+  for i in 0..<N:
+    arr[i] = s ** pointsJens[i] 
+    s *= etaInv
+
+  # do an IFFT
+  return inverseGroupFFT( arr , D )
+
+proc testJordiConversion*(N: int, tau: F, delta: F): bool =
+
+  let D  : Domain = createDomain(   N )
+  let D2 : Domain = createDomain( 2*N )
+
+  let deltaInv : F  = invFr(delta)
+  let ztauG1   : G1 = (smallPowFr(tau,N) - oneFr) ** gen1     # (tau^N - 1)
+
+  #
+  # in the original paper, these are the curve points
+  #   [ delta^-1 * tau^i * Z(tau) ] 
+  #
+  let pointsJens: seq[G1] = collect( newSeq , (for i in 0..<N: 
+        (deltaInv * smallPowFr(tau,i)) ** ztauG1 ))
+
+  #
+  # in the Snarkjs implementation, these are the curve points
+  #   [ delta^-1 * L_{2i+1} (tau) ]
+  # where L_k are the Lagrange polynomials on the refined domain
+  #
+  let pointsJordi: seq[G1] = collect( newSeq , (for i in 0..<N: 
+        (deltaInv * evalLagrangePolyAt(D2, 2*i+1, tau)) ** gen1 ))
+
+  let jensFromJordi = convertPointsFromJordi( D , pointsJordi )
+  let jordiFromJens = convertPointsToJordi(   D , pointsJens  )
+  
+  let ok1 = isEqualG1Seq( pointsJens  , jensFromJordi )
+  let ok2 = isEqualG1Seq( pointsJordi , jordiFromJens )
+
+  # echo "fromJordi = " & ($ok1)
+  # echo "toJordi   = " & ($ok2)
+
+  return (ok1 and ok2)
+
+#-------------------------------------------------------------------------------
+# for testing purposes, computing the basics from tau
+
+func computePowersOfScalar*(N: int, tau: F): seq[F] =
+  var ts: seq[F] = newSeq[F]( N )
+  var s: F = oneFr
+  for i in 0..<N:
+    ts[i] = s
+    s *= tau
+  return ts
+
+func computePowersOfTau*(N: int, tau: F): PowersOfTau =
+  var hs: seq[G1] = newSeq[G1]( N )
+  var s: F = oneFr
+  for i in 0..<N:
+    hs[i] = s ** gen1
+    s *= tau
+  return PowersOfTau(elements: hs)
+
+# direct computation
+func computeLagrangeOfTauV1*(D: Domain, tau: F): LagrangeOfTau =
+  let N  = D.domainSize
+  var hs: seq[G1] = newSeq[G1]( N )
+  for k in 0..<N:
+    hs[k] = evalLagrangePolyAt(D, k, tau) ** gen1
+  return LagrangeOfTau(domain: D, elements:hs)
+
+# via Fourier transform
+func computeLagrangeOfTauV2*(D: Domain, tau: F): LagrangeOfTau =
+  let N  = D.domainSize
+  let ts = computePowersOfScalar(N, tau) 
+  let ls = inverseNTT(ts, D)
+  var hs: seq[G1] = newSeq[G1]( N )
+  for i in 0..<N:
+    hs[i] = ls[i] ** gen1
+  return LagrangeOfTau(domain: D, elements:hs)
+
+#-------------------------------------------------------------------------------

groth16/math/convolution.nim

@@ -0,0 +1,70 @@
+
+import constantine/math/arithmetic
+import constantine/named/properties_fields
+
+import groth16/bn128
+import groth16/bn128/curves
+import groth16/bn128/arrays
+
+import groth16/math/domain
+import groth16/math/ntt
+import groth16/math/group_fft
+
+#-------------------------------------------------------------------------------
+
+type F = Fr[BN254_Snarks]
+
+func fieldConvolution*(xs: seq[F], ys: seq[F]): seq[F] = 
+  let N = xs.len
+  assert( N == ys.len )
+  let D = createDomain( N )
+
+  let us = forwardNTT(xs, D) 
+  let vs = forwardNTT(ys, D)
+  let ws = pointwiseProdFr(us, vs)
+
+  return inverseNTT(ws, D)
+
+#-------------------------------------------------------------------------------
+
+func groupConvolution*(xs: seq[F], gs: seq[G1]): seq[G1] =
+  let N = xs.len
+  assert( N == gs.len )
+  let D = createDomain( N )
+
+  let us =      forwardNTT(xs, D) 
+  let hs = forwardGroupFFT(gs, D)
+  let rs = pointwiseScaleG1(us, hs)
+
+  return inverseGroupFFT(rs, D)
+
+#-------------------------------------------------------------------------------
+# naive reference implementations for testing
+
+func naiveFieldConvolution*(xs: seq[F], ys: seq[F]): seq[F] = 
+  let N = xs.len
+  assert( N == ys.len )
+
+  var zs: seq[F] = newSeq[F]( N )
+  for k in 0..<N:
+    var acc: F = zeroFr 
+    for i in 0..<N:
+      acc += xs[ (k-i+N) mod N ] * ys[ i ]
+    zs[k] = acc
+
+  return zs  
+
+func naiveGroupConvolution*(xs: seq[F], gs: seq[G1]): seq[G1] = 
+  let N = xs.len
+  assert( N == gs.len )
+
+  var hs: seq[G1] = newSeq[G1]( N )
+  for k in 0..<N:
+    var acc: G1 = infG1 
+    for i in 0..<N:
+      acc += xs[ (k-i+N) mod N ] ** gs[ i ]
+    hs[k] = acc
+
+  return hs
+
+#-------------------------------------------------------------------------------

groth16/partial/README.md

@@ -0,0 +1,16 @@
+
+Partial proofs
+--------------
+
+If a large part of your witness is constant (or changing infrequently),
+you can precalculate a significant portion of the proof.
+
+See [this forum post](https://forum.research.logos.co/t/speeding-up-rln-proofs/663) 
+for more details.
+
+This is essentially free (the overhead is minimal, so already when doing 2 proofs
+with same shared part of the witness, it's worth do it).
+
+A more advanced version, but also with a more heave trade-off, is developed in
+[`groth16/dynamic`](../dynamic/).
+

groth16/partial/precalc.nim

@@ -6,7 +6,7 @@
 # between proofs - then such precalculation can result in a potentially big speedup.
 #
 # an example use is RLN proofs, where the circuit is dominated by the Merkle 
-#vinclusion proof check, which is expected to change relatively rarely.
+# inclusion proof check, which is expected to change relatively rarely.
 #
 
 {.push raises:[].}

groth16/partial/types.nim

@@ -30,3 +30,12 @@ type
 
 func makePartialWitness*(vals: seq[Option[Fr[BN254_Snarks]]]): PartialWitness =
   return PartialWitness(values: vals)
+
+func partialWitnessMask*(pw: PartialWitness): seq[bool] = 
+  let N = pw.values.len
+  var bs: seq[bool] = newSeq[bool]( N )
+  for i in 0..<N:
+    bs[i] = isSome(pw.values[i])
+  return bs 
+
+#-------------------------------------------------------------------------------

tests/groth16/testFFT.nim

@@ -16,6 +16,7 @@ import groth16/bn128/arrays
 import groth16/math/domain
 import groth16/math/ntt
 import groth16/math/group_fft
+import groth16/math/convolution
 
 #-------------------------------------------------------------------------------
 
@@ -69,4 +70,18 @@ suite "group FFT checks (for the group G1)":
 
 #-------------------------------------------------------------------------------
   
+suite "convolution checks":
 
+  let N = 64
+
+  let xs = randFrSeq(N)
+  let ys = randFrSeq(N)
+  let gs = randG1Seq(N)
+
+  test "FFT field convolutions vs. slow reference":
+    check isEqualFrSeq( fieldConvolution(xs,ys) , naiveFieldConvolution(xs,ys) )
+
+  test "FFT group convolutions vs. slow reference":
+    check isEqualG1Seq( groupConvolution(xs,gs) , naiveGroupConvolution(xs,gs) )
+
+#-------------------------------------------------------------------------------

tests/groth16/testLagrangeTau.nim

@@ -0,0 +1,56 @@
+
+{.used.} 
+
+import std/unittest
+
+# import constantine/math/arithmetic
+# import constantine/named/properties_fields
+
+#import groth16/bn128/fields
+#import groth16/bn128/curves
+import groth16/bn128/arrays
+import groth16/bn128/rnd
+# import groth16/bn128/debug
+
+import groth16/math/domain
+import groth16/math/convert
+#import groth16/math/ntt
+#import groth16/math/group_fft
+
+import groth16/dynamic/types
+# import groth16/dynamic/setup
+
+#-------------------------------------------------------------------------------
+
+suite "Lagrange-Tau":
+
+  let N   : int    = 128
+  let D   : Domain = createDomain( N )
+  let tau : F      = randFr()
+
+  let powers_of_tau = computePowersOfTau(N , tau)
+  let lagrange_v0   = powersToLagrange( powers_of_tau )
+  let lagrange_v1   = computeLagrangeOfTauV1(D , tau )
+  let lagrange_v2   = computeLagrangeOfTauV2(D , tau )
+
+  test "converting from powers of tau == computing directly via evaluation formula":
+    check isEqualG1Seq( lagrange_v0.elements, lagrange_v1.elements )
+
+  test "converting from powers of tau == computing via scalar FFT":
+    check isEqualG1Seq( lagrange_v0.elements, lagrange_v2.elements )
+
+  test "computing directly via evaluation formulaa == computing via scalar FFT":
+    check isEqualG1Seq( lagrange_v1.elements, lagrange_v2.elements )
+
+#-------------------------------------------------------------------------------
+
+suite "to/from Jordi":
+
+  let N = 128
+  let tau   = randFr()
+  let delta = randFr()
+
+  test "to/from Jordi":  
+    check testJordiConversion(N, tau, delta)
+
+#-------------------------------------------------------------------------------

tests/groth16/testPtCompression.nim

@@ -11,7 +11,6 @@ import constantine/named/properties_fields
 import groth16/bn128/fields
 import groth16/bn128/curves
 import groth16/bn128/rnd
-
 # import groth16/bn128/debug
 
 #-------------------------------------------------------------------------------

tests/test.nim

@@ -1,9 +1,10 @@
 
 import ./groth16/testField
+import ./groth16/testCurve
 import ./groth16/testPoly
 import ./groth16/testFFT
+import ./groth16/testLagrangeTau
 import ./groth16/testPtCompression
-import ./groth16/testCurve
 import ./groth16/testProver
 import ./groth16/testMultithreading