Conversation
|
I just had another look at SnarkJS to see if I can figure out the double Montgomery encoding situation. As it turns out at least I could identify the SnarkJS code that writes the coefficients as doubly Montgomery encoded. This happens here: const nR2 = curve.Fr.mul(n, R2r);
curve.Fr.toRprLE(buffCoeff, 12, nR2);where const R2r = curve.Fr.e(Scalar.mod(Scalar.mul(Rr,Rr), primeR));Which at the very least finally ends the question of whether I'm imagining things or not. I still have no clue why this is done to be honest. The commits touching these lines are not very illuminating either. |
2be3509 to
47ba6f9
Compare
Otherwise when performing operations on EC points in Jacobian, depending on import these can be missing.
Previously these were generated from the wrong branch of the
`{.booldefine.}` variables we had in the code when the PR was still a
draft. Instead of parsing the data as doubly Montgomery encoded as
they actually are, we parsed them as regular Montgomery encoded
resulting in wrong test vectors.
29f1a4e to
8efc5bb
Compare
|
I guess the code at the moment only works on the 64-bit backend. That rings a bell, it was probably due to the parsing of the binary files. Need to check. |
mratsim
left a comment
There was a problem hiding this comment.
LGTM.
I have added several notes for further refactoring and will follow-up with an issue once that is merged. I think the PR is too big at the moment to add further non-critical fix.
Critical if possible would be to fix the 32-bit tests.
The FFT part is not worth changing at the moment but is important to merge now because the following future work use FFT:
- FRI for zkVM acceleration
- Whir for Lean Ethereum
- PeerDAS. Ethereum 2D KZG commitments with Error Correcting Codes
| result &= "\n" & sp & ")" | ||
|
|
||
| func toDecimal*[EC: EC_ShortW_Prj or EC_ShortW_Jac or EC_ShortW_Aff or EC_ShortW_JacExt](P: EC, indent: static int = 0): string = | ||
| ## Stringify an elliptic curve point to Hex |
There was a problem hiding this comment.
| ## Stringify an elliptic curve point to Hex | |
| ## Stringify an elliptic curve point to decimal |
| ## Note. Leading zeros are not removed. | ||
| ## Output as decimal. | ||
| ## | ||
| ## WARNING: NOT constant time! |
There was a problem hiding this comment.
Actually it is constant-time see low-level impl:
constantine/constantine/math/io/io_bigints.nim
Lines 371 to 385 in b3f4ebd
The mid-level impl is wrong:
constantine/constantine/math/io/io_fields.nim
Lines 110 to 117 in b3f4ebd
But it's fine, afaik, it's only used for debugging purposes.
| accum.add toDecimal(f) | ||
|
|
||
| func appendDecimal*(accum: var string, f: ExtensionField, indent = 0, order: static Endianness = bigEndian) = | ||
| ## Stringify a tower field element to hex. |
There was a problem hiding this comment.
| ## Stringify a tower field element to hex. | |
| ## Stringify a tower field element to decimal. |
| for i in 0 ..< half: | ||
| # FFT Butterfly | ||
| y_times_root .scalarMul_vartime(output[i+half], rootsOfUnity[i]) | ||
| y_times_root .scalarMul_vartime(rootsOfUnity[i], output[i+half]) |
There was a problem hiding this comment.
Curious the PR that introduced scalarMul_vartime also rewrote this part of the FFT and the argument order should be scalarMul_vartime(EC, scalar).
There was a problem hiding this comment.
Ah my mistake, the argument order when out-of-place is correct. (EC out, scalar, EC in)
| FFTS_SizeNotPowerOfTwo = "Input must be of a power of 2 length" | ||
|
|
||
| FFT_Descriptor*[F] = object # `F` is either `Fp[Name]` or `Fr[Name]` | ||
| ## Metadata for FFT on Elliptic Curve |
There was a problem hiding this comment.
on fields.
Note that there was an implementation here: https://github.com/mratsim/constantine/blob/v0.2.0/research/kzg/fft_fr.nim
| beta2*: EC_ShortW_Aff[Fp2[Name], G2] ## g₂^β | ||
| gamma2*: EC_ShortW_Aff[Fp2[Name], G2] ## g₂^γ | ||
| delta2*: EC_ShortW_Aff[Fp2[Name], G2] ## g₂^δ | ||
| gamma_abc*: seq[EC_ShortW_Aff[Fp[Name], G1]] ## == `zkey.ic` == `[g₁^{ \frac{β·A_i(τ) + α·B_i(τ) + C_0(τ)}{ γ }}]` |
There was a problem hiding this comment.
TODO: we will probably want to remove all seqs here in the future to unify with GPU impl.
Also we just don't use seq in Constantine except at the very edge for easy integration with Nim
| r1cs: r1cs | ||
| ) | ||
|
|
||
| proc calcAp[Name: static Algebra](ctx: Groth16Prover[Name], wt: seq[Fr[Name]]): EC_ShortW_Jac[Fp[Name], G1] {.noinit.} = |
There was a problem hiding this comment.
All non-var parameters should be openArray instead of seq
| result = B1_p | ||
|
|
||
|
|
||
| proc buildABC[Name: static Algebra](ctx: Groth16Prover[Name], wt: seq[Fr[Name]]): tuple[A, B, C: seq[Fr[Name]]] = |
There was a problem hiding this comment.
TODO: remove seq creation and prefer the caller to pass allocated buffers if possible
| proc transform[Name: static Algebra](args: seq[Fr[Name]], inc: Fr[Name]): seq[Fr[Name]] = | ||
| ## Applies (multiplies) increasing powers of `inc` to each element | ||
| ## of `args`, i.e. | ||
| ## | ||
| ## `{ a[0], a[1]·inc, a[2]·inc², a[3]·inc³, ... }`. | ||
| ## | ||
| ## In our case `inc` is usually a root of unity of the power given by | ||
| ## `log2( FFT order ) + 1`. | ||
| result = newSeq[Fr[Name]](args.len) | ||
| var cur = Fr[Name].fromUint(1.uint64) | ||
| for i in 0 ..< args.len: | ||
| result[i].prod(args[i], cur) | ||
| cur *= inc |
There was a problem hiding this comment.
Note: This patterns also exist in KZG and IPA and other commitments:
See
constantine/tests/t_ethereum_verkle_ipa_test_helper.nim
Lines 36 to 52 in b3f4ebd
So there is a computePowers proc, though it does allocate more.
constantine/constantine/math/arithmetic/finite_fields.nim
Lines 841 to 854 in b3f4ebd
| pairing(proof.C, vk.delta2) | ||
|
|
||
| # 4. Check pairing equality | ||
| result = bool(lhs == rhs) |
There was a problem hiding this comment.
TODO: This can be rewritten by negating either proof.A or proof.B then use multipairing and check that the result is the pairing group neutral element (i.e. 1 in GT/Fp12)
Greptile SummaryThis PR adds a full Groth16 prover and verifier to Constantine, including binary parsers for SnarkJS-produced Key issues found:
Confidence Score: 2/5
Important Files Changed
Sequence DiagramsequenceDiagram
participant User
participant Groth16Prover
participant BinaryParsers
participant FFT as FFT (fft_fields)
participant MSM as MSM (ec_multi_scalar_mul)
participant Pairing as Pairings
User->>BinaryParsers: parseZkeyFile() → ZkeyBin.toZkey[Name]()
User->>BinaryParsers: parseWtnsFile() → WtnsBin.toWtns[Name]()
User->>BinaryParsers: parseR1csFile() → R1csBin.toR1CS()
User->>Groth16Prover: init(zkey, wtns, r1cs)
User->>Groth16Prover: prove() [samples r, s]
Groth16Prover->>MSM: calcAp: multiScalarMul(witnesses, A[])
Groth16Prover->>MSM: calcBp: multiScalarMul(witnesses, B2[])
Groth16Prover->>MSM: calcB1: multiScalarMul(witnesses, B1[])
Groth16Prover->>Groth16Prover: buildABC(coeffs, witnesses)
Groth16Prover->>FFT: itf(A_poly): ifft → transform → fft
Groth16Prover->>FFT: itf(B_poly): ifft → transform → fft
Groth16Prover->>FFT: itf(C_poly): ifft → transform → fft
Groth16Prover->>MSM: calcCp: multiScalarMul(priv, C[]) + multiScalarMul(jabc, H[])
Groth16Prover-->>User: Groth16Proof(A, B, C)
User->>Pairing: verify(vk, proof, publicInputs)
Pairing->>MSM: multiScalarMul(publicInputs, IC[])
Pairing->>Pairing: e(A,B) == e(α,β)·e(I,γ)·e(C,δ)
Pairing-->>User: bool
|
| proc prove*[Name: static Algebra](ctx: Groth16Prover[Name]): Groth16Proof {.noinit.} = | ||
| ## Generate a proof given the Groth16 prover context data. | ||
| ## | ||
| ## This implies calculating the proof elements `π = (g₁^A, g₁^C, g₂^B)` | ||
| ## | ||
| ## See `calcAp`, `calcBp` and `calcCp` on how these elements are computed. | ||
| # 1. Sample the random field elements `r` and `s` for the proof | ||
| ctx.r = randomFieldElement(Fr[Name]) | ||
| ctx.s = randomFieldElement(Fr[Name]) | ||
| # 2. get the witness data needed for all proof elements | ||
| let wt = ctx.wtns.witnesses | ||
| # 3. compute the individual proof elements | ||
| let A_p = ctx.calcAp(wt) | ||
| let B2_p = ctx.calcBp(wt) | ||
| let B1_p = ctx.calcB1(wt) | ||
| let C_p = ctx.calcCp(A_p, B1_p, wt) | ||
|
|
||
| result = Groth16Proof(A: A_p.getAffine(), B: B2_p.getAffine(), C: C_p.getAffine()) |
There was a problem hiding this comment.
Missing generic parameter in prove* return type and constructor
The return type is declared as Groth16Proof but Groth16Proof is a generic type that requires the Name parameter. Without it the declaration is incomplete and would fail to compile. The same issue appears in the result constructor on line 255.
Additionally, ctx is an immutable parameter (no var) yet the body assigns to ctx.r and ctx.s on lines 245–246 — Nim does not allow mutating a non-var parameter, which is another compile error.
| proc prove*[Name: static Algebra](ctx: Groth16Prover[Name]): Groth16Proof {.noinit.} = | |
| ## Generate a proof given the Groth16 prover context data. | |
| ## | |
| ## This implies calculating the proof elements `π = (g₁^A, g₁^C, g₂^B)` | |
| ## | |
| ## See `calcAp`, `calcBp` and `calcCp` on how these elements are computed. | |
| # 1. Sample the random field elements `r` and `s` for the proof | |
| ctx.r = randomFieldElement(Fr[Name]) | |
| ctx.s = randomFieldElement(Fr[Name]) | |
| # 2. get the witness data needed for all proof elements | |
| let wt = ctx.wtns.witnesses | |
| # 3. compute the individual proof elements | |
| let A_p = ctx.calcAp(wt) | |
| let B2_p = ctx.calcBp(wt) | |
| let B1_p = ctx.calcB1(wt) | |
| let C_p = ctx.calcCp(A_p, B1_p, wt) | |
| result = Groth16Proof(A: A_p.getAffine(), B: B2_p.getAffine(), C: C_p.getAffine()) | |
| proc prove*[Name: static Algebra](ctx: var Groth16Prover[Name]): Groth16Proof[Name] {.noinit.} = |
And on line 255:
| proc prove*[Name: static Algebra](ctx: Groth16Prover[Name]): Groth16Proof {.noinit.} = | |
| ## Generate a proof given the Groth16 prover context data. | |
| ## | |
| ## This implies calculating the proof elements `π = (g₁^A, g₁^C, g₂^B)` | |
| ## | |
| ## See `calcAp`, `calcBp` and `calcCp` on how these elements are computed. | |
| # 1. Sample the random field elements `r` and `s` for the proof | |
| ctx.r = randomFieldElement(Fr[Name]) | |
| ctx.s = randomFieldElement(Fr[Name]) | |
| # 2. get the witness data needed for all proof elements | |
| let wt = ctx.wtns.witnesses | |
| # 3. compute the individual proof elements | |
| let A_p = ctx.calcAp(wt) | |
| let B2_p = ctx.calcBp(wt) | |
| let B1_p = ctx.calcB1(wt) | |
| let C_p = ctx.calcCp(A_p, B1_p, wt) | |
| result = Groth16Proof(A: A_p.getAffine(), B: B2_p.getAffine(), C: C_p.getAffine()) | |
| result = Groth16Proof[Name](A: A_p.getAffine(), B: B2_p.getAffine(), C: C_p.getAffine()) |
| let fftD = FFTDescriptor[Fr[Name]].init(abc[0].len * 2) | ||
|
|
||
| let A = itf(abc[0]) | ||
| let B = itf(abc[1]) | ||
| let C = itf(abc[2]) |
There was a problem hiding this comment.
fftD is allocated but never used or freed — memory leak
FFTDescriptor allocates heap memory for rootsOfUnity (via allocHeapArrayAligned). The descriptor fftD is initialised on line 194 but then the actual FFT work is performed by itf, which internally calls the convenience wrappers ifft_vartime/fft_vartime — each of which creates, uses, and deletes its own descriptor via defer: fftDesc.delete().
As a result fftD is:
- Never used — the
itfcalls ignore it entirely. - Never freed —
fftD.delete()is never called, leaking the allocated memory every timecalcCpruns.
The line should be removed (the FFT descriptors are managed inside itf already).
| let fftD = FFTDescriptor[Fr[Name]].init(abc[0].len * 2) | |
| let A = itf(abc[0]) | |
| let B = itf(abc[1]) | |
| let C = itf(abc[2]) |
| proc calcB1[Name: static Algebra](ctx: Groth16Prover[Name], wt: seq[Fr[Name]]): EC_ShortW_Jac[Fp[Name], G1] {.noinit.} = | ||
| let g16h = ctx.zkey.g16h | ||
|
|
||
| let beta1 = g16h.beta1 | ||
| let delta1 = g16h.delta1 | ||
| result = beta1.getJacobian + ctx.s * delta1 | ||
|
|
||
| # Get the B1 data | ||
| let Bs = ctx.zkey.B1 | ||
| doAssert Bs.len == wt.len | ||
| var B1_p {.noinit.}: EC_ShortW_Jac[Fp[Name], G1] | ||
| B1_p.multiScalarMul_vartime(wt, Bs) | ||
| # Add the independent terms, `β₁ + [s] · δ₁` | ||
| B1_p += beta1.getJacobian + ctx.s * delta1 | ||
|
|
||
| result = B1_p |
There was a problem hiding this comment.
First result assignment in calcB1 is dead code
Line 104 assigns result = beta1.getJacobian + ctx.s * delta1. However, result is unconditionally overwritten on line 114 (result = B1_p), so the first assignment is wasted computation. The correct formula is already applied via B1_p += beta1.getJacobian + ctx.s * delta1 on line 112, so the first line should be removed to avoid confusion.
| proc calcB1[Name: static Algebra](ctx: Groth16Prover[Name], wt: seq[Fr[Name]]): EC_ShortW_Jac[Fp[Name], G1] {.noinit.} = | |
| let g16h = ctx.zkey.g16h | |
| let beta1 = g16h.beta1 | |
| let delta1 = g16h.delta1 | |
| result = beta1.getJacobian + ctx.s * delta1 | |
| # Get the B1 data | |
| let Bs = ctx.zkey.B1 | |
| doAssert Bs.len == wt.len | |
| var B1_p {.noinit.}: EC_ShortW_Jac[Fp[Name], G1] | |
| B1_p.multiScalarMul_vartime(wt, Bs) | |
| # Add the independent terms, `β₁ + [s] · δ₁` | |
| B1_p += beta1.getJacobian + ctx.s * delta1 | |
| result = B1_p | |
| proc calcB1[Name: static Algebra](ctx: Groth16Prover[Name], wt: seq[Fr[Name]]): EC_ShortW_Jac[Fp[Name], G1] {.noinit.} = | |
| let g16h = ctx.zkey.g16h | |
| let beta1 = g16h.beta1 | |
| let delta1 = g16h.delta1 | |
| # Get the B1 data | |
| let Bs = ctx.zkey.B1 | |
| doAssert Bs.len == wt.len | |
| var B1_p {.noinit.}: EC_ShortW_Jac[Fp[Name], G1] | |
| B1_p.multiScalarMul_vartime(wt, Bs) | |
| # Add the independent terms, `β₁ + [s] · δ₁` | |
| B1_p += beta1.getJacobian + ctx.s * delta1 | |
| result = B1_p |
| c = coeffs.cs[i].section | ||
| s = coeffs.cs[i].index | ||
| coef = coeffs.cs[i].value | ||
| assert s.int < wt.len |
There was a problem hiding this comment.
assert is silently disabled in release builds
assert in Nim is equivalent to the C assert() macro — it is compiled away when the program is built with -d:release or -d:danger. This bounds check would then be silently skipped in production, risking an out-of-bounds access. Use doAssert instead, which is unconditional.
| assert s.int < wt.len | |
| doAssert s.int < wt.len |
|
|
||
| proc toEcG2*[Name: static Algebra](s: seq[byte]): EC_ShortW_Aff[Fp2[Name], G2] = | ||
| let x = toFp2[Name](s[0 .. 63]) | ||
| let y = toFp2[Name](s[64 .. ^1]) | ||
| result.x = x | ||
| result.y = y | ||
| if not bool(result.isNeutral()): |
There was a problem hiding this comment.
Hardcoded 32-byte slice sizes break for other field sizes
toEcG1 slices the input as s[0 .. 31] and s[32 .. ^1], hardcoding the assumption that each field element is exactly 32 bytes. This is correct for BN254 (256-bit prime → 32 bytes) but silently produces wrong results for any curve whose base field element has a different byte length (e.g., BLS12-381 with 48-byte Fp elements). The same hardcoding exists in toFp2 (x[0..31], x[32..63]) and toEcG2 (s[0..63], s[64..^1]).
Consider deriving the slice boundaries from the curve's actual field size, e.g. Fp[Name].bits().ceilDiv_vartime(8), instead of the magic constant 32.
Edit: The PR previously contained a long description of my debugging of the Groth16 prover implementation. For anyone interested, it can still be found here.
This PR implements a Groth16 prover and verifier, which requires SnarkJS produced
.wtns,.zkeyand.r1csbinary files.