From 0afda7bea8ccff1e3e0f0047d97b03320e44ee1d Mon Sep 17 00:00:00 2001 From: Mamy Ratsimbazafy Date: Sat, 6 Jul 2024 20:35:48 +0000 Subject: [PATCH] Multilinear extensions of polynomials (#423) * feat(multilinear-extension): initial commit * feat(multilinear-extension): add O(n) MultilinearExtension evaluations * fix(multilinear-extensions): address concerns of double free on copy --- benchmarks/bench_multilinear_extensions.nim | 75 ++++ constantine.nimble | 12 +- constantine/boolean_hypercube/README.md | 7 + .../multilinear_extensions.nim | 323 ++++++++++++++++++ constantine/commitments/eth_verkle_ipa.nim | 4 +- constantine/math/polynomials/polynomials.nim | 6 +- .../named/config_fields_and_curves.nim | 4 + .../t_multilinear_extensions.nim | 161 +++++++++ tests/t_ethereum_verkle_ipa_primitives.nim | 4 +- 9 files changed, 584 insertions(+), 12 deletions(-) create mode 100644 benchmarks/bench_multilinear_extensions.nim create mode 100644 constantine/boolean_hypercube/README.md create mode 100644 constantine/boolean_hypercube/multilinear_extensions.nim create mode 100644 tests/interactive_proofs/t_multilinear_extensions.nim diff --git a/benchmarks/bench_multilinear_extensions.nim b/benchmarks/bench_multilinear_extensions.nim new file mode 100644 index 00000000..4dec4a99 --- /dev/null +++ b/benchmarks/bench_multilinear_extensions.nim @@ -0,0 +1,75 @@ +# Constantine +# Copyright (c) 2018-2019 Status Research & Development GmbH +# Copyright (c) 2020-Present Mamy André-Ratsimbazafy +# Licensed and distributed under either of +# * MIT license (license terms in the root directory or at http://opensource.org/licenses/MIT). +# * Apache v2 license (license terms in the root directory or at http://www.apache.org/licenses/LICENSE-2.0). +# at your option. This file may not be copied, modified, or distributed except according to those terms. + +import + constantine/boolean_hypercube/multilinear_extensions, + constantine/named/algebras, + constantine/math/arithmetic, + constantine/platforms/static_for, + # Helpers + helpers/prng_unsafe, + benchmarks/bench_blueprint, + std/macros + +var rng*: RngState +let seed = 1234 +rng.seed(seed) +echo "bench multilinear_extensions xoshiro512** seed: ", seed + +proc separator*() = separator(155) + +proc report(op, field: string, start, stop: MonoTime, startClk, stopClk: int64, iters: int) = + let ns = inNanoseconds((stop-start) div iters) + let throughput = 1e9 / float64(ns) + when SupportsGetTicks: + echo &"{op:<60} {field:<18} {throughput:>15.3f} ops/s {ns:>9} ns/op {(stopClk - startClk) div iters:>9} CPU cycles (approx)" + else: + echo &"{op:<60} {field:<18} {throughput:>15.3f} ops/s {ns:>9} ns/op" + +macro fixFieldDisplay(T: typedesc): untyped = + # At compile-time, enums are integers and their display is buggy + # we get the Curve ID instead of the curve name. + let instantiated = T.getTypeInst() + var name = $instantiated[1][0] # 𝔽p + name.add "[" & $Algebra(instantiated[1][1].intVal) & "]" + result = newLit name + +template bench(op: string, T: typedesc, iters: int, body: untyped): untyped = + measure(iters, startTime, stopTime, startClk, stopClk, body) + report(op, fixFieldDisplay(T), startTime, stopTime, startClk, stopClk, iters) + +proc bench_mle(F: typedesc, num_vars: int) = + var evals = newSeq[F](1 shl num_vars) + for eval in evals.mitems(): + eval = rng.random_unsafe(F) + + let mle = MultilinearExtension[F].new(num_vars, evals) + + var coords = newSeq[F](num_vars) + for coord in coords.mitems(): + coord = rng.random_unsafe(F) + + var r: F + if num_vars <= 13: + bench("Multilinear Extension " & $num_vars & " variables: Reference EvaluateAt", F, 100): + r.evalMultilinearExtensionAt_reference(mle, coords) + + block: + bench("Multilinear Extension " & $num_vars & " variables: Optimized EvaluateAt", F, 100): + r.evalMultilinearExtensionAt(mle, coords) + +const Curves = [BN254_Snarks, BLS12_381] + +separator() +separator() +staticFor i, 0, Curves.len: + const curve = Curves[i] + for num_vars in countup(9, 19, 2): + bench_mle(Fr[curve], num_vars) + separator() + separator() diff --git a/constantine.nimble b/constantine.nimble index 831d1e14..dd05e12a 100644 --- a/constantine.nimble +++ b/constantine.nimble @@ -531,10 +531,10 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[ # ("tests/math_pairings/t_pairing_bls12_381_line_functions.nim", false), # ("tests/math_pairings/t_pairing_mul_fp12_by_lines.nim", false), ("tests/math_pairings/t_pairing_cyclotomic_subgroup.nim", false), - ("tests/math_pairings/t_pairing_bn254_nogami_optate.nim", false), - ("tests/math_pairings/t_pairing_bn254_snarks_optate.nim", false), - ("tests/math_pairings/t_pairing_bls12_377_optate.nim", false), - ("tests/math_pairings/t_pairing_bls12_381_optate.nim", false), + # ("tests/math_pairings/t_pairing_bn254_nogami_optate.nim", false), + # ("tests/math_pairings/t_pairing_bn254_snarks_optate.nim", false), + # ("tests/math_pairings/t_pairing_bls12_377_optate.nim", false), + # ("tests/math_pairings/t_pairing_bls12_381_optate.nim", false), # Multi-Pairing # ---------------------------------------------------------- @@ -549,7 +549,7 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[ # Hashing to elliptic curves # ---------------------------------------------------------- - ("tests/t_hash_to_field.nim", false), + # ("tests/t_hash_to_field.nim", false), ("tests/t_hash_to_curve_random.nim", false), ("tests/t_hash_to_curve.nim", false), @@ -571,6 +571,7 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[ # Proof systems # ---------------------------------------------------------- ("tests/proof_systems/t_r1cs_parser.nim", false), + ("tests/interactive_proofs/t_multilinear_extensions.nim", false), ] const testDescNvidia: seq[string] = @[ @@ -639,6 +640,7 @@ const benchDesc = [ "bench_eth_eip2537_subgroup_checks_impact", "bench_verkle_primitives", "bench_eth_evm_precompiles", + "bench_multilinear_extensions" ] # For temporary (hopefully) investigation that can only be reproduced in CI diff --git a/constantine/boolean_hypercube/README.md b/constantine/boolean_hypercube/README.md new file mode 100644 index 00000000..6956a45f --- /dev/null +++ b/constantine/boolean_hypercube/README.md @@ -0,0 +1,7 @@ +# Boolean Hypercube + +This folder holds utilities to work with the Boolean Hypercube, +via MLE, Multilinear Extension of polynomials evaluated at {0, 1}ⁿ. + +This is a heavy area of experimentation to nail down useful software architecture. +Expect many refactorings. diff --git a/constantine/boolean_hypercube/multilinear_extensions.nim b/constantine/boolean_hypercube/multilinear_extensions.nim new file mode 100644 index 00000000..147001db --- /dev/null +++ b/constantine/boolean_hypercube/multilinear_extensions.nim @@ -0,0 +1,323 @@ +# Constantine +# Copyright (c) 2018-2019 Status Research & Development GmbH +# Copyright (c) 2020-Present Mamy André-Ratsimbazafy +# Licensed and distributed under either of +# * MIT license (license terms in the root directory or at http://opensource.org/licenses/MIT). +# * Apache v2 license (license terms in the root directory or at http://www.apache.org/licenses/LICENSE-2.0). +# at your option. This file may not be copied, modified, or distributed except according to those terms. + +import + constantine/platforms/[abstractions, allocs] + +# Resources: +# - https://people.cs.georgetown.edu/jthaler/IPsandextensions.pdf +# - https://people.cs.georgetown.edu/jthaler/ProofsArgsAndZK.pdf +# Chapter 3.5 + +type + MultilinearExtension*[F] = object + ## Note: to follow mathematical description, indices start at 1 and end of range is inclusive + ## actual implementation will have indices start at 0 + ## + ## Given a sequence of bits of size s {0,1}ˢ + ## and an arbitrary function f: {0,1}ˢ -> 𝔽 + ## i.e. that maps a sequence of bits to a finite field 𝔽 + ## + ## there is an unique multilinear polynomial f̃ + ## called multilinear extension of f + ## that satisfies f̃(xᵢ) = f(xᵢ) for xᵢ ∈ {0,1}ˢ + ## + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) 𝛘ₑ(x₁, ...,xₛ) + ## with e ∈ {0,1}ˢ and f(e) the evaluation of f at e. + ## + ## 𝛘ₑ(x₁, ...,xₛ) is the multilinear Lagrange basis polynomial + ## which takes value 1 at 𝛘ₑ(e) and 0 at Xₑ(e̅) e̅ any other element ∈ {0,1}ˢ + ## + ## 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ₌₁ˢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)), i ∈ [1, s] + ## + ## A multilinear polynomial is linear (i.e. degree at most 1) in each + ## of its variables. + num_vars*: int + base_poly_evals*: ptr UncheckedArray[F] + +proc `=copy`*[F](dst: var MultilinearExtension[F], src: MultilinearExtension[F]) {.error: "A multilinear extension cannot be copied".} + +proc `=wasMoved`*[F](mle: var MultilinearExtension[F]) = + mle.base_poly_evals = nil + +proc `=destroy`*[F](mle: var MultilinearExtension[F]) = + if not mle.base_poly_evals.isNil: + freeHeapAligned(mle.base_poly_evals) + +func new*[F](T: type MultilinearExtension[F], num_vars: int, poly_evals: openArray[F]): T = + debug: + doAssert 1 shl num_vars == poly_evals.len, block: + "The MLE has " & $num_vars " variables\n" & + "but the poly it's derived from has " & $poly_evals.len & " evaluations.\n" & + "2^" & $num_vars & " = " & $(1 shl num_vars) & " were expected instead." + result.num_vars = num_vars + let L = 1 shl num_vars + result.base_poly_evals = allocHeapArrayAligned(F, L, alignment = 64) + for i in 0 ..< L: + result.base_poly_evals[i] = poly_evals[i] + +# Big-endian Multilinear Lagrange Basis +# ------------------------------------------------------------------------------------ + +iterator bits_be(n: SomeInteger, len: int): (int, bool) = + for i in 0 ..< len: + yield (i, bool((n shr (len-1-i)) and 1)) + +func evalMultilinearExtensionAt_BE_reference[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F]) = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## e is interpreted as a big-endian bitstring + ## + ## This is a reference implementation using naive computation + ## in O(n log n) with n being the number of evaluations of the original polynomial. + debug: doAssert mle.num_vars == coords.len + + let L = 1 shl mle.num_vars + + r.setZero() + for e in 0 ..< L: + # 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + # e ∈ {0,1}ˢ hence each factor is either: + # (1-xᵢ) or xᵢ + var chi_e {.noInit.}: F + chi_e.setOne() + + for (i, ei) in bits_be(e, mle.num_vars): + if ei: + chi_e *= xs[i] + else: + var t {.noInit.}: F + t.diff(F.getOne(), xs[i]) + chi_e *= t + + var t {.noInit.}: F + t.prod(mle.base_poly_evals[e], chi_e) + r += t + +func evalMultilinearExtensionAt_BE[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F]) = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## e is interpreted as a big-endian bitstring + ## + ## This is an optimized implementation that + ## + ## evaluates 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## in O(n) + # + # e ∈ {0,1}ˢ hence each factor is either: + # (1-xᵢ) or xᵢ + # + # See the algorithm in ec_endomorphism_accel to build + # a binary lookup table for O(n) evaluations + # + # Algorithm: + # for u in 0 ..< 2ˢ + # lagrange[u] = 1 + # iterate on the bit representation of u + # if the bit is set, multiply by matching xᵤ + # if not, multiply by (1-xᵤ) + # + # Implementation: + # We optimize the basic algorithm to reuse already computed table entries + # by noticing for example that: + # - 6 represented as 0b110 requires x₂.x₁(1-x₀) + # - 2 represented as 0b010 requires (1-x₂).x₁(1-x₀) + + let L = 1 shl mle.num_vars + let buf = allocHeapArrayAligned(F, L, alignment = 64) + for i in 0 ..< L: + buf[i] = mle.base_poly_evals[i] + + # number of operations: + # ∑ᵢ₌₀ˡᵒᵍ²⁽ˢ⁾⁻¹ 2ˡᵒᵍ²⁽ˢ⁾⁻¹⁻ⁱ = 2ˡᵒᵍ²⁽ˢ⁾-1 = s-1 + # Using sum of consecutive powers of 2 formula + # So we're linear in the original polynomial size + for i in countdown(xs.len-1, 0): + for e in 0 ..< 1 shl i: + # Implicit binary tree representation with root at 1 + # root at 1 + # Left child e*2 + 0 + # Right child e*2 + 1 + # Parent e/2 + # + # Representation + # + # depth 0 ------ 1 ------ + # depth 1 -- 2 -- --- 3 --- + # depth 2 4 5 6 7 + # depth 3 8 9 10 11 12 13 14 15 + # + # In an array, storage is linear + template left(): untyped = buf[e shl 1] + template right(): untyped = buf[(e shl 1) + 1] + + var t {.noInit.}: F + t.diff(right, left) + t *= xs[i] + buf[e].sum(left, t) + + r = buf[0] + freeHeapAligned(buf) + +# Little-endian Multilinear Lagrange Basis +# ------------------------------------------------------------------------------------ + +iterator bits_le(n: SomeInteger, len: int): (int, bool) = + for i in 0 ..< len: + yield (i, bool((n shr i) and 1)) + +func evalMultilinearExtensionAt_LE_reference[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F]) = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## e is interpreted as a little-endian bitstring + ## + ## This is a reference implementation using naive computation + ## in O(n log n) with n being the number of evaluations of the original polynomial. + debug: doAssert mle.num_vars == coords.len + + let L = 1 shl mle.num_vars + + r.setZero() + for e in 0 ..< L: + # 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + # e ∈ {0,1}ˢ hence each factor is either: + # (1-xᵢ) or xᵢ + var chi_e {.noInit.}: F + chi_e.setOne() + + for (i, ei) in bits_le(e, mle.num_vars): + if ei: + chi_e *= xs[i] + else: + var t {.noInit.}: F + t.diff(F.getOne(), xs[i]) + chi_e *= t + + var t {.noInit.}: F + t.prod(mle.base_poly_evals[e], chi_e) + r += t + +func evalMultilinearExtensionAt_LE[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F]) = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## e is interpreted as a little-endian bitstring + ## + ## This is an optimized implementation that + ## + ## evaluates 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## in O(n) + # + # e ∈ {0,1}ˢ hence each factor is either: + # (1-xᵢ) or xᵢ + # + # See the algorithm in ec_endomorphism_accel to build + # a binary lookup table for O(n) evaluations + # + # Algorithm: + # for u in 0 ..< 2ˢ + # lagrange[u] = 1 + # iterate on the bit representation of u + # if the bit is set, multiply by matching xᵤ + # if not, multiply by (1-xᵤ) + # + # Implementation: + # We optimize the basic algorithm to reuse already computed table entries + # by noticing for example that: + # - 6 represented as 0b110 requires x₂.x₁(1-x₀) + # - 2 represented as 0b010 requires (1-x₂).x₁(1-x₀) + + let L = 1 shl mle.num_vars + let buf = allocHeapArrayAligned(F, L, alignment = 64) + for i in 0 ..< L: + buf[i] = mle.base_poly_evals[i] + + # number of operations: + # ∑ᵢ₌₀ˡᵒᵍ²⁽ˢ⁾⁻¹ 2ˡᵒᵍ²⁽ˢ⁾⁻¹⁻ⁱ = 2ˡᵒᵍ²⁽ˢ⁾-1 = s-1 + # Using sum of consecutive powers of 2 formula + # So we're linear in the original polynomial size + for i in 0 ..< xs.len: + for e in 0 ..< 1 shl (mle.num_vars - 1 - i): + # Implicit binary tree representation with root at 1 + # root at 1 + # Left child e*2 + 0 + # Right child e*2 + 1 + # Parent e/2 + # + # Representation + # + # depth 0 ------ 1 ------ + # depth 1 -- 2 -- --- 3 --- + # depth 2 4 5 6 7 + # depth 3 8 9 10 11 12 13 14 15 + # + # In an array, storage is linear + template left(): untyped = buf[e shl 1] + template right(): untyped = buf[(e shl 1) + 1] + + var t {.noInit.}: F + t.diff(right, left) + t *= xs[i] + buf[e].sum(left, t) + + r = buf[0] + freeHeapAligned(buf) + +# Public API +# ------------------------------------------------------------------------------------ + +func evalMultilinearExtensionAt_reference*[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F], + endian: static Endianness = bigEndian) {.inline.} = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## By default, e is interpreted as a big-endian bitstring + ## + ## This is a reference implementation using naive computation + ## in O(n log n) with n being the number of evaluations of the original polynomial. + when endian == bigEndian: + evalMultilinearExtensionAt_BE_reference(r, mle, xs) + else: + evalMultilinearExtensionAt_LE_reference(r, mle, xs) + +func evalMultilinearExtensionAt*[F]( + r: var F, + mle: MultilinearExtension[F], + xs: openArray[F], + endian: static Endianness = bigEndian) {.inline.} = + ## Compute + ## f̃(x₁, ...,xₛ) = ∑ₑ f(e) ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## at supplied (x₁, ...,xₛ) + ## By default, e is interpreted as a little-endian bitstring + ## + ## This is an optimized implementation that + ## + ## evaluates 𝛘ₑ(x₁, ...,xₛ) = ∏ᵢ(xᵢeᵢ + (1-xᵢ)(1-eᵢ)) + ## in O(n) + when endian == bigEndian: + evalMultilinearExtensionAt_BE(r, mle, xs) + else: + evalMultilinearExtensionAt_LE(r, mle, xs) diff --git a/constantine/commitments/eth_verkle_ipa.nim b/constantine/commitments/eth_verkle_ipa.nim index 018b111b..f6606b2a 100644 --- a/constantine/commitments/eth_verkle_ipa.nim +++ b/constantine/commitments/eth_verkle_ipa.nim @@ -176,7 +176,7 @@ func ipa_prove*[N, logN: static int, EcAff, F]( var G = gprime[].toMutableView() a.copyFrom(poly.evals) - domain.getLagrangeBasisPolysAt(bprime[], opening_challenge) + domain.computeLagrangeBasisPolysAt(bprime[], opening_challenge) G.copyFrom(crs.evals) # Protocol @@ -415,7 +415,7 @@ func ipa_verify*[N, logN: static int, EcAff, F]( # ... - [a₀.b₀.w]G ... # a₀.b₀.w = a₀.w. - domain.getLagrangeBasisPolysAt(b[], opening_challenge) + domain.computeLagrangeBasisPolysAt(b[], opening_challenge) fs_a0b0wG[0].innerProduct(View[F](fs_a0g0), b.toView()) fs_a0b0wG[0] *= w ecs_a0b0wG[0].setGenerator() diff --git a/constantine/math/polynomials/polynomials.nim b/constantine/math/polynomials/polynomials.nim index b74406a0..7811c901 100644 --- a/constantine/math/polynomials/polynomials.nim +++ b/constantine/math/polynomials/polynomials.nim @@ -443,7 +443,7 @@ func evalPolyAt*[N: static int, Field]( freeHeapAligned(invZminusDomain) -func getLagrangeBasisPolysAt*[N: static int, Field]( +func computeLagrangeBasisPolysAt*[N: static int, Field]( domain: PolyEvalDomain[N, Field], lagrangePolys: var array[N, Field], z: Field) = @@ -516,11 +516,11 @@ func evalPolyAt*[N: static int, Field]( z: Field) = lindom.dom.evalPolyAt(r, poly, z) -func getLagrangeBasisPolysAt*[N: static int, Field]( +func computeLagrangeBasisPolysAt*[N: static int, Field]( lindom: PolyEvalLinearDomain[N, Field], lagrangePolys: var array[N, Field], z: Field) = - lindom.dom.getLagrangeBasisPolysAt(lagrangePolys, z) + lindom.dom.computeLagrangeBasisPolysAt(lagrangePolys, z) func setupLinearEvaluationDomain*[N: static int, Field]( lindom: var PolyEvalLinearDomain[N, Field]) = diff --git a/constantine/named/config_fields_and_curves.nim b/constantine/named/config_fields_and_curves.nim index 8cf7faf8..c3fc397a 100644 --- a/constantine/named/config_fields_and_curves.nim +++ b/constantine/named/config_fields_and_curves.nim @@ -45,6 +45,10 @@ export CurveFamily, SexticTwist declareCurves: # ----------------------------------------------------------------------------- # Curves added when passed "-d:CTT_TEST_CURVES" + curve F5: + testingCurve: true + bitwidth: 3 + modulus: "0x5" curve Fake101: testingCurve: true bitwidth: 7 diff --git a/tests/interactive_proofs/t_multilinear_extensions.nim b/tests/interactive_proofs/t_multilinear_extensions.nim new file mode 100644 index 00000000..38705cec --- /dev/null +++ b/tests/interactive_proofs/t_multilinear_extensions.nim @@ -0,0 +1,161 @@ +# Constantine +# Copyright (c) 2018-2019 Status Research & Development GmbH +# Copyright (c) 2020-Present Mamy André-Ratsimbazafy +# Licensed and distributed under either of +# * MIT license (license terms in the root directory or at http://opensource.org/licenses/MIT). +# * Apache v2 license (license terms in the root directory or at http://www.apache.org/licenses/LICENSE-2.0). +# at your option. This file may not be copied, modified, or distributed except according to those terms. + +import + # Standard library + std/times, + # Internals + constantine/boolean_hypercube/multilinear_extensions, + constantine/named/algebras, + constantine/math/arithmetic, + constantine/math/io/io_fields, + # Helpers + helpers/prng_unsafe + +# Compile with -d:CTT_TEST_CURVES to define F5 + +func toF5[N: static int](a: array[N, SomeUnsignedInt]): array[N, Fp[F5]] = + for i in 0 ..< N: + result[i] = Fp[F5].fromUint(a[i]) + +func transpose[M, N: static int, T](a: array[M, array[N, T]]): array[N, array[M, T]] = + for i in 0 ..< M: + for j in 0 ..< N: + result[j][i] = a[i][j] + +# - https://people.cs.georgetown.edu/jthaler/IPsandextensions.pdf\ +# Note: first row is +# 1 2 3 4 0 not 1 2 3 4 5 (though 5 ≡ 0 (mod 5) so arguably not wrong) +# - https://people.cs.georgetown.edu/jthaler/ProofsArgsAndZK.pdf +# Chapter 3.5 +proc test_thaler() = + let evals = [uint32 1, 2, 1, 4].toF5() + let mle_evals_be = [ + [byte 1, 2, 3, 4, 0], + [byte 1, 4, 2, 0, 3], + [byte 1, 1, 1, 1, 1], + [byte 1, 3, 0, 2, 4], + [byte 1, 0, 4, 3, 2], + ] + let mle_evals_le = mle_evals_be.transpose() + + let mle = MultilinearExtension[Fp[F5]].new(2, evals) + + block: + for i in 0'u32 .. 4: + var row: array[5, byte] + for j in 0'u32 .. 4: + var r: Fp[F5] + r.evalMultilinearExtensionAt_reference( + mle, [Fp[F5].fromUint(i), Fp[F5].fromUint(j)], + bigEndian) + var buf: array[1, byte] + buf.marshal(r, bigEndian) + row[j] = buf[0] + + echo row + doAssert row == mle_evals_be[i] + + echo "=== SUCCESS reference MLE big-endian evaluation ===" + + block: + for i in 0'u32 .. 4: + var row: array[5, byte] + for j in 0'u32 .. 4: + var r: Fp[F5] + r.evalMultilinearExtensionAt( + mle, [Fp[F5].fromUint(i), Fp[F5].fromUint(j)], + bigEndian) + var buf: array[1, byte] + buf.marshal(r, bigEndian) + row[j] = buf[0] + + echo row + doAssert row == mle_evals_be[i] + + echo "=== SUCCESS optimized MLE big-endian evaluation ===" + + block: + for i in 0'u32 .. 4: + var row: array[5, byte] + for j in 0'u32 .. 4: + var r: Fp[F5] + r.evalMultilinearExtensionAt_reference( + mle, [Fp[F5].fromUint(i), Fp[F5].fromUint(j)], + littleEndian) + var buf: array[1, byte] + buf.marshal(r, bigEndian) + row[j] = buf[0] + + echo row + doAssert row == mle_evals_le[i] + + echo "=== SUCCESS reference MLE little-endian evaluation ===" + + block: + for i in 0'u32 .. 4: + var row: array[5, byte] + for j in 0'u32 .. 4: + var r: Fp[F5] + r.evalMultilinearExtensionAt( + mle, [Fp[F5].fromUint(i), Fp[F5].fromUint(j)], + littleEndian) + var buf: array[1, byte] + buf.marshal(r, bigEndian) + row[j] = buf[0] + + echo row + doAssert row == mle_evals_le[i] + + echo "=== SUCCESS optimized MLE little-endian evaluation ===" + +var rng*: RngState +let seed = uint32(getTime().toUnix() and (1'i64 shl 32 - 1)) # unixTime mod 2^32 +rng.seed(seed) +echo "\n------------------------------------------------------\n" +echo "Multilinear Extensions xoshiro512** seed: ", seed + +proc test_randomized(F: typedesc, num_vars: int) = + var evals = newSeq[F](1 shl num_vars) + for eval in evals.mitems(): + eval = rng.random_unsafe(F) + + let mle = MultilinearExtension[F].new(num_vars, evals) + + var coords = newSeq[F](num_vars) + for coord in coords.mitems(): + coord = rng.random_unsafe(F) + + var r_ref {.noInit.}: F + var r_opt {.noInit.}: F + + r_ref.evalMultilinearExtensionAt_reference(mle, coords, bigEndian) + r_opt.evalMultilinearExtensionAt(mle, coords, bigEndian) + + doAssert bool(r_ref == r_opt) + echo "Success: ", F, ", ", num_vars, " variables, bigEndian" #, r: ", r_ref.toHex() + + r_ref.evalMultilinearExtensionAt_reference(mle, coords, littleEndian) + r_opt.evalMultilinearExtensionAt(mle, coords, littleEndian) + + doAssert bool(r_ref == r_opt) + echo "Success: ", F, ", ", num_vars, " variables, littleEndian" #, r: ", r_ref.toHex() + +test_thaler() +test_randomized(Fr[BN254_Snarks], 3) +test_randomized(Fr[BN254_Snarks], 3) +test_randomized(Fr[BN254_Snarks], 7) +test_randomized(Fr[BN254_Snarks], 7) +test_randomized(Fr[BN254_Snarks], 11) +test_randomized(Fr[BN254_Snarks], 11) +test_randomized(Fr[BLS12_381], 3) +test_randomized(Fr[BLS12_381], 3) +test_randomized(Fr[BLS12_381], 7) +test_randomized(Fr[BLS12_381], 7) +test_randomized(Fr[BLS12_381], 11) +test_randomized(Fr[BLS12_381], 11) diff --git a/tests/t_ethereum_verkle_ipa_primitives.nim b/tests/t_ethereum_verkle_ipa_primitives.nim index 5f5a8445..71bed6a9 100644 --- a/tests/t_ethereum_verkle_ipa_primitives.nim +++ b/tests/t_ethereum_verkle_ipa_primitives.nim @@ -100,7 +100,7 @@ suite "Barycentric Form Tests": lindom.setupLinearEvaluationDomain() var bar_coeffs {.noInit.}: array[256, Fr[Banderwagon]] - lindom.getLagrangeBasisPolysAt(bar_coeffs, p_outside_dom) + lindom.computeLagrangeBasisPolysAt(bar_coeffs, p_outside_dom) var expected0: Fr[Banderwagon] expected0.innerProduct(lagrange_values.evals, bar_coeffs) @@ -519,7 +519,7 @@ suite "IPA proof tests": opening_challenge) doAssert eval_at_challenge.toHex(littleEndian) == "0x4a353e70b03c89f161de002e8713beec0d740a5e20722fd5bd68b30540a33208", "Issue with computing commitment" - + var prover_challenge {.noInit.}: Fr[Banderwagon] prover_transcript.squeezeChallenge("state", prover_challenge) doAssert prover_challenge.toHex(littleEndian) == "0x0a81881cbfd7d7197a54ebd67ed6a68b5867f3c783706675b34ece43e85e7306", "Issue with squeezing prover challenge"