diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 6122022..64fd5c9 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -55,6 +55,11 @@ jobs: name: '${{ matrix.target.os }}-${{ matrix.target.cpu }}-${{ matrix.target.TEST_LANG }} (${{ matrix.branch }})' runs-on: ${{ matrix.builder }} steps: + - name: Cancel Previous Runs + uses: styfle/cancel-workflow-action@0.5.0 + with: + access_token: ${{ github.token }} + - name: Checkout constantine uses: actions/checkout@v2 with: @@ -197,13 +202,20 @@ jobs: run: | nimble refresh nimble install -y gmp stew - - name: Run Constantine tests (with GMP) + - name: Run Constantine tests (with Assembler & with GMP) if: (runner.os == 'Linux' || runner.os == 'macOS') && matrix.target.cpu == 'amd64' shell: bash run: | export UCPU="$cpu" cd constantine nimble test_parallel + - name: Run Constantine tests (no Assembler & with GMP) + if: (runner.os == 'Linux' || runner.os == 'macOS') && matrix.target.cpu == 'amd64' + shell: bash + run: | + export UCPU="$cpu" + cd constantine + nimble test_parallel_no_assembler - name: Run Constantine tests (without GMP) if: runner.os == 'Linux' && matrix.target.cpu == 'i386' shell: bash @@ -211,3 +223,10 @@ jobs: export UCPU="$cpu" cd constantine nimble test_parallel_no_gmp + - name: Run Constantine tests (without Assembler or GMP) + if: runner.os == 'Linux' && matrix.target.cpu == 'i386' + shell: bash + run: | + export UCPU="$cpu" + cd constantine + nimble test_parallel_no_gmp_no_assembler diff --git a/azure-pipelines.yml b/azure-pipelines.yml index 66646de..fa721cf 100644 --- a/azure-pipelines.yml +++ b/azure-pipelines.yml @@ -240,5 +240,5 @@ steps: echo "PATH=${PATH}" export ucpu=${UCPU} nimble test_no_gmp - displayName: 'Testing the package (without Assembler or GMP)' + displayName: 'Testing the package (without GMP and ASM on Windows)' condition: eq(variables['Agent.OS'], 'Windows_NT') diff --git a/benchmarks/bench_ec_g1.nim b/benchmarks/bench_ec_g1.nim index e3454e0..7684fdb 100644 --- a/benchmarks/bench_ec_g1.nim +++ b/benchmarks/bench_ec_g1.nim @@ -26,8 +26,8 @@ import # ############################################################ -const Iters = 1_000_000 -const MulIters = 1000 +const Iters = 10_000 +const MulIters = 100 const AvailableCurves = [ # P224, # BN254_Nogami, diff --git a/benchmarks/bench_ec_g2.nim b/benchmarks/bench_ec_g2.nim index 37f5f41..c37f299 100644 --- a/benchmarks/bench_ec_g2.nim +++ b/benchmarks/bench_ec_g2.nim @@ -27,7 +27,7 @@ import # ############################################################ -const Iters = 500_000 +const Iters = 10_000 const MulIters = 500 const AvailableCurves = [ # P224, diff --git a/benchmarks/bench_fp6.nim b/benchmarks/bench_fp6.nim index 717e1f6..e843bc9 100644 --- a/benchmarks/bench_fp6.nim +++ b/benchmarks/bench_fp6.nim @@ -23,7 +23,7 @@ import # ############################################################ -const Iters = 1_000_000 +const Iters = 100_000 const InvIters = 1000 const AvailableCurves = [ # Pairing-Friendly curves diff --git a/constantine.nimble b/constantine.nimble index a2d74f6..4c59b62 100644 --- a/constantine.nimble +++ b/constantine.nimble @@ -46,6 +46,9 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[ ("tests/t_fp12_bn254_snarks.nim", false), ("tests/t_fp12_bls12_377.nim", false), ("tests/t_fp12_bls12_381.nim", false), + ("tests/t_fp12_exponentiation.nim", false), + + ("tests/t_fp4_frobenius.nim", false), # Elliptic curve arithmetic G1 ("tests/t_ec_wstrass_prj_g1_add_double.nim", false), ("tests/t_ec_wstrass_prj_g1_mul_sanity.nim", false), @@ -66,7 +69,12 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[ ("tests/t_ec_sage_bn254.nim", false), ("tests/t_ec_sage_bls12_381.nim", false), # Edge cases highlighted by past bugs - ("tests/t_ec_wstrass_prj_edge_cases.nim", false) + ("tests/t_ec_wstrass_prj_edge_cases.nim", false), + # Pairing + ("tests/t_pairing_fp12_sparse.nim", false), + ("tests/t_pairing_bn254_nogami_optate.nim", false), + ("tests/t_pairing_bn254_snarks_optate.nim", false), + ("tests/t_pairing_bls12_381_optate.nim", false) ] # For temporary (hopefully) investigation that can only be reproduced in CI @@ -139,8 +147,6 @@ task test, "Run all tests": else: test "-d:Constantine32", td.path - # Benchmarks compile and run - # ignore Windows 32-bit for the moment # Ensure benchmarks stay relevant. Ignore Windows 32-bit at the moment if not defined(windows) or not (existsEnv"UCPU" or getEnv"UCPU" == "i686"): runBench("bench_fp") @@ -167,9 +173,6 @@ task test_no_gmp, "Run tests that don't require GMP": else: test "-d:Constantine32", td.path - - # Benchmarks compile and run - # ignore Windows 32-bit for the moment # Ensure benchmarks stay relevant. Ignore Windows 32-bit at the moment if not defined(windows) or not (existsEnv"UCPU" or getEnv"UCPU" == "i686"): runBench("bench_fp") @@ -298,6 +301,47 @@ task test_parallel_no_gmp, "Run all tests in parallel (via GNU parallel)": runBench("bench_ec_g1") runBench("bench_ec_g2") +task test_parallel_no_gmp_no_assembler, "Run all tests in parallel (via GNU parallel)": + # -d:testingCurves is configured in a *.nim.cfg for convenience + let cmdFile = true # open(buildParallel, mode = fmWrite) # Nimscript doesn't support IO :/ + exec "> " & buildParallel + + for td in testDesc: + if not td.useGMP: + if td.path in useDebug: + test "-d:debugConstantine -d:ConstantineASM=false", td.path, cmdFile + else: + test "-d:ConstantineASM=false", td.path, cmdFile + + # cmdFile.close() + # Execute everything in parallel with GNU parallel + exec "parallel --keep-order --group < " & buildParallel + + exec "> " & buildParallel + if sizeof(int) == 8: # 32-bit tests on 64-bit arch + for td in testDesc: + if not td.useGMP: + if td.path in useDebug: + test "-d:Constantine32 -d:debugConstantine", td.path, cmdFile + else: + test "-d:Constantine32", td.path, cmdFile + # cmdFile.close() + # Execute everything in parallel with GNU parallel + exec "parallel --keep-order --group < " & buildParallel + + # Now run the benchmarks + # + # Benchmarks compile and run + # ignore Windows 32-bit for the moment + # Ensure benchmarks stay relevant. Ignore Windows 32-bit at the moment + if not defined(windows) or not (existsEnv"UCPU" or getEnv"UCPU" == "i686"): + runBench("bench_fp") + runBench("bench_fp2") + runBench("bench_fp6") + runBench("bench_fp12") + runBench("bench_ec_g1") + runBench("bench_ec_g2") + task bench_fp, "Run benchmark 𝔽p with your default compiler": runBench("bench_fp") diff --git a/constantine/arithmetic/bigints.nim b/constantine/arithmetic/bigints.nim index 0f3c9aa..b5f5018 100644 --- a/constantine/arithmetic/bigints.nim +++ b/constantine/arithmetic/bigints.nim @@ -318,6 +318,87 @@ func bit0*(a: BigInt): Ct[uint8] = ## Access the least significant bit ct(a.limbs[0] and SecretWord(1), uint8) +# Multiplication by small cosntants +# ------------------------------------------------------------ + +func `*=`*(a: var BigInt, b: static int) {.inline.} = + ## Multiplication by a small integer known at compile-time + # Implementation: + # + # we hardcode addition chains for small integer + const negate = b < 0 + const b = if negate: -b + else: b + when negate: + a.neg(a) + when b == 0: + a.setZero() + elif b == 1: + return + elif b == 2: + discard a.double() + elif b == 3: + let t1 = a + discard a.double() + a += t1 + elif b == 4: + discard a.double() + discard a.double() + elif b == 5: + let t1 = a + discard a.double() + discard a.double() + a += t1 + elif b == 6: + discard a.double() + let t2 = a + discard a.double() # 4 + a += t2 + elif b == 7: + let t1 = a + discard a.double() + let t2 = a + discard a.double() # 4 + a += t2 + a += t1 + elif b == 8: + discard a.double() + discard a.double() + discard a.double() + elif b == 9: + let t1 = a + discard a.double() + discard a.double() + discard a.double() # 8 + a += t1 + elif b == 10: + discard a.double() + let t2 = a + discard a.double() + discard a.double() # 8 + a += t2 + elif b == 11: + let t1 = a + discard a.double() + let t2 = a + discard a.double() + discard a.double() # 8 + a += t2 + a += t1 + elif b == 12: + discard a.double() + discard a.double() # 4 + let t4 = a + discard a.double() # 8 + a += t4 + else: + {.error: "Multiplication by this small int not implemented".} + +func `*`*(b: static int, a: BigInt): BigInt {.noinit, inline.} = + ## Multiplication by a small integer known at compile-time + result = a + result *= b + # ############################################################ # # Modular BigInt @@ -417,58 +498,6 @@ func montySquare*(r: var BigInt, a, M: BigInt, negInvModWord: static BaseType, c ## to avoid duplicating with Nim zero-init policy montySquare(r.limbs, a.limbs, M.limbs, negInvModWord, canUseNoCarryMontyMul) -func montyPow*[mBits, eBits: static int]( - a: var BigInt[mBits], exponent: BigInt[eBits], - M, one: BigInt[mBits], negInvModWord: static BaseType, windowSize: static int, - canUseNoCarryMontyMul, canUseNoCarryMontySquare: static bool - ) = - ## Compute a <- a^exponent (mod M) - ## ``a`` in the Montgomery domain - ## ``exponent`` is any BigInt, in the canonical domain - ## - ## This uses fixed window optimization - ## A window size in the range [1, 5] must be chosen - ## - ## This is constant-time: the window optimization does - ## not reveal the exponent bits or hamming weight - mixin exportRawUint # exported in io_bigints which depends on this module ... - - var expBE {.noInit.}: array[(ebits + 7) div 8, byte] - expBE.exportRawUint(exponent, bigEndian) - - const scratchLen = if windowSize == 1: 2 - else: (1 shl windowSize) + 1 - var scratchSpace {.noInit.}: array[scratchLen, Limbs[mBits.wordsRequired]] - montyPow(a.limbs, expBE, M.limbs, one.limbs, negInvModWord, scratchSpace, canUseNoCarryMontyMul, canUseNoCarryMontySquare) - -func montyPowUnsafeExponent*[mBits, eBits: static int]( - a: var BigInt[mBits], exponent: BigInt[eBits], - M, one: BigInt[mBits], negInvModWord: static BaseType, windowSize: static int, - canUseNoCarryMontyMul, canUseNoCarryMontySquare: static bool - ) = - ## Compute a <- a^exponent (mod M) - ## ``a`` in the Montgomery domain - ## ``exponent`` is any BigInt, in the canonical domain - ## - ## Warning ⚠️ : - ## This is an optimization for public exponent - ## Otherwise bits of the exponent can be retrieved with: - ## - memory access analysis - ## - power analysis - ## - timing analysis - ## - ## This uses fixed window optimization - ## A window size in the range [1, 5] must be chosen - mixin exportRawUint # exported in io_bigints which depends on this module ... - - var expBE {.noInit.}: array[(ebits + 7) div 8, byte] - expBE.exportRawUint(exponent, bigEndian) - - const scratchLen = if windowSize == 1: 2 - else: (1 shl windowSize) + 1 - var scratchSpace {.noInit.}: array[scratchLen, Limbs[mBits.wordsRequired]] - montyPowUnsafeExponent(a.limbs, expBE, M.limbs, one.limbs, negInvModWord, scratchSpace, canUseNoCarryMontyMul, canUseNoCarryMontySquare) - func montyPow*[mBits: static int]( a: var BigInt[mBits], exponent: openarray[byte], M, one: BigInt[mBits], negInvModWord: static BaseType, windowSize: static int, @@ -517,5 +546,51 @@ func montyPowUnsafeExponent*[mBits: static int]( var scratchSpace {.noInit.}: array[scratchLen, Limbs[mBits.wordsRequired]] montyPowUnsafeExponent(a.limbs, exponent, M.limbs, one.limbs, negInvModWord, scratchSpace, canUseNoCarryMontyMul, canUseNoCarryMontySquare) +func montyPow*[mBits, eBits: static int]( + a: var BigInt[mBits], exponent: BigInt[eBits], + M, one: BigInt[mBits], negInvModWord: static BaseType, windowSize: static int, + canUseNoCarryMontyMul, canUseNoCarryMontySquare: static bool + ) = + ## Compute a <- a^exponent (mod M) + ## ``a`` in the Montgomery domain + ## ``exponent`` is any BigInt, in the canonical domain + ## + ## This uses fixed window optimization + ## A window size in the range [1, 5] must be chosen + ## + ## This is constant-time: the window optimization does + ## not reveal the exponent bits or hamming weight + mixin exportRawUint # exported in io_bigints which depends on this module ... + + var expBE {.noInit.}: array[(ebits + 7) div 8, byte] + expBE.exportRawUint(exponent, bigEndian) + + montyPow(a, expBE, M, one, negInvModWord, windowSize, canUseNoCarryMontyMul, canUseNoCarryMontySquare) + +func montyPowUnsafeExponent*[mBits, eBits: static int]( + a: var BigInt[mBits], exponent: BigInt[eBits], + M, one: BigInt[mBits], negInvModWord: static BaseType, windowSize: static int, + canUseNoCarryMontyMul, canUseNoCarryMontySquare: static bool + ) = + ## Compute a <- a^exponent (mod M) + ## ``a`` in the Montgomery domain + ## ``exponent`` is any BigInt, in the canonical domain + ## + ## Warning ⚠️ : + ## This is an optimization for public exponent + ## Otherwise bits of the exponent can be retrieved with: + ## - memory access analysis + ## - power analysis + ## - timing analysis + ## + ## This uses fixed window optimization + ## A window size in the range [1, 5] must be chosen + mixin exportRawUint # exported in io_bigints which depends on this module ... + + var expBE {.noInit.}: array[(ebits + 7) div 8, byte] + expBE.exportRawUint(exponent, bigEndian) + + montyPowUnsafeExponent(a, expBE, M, one, negInvModWord, windowSize, canUseNoCarryMontyMul, canUseNoCarryMontySquare) + {.pop.} # inline {.pop.} # raises no exceptions diff --git a/constantine/arithmetic/finite_fields.nim b/constantine/arithmetic/finite_fields.nim index e4b1983..0bae46c 100644 --- a/constantine/arithmetic/finite_fields.nim +++ b/constantine/arithmetic/finite_fields.nim @@ -225,6 +225,10 @@ func neg*(r: var Fp, a: Fp) {.inline.} = else: discard r.mres.diff(Fp.C.Mod, a.mres) +func neg*(a: var Fp) {.inline.} = + ## Negate modulo p + a.neg(a) + func div2*(a: var Fp) {.inline.} = ## Modular division by 2 a.mres.div2_modular(Fp.C.getPrimePlus1div2()) diff --git a/constantine/config/curves.nim b/constantine/config/curves.nim index a2debf6..83629d1 100644 --- a/constantine/config/curves.nim +++ b/constantine/config/curves.nim @@ -94,6 +94,7 @@ macro get_CNR_Fp2*(C: static Curve): untyped = result = bindSym($C & "_nonresidue_cube_fp2") macro getSexticTwist*(C: static Curve): untyped = + ## Returns if D-Twist or M-Twist result = bindSym($C & "_sexticTwist") macro get_SNR_Fp2*(C: static Curve): untyped = diff --git a/constantine/config/curves_declaration.nim b/constantine/config/curves_declaration.nim index 3efdfc7..fda3bc4 100644 --- a/constantine/config/curves_declaration.nim +++ b/constantine/config/curves_declaration.nim @@ -85,6 +85,19 @@ declareCurves: family: BarretoNaehrig # Equation: Y^2 = X^3 + 2 # u: -(2^62 + 2^55 + 1) + + order: "0x2523648240000001ba344d8000000007ff9f800000000010a10000000000000d" + orderBitwidth: 254 + cofactor: 1 + eq_form: ShortWeierstrass + coef_a: 0 + coef_b: 2 + nonresidue_quad_fp: -1 # -1 is not a square in 𝔽p + nonresidue_cube_fp2: (1, 1) # 1+𝑖 1+𝑖 is not a cube in 𝔽p² + + sexticTwist: D_Twist + sexticNonResidue_fp2: (1, 1) # 1+𝑖 + curve BN254_Snarks: # Zero-Knowledge proofs curve (SNARKS, STARKS, Ethereum) bitwidth: 254 modulus: "0x30644e72e131a029b85045b68181585d97816a916871ca8d3c208c16d87cfd47" diff --git a/constantine/elliptic/ec_weierstrass_affine.nim b/constantine/elliptic/ec_weierstrass_affine.nim index a25c167..0e611bf 100644 --- a/constantine/elliptic/ec_weierstrass_affine.nim +++ b/constantine/elliptic/ec_weierstrass_affine.nim @@ -13,6 +13,20 @@ import ../towers, ../io/io_bigints +# ############################################################ +# +# Elliptic Curve in Short Weierstrass form +# with Projective Coordinates +# +# ############################################################ + +type ECP_SWei_Aff*[F] = object + ## Elliptic curve point for a curve in Short Weierstrass form + ## y² = x³ + a x + b + ## + ## over a field F + x*, y*: F + func curve_eq_rhs*[F](y2: var F, x: F) = ## Compute the curve equation right-hand-side from field element `x` ## i.e. `y²` in `y² = x³ + a x + b` @@ -65,3 +79,35 @@ func isOnCurve*[F](x, y: F): SecretBool = rhs.curve_eq_rhs(x) return y2 == rhs + +func trySetFromCoordX*[F](P: var ECP_SWei_Aff[F], x: F): SecretBool = + ## Try to create a point the elliptic curve + ## y² = x³ + a x + b (affine coordinate) + ## + ## The `Z` coordinates is set to 1 + ## + ## return true and update `P` if `x` leads to a valid point + ## return false otherwise, in that case `P` is undefined. + ## + ## Note: Dedicated robust procedures for hashing-to-curve + ## will be provided, this is intended for testing purposes. + P.y.curve_eq_rhs(x) + # TODO: supports non p ≡ 3 (mod 4) modulus like BLS12-377 + result = sqrt_if_square(P.y) + +func neg*(P: var ECP_SWei_Aff, Q: ECP_SWei_Aff) = + ## Negate ``P`` + P.x = Q.x + P.y.neg(Q.y) + +func neg*(P: var ECP_SWei_Aff) = + ## Negate ``P`` + P.y.neg() + +func cneg*(P: var ECP_SWei_Aff, ctl: CTBool) = + ## Conditional negation. + ## Negate if ``ctl`` is true + var Q{.noInit.}: typeof(P) + Q.x = P.x + Q.y.neg(P.y) + P.ccopy(Q, ctl) diff --git a/constantine/elliptic/ec_weierstrass_projective.nim b/constantine/elliptic/ec_weierstrass_projective.nim index 97ea1f9..697ffe1 100644 --- a/constantine/elliptic/ec_weierstrass_projective.nim +++ b/constantine/elliptic/ec_weierstrass_projective.nim @@ -15,7 +15,7 @@ import # ############################################################ # -# Elliptic Curve in Weierstrass form +# Elliptic Curve in Short Weierstrass form # with Projective Coordinates # # ############################################################ @@ -104,9 +104,15 @@ func trySetFromCoordX*[F](P: var ECP_SWei_Proj[F], x: F): SecretBool = P.x = x P.z.setOne() +func neg*(P: var ECP_SWei_Proj, Q: ECP_SWei_Proj) = + ## Negate ``P`` + P.x = Q.x + P.y.neg(Q.y) + P.z = Q.z + func neg*(P: var ECP_SWei_Proj) = ## Negate ``P`` - P.y.neg(P.y) + P.y.neg() func cneg*(P: var ECP_SWei_Proj, ctl: CTBool) = ## Conditional negation. @@ -309,13 +315,17 @@ func diff*[F](r: var ECP_SWei_Proj[F], nQ.neg() r.sum(P, nQ) -func affineFromProjective*[F](aff: var ECP_SWei_Proj[F], proj: ECP_SWei_Proj) = +func affineFromProjective*[F](aff: var ECP_SWei_Aff[F], proj: ECP_SWei_Proj) = # TODO: for now we reuse projective coordinate backend # TODO: simultaneous inversion var invZ {.noInit.}: F invZ.inv(proj.z) - aff.z.setOne() aff.x.prod(proj.x, invZ) aff.y.prod(proj.y, invZ) + +func projectiveFromAffine*[F](proj: var ECP_SWei_Proj, aff: ECP_SWei_Aff[F]) {.inline.} = + proj.x = aff.x + proj.y = aff.y + proj.z.setOne() diff --git a/constantine/hash_to_curve/cofactors.nim b/constantine/hash_to_curve/cofactors.nim new file mode 100644 index 0000000..174c656 --- /dev/null +++ b/constantine/hash_to_curve/cofactors.nim @@ -0,0 +1,69 @@ +# 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/[tables, unittest, times], + # Internals + ../config/common, + ../arithmetic, + ../primitives, + ../towers, + ../config/curves, + ../io/io_bigints, + ../elliptic/[ec_weierstrass_projective, ec_scalar_mul] + +# ############################################################ +# +# Clear Cofactor - Naive +# +# ############################################################ + +const Cofactor_Eff_BN254_Nogami_G1 = BigInt[1].fromHex"0x1" +const Cofactor_Eff_BN254_Nogami_G2 = BigInt[254].fromHex"0x2523648240000001ba344d8000000008c2a2800000000016ad00000000000019" + ## G2.order // r + +const Cofactor_Eff_BN254_Snarks_G1 = BigInt[1].fromHex"0x1" +const Cofactor_Eff_BN254_Snarks_G2 = BigInt[254].fromHex"0x30644e72e131a029b85045b68181585e06ceecda572a2489345f2299c0f9fa8d" + ## G2.order // r + +# https://tools.ietf.org/html/draft-irtf-cfrg-hash-to-curve-09#section-8.8 +const Cofactor_Eff_BLS12_381_G1 = BigInt[64].fromHex"0xd201000000010001" + ## P -> (1 - x) P +const Cofactor_Eff_BLS12_381_G2 = BigInt[636].fromHex"0xbc69f08f2ee75b3584c6a0ea91b352888e2a8e9145ad7689986ff031508ffe1329c2f178731db956d82bf015d1212b02ec0ec69d7477c1ae954cbc06689f6a359894c0adebbf6b4e8020005aaa95551" + ## P -> (x^2 - x - 1) P + (x - 1) psi(P) + psi(psi(2P)) + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp[BN254_Nogami]]) {.inline.} = + ## Clear the cofactor of BN254_Nogami G1 + ## BN curve have a G1 cofactor of 1 so this is a no-op + discard + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp2[BN254_Nogami]]) {.inline.} = + ## Clear the cofactor of BN254_Snarks G2 + # Endomorphism acceleration cannot be used if cofactor is not cleared + P.scalarMulGeneric(Cofactor_Eff_BN254_Nogami_G2) + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp[BN254_Snarks]]) {.inline.} = + ## Clear the cofactor of BN254_Snarks G1 + ## BN curve have a G1 cofactor of 1 so this is a no-op + discard + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp2[BN254_Snarks]]) {.inline.} = + ## Clear the cofactor of BN254_Snarks G2 + # Endomorphism acceleration cannot be used if cofactor is not cleared + P.scalarMulGeneric(Cofactor_Eff_BN254_Snarks_G2) + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp[BLS12_381]]) {.inline.} = + ## Clear the cofactor of BLS12_381 G1 + ## BN curve have a G1 cofactor of 1 so this is a no-op + P.scalarMulGeneric(Cofactor_Eff_BLS12_381_G1) + +func clearCofactorReference*(P: var ECP_SWei_Proj[Fp2[BLS12_381]]) {.inline.} = + ## Clear the cofactor of BLS12_381 G2 + # Endomorphism acceleration cannot be used if cofactor is not cleared + P.scalarMulGeneric(Cofactor_Eff_BLS12_381_G2) diff --git a/constantine/io/io_ec.nim b/constantine/io/io_ec.nim index c548ef2..9305988 100644 --- a/constantine/io/io_ec.nim +++ b/constantine/io/io_ec.nim @@ -40,7 +40,7 @@ func toHex*(P: ECP_SWei_Proj): string = ## ## This proc output may change format in the future - var aff {.noInit.}: typeof(P) + var aff {.noInit.}: ECP_SWei_Aff[ECP_SWei_Proj.F] aff.affineFromProjective(P) result = "ECP[" & $aff.F & "](\n x: " @@ -70,3 +70,23 @@ func fromHex*(dst: var ECP_SWei_Proj, x0, x1, y0, y1: string): bool {.raises: [V dst.y.fromHex(y0, y1) dst.z.setOne() return bool(isOnCurve(dst.x, dst.y)) + +func fromHex*(dst: var ECP_SWei_Aff, x, y: string): bool {.raises: [ValueError].}= + ## Convert hex strings to a G1 curve point + ## Returns `false` + ## if there is no point with coordinates (`x`, `y`) on the curve + ## In that case, `dst` content is undefined. + static: doAssert dst.F is Fp, "dst must be on G1, an elliptic curve over 𝔽p" + dst.x.fromHex(x) + dst.y.fromHex(y) + return bool(isOnCurve(dst.x, dst.y)) + +func fromHex*(dst: var ECP_SWei_Aff, x0, x1, y0, y1: string): bool {.raises: [ValueError].}= + ## Convert hex strings to a G2 curve point + ## Returns `false` + ## if there is no point with coordinates (`x`, `y`) on the curve + ## In that case, `dst` content is undefined. + static: doAssert dst.F is Fp2, "dst must be on G2, an elliptic curve over 𝔽p2" + dst.x.fromHex(x0, x1) + dst.y.fromHex(y0, y1) + return bool(isOnCurve(dst.x, dst.y)) diff --git a/constantine/io/io_towers.nim b/constantine/io/io_towers.nim index 27e0334..41ab40c 100644 --- a/constantine/io/io_towers.nim +++ b/constantine/io/io_towers.nim @@ -26,7 +26,7 @@ import # # ############################################################ -func appendHex*(accum: var string, f: Fp2 or Fp6 or Fp12, order: static Endianness = bigEndian) = +func appendHex*(accum: var string, f: Fp2 or Fp4 or Fp6 or Fp12, order: static Endianness = bigEndian) = ## Hex accumulator accum.add static($f.typeof.genericHead() & '(') for fieldName, fieldValue in fieldPairs(f): @@ -36,7 +36,7 @@ func appendHex*(accum: var string, f: Fp2 or Fp6 or Fp12, order: static Endianne accum.appendHex(fieldValue, order) accum.add ")" -func toHex*(f: Fp2 or Fp6 or Fp12, order: static Endianness = bigEndian): string = +func toHex*(f: Fp2 or Fp4 or Fp6 or Fp12, order: static Endianness = bigEndian): string = ## Stringify a tower field element to hex. ## Note. Leading zeros are not removed. ## Result is prefixed with 0x @@ -59,3 +59,49 @@ func fromHex*(T: typedesc[Fp2], c0, c1: string): T {.raises: [ValueError].}= ## with dst = c0 + β * c1 ## β is the quadratic non-residue chosen to construct 𝔽p2 result.fromHex(c0, c1) + +func fromHex*(dst: var Fp4, + c0, c1, c2, c3: string) {.raises: [ValueError].}= + ## Convert 4 coordinates to an element of 𝔽p4 + dst.c0.fromHex(c0, c1) + dst.c1.fromHex(c2, c3) + +func fromHex*(T: typedesc[Fp4], + c0, c1, c2: string, + c3, c4, c5: string): T {.raises: [ValueError].}= + ## Convert 4 coordinates to an element of 𝔽p4 + result.fromHex(c0, c1, c2, c3) + +func fromHex*(dst: var Fp6, + c0, c1, c2: string, + c3, c4, c5: string) {.raises: [ValueError].}= + ## Convert 6 coordinates to an element of 𝔽p6 + dst.c0.fromHex(c0, c1) + dst.c1.fromHex(c2, c3) + dst.c2.fromHex(c4, c5) + +func fromHex*(T: typedesc[Fp6], + c0, c1, c2: string, + c3, c4, c5: string): T {.raises: [ValueError].}= + ## Convert 6 coordinates to an element of 𝔽p6 + result.fromHex(c0, c1, c2, c3, c4, c5) + +func fromHex*(dst: var Fp12, + c0, c1, c2, c3: string, + c4, c5, c6, c7: string, + c8, c9, c10, c11: string) {.raises: [ValueError].}= + ## Convert 12 coordinates to an element of 𝔽p12 + when Fp12.c0 is Fp6: + dst.c0.fromHex(c0, c1, c2, c3, c4, c5) + dst.c1.fromHex(c6, c7, c8, c9, c10, c11) + else: + dst.c0.fromHex(c0, c1, c2, c3) + dst.c1.fromHex(c4, c5, c6, c7) + dst.c2.fromHex(c8, c9, c10, c11) + +func fromHex*(T: typedesc[Fp12], + c0, c1, c2, c3: string, + c4, c5, c6, c7: string, + c8, c9, c10, c11: string): T {.raises: [ValueError].}= + ## Convert 12 coordinates to an element of 𝔽p12 + result.fromHex(c0, c1, c2, c3, c4, c5, c6, c7, c8, c9, c10, c11) diff --git a/constantine/isogeny/frobenius.nim b/constantine/isogeny/frobenius.nim index 6ea3164..f9d36b7 100644 --- a/constantine/isogeny/frobenius.nim +++ b/constantine/isogeny/frobenius.nim @@ -12,7 +12,7 @@ import ../io/io_towers, ../towers, ../arithmetic -# Frobenius automorphism +# Frobenius Map # ------------------------------------------------------------ # # https://en.wikipedia.org/wiki/Frobenius_endomorphism @@ -35,91 +35,13 @@ import # whether u = √-1 = i # or √-2 or √-5 -func frobenius*(r: var Fp2, a: Fp2, k: static int = 1) {.inline.} = +func frobenius_map*(r: var Fp2, a: Fp2, k: static int = 1) {.inline.} = ## Computes a^(p^k) ## The p-power frobenius automorphism on 𝔽p2 - r.c0 = a.c0 when (k and 1) == 1: - r.c1.neg(a.c1) + r.conj(a) else: - r.c1 = a.c1 - -# Frobenius endomorphism -# ------------------------------------------------------------ -# TODO: generate those constants via Sage in a Json file -# and parse at compile-time - -# Constants: -# Assuming embedding degree of 12 and a sextic twist -# with SNR the sextic non-residue -# -# BN254_Snarks is a D-Twist: SNR^((p-1)/6) -const FrobConst_BN254_Snarks_psi1_coef1 = Fp2[BN254_Snarks].fromHex( - "0x1284b71c2865a7dfe8b99fdd76e68b605c521e08292f2176d60b35dadcc9e470", - "0x246996f3b4fae7e6a6327cfe12150b8e747992778eeec7e5ca5cf05f80f362ac" -) -# SNR^((p-1)/3) -const FrobConst_BN254_Snarks_psi1_coef2 = Fp2[BN254_Snarks].fromHex( - "0x2fb347984f7911f74c0bec3cf559b143b78cc310c2c3330c99e39557176f553d", - "0x16c9e55061ebae204ba4cc8bd75a079432ae2a1d0b7c9dce1665d51c640fcba2" -) -# SNR^((p-1)/2) -const FrobConst_BN254_Snarks_psi1_coef3 = Fp2[BN254_Snarks].fromHex( - "0x63cf305489af5dcdc5ec698b6e2f9b9dbaae0eda9c95998dc54014671a0135a", - "0x7c03cbcac41049a0704b5a7ec796f2b21807dc98fa25bd282d37f632623b0e3" -) -# norm(SNR)^((p-1)/3) -const FrobConst_BN254_Snarks_psi2_coef2 = Fp2[BN254_Snarks].fromHex( - "0x30644e72e131a0295e6dd9e7e0acccb0c28f069fbb966e3de4bd44e5607cfd48", - "0x0" -) - -# BLS12_377 is a D-Twist: SNR^((p-1)/6) -const FrobConst_BLS12_377_psi1_coef1 = Fp2[BLS12_377].fromHex( - "0x9a9975399c019633c1e30682567f915c8a45e0f94ebc8ec681bf34a3aa559db57668e558eb0188e938a9d1104f2031", - "0x0" -) -# SNR^((p-1)/3) -const FrobConst_BLS12_377_psi1_coef2 = Fp2[BLS12_377].fromHex( - "0x9b3af05dd14f6ec619aaf7d34594aabc5ed1347970dec00452217cc900000008508c00000000002", - "0x0" -) -# SNR^((p-1)/2) -const FrobConst_BLS12_377_psi1_coef3 = Fp2[BLS12_377].fromHex( - "0x1680a40796537cac0c534db1a79beb1400398f50ad1dec1bce649cf436b0f6299588459bff27d8e6e76d5ecf1391c63", - "0x0" -) -# norm(SNR)^((p-1)/3) -const FrobConst_BLS12_377_psi2_coef2 = Fp2[BLS12_377].fromHex( - "0x9b3af05dd14f6ec619aaf7d34594aabc5ed1347970dec00452217cc900000008508c00000000001", - "0x0" -) - -# BLS12_381 is a M-twist: (1/SNR)^((p-1)/6) -const FrobConst_BLS12_381_psi1_coef1 = Fp2[BLS12_381].fromHex( - "0x5b2cfd9013a5fd8df47fa6b48b1e045f39816240c0b8fee8beadf4d8e9c0566c63a3e6e257f87329b18fae980078116", - "0x5b2cfd9013a5fd8df47fa6b48b1e045f39816240c0b8fee8beadf4d8e9c0566c63a3e6e257f87329b18fae980078116" -) -# (1/SNR)^((p-1)/3) -const FrobConst_BLS12_381_psi1_coef2 = Fp2[BLS12_381].fromHex( - "0x0", - "0x1a0111ea397fe699ec02408663d4de85aa0d857d89759ad4897d29650fb85f9b409427eb4f49fffd8bfd00000000aaad" -) -# (1/SNR)^((p-1)/2) -const FrobConst_BLS12_381_psi1_coef3 = Fp2[BLS12_381].fromHex( - "0x135203e60180a68ee2e9c448d77a2cd91c3dedd930b1cf60ef396489f61eb45e304466cf3e67fa0af1ee7b04121bdea2", - "0x6af0e0437ff400b6831e36d6bd17ffe48395dabc2d3435e77f76e17009241c5ee67992f72ec05f4c81084fbede3cc09" -) -# norm(SNR)^((p-1)/3) -const FrobConst_BLS12_381_psi2_coef2 = Fp2[BLS12_381].fromHex( - "0x1a0111ea397fe699ec02408663d4de85aa0d857d89759ad4897d29650fb85f9b409427eb4f49fffd8bfd00000000aaac", - "0x0" -) - -{.experimental: "dynamicBindSym".} - -macro frobConst(C: static Curve, psipow, coefpow: static int): untyped = - return bindSym("FrobConst_" & $C & "_psi" & $psipow & "_coef" & $coefpow) + r = a template mulCheckSparse[Fp2](a: var Fp2, b: Fp2) = when b.c0.isZero().bool: @@ -129,11 +51,129 @@ template mulCheckSparse[Fp2](a: var Fp2, b: Fp2) = else: a *= b +# Frobenius map - on extension fields +# ----------------------------------------------------------------- + +{.experimental: "dynamicBindSym".} + +macro frobMapConst(C: static Curve, coefpow, frobpow: static int): untyped = + return bindSym("FrobMapConst_" & $C & "_coef" & $coefpow & "_pow" & $frobpow) + +# SNR^((p-1)/6)^3 +const FrobMapConst_BLS12_381_coef3_pow1 = Fp2[BLS12_381].fromHex( + "0x6af0e0437ff400b6831e36d6bd17ffe48395dabc2d3435e77f76e17009241c5ee67992f72ec05f4c81084fbede3cc09", + "0x6af0e0437ff400b6831e36d6bd17ffe48395dabc2d3435e77f76e17009241c5ee67992f72ec05f4c81084fbede3cc09" + ) + +func frobenius_map*(r: var Fp4, a: Fp4, k: static int = 1) {.inline.} = + ## Computes a^(p^k) + ## The p-power frobenius automorphism on 𝔽p4 + r.c0.frobenius_map(a.c0, k) + r.c1.frobenius_map(a.c1, k) + + r.c1.mulCheckSparse frobMapConst(Fp4.C, 3, k) + +# ψ (Psi) - Untwist-Frobenius-Twist Endomorphisms on twisted curves +# ----------------------------------------------------------------- +# TODO: generate those constants via Sage in a Json file +# and parse at compile-time + +# Constants: +# Assuming embedding degree of 12 and a sextic twist +# with SNR the sextic non-residue +# +# BN254_Snarks is a D-Twist: SNR^((p-1)/6) +const FrobPsiConst_BN254_Snarks_psi1_coef1 = Fp2[BN254_Snarks].fromHex( + "0x1284b71c2865a7dfe8b99fdd76e68b605c521e08292f2176d60b35dadcc9e470", + "0x246996f3b4fae7e6a6327cfe12150b8e747992778eeec7e5ca5cf05f80f362ac" +) +# SNR^((p-1)/3) +const FrobPsiConst_BN254_Snarks_psi1_coef2 = Fp2[BN254_Snarks].fromHex( + "0x2fb347984f7911f74c0bec3cf559b143b78cc310c2c3330c99e39557176f553d", + "0x16c9e55061ebae204ba4cc8bd75a079432ae2a1d0b7c9dce1665d51c640fcba2" +) +# SNR^((p-1)/2) +const FrobPsiConst_BN254_Snarks_psi1_coef3 = Fp2[BN254_Snarks].fromHex( + "0x63cf305489af5dcdc5ec698b6e2f9b9dbaae0eda9c95998dc54014671a0135a", + "0x7c03cbcac41049a0704b5a7ec796f2b21807dc98fa25bd282d37f632623b0e3" +) +# norm(SNR)^((p-1)/3) +const FrobPsiConst_BN254_Snarks_psi2_coef2 = Fp2[BN254_Snarks].fromHex( + "0x30644e72e131a0295e6dd9e7e0acccb0c28f069fbb966e3de4bd44e5607cfd48", + "0x0" +) + +# BN254_Nogami is a D-Twist: SNR^((p-1)/6) +const FrobPsiConst_BN254_Nogami_psi1_coef1 = Fp2[BN254_Nogami].fromHex( + "0x1b377619212e7c8cb6499b50a846953f850974924d3f77c2e17de6c06f2a6de9", + "0x9ebee691ed1837503eab22f57b96ac8dc178b6db2c08850c582193f90d5922a" +) +# SNR^((p-1)/3) +const FrobPsiConst_BN254_Nogami_psi1_coef2 = Fp2[BN254_Nogami].fromHex( + "0x0", + "0x25236482400000017080eb4000000006181800000000000cd98000000000000b" +) +# SNR^((p-1)/2) +const FrobPsiConst_BN254_Nogami_psi1_coef3 = Fp2[BN254_Nogami].fromHex( + "0x23dfc9d1a39f4db8c69b87a8848aa075a7333a0e62d78cbf4b1b8eeae58b81c5", + "0x23dfc9d1a39f4db8c69b87a8848aa075a7333a0e62d78cbf4b1b8eeae58b81c5" +) +# norm(SNR)^((p-1)/3) +const FrobPsiConst_BN254_Nogami_psi2_coef2 = Fp2[BN254_Nogami].fromHex( + "0x49b36240000000024909000000000006cd80000000000007", + "0x0" +) + +# BLS12_377 is a D-Twist: SNR^((p-1)/6) +const FrobPsiConst_BLS12_377_psi1_coef1 = Fp2[BLS12_377].fromHex( + "0x9a9975399c019633c1e30682567f915c8a45e0f94ebc8ec681bf34a3aa559db57668e558eb0188e938a9d1104f2031", + "0x0" +) +# SNR^((p-1)/3) +const FrobPsiConst_BLS12_377_psi1_coef2 = Fp2[BLS12_377].fromHex( + "0x9b3af05dd14f6ec619aaf7d34594aabc5ed1347970dec00452217cc900000008508c00000000002", + "0x0" +) +# SNR^((p-1)/2) +const FrobPsiConst_BLS12_377_psi1_coef3 = Fp2[BLS12_377].fromHex( + "0x1680a40796537cac0c534db1a79beb1400398f50ad1dec1bce649cf436b0f6299588459bff27d8e6e76d5ecf1391c63", + "0x0" +) +# norm(SNR)^((p-1)/3) +const FrobPsiConst_BLS12_377_psi2_coef2 = Fp2[BLS12_377].fromHex( + "0x9b3af05dd14f6ec619aaf7d34594aabc5ed1347970dec00452217cc900000008508c00000000001", + "0x0" +) + +# BLS12_381 is a M-twist: (1/SNR)^((p-1)/6) +const FrobPsiConst_BLS12_381_psi1_coef1 = Fp2[BLS12_381].fromHex( + "0x5b2cfd9013a5fd8df47fa6b48b1e045f39816240c0b8fee8beadf4d8e9c0566c63a3e6e257f87329b18fae980078116", + "0x5b2cfd9013a5fd8df47fa6b48b1e045f39816240c0b8fee8beadf4d8e9c0566c63a3e6e257f87329b18fae980078116" +) +# (1/SNR)^((p-1)/3) +const FrobPsiConst_BLS12_381_psi1_coef2 = Fp2[BLS12_381].fromHex( + "0x0", + "0x1a0111ea397fe699ec02408663d4de85aa0d857d89759ad4897d29650fb85f9b409427eb4f49fffd8bfd00000000aaad" +) +# (1/SNR)^((p-1)/2) +const FrobPsiConst_BLS12_381_psi1_coef3 = Fp2[BLS12_381].fromHex( + "0x135203e60180a68ee2e9c448d77a2cd91c3dedd930b1cf60ef396489f61eb45e304466cf3e67fa0af1ee7b04121bdea2", + "0x6af0e0437ff400b6831e36d6bd17ffe48395dabc2d3435e77f76e17009241c5ee67992f72ec05f4c81084fbede3cc09" +) +# norm(SNR)^((p-1)/3) +const FrobPsiConst_BLS12_381_psi2_coef2 = Fp2[BLS12_381].fromHex( + "0x1a0111ea397fe699ec02408663d4de85aa0d857d89759ad4897d29650fb85f9b409427eb4f49fffd8bfd00000000aaac", + "0x0" +) + +macro frobPsiConst(C: static Curve, psipow, coefpow: static int): untyped = + return bindSym("FrobPsiConst_" & $C & "_psi" & $psipow & "_coef" & $coefpow) + func frobenius_psi*[PointG2](r: var PointG2, P: PointG2) = ## "Untwist-Frobenius-Twist" endomorphism ## r = ψ(P) for coordR, coordP in fields(r, P): - coordR.frobenius(coordP, 1) + coordR.frobenius_map(coordP, 1) # With ξ (xi) the sextic non-residue # c = ξ^((p-1)/6) for D-Twist @@ -142,14 +182,14 @@ func frobenius_psi*[PointG2](r: var PointG2, P: PointG2) = # c1_2 = c² # c1_3 = c³ - r.x.mulCheckSparse frobConst(PointG2.F.C, psipow=1, coefpow=2) - r.y.mulCheckSparse frobConst(PointG2.F.C, psipow=1, coefpow=3) + r.x.mulCheckSparse frobPsiConst(PointG2.F.C, psipow=1, coefpow=2) + r.y.mulCheckSparse frobPsiConst(PointG2.F.C, psipow=1, coefpow=3) func frobenius_psi2*[PointG2](r: var PointG2, P: PointG2) = ## "Untwist-Frobenius-Twist" endomorphism applied twice ## r = ψ(ψ(P)) for coordR, coordP in fields(r, P): - coordR.frobenius(coordP, 2) + coordR.frobenius_map(coordP, 2) # With ξ (xi) the sextic non-residue # c = ξ for D-Twist @@ -181,5 +221,5 @@ func frobenius_psi2*[PointG2](r: var PointG2, P: PointG2) = # c2_3 ≡ -1 (mod p²) # QED - r.x.mulCheckSparse frobConst(PointG2.F.C, psipow=2, coefpow=2) + r.x.mulCheckSparse frobPsiConst(PointG2.F.C, psipow=2, coefpow=2) r.y.neg(r.y) diff --git a/constantine/pairing/README.md b/constantine/pairing/README.md index 486c87d..be35783 100644 --- a/constantine/pairing/README.md +++ b/constantine/pairing/README.md @@ -15,10 +15,28 @@ Ben Lynn, 2007\ https://crypto.stanford.edu/pbc/thesis.pdf +- Faster Pairing Computations on Curves with High-Degree Twists + Craig Costello, Tanja Lange, and Michael Naehrig, 2009 + https://eprint.iacr.org/2009/615.pdf + +- High-Speed Software Implementation of the Optimal Ate Pairing over Barreto-Naehrig Curves\ + Jean-Luc Beuchat and Jorge Enrique González Díaz and Shigeo Mitsunari and Eiji Okamoto and Francisco Rodríguez-Henríquez and Tadanori Teruya, 2010\ + https://eprint.iacr.org/2010/354 + - An Analysis of Affine Coordinates for Pairing Computation\ Kristin Lauter, Peter L. Montgomery, and Michael Naehrig, 2010\ https://eprint.iacr.org/2010/363.pdf +- Faster Explicit Formulas for Computing Pairings over Ordinary Curves\ + Diego F. Aranha and Koray Karabina and Patrick Longa and Catherine H. Gebotys and Julio López, 2010\ + https://eprint.iacr.org/2010/526.pdf\ + https://www.iacr.org/archive/eurocrypt2011/66320047/66320047.pdf + +- Avoiding Full Extension Field Arithmetic in Pairing Computations + Craig Costello, Colin Boyd, + Juan Manuel Gonzalez Nieto, and Kenneth Koon-Ho Wong, 2010 + https://eprint.iacr.org/2010/104.pdf + - Pairings for beginners\ Craig Costello, 2012 (?)\ http://www.craigcostello.com.au/pairings/PairingsForBeginners.pdf @@ -28,6 +46,16 @@ Craig Costello, 2012\ https://eprints.qut.edu.au/61037/1/Craig_Costello_Thesis.pdf +- Efficient Implementation of Bilinear Pairings on ARM Processors + Gurleen Grewal, Reza Azarderakhsh, + Patrick Longa, Shi Hu, and David Jao, 2012 + https://eprint.iacr.org/2012/408.pdf + +- The Realm of the Pairings\ + Diego F. Aranha and Paulo S. L. M. Barreto and Patrick Longa and Jefferson E. Ricardini\ + https://eprint.iacr.org/2013/722.pdf\ + http://sac2013.irmacs.sfu.ca/slides/s1.pdf + - Efficient Implementations of Pairing-Based Cryptography on Embedded Systems\ Master Thesis\ Rajeev Verma, 2015\ @@ -38,10 +66,21 @@ Razvan Barbulescu, Nadia El Mrabet, and Loubna Ghammam, 2019\ https://hal.archives-ouvertes.fr/hal-02129868/file/2019-485.pdf -- A short-list of pairing-friendly curves resistantto Special TNFS at the 128-bit security level\ - Aurore Guillevic\ +- Pairing Implementation Revisited\ + Mike Scott, 2019\ + https://eprint.iacr.org/2019/077.pdf + +- A short-list of pairing-friendly curves resistant to Special TNFS at the 128-bit security level\ + Aurore Guillevic, 2019\ https://eprint.iacr.org/2019/1371.pdf +- Efficient Final Exponentiation + via Cyclotomic Structure for Pairings + over Families of Elliptic Curves + Daiki Hayashida and Kenichiro Hayasaka + and Tadanori Teruya, 2020 + https://eprint.iacr.org/2020/875.pdf + ### Presentations - Introduction to pairings\ diff --git a/constantine/pairing/gt_fp12.nim b/constantine/pairing/gt_fp12.nim new file mode 100644 index 0000000..1fb9c97 --- /dev/null +++ b/constantine/pairing/gt_fp12.nim @@ -0,0 +1,168 @@ +# 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 + ../primitives, + ../config/[common, curves], + ../arithmetic, + ../towers, + ./lines_projective + +# ############################################################ +# +# Sparse Multiplication +# and cyclotomic squaring +# for elements of Gₜ = E(Fp¹²) +# +# ############################################################ + +# - Pairing Implementation Revisited +# Michael Scott, 2019 +# https://eprint.iacr.org/2019/077 +# +# - Efficient Implementation of Bilinear Pairings on ARM Processors +# Gurleen Grewal, Reza Azarderakhsh, +# Patrick Longa, Shi Hu, and David Jao, 2012 +# https://eprint.iacr.org/2012/408.pdf +# +# - High-Speed Software Implementation of the Optimal Ate Pairing over Barreto-Naehrig Curves\ +# Jean-Luc Beuchat and Jorge Enrique González Díaz and Shigeo Mitsunari and Eiji Okamoto and Francisco Rodríguez-Henríquez and Tadanori Teruya, 2010\ +# https://eprint.iacr.org/2010/354 +# +# - Faster Pairing Computations on Curves with High-Degree Twists +# Craig Costello, Tanja Lange, and Michael Naehrig, 2009 +# https://eprint.iacr.org/2009/615.pdf + +# TODO: we assume an embedding degree k of 12 and a sextic twist. +# -> Generalize to KSS (k=18), BLS24 and BLS48 curves +# +# TODO: we assume a 2-3-2 towering scheme +# +# TODO: merge that in the quadratic/cubic files + +# 𝔽p12 - Sparse functions +# ---------------------------------------------------------------- + +func mul_sparse_by_0y0*[C: static Curve](r: var Fp6[C], a: Fp6[C], b: Fp2[C]) = + ## Sparse multiplication of an Fp6 element + ## with coordinates (a₀, a₁, a₂) by (0, b₁, 0) + # TODO: make generic and move to tower_field_extensions + + # v0 = a0 b0 = 0 + # v1 = a1 b1 + # v2 = a2 b2 = 0 + # + # r0 = ξ ((a1 + a2) * (b1 + b2) - v1 - v2) + v0 + # = ξ (a1 b1 + a2 b1 - v1) + # = ξ a2 b1 + # r1 = (a0 + a1) * (b0 + b1) - v0 - v1 + β v2 + # = a0 b1 + a1 b1 - v1 + # = a0 b1 + # r2 = (a0 + a2) * (b0 + b2) - v0 - v2 + v1 + # = v1 + # = a1 b1 + + r.c0.prod(a.c2, b) + r.c0 *= ξ + r.c1.prod(a.c0, b) + r.c2.prod(a.c1, b) + +func mul_by_line_xy0*[C: static Curve, twist: static SexticTwist]( + r: var Fp6[C], + a: Fp6[C], + b: Line[Fp2[C], twist]) = + ## Sparse multiplication of an 𝔽p6 + ## with coordinates (a₀, a₁, a₂) by a line (x, y, 0) + ## The z coordinates in the line will be ignored. + ## `r` and `a` must not alias + var + v0 {.noInit.}: Fp2[C] + v1 {.noInit.}: Fp2[C] + + v0.prod(a.c0, b.x) + v1.prod(a.c1, b.y) + r.c0.prod(a.c2, b.y) + r.c0 *= ξ + r.c0 += v0 + + r.c1.sum(a.c0, a.c1) # Error when r and a alias as r.c0 was updated + r.c2.sum(b.x, b.y) + r.c1 *= r.c2 + r.c1 -= v0 + r.c1 -= v1 + + r.c2.prod(a.c2, b.x) + r.c2 += v1 + +func mul_sparse_by_line_xy00z0*[C: static Curve, Tw: static SexticTwist]( + f: var Fp12[C], l: Line[Fp2[C], Tw]) = + ## Sparse multiplication of an 𝔽p12 element + ## by a sparse 𝔽p12 element coming from an D-Twist line function. + ## The sparse element is represented by a packed Line type + ## with coordinate (x,y,z) matching 𝔽p12 coordinates xy00z0 (TODO: verify this) + + static: doAssert f.c0.typeof is Fp6, "This assumes 𝔽p12 as a quadratic extension of 𝔽p6" + + var + v0 {.noInit.}: Fp6[C] + v1 {.noInit.}: Fp6[C] + v2 {.noInit.}: Line[Fp2[C], Tw] + v3 {.noInit.}: Fp6[C] + + v0.mul_by_line_xy0(f.c0, l) + v1.mul_sparse_by_0y0(f.c1, l.z) + + v2.x = l.x + v2.y.sum(l.y, l.z) + f.c1 += f.c0 + v3.mul_by_line_xy0(f.c1, v2) + v3 -= v0 + v3 -= v1 + f.c1 = v3 + + v3.c0 = ξ * v1.c2 + v3.c0 += v0.c0 + v3.c1.sum(v0.c1, v1.c0) + v3.c2.sum(v0.c2, v1.c1) + f.c0 = v3 + +func mul_sparse_by_line_xyz000*[C: static Curve, Tw: static SexticTwist]( + f: var Fp12[C], l: Line[Fp2[C], Tw]) = + ## Sparse multiplication of an 𝔽p12 element + ## by a sparse 𝔽p12 element coming from an D-Twist line function. + ## The sparse element is represented by a packed Line type + ## with coordinates (x,y,z) matching 𝔽p12 coordinates xyz000 + + static: doAssert f.c0.typeof is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4" + + var v: Fp12[C] + v.c0.c0 = l.x + v.c0.c1 = l.y + v.c1.c0 = l.z + + f *= v + +func mul_sparse_by_line_xy000z*[C: static Curve, Tw: static SexticTwist]( + f: var Fp12[C], l: Line[Fp2[C], Tw]) = + + static: doAssert f.c0.typeof is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4" + + var v: Fp12[C] + v.c0.c0 = l.x + v.c0.c1 = l.y + v.c2.c1 = l.z + + f *= v + +# Gₜ = 𝔽p12 - Cyclotomic functions +# ---------------------------------------------------------------- +# A cyclotomic group is a subgroup of Fp^n defined by +# +# GΦₙ(p) = {α ∈ Fpⁿ : α^Φₙ(p) = 1} +# +# The result of any pairing is in a cyclotomic subgroup diff --git a/constantine/pairing/lines_projective.nim b/constantine/pairing/lines_projective.nim new file mode 100644 index 0000000..9e115d2 --- /dev/null +++ b/constantine/pairing/lines_projective.nim @@ -0,0 +1,234 @@ +# 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 + std/typetraits, + ../primitives, + ../config/[common, curves], + ../arithmetic, + ../towers, + ../elliptic/[ + ec_weierstrass_affine, + ec_weierstrass_projective + ], + ../io/io_towers + +# ############################################################ +# +# Miller Loop's Line Evaluation +# with projective coordinates +# +# ############################################################ +# +# - Pairing Implementation Revisited +# Michael Scott, 2019 +# https://eprint.iacr.org/2019/077 +# +# - The Realm of the Pairings +# Diego F. Aranha and Paulo S. L. M. Barreto +# and Patrick Longa and Jefferson E. Ricardini, 2013 +# https://eprint.iacr.org/2013/722.pdf +# http://sac2013.irmacs.sfu.ca/slides/s1.pdf + +# +# TODO: Implement fused line doubling and addition +# from Costello2009 or Aranha2010 +# We don't need the complete formulae in the Miller Loop + +type + Line*[F; twist: static SexticTwist] = object + ## Packed line representation over a E'(Fp^k/d) + ## with k the embedding degree and d the twist degree + ## i.e. for a curve with embedding degree 12 and sextic twist + ## F is Fp2 + ## + ## Assuming a Sextic Twist + ## + ## Out of 6 Fp2 coordinates, 3 are 0 and + ## the non-zero coordinates depend on the twist kind. + ## + ## For a D-twist, + ## (x, y, z) corresponds to an sparse element of Fp12 + ## with Fp2 coordinates: xy00z0 + ## For a M-Twist + ## (x, y, z) corresponds to an sparse element of Fp12 + ## with Fp2 coordinates: xyz000 + x*, y*, z*: F + +func toHex*(line: Line, order: static Endianness = bigEndian): string = + result = static($line.typeof.genericHead() & '(') + for fieldName, fieldValue in fieldPairs(line): + when fieldName != "x": + result.add ", " + result.add fieldName & ": " + result.appendHex(fieldValue, order) + result.add ")" + +# Line evaluation +# -------------------------------------------------- + +func `*=`(a: var Fp2, b: Fp) = + ## Multiply an element of Fp2 by an element of Fp + # TODO: make generic and move to tower_field_extensions + a.c0 *= b + a.c1 *= b + +func line_update(line: var Line, P: ECP_SWei_Aff) = + ## Update the line evaluation with P + ## after addition or doubling + ## P in G1 + line.x *= P.y + line.z *= P.x + +func line_eval_double*(line: var Line, T: ECP_SWei_Proj) = + ## Evaluate the line function for doubling + ## i.e. the tangent at T + ## + ## With T in homogenous projective coordinates (X, Y, Z) + ## And ξ the sextic non residue to construct 𝔽p4 / 𝔽p6 / 𝔽p12 + ## + ## M-Twist: + ## A = -2ξ Y.Z + ## B = 3bξ Z² - Y² + ## C = 3 X² + ## + ## D-Twist are scaled by ξ to avoid dividing by ξ: + ## A = -2ξ Y.Z + ## B = 3b Z² - ξY² + ## C = 3ξ X² + ## + ## Instead of + ## - equation 10 from The Real of pairing, Aranha et al, 2013 + ## - or chapter 3 from pairing Implementation Revisited, Scott 2019 + ## A = -2 Y.Z + ## B = 3b/ξ Z² - Y² + ## C = 3 X² + ## + ## A constant factor will be wiped by the final exponentiation + ## as for all non-zero α ∈ GF(pᵐ) + ## with + ## - p odd prime + ## - and gcd(α,pᵐ) = 1 (i.e. the extension field pᵐ is using irreducible polynomials) + ## + ## Little Fermat holds and we have + ## α^(pᵐ - 1) ≡ 1 (mod pᵐ) + ## + ## The final exponent is of the form + ## (pᵏ-1)/r + ## + ## A constant factor on twisted coordinates pᵏᐟᵈ + ## is a constant factor on pᵏ with d the twisting degree + ## and so will be elminated. QED. + var v {.noInit.}: Line.F + const b3 = 3 * ECP_SWei_Proj.F.C.getCoefB() + + template A: untyped = line.x + template B: untyped = line.y + template C: untyped = line.z + + A = T.y # A = Y + v = T.y # v = Y + B = T.z # B = Z + C = T.x # C = X + + A *= B # A = Y.Z + C.square() # C = X² + v.square() # v = Y² + B.square() # B = Z² + + A.double() # A = 2 Y.Z + A.neg() # A = -2 Y.Z + A *= SexticNonResidue # A = -2 ξ Y.Z + + B *= b3 # B = 3b Z² + C *= 3 # C = 3X² + when ECP_SWei_Proj.F.C.getSexticTwist() == M_Twist: + B *= SexticNonResidue # B = 3b' Z² = 3bξ Z² + elif ECP_SWei_Proj.F.C.getSexticTwist() == D_Twist: + v *= SexticNonResidue # v = ξ Y² + C *= SexticNonResidue # C = 3ξ X² + else: + {.error: "unreachable".} + + B -= v # B = 3bξ Z² - Y² (M-twist) + # B = 3b Z² - ξ Y² (D-twist) + +func line_eval_add*(line: var Line, T: ECP_SWei_Proj, Q: ECP_SWei_Aff) = + ## Evaluate the line function for addition + ## i.e. the line between T and Q + ## + ## With T in homogenous projective coordinates (X, Y, Z) + ## And ξ the sextic non residue to construct 𝔽p4 / 𝔽p6 / 𝔽p12 + ## + ## M-Twist: + ## A = ξ (X₁ - Z₁X₂) + ## B = (Y₁ - Z₁Y₂) X₂ - (X₁ - Z₁X₂) Y₂ + ## C = - (Y₁ - Z₁Y₂) + ## + ## D-Twist: + ## A = X₁ - Z₁X₂ + ## B = (Y₁ - Z₁Y₂) X₂ - (X₁ - Z₁X₂) Y₂ + ## C = - (Y₁ - Z₁Y₂) + ## + ## Note: There is no need for complete formula as + ## we have T ∉ [Q, -Q] in the Miller loop doubling-and-add + ## i.e. the line cannot be vertical + var v {.noInit.}: Line.F + + template A: untyped = line.x + template B: untyped = line.y + template C: untyped = line.z + + A = T.x # A = X₁ + v = T.z # v = Z₁ + B = T.z # B = Z₁ + C = T.y # C = Y₁ + + v *= Q.y # v = Z₁Y₂ + B *= Q.x # B = Z₁X₂ + + A -= B # A = X₁-Z₁X₂ + C -= v # C = Y₁-Z₁Y₂ + + v = A # v = X₁-Z₁X₂ + when ECP_SWei_Proj.F.C.getSexticTwist() == M_Twist: + A *= SexticNonResidue # A = ξ (X₁ - Z₁X₂) + + v *= Q.y # v = (X₁-Z₁X₂) Y₂ + B = C # B = Y₁-Z₁Y₂ + B *= Q.x # B = (Y₁-Z₁Y₂) X₂ + B -= v # B = (Y₁-Z₁Y₂) X₂ - (X₁-Z₁X₂) Y₂ + + C.neg() # C = -(Y₁-Z₁Y₂) + +func line_double*(line: var Line, T: var ECP_SWei_Proj, P: ECP_SWei_Aff) = + ## Doubling step of the Miller loop + ## T in G2, P in G1 + ## + ## Compute lt,t(P) + ## + # TODO fused line doubling from Costello 2009, Grewal 2012, Aranha 2013 + line_eval_double(line, T) + line.line_update(P) + T.double() + +func line_add*[C]( + line: var Line, + T: var ECP_SWei_Proj[Fp2[C]], + Q: ECP_SWei_Aff[Fp2[C]], P: ECP_SWei_Aff[Fp[C]]) = + ## Addition step of the Miller loop + ## T and Q in G2, P in G1 + ## + ## Compute lt,q(P) + # TODO fused line addition from Costello 2009, Grewal 2012, Aranha 2013 + line_eval_add(line, T, Q) + line.line_update(P) + # TODO: mixed addition + var QProj {.noInit.}: ECP_SWei_Proj[Fp2[C]] + QProj.projectiveFromAffine(Q) + T += QProj diff --git a/constantine/pairing/pairing_bls12.nim b/constantine/pairing/pairing_bls12.nim new file mode 100644 index 0000000..12a8079 --- /dev/null +++ b/constantine/pairing/pairing_bls12.nim @@ -0,0 +1,152 @@ +# 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 + ../primitives, + ../config/[common, curves], + ../arithmetic, + ../towers, + ../io/io_bigints, + ../elliptic/[ + ec_weierstrass_affine, + ec_weierstrass_projective + ], + ./lines_projective, + ./gt_fp12 + +# ############################################################ +# +# Optimal ATE pairing for +# BLS12 curves +# +# ############################################################ + +# - Efficient Final Exponentiation +# via Cyclotomic Structure for Pairings +# over Families of Elliptic Curves +# Daiki Hayashida and Kenichiro Hayasaka +# and Tadanori Teruya, 2020 +# https://eprint.iacr.org/2020/875.pdf +# +# - Faster Pairing Computations on Curves with High-Degree Twists +# Craig Costello, Tanja Lange, and Michael Naehrig, 2009 +# https://eprint.iacr.org/2009/615.pdf + +# TODO: implement quadruple-and-add and octuple-and-add +# from Costello2009 to trade multiplications in Fpᵏ +# for multiplications in Fp + +# TODO: should be part of curve parameters +const BLS12_381_param = block: + # BLS Miller loop is parametrized by u + BigInt[64+2].fromHex("0xd201000000010000") # +2 so that we can take *3 and NAF encode it + +const BLS12_381_param_isNeg = true + +const BLS12_381_finalexponent = block: + # (p^12 - 1) / r + # BigInt[4314].fromHex"0x2ee1db5dcc825b7e1bda9c0496a1c0a89ee0193d4977b3f7d4507d07363baa13f8d14a917848517badc3a43d1073776ab353f2c30698e8cc7deada9c0aadff5e9cfee9a074e43b9a660835cc872ee83ff3a0f0f1c0ad0d6106feaf4e347aa68ad49466fa927e7bb9375331807a0dce2630d9aa4b113f414386b0e8819328148978e2b0dd39099b86e1ab656d2670d93e4d7acdd350da5359bc73ab61a0c5bf24c374693c49f570bcd2b01f3077ffb10bf24dde41064837f27611212596bc293c8d4c01f25118790f4684d0b9c40a68eb74bb22a40ee7169cdc1041296532fef459f12438dfc8e2886ef965e61a474c5c85b0129127a1b5ad0463434724538411d1676a53b5a62eb34c05739334f46c02c3f0bd0c55d3109cd15948d0a1fad20044ce6ad4c6bec3ec03ef19592004cedd556952c6d8823b19dadd7c2498345c6e5308f1c511291097db60b1749bf9b71a9f9e0100418a3ef0bc627751bbd81367066bca6a4c1b6dcfc5cceb73fc56947a403577dfa9e13c24ea820b09c1d9f7c31759c3635de3f7a3639991708e88adce88177456c49637fd7961be1a4c7e79fb02faa732e2f3ec2bea83d196283313492caa9d4aff1c910e9622d2a73f62537f2701aaef6539314043f7bbce5b78c7869aeb2181a67e49eeed2161daf3f881bd88592d767f67c4717489119226c2f011d4cab803e9d71650a6f80698e2f8491d12191a04406fbc8fbd5f48925f98630e68bfb24c0bcb9b55df57510" + # (p^12 - 1) / r * 3 + BigInt[4316].fromHex"0x8ca592196587127a538fd40dc3e541f9dca04bb7dc671be77cf17715a2b2fe3bea73dfb468d8f473094aecb7315a664019fbd84913caba6579c08fd42009fe1bd6fcbce15eacb2cf3218a165958cb8bfdae2d2d54207282314fc0dea9d6ff3a07dbd34efb77b732ba5f994816e296a72928cfee133bdc3ca9412b984b9783d9c6aa81297ab1cd294a502304773528bbae8706979f28efa0d355b0224e2513d6e4a5d3bb4dde0523678105d9167ff1323d6e99ac312d8a7d762336370c4347bb5a7e405d6f3496b2dd38e722d4c1f3ac25e3167ec2cb543d69430c37c2f98fcdd0dd36caa9f5aa7994cec31b24ed5e515911037b376e521070d29c9d56cfa8c3574363efb20f28c19e4105ab99edd44084bd23725017931d6740bda71e5f07600ce6b407e543c4bc40bcd4c0b600e6c98003bf8548986b14d9098746dc89d154af91ad54f337b31c79222145dd3ed254fdeda0300c49ebcd2352765f533883a3513435f3ee452496f5166c25bf503bd6ec0a0679efda3b46ebf86211d458de749460d4a2a19abe6ea2accb451ab9a096b98465d044dc2a7f86c253a4ee57b6df108eff598a8dbc483bf8b74c2789939db85ffd7e0fd55b32bc26877f5be26fa7d750500ce2fab93c0cbe7336b126a5693d0c16484f37addccc7642590dbe98538990b88637e374d545d9b34b67448d0357e60280bbd8542f1f4e813caa8e8db57364b4e0cc14f35af381dd9b71ec9292b3a3f16e42362d2019e05f30" + +func millerLoopGenericBLS12[C: static Curve]( + f: var Fp12[C], + P: ECP_SWei_Aff[Fp[C]], + Q: ECP_SWei_Aff[Fp2[C]] + ) = + ## Generic Miller Loop for BLS12 curve + ## Computes f{u,Q}(P) with u the BLS curve parameter + # TODO: retrieve the curve parameter from the curve declaration + + # Boundary cases + # Loop start + # The litterature starts from both L-1 or L-2: + # L-1: + # - Scott2019, Pairing Implementation Revisited, Algorithm 1 + # - Aranha2010, Faster Explicit Formulas ..., Algorithm 1 + # L-2 + # - Beuchat2010, High-Speed Software Implementation ..., Algorithm 1 + # - Aranha2013, The Realm of The Pairings, Algorithm 1 + # - Costello, Thesis, Algorithm 2.1 + # - Costello2012, Pairings for Beginners, Algorithm 5.1 + # + # Even the guide to pairing based cryptography has both + # Chapter 3: L-1 (Algorithm 3.1) + # Chapter 11: L-2 (Algorithm 11.1) but it explains why L-2 (unrolling) + # Loop end + # - Some implementation, for example Beuchat2010 or the Guide to Pairing-Based Cryptography + # have extra line additions after the main loop, + # this is needed for BN curves. + # - With r the order of G1 / G2 / GT, + # we have [r]T = Inf + # Hence, [r-1]T = -T + # so either we use complete addition + # or we special case line addition of T and -T (it's a vertical line) + # or we ensure the loop is done for a number of iterations strictly less + # than the curve order which is the case for BLS12 curves + + static: doAssert C == BLS12_381, "Only BLS12-381 is supported at the moment" + + var + T {.noInit.}: ECP_SWei_Proj[Fp2[C]] + line {.noInit.}: Line[Fp2[C], C.getSexticTwist()] + nQ{.noInit.}: typeof(Q) + + T.projectiveFromAffine(Q) + nQ.neg(Q) + f.setOne() + + template u: untyped = BLS12_381_param + let u3 = 3*BLS12_381_param + for i in countdown(u3.bits - 2, 1): + f.square() + line.line_double(T, P) + f.mul_sparse_by_line_xy000z(line) + + let naf = u3.bit(i).int8 - u.bit(i).int8 # This can throw exception + if naf == 1: + line.line_add(T, Q, P) + f.mul_sparse_by_line_xy000z(line) + elif naf == -1: + line.line_add(T, nQ, P) + f.mul_sparse_by_line_xy000z(line) + + when BLS12_381_param_isNeg: # TODO generic + # In GT, x^-1 == conjugate(x) + # Remark 7.1, chapter 7.1.1 of Guide to Pairing-Based Cryptography, El Mrabet, 2017 + f.conj() + +func finalExpGeneric[C: static Curve](f: var Fp12[C]) = + ## A generic and slow implementation of final exponentiation + ## for sanity checks purposes. + static: doAssert C == BLS12_381, "Only BLS12-381 is supported at the moment" + f.powUnsafeExponent(BLS12_381_finalexponent, window = 3) + +func pairing_bls12_reference*[C](gt: var Fp12[C], P: ECP_SWei_Proj[Fp[C]], Q: ECP_SWei_Proj[Fp2[C]]) = + ## Compute the optimal Ate Pairing for BLS12 curves + ## Input: P ∈ G1, Q ∈ G2 + ## Output: e(P, Q) ∈ Gt + ## + ## Reference implementation + var Paff {.noInit.}: ECP_SWei_Aff[Fp[C]] + var Qaff {.noInit.}: ECP_SWei_Aff[Fp2[C]] + Paff.affineFromProjective(P) + Qaff.affineFromProjective(Q) + gt.millerLoopGenericBLS12(Paff, Qaff) + gt.finalExpGeneric() + +func finalExpEasy[C: static Curve](f: var Fp12[C]) = + ## Easy part of the final exponentiation + ## We need to clear the GT cofactor to obtain + ## an unique GT representation + ## (reminder, GT is a multiplicative group hence we exponentiate by the cofactor) + ## + ## With an embedding degree of 12, the easy part of final exponentiation is + ## + ## f^(p⁶−1)(p²+1) + discard diff --git a/constantine/pairing/pairing_bn.nim b/constantine/pairing/pairing_bn.nim new file mode 100644 index 0000000..5e928a2 --- /dev/null +++ b/constantine/pairing/pairing_bn.nim @@ -0,0 +1,178 @@ +# 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 + std/macros, + ../primitives, + ../config/[common, curves], + ../arithmetic, + ../towers, + ../io/io_bigints, + ../elliptic/[ + ec_weierstrass_affine, + ec_weierstrass_projective + ], + ./lines_projective, + ./gt_fp12, + ../isogeny/frobenius + +# ############################################################ +# +# Optimal ATE pairing for +# BN curves +# +# ############################################################ + +# - Efficient Final Exponentiation +# via Cyclotomic Structure for Pairings +# over Families of Elliptic Curves +# Daiki Hayashida and Kenichiro Hayasaka +# and Tadanori Teruya, 2020 +# https://eprint.iacr.org/2020/875.pdf +# +# - Faster Pairing Computations on Curves with High-Degree Twists +# Craig Costello, Tanja Lange, and Michael Naehrig, 2009 +# https://eprint.iacr.org/2009/615.pdf + +# TODO: implement quadruple-and-add and octuple-and-add +# from Costello2009 to trade multiplications in Fpᵏ +# for multiplications in Fp + +# TODO: should be part of curve parameters +const BN254_Snarks_ate_param = block: + # BN Miller loop is parametrized by 6u+2 + BigInt[67].fromHex"0x19d797039be763ba8" + +const BN254_Snarks_ate_param_isNeg = false + +const BN254_Snarks_finalexponent = block: + # (p^12 - 1) / r + BigInt[2790].fromHex"0x2f4b6dc97020fddadf107d20bc842d43bf6369b1ff6a1c71015f3f7be2e1e30a73bb94fec0daf15466b2383a5d3ec3d15ad524d8f70c54efee1bd8c3b21377e563a09a1b705887e72eceaddea3790364a61f676baaf977870e88d5c6c8fef0781361e443ae77f5b63a2a2264487f2940a8b1ddb3d15062cd0fb2015dfc6668449aed3cc48a82d0d602d268c7daab6a41294c0cc4ebe5664568dfc50e1648a45a4a1e3a5195846a3ed011a337a02088ec80e0ebae8755cfe107acf3aafb40494e406f804216bb10cf430b0f37856b42db8dc5514724ee93dfb10826f0dd4a0364b9580291d2cd65664814fde37ca80bb4ea44eacc5e641bbadf423f9a2cbf813b8d145da90029baee7ddadda71c7f3811c4105262945bba1668c3be69a3c230974d83561841d766f9c9d570bb7fbe04c7e8a6c3c760c0de81def35692da361102b6b9b2b918837fa97896e84abb40a4efb7e54523a486964b64ca86f120" + +const BN254_Nogami_ate_param = block: + # BN Miller loop is parametrized by 6u+2 + BigInt[67].fromHex"0x18300000000000004" # 65+2 bit for NAF x3 encoding + +const BN254_Nogami_ate_param_isNeg = true + +const BN254_Nogami_finalexponent = block: + # (p^12 - 1) / r + BigInt[2786].fromHex"0x2928fbb36b391596ee3fe4cbe857330da83e46fedf04d235a4a8daf5ff9f6eabcb4e3f20aa06f0a0d96b24f9af0cbbce750d61627dcbf5fec9139b8f1c46c86b49b4f8a202af26e4504f2c0f56570e9bd5b94c403f385d1908556486e24b396ddc2cdf13d06542f84fe8e82ccbad7b7423fc1ef4e8cc73d605e3e867c0a75f45ea7f6356d9846ce35d5a34f30396938818ad41914b97b99c289a7259b5d2e09477a77bd3c409b19f19e893f8ade90b0aed1b5fc8a07a3cebb41d4e9eee96b21a832ddb1e93e113edfb704fa532848c18593cd0ee90444a1b3499a800177ea38bdec62ec5191f2b6bbee449722f98d2173ad33077545c2ad10347e125a56fb40f086e9a4e62ad336a72c8b202ac3c1473d73b93d93dc0795ca0ca39226e7b4c1bb92f99248ec0806e0ad70744e9f2238736790f5185ea4c70808442a7d530c6ccd56b55a6973867ec6c73599bbd020bbe105da9c6b5c009ad8946cd6f0" + +{.experimental: "dynamicBindSym".} + +macro get(C: static Curve, value: untyped): untyped = + return bindSym($C & "_" & $value) + +func millerLoopGenericBN[C: static Curve]( + f: var Fp12[C], + P: ECP_SWei_Aff[Fp[C]], + Q: ECP_SWei_Aff[Fp2[C]] + ) = + ## Generic Miller Loop for BN curves + ## Computes f{6u+2,Q}(P) with u the BN curve parameter + # TODO: retrieve the curve parameter from the curve declaration + + # TODO - boundary cases + # Loop start + # The literatture starts from both L-1 or L-2: + # L-1: + # - Scott2019, Pairing Implementation Revisited, Algorithm 1 + # - Aranha2010, Faster Explicit Formulas ..., Algorithm 1 + # L-2 + # - Beuchat2010, High-Speed Software Implementation ..., Algorithm 1 + # - Aranha2013, The Realm of The Pairings, Algorithm 1 + # - Costello, Thesis, Algorithm 2.1 + # - Costello2012, Pairings for Beginners, Algorithm 5.1 + # + # Even the guide to pairing based cryptography has both + # Chapter 3: L-1 (Algorithm 3.1) + # Chapter 11: L-2 (Algorithm 11.1) but it explains why L-2 (unrolling) + # Loop end + # - Some implementation, for example Beuchat2010 or the Guide to Pairing-Based Cryptography + # have an extra line addition after the main loop, this seems related to + # the NAF recoding and not Miller Loop + # - With r the order of G1 / G2 / GT, + # we have [r]T = Inf + # Hence, [r-1]T = -T + # so either we use complete addition + # or we special case line addition of T and -T (it's a vertical line) + # or we ensure the loop is done for a number of iterations strictly less + # than the curve order which is the case for BN curves + + var + T {.noInit.}: ECP_SWei_Proj[Fp2[C]] + line {.noInit.}: Line[Fp2[C], C.getSexticTwist()] + nQ{.noInit.}: typeof(Q) + + T.projectiveFromAffine(Q) + nQ.neg(Q) + f.setOne() + + template u: untyped = C.get(ate_param) + let u3 = 3*C.get(ate_param) + for i in countdown(u3.bits - 2, 1): + f.square() + line.line_double(T, P) + f.mul_sparse_by_line_xyz000(line) + + let naf = u3.bit(i).int8 - u.bit(i).int8 # This can throw exception + if naf == 1: + line.line_add(T, Q, P) + f.mul_sparse_by_line_xyz000(line) + elif naf == -1: + line.line_add(T, nQ, P) + f.mul_sparse_by_line_xyz000(line) + + when C.get(ate_param_isNeg): # TODO generic + # In GT, x^-1 == conjugate(x) + # Remark 7.1, chapter 7.1.1 of Guide to Pairing-Based Cryptography, El Mrabet, 2017 + f.conj() + + # Ate pairing for BN curves need adjustment after Miller loop + when C.get(ate_param_isNeg): + T.neg() + var V {.noInit.}: typeof(Q) + + V.frobenius_psi(Q) + line.line_add(T, V, P) + f.mul_sparse_by_line_xyz000(line) + + V.frobenius_psi2(Q) + V.neg() + line.line_add(T, V, P) + f.mul_sparse_by_line_xyz000(line) + +func finalExpGeneric[C: static Curve](f: var Fp12[C]) = + ## A generic and slow implementation of final exponentiation + ## for sanity checks purposes. + f.powUnsafeExponent(C.get(finalexponent), window = 3) + +func pairing_bn_reference*[C](gt: var Fp12[C], P: ECP_SWei_Proj[Fp[C]], Q: ECP_SWei_Proj[Fp2[C]]) = + ## Compute the optimal Ate Pairing for BN curves + ## Input: P ∈ G1, Q ∈ G2 + ## Output: e(P, Q) ∈ Gt + ## + ## Reference implementation + var Paff {.noInit.}: ECP_SWei_Aff[Fp[C]] + var Qaff {.noInit.}: ECP_SWei_Aff[Fp2[C]] + Paff.affineFromProjective(P) + Qaff.affineFromProjective(Q) + gt.millerLoopGenericBN(Paff, Qaff) + gt.finalExpGeneric() + +func finalExpEasy[C: static Curve](f: var Fp12[C]) = + ## Easy part of the final exponentiation + ## We need to clear the GT cofactor to obtain + ## an unique GT representation + ## (reminder, GT is a multiplicative group hence we exponentiate by the cofactor) + ## + ## With an embedding degree of 12, the easy part of final exponentiation is + ## + ## f^(p⁶−1)(p²+1) + discard diff --git a/constantine/tower_field_extensions/cubic_extensions.nim b/constantine/tower_field_extensions/cubic_extensions.nim index 0681420..c61bda1 100644 --- a/constantine/tower_field_extensions/cubic_extensions.nim +++ b/constantine/tower_field_extensions/cubic_extensions.nim @@ -190,9 +190,31 @@ func inv*(r: var CubicExt, a: CubicExt) = r.c1.prod(v1, v3) r.c2.prod(v2, v3) +func inv*(a: var CubicExt) = + ## Compute the multiplicative inverse of ``a`` + ## + ## The inverse of 0 is 0. + ## Incidentally this avoids extra check + ## to convert Jacobian and Projective coordinates + ## to affine for elliptic curve + let t = a + a.inv(t) + func `*=`*(a: var CubicExt, b: CubicExt) {.inline.} = ## In-place multiplication # On higher extension field like 𝔽p6, # if `prod` is called on shared in and out buffer, the result is wrong let t = a a.prod(t, b) + +func conj*(a: var CubicExt) {.inline.} = + ## Computes the conjugate in-place + mixin conj, conjneg + a.c0.conj() + a.c1.conjneg() + a.c2.conj() + +func square*(a: var CubicExt) {.inline.} = + ## In-place squaring + let t = a + a.square(t) diff --git a/constantine/tower_field_extensions/exponentiations.nim b/constantine/tower_field_extensions/exponentiations.nim index 955af3d..ec48454 100644 --- a/constantine/tower_field_extensions/exponentiations.nim +++ b/constantine/tower_field_extensions/exponentiations.nim @@ -11,6 +11,7 @@ import ../arithmetic, ../config/[common, curves], ../primitives, + ../io/io_bigints, ./tower_common, ./quadratic_extensions, ./cubic_extensions @@ -60,7 +61,7 @@ func powPrologue[F](a: var F, scratchspace: var openarray[F]): uint = else: scratchspace[2] = a for k in 2 ..< 1 shl result: - scratchspace[k+1] + scratchspace[k+1].prod(scratchspace[k], a) a.setOne() func powSquarings[F]( @@ -69,7 +70,7 @@ func powSquarings[F]( tmp: var F, window: uint, acc, acc_len: var uint, - e: var uint + e: var int ): tuple[k, bits: uint] {.inline.}= ## Squaring step of exponentiation by squaring ## Get the next k bits in range [1, window) @@ -109,10 +110,10 @@ func powSquarings[F]( return (k, bits) -func powUnsafeExponent( - a: var ExtensionField, +func powUnsafeExponent[F]( + a: var F, exponent: openArray[byte], - scratchspace: var openArray[byte] + scratchspace: var openArray[F] ) = ## Extension field exponentiation r = a^exponent (mod p^m) ## @@ -145,6 +146,55 @@ func powUnsafeExponent( scratchspace[0].prod(a, scratchspace[1]) a = scratchspace[0] +func powUnsafeExponent*[F]( + a: var F, + exponent: openArray[byte], + window: static int + ) = + ## Extension field exponentiation r = a^exponent (mod p^m) + ## exponent is an big integer in canonical octet-string format + ## + ## Window is used for window optimization. + ## 2^window field elements are allocated for scratchspace. + ## + ## - On Fp2, with a 256-bit base field, a window of size 5 requires + ## 2*256*2^5 = 16KiB + ## + ## Warning ⚠️ : + ## This is an optimization for public exponent + ## Otherwise bits of the exponent can be retrieved with: + ## - memory access analysis + ## - power analysis + ## - timing analysis + const scratchLen = if window == 1: 2 + else: (1 shl window) + 1 + var scratchSpace {.noInit.}: array[scratchLen, typeof(a)] + a.powUnsafeExponent(exponent, scratchspace) + +func powUnsafeExponent*[F; bits: static int]( + a: var F, + exponent: BigInt[bits], + window: static int + ) = + ## Extension field exponentiation r = a^exponent (mod p^m) + ## exponent is an big integer in canonical octet-string format + ## + ## Window is used for window optimization. + ## 2^window field elements are allocated for scratchspace. + ## + ## - On Fp2, with a 256-bit base field, a window of size 5 requires + ## 2*256*2^5 = 16KiB + ## + ## Warning ⚠️ : + ## This is an optimization for public exponent + ## Otherwise bits of the exponent can be retrieved with: + ## - memory access analysis + ## - power analysis + ## - timing analysis + var expBE {.noInit.}: array[(bits + 7) div 8, byte] + expBE.exportRawUint(exponent, bigEndian) + a.powUnsafeExponent(expBE, window) + # Square root # ----------------------------------------------------------- # diff --git a/constantine/tower_field_extensions/quadratic_extensions.nim b/constantine/tower_field_extensions/quadratic_extensions.nim index 6a5e65f..d004ba4 100644 --- a/constantine/tower_field_extensions/quadratic_extensions.nim +++ b/constantine/tower_field_extensions/quadratic_extensions.nim @@ -246,6 +246,19 @@ func mul_sparse_generic_by_x0(r: var QuadraticExt, a, sparseB: QuadraticExt) = # Exported symbols # ------------------------------------------------------------------- +func conj*(a: var QuadraticExt) {.inline.} = + ## Computes the conjugate in-place + a.c1.neg() + +func conj*(r: var QuadraticExt, a: QuadraticExt) {.inline.} = + ## Computes the conjugate out-of-place + r.c0 = a.c0 + r.c1.neg(a.c1) + +func conjneg*(a: var QuadraticExt) {.inline.} = + ## Computes the negated conjugate in-place + a.c0.neg() + func square*(r: var QuadraticExt, a: QuadraticExt) {.inline.} = mixin fromComplexExtension when r.fromComplexExtension(): @@ -294,6 +307,15 @@ func inv*(r: var QuadraticExt, a: QuadraticExt) = v0.neg(v1) # v0 = -1 / (a0² - β a1²) r.c1.prod(a.c1, v0) # r1 = -a1 / (a0² - β a1²) +func inv*(a: var QuadraticExt) = + ## Compute the multiplicative inverse of ``a`` + ## + ## The inverse of 0 is 0. + ## Incidentally this avoids extra check + ## to convert Jacobian and Projective coordinates + ## to affine for elliptic curve + a.inv(a) + func `*=`*(a: var QuadraticExt, b: QuadraticExt) {.inline.} = ## In-place multiplication # On higher extension field like 𝔽p12, diff --git a/constantine/tower_field_extensions/tower_common.nim b/constantine/tower_field_extensions/tower_common.nim index c05f01d..50fed56 100644 --- a/constantine/tower_field_extensions/tower_common.nim +++ b/constantine/tower_field_extensions/tower_common.nim @@ -107,6 +107,11 @@ func neg*(r: var ExtensionField, a: ExtensionField) = for fR, fA in fields(r, a): fR.neg(fA) +func neg*(a: var ExtensionField) = + ## Field out-of-place negation + for fA in fields(a): + fA.neg() + func `+=`*(a: var ExtensionField, b: ExtensionField) = ## Addition in the extension field for fA, fB in fields(a, b): diff --git a/constantine/towers.nim b/constantine/towers.nim index 3e14ec8..5d77479 100644 --- a/constantine/towers.nim +++ b/constantine/towers.nim @@ -121,18 +121,22 @@ func `/=`*(a: var Fp2, _: typedesc[SexticNonResidue]) {.inline.} = # ---------------------------------------------------------------- type + Fp4*[C: static Curve] = object + c0*, c1*: Fp2[C] + Fp6*[C: static Curve] = object c0*, c1*, c2*: Fp2[C] - ξ = NonResidue - # We call the non-residue ξ on 𝔽p6 to avoid confusion between non-residue + ξ* = NonResidue + # We call the non-residue ξ on 𝔽p4/𝔽p6 to avoid confusion + # between non-residue # of different tower level func `*`*(_: typedesc[ξ], a: Fp2): Fp2 {.inline, noInit.} = - ## Multiply an element of 𝔽p2 by the cubic non-residue - ## chosen to construct 𝔽p6 + ## Multiply an element of 𝔽p2 by the quadratic and cubic non-residue + ## chosen to construct 𝔽p4/𝔽p6 # Yet another const tuple unpacking bug - const u = Fp2.C.get_CNR_Fp2()[0] # Cubic non-residue to construct 𝔽p6 + const u = Fp2.C.get_CNR_Fp2()[0] # Quadratic & Cubic non-residue to construct 𝔽p4/𝔽p6 const v = Fp2.C.get_CNR_Fp2()[1] const Beta = Fp2.C.get_QNR_Fp() # Quadratic non-residue to construct 𝔽p2 # ξ = u + v x @@ -176,21 +180,20 @@ func `*=`*(a: var Fp2, _: typedesc[ξ]) {.inline.} = type Fp12*[C: static Curve] = object - c0*, c1*: Fp6[C] + c0*, c1*, c2*: Fp4[C] γ = NonResidue # We call the non-residue γ (Gamma) on 𝔽p6 to avoid confusion between non-residue # of different tower level -func `*`*(_: typedesc[γ], a: Fp6): Fp6 {.noInit, inline.} = +func `*`*(_: typedesc[γ], a: Fp4): Fp4 {.noInit, inline.} = ## Multiply an element of 𝔽p6 by the cubic non-residue ## chosen to construct 𝔽p12 ## For all curves γ = v with v the factor for 𝔽p6 coordinate ## and v³ = ξ ## (c0 + c1 v + c2 v²) v => ξ c2 + c0 v + c1 v² - result.c0 = ξ * a.c2 + result.c0 = ξ * a.c1 result.c1 = a.c0 - result.c2 = a.c1 -func `*=`*(a: var Fp6, _: typedesc[γ]) {.inline.} = +func `*=`*(a: var Fp4, _: typedesc[γ]) {.inline.} = a = γ * a diff --git a/helpers/prng_unsafe.nim b/helpers/prng_unsafe.nim index 6d98a19..b2756b5 100644 --- a/helpers/prng_unsafe.nim +++ b/helpers/prng_unsafe.nim @@ -226,7 +226,7 @@ func random_long01Seq[T](rng: var RngState, a: var T, C: static Curve) = # Elliptic curves # ------------------------------------------------------------ -func random_unsafe[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = +func random_unsafe[F](rng: var RngState, a: var (ECP_SWei_Proj[F] or ECP_SWei_Aff[F])) = ## Initialize a random curve point with Z coordinate == 1 ## Unsafe: for testing and benchmarking purposes only var fieldElem {.noInit.}: F @@ -251,7 +251,7 @@ func random_unsafe_with_randZ[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = rng.random_unsafe(fieldElem, F.C) success = trySetFromCoordsXandZ(a, fieldElem, Z) -func random_highHammingWeight[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = +func random_highHammingWeight[F](rng: var RngState, a: var (ECP_SWei_Proj[F] or ECP_SWei_Aff[F])) = ## Initialize a random curve point with Z coordinate == 1 ## This will be generated with a biaised RNG with high Hamming Weight ## to trigger carry bugs @@ -264,7 +264,7 @@ func random_highHammingWeight[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = rng.random_highHammingWeight(fieldElem, F.C) success = trySetFromCoordX(a, fieldElem) -func random_highHammingWeight_with_randZ[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = +func random_highHammingWeight_with_randZ[F](rng: var RngState, a: var (ECP_SWei_Proj[F] or ECP_SWei_Aff[F])) = ## Initialize a random curve point with Z coordinate == 1 ## This will be generated with a biaised RNG with high Hamming Weight ## to trigger carry bugs @@ -278,7 +278,7 @@ func random_highHammingWeight_with_randZ[F](rng: var RngState, a: var ECP_SWei_P rng.random_highHammingWeight(fieldElem, F.C) success = trySetFromCoordsXandZ(a, fieldElem, Z) -func random_long01Seq[F](rng: var RngState, a: var ECP_SWei_Proj[F]) = +func random_long01Seq[F](rng: var RngState, a: var (ECP_SWei_Proj[F] or ECP_SWei_Aff[F])) = ## Initialize a random curve point with Z coordinate == 1 ## This will be generated with a biaised RNG ## that produces long bitstrings of 0 and 1 @@ -313,7 +313,7 @@ func random_long01Seq_with_randZ[F](rng: var RngState, a: var ECP_SWei_Proj[F]) func random_unsafe*(rng: var RngState, T: typedesc): T = ## Create a random Field or Extension Field or Curve Element ## Unsafe: for testing and benchmarking purposes only - when T is ECP_SWei_Proj: + when T is (ECP_SWei_Proj or ECP_SWei_Aff): rng.random_unsafe(result) elif T is SomeNumber: cast[T](rng.next()) # TODO: Rely on casting integer actually converting in C (i.e. uint64->uint32 is valid) @@ -330,7 +330,7 @@ func random_unsafe_with_randZ*(rng: var RngState, T: typedesc[ECP_SWei_Proj]): T func random_highHammingWeight*(rng: var RngState, T: typedesc): T = ## Create a random Field or Extension Field or Curve Element ## Skewed towards high Hamming Weight - when T is ECP_SWei_Proj: + when T is (ECP_SWei_Proj or ECP_SWei_Aff): rng.random_highHammingWeight(result) elif T is SomeNumber: cast[T](rng.next()) # TODO: Rely on casting integer actually converting in C (i.e. uint64->uint32 is valid) @@ -347,7 +347,7 @@ func random_highHammingWeight_with_randZ*(rng: var RngState, T: typedesc[ECP_SWe func random_long01Seq*(rng: var RngState, T: typedesc): T = ## Create a random Field or Extension Field or Curve Element ## Skewed towards long bitstrings of 0 or 1 - when T is ECP_SWei_Proj: + when T is (ECP_SWei_Proj or ECP_SWei_Aff): rng.random_long01Seq(result) elif T is SomeNumber: cast[T](rng.next()) # TODO: Rely on casting integer actually converting in C (i.e. uint64->uint32 is valid) diff --git a/helpers/static_for.nim b/helpers/static_for.nim index ffedcad..f78305b 100644 --- a/helpers/static_for.nim +++ b/helpers/static_for.nim @@ -34,6 +34,14 @@ macro staticFor*(ident: untyped{nkIdent}, choices: typed, body: untyped): untype ## staticFor(curve, TestCurves): ## body ## and unroll the body for each curve in TestCurves + + let choices = if choices.kind == nnkSym: + # Unpack symbol + choices.getImpl() + else: + choices.expectKind(nnkBracket) + choices + result = newStmtList() for choice in choices: result.add nnkBlockStmt.newTree( diff --git a/sage/frobenius_bls12_381.sage b/sage/frobenius_bls12_381.sage index 178fac7..3da1701 100644 --- a/sage/frobenius_bls12_381.sage +++ b/sage/frobenius_bls12_381.sage @@ -50,7 +50,15 @@ def fp2_to_hex(a): v = vector(a) return Integer(v[0]).hex() + ' + β * ' + Integer(v[1]).hex() -# Frobenius constants (D type: use SNR, M type use 1/SNR) +# Frobenius map constants (D type: use SNR, M type use 1/SNR) +print('\nFrobenius extension field constants') +FrobConst_map = SNR^((p-1)/6) +print('FrobConst_map : ' + fp2_to_hex(FrobConst_map)) +FrobConst_map_fp4 = FrobConst_map^(12//4) +print('FrobConst_map_fp4_y : ' + fp2_to_hex(FrobConst_map_fp4)) + +# Frobenius psi constants (D type: use SNR, M type use 1/SNR) +print('\nψ (Psi) - Untwist-Frobenius-Twist constants') print('1/sextic_non_residue: ' + fp2_to_hex(1/SNR)) FrobConst_psi = (1/SNR)^((p-1)/6) FrobConst_psi_2 = FrobConst_psi * FrobConst_psi diff --git a/sage/frobenius_bn254_nogami.sage b/sage/frobenius_bn254_nogami.sage index 2f138e1..97c0132 100644 --- a/sage/frobenius_bn254_nogami.sage +++ b/sage/frobenius_bn254_nogami.sage @@ -8,7 +8,7 @@ # ############################################################ # -# BLS12-381 +# BN254-Nogami # Frobenius Endomorphism # Untwist-Frobenius-Twist isogeny # diff --git a/sage/frobenius_bn254_snarks.sage b/sage/frobenius_bn254_snarks.sage index a63d69d..3eaa34d 100644 --- a/sage/frobenius_bn254_snarks.sage +++ b/sage/frobenius_bn254_snarks.sage @@ -8,7 +8,7 @@ # ############################################################ # -# BLS12-381 +# BN254-Snarks # Frobenius Endomorphism # Untwist-Frobenius-Twist isogeny # diff --git a/tests/t_bigints_mod_vs_gmp.nim b/tests/t_bigints_mod_vs_gmp.nim index 2f96cd6..6a9c654 100644 --- a/tests/t_bigints_mod_vs_gmp.nim +++ b/tests/t_bigints_mod_vs_gmp.nim @@ -92,7 +92,7 @@ proc main() = mpz_init(m) mpz_init(r) - testRandomModSizes(24, aBits, mBits): + testRandomModSizes(12, aBits, mBits): # echo "--------------------------------------------------------------------------------" echo "Testing: random dividend (" & align($aBits, 4) & "-bit) -- random modulus (" & align($mBits, 4) & "-bit)" diff --git a/tests/t_bigints_mul_high_words_vs_gmp.nim b/tests/t_bigints_mul_high_words_vs_gmp.nim index c2daafd..b9a7418 100644 --- a/tests/t_bigints_mul_high_words_vs_gmp.nim +++ b/tests/t_bigints_mul_high_words_vs_gmp.nim @@ -67,7 +67,7 @@ proc main() = mpz_init(a) mpz_init(b) - testRandomModSizes(24, rBits, aBits, bBits, wordsStartIndex): + testRandomModSizes(12, rBits, aBits, bBits, wordsStartIndex): # echo "--------------------------------------------------------------------------------" echo "Testing: random mul_high_words r (", align($rBits, 4), "-bit, keeping from ", wordsStartIndex, diff --git a/tests/t_bigints_mul_vs_gmp.nim b/tests/t_bigints_mul_vs_gmp.nim index db960cb..cb4d3e2 100644 --- a/tests/t_bigints_mul_vs_gmp.nim +++ b/tests/t_bigints_mul_vs_gmp.nim @@ -64,7 +64,7 @@ proc main() = mpz_init(a) mpz_init(b) - testRandomModSizes(24, rBits, aBits, bBits): + testRandomModSizes(12, rBits, aBits, bBits): # echo "--------------------------------------------------------------------------------" echo "Testing: random mul r (", align($rBits, 4), "-bit) <- a (", align($aBits, 4), "-bit) * b (", align($bBits, 4), "-bit) (full mul bits: ", align($(aBits+bBits), 4), "), r large enough? ", rBits >= aBits+bBits diff --git a/tests/t_ec_template.nim b/tests/t_ec_template.nim index 33d4883..0840475 100644 --- a/tests/t_ec_template.nim +++ b/tests/t_ec_template.nim @@ -31,21 +31,21 @@ type HighHammingWeight Long01Sequence -func random_point*(rng: var RngState, F: typedesc, randZ: static bool, gen: static RandomGen): F {.inline, noInit.} = - when not randZ: - when gen == Uniform: - result = rng.random_unsafe(F) +func random_point*(rng: var RngState, EC: typedesc, randZ: bool, gen: RandomGen): EC {.noInit.} = + if not randZ: + if gen == Uniform: + result = rng.random_unsafe(EC) elif gen == HighHammingWeight: - result = rng.random_highHammingWeight(F) + result = rng.random_highHammingWeight(EC) else: - result = rng.random_long01Seq(F) + result = rng.random_long01Seq(EC) else: - when gen == Uniform: - result = rng.random_unsafe_with_randZ(F) + if gen == Uniform: + result = rng.random_unsafe_with_randZ(EC) elif gen == HighHammingWeight: - result = rng.random_highHammingWeight_with_randZ(F) + result = rng.random_highHammingWeight_with_randZ(EC) else: - result = rng.random_long01Seq_with_randZ(F) + result = rng.random_long01Seq_with_randZ(EC) proc run_EC_addition_tests*( ec: typedesc, @@ -69,7 +69,7 @@ proc run_EC_addition_tests*( suite testSuiteDesc & " - " & $ec & " - [" & $WordBitwidth & "-bit mode]": test "The infinity point is the neutral element w.r.t. to EC " & G1_or_G2 & " addition": - proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, randZ: bool, gen: RandomGen) = var inf {.noInit.}: EC inf.setInf() check: bool inf.isInf() @@ -92,7 +92,7 @@ proc run_EC_addition_tests*( test(ec, randZ = true, gen = Long01Sequence) test "Adding opposites gives an infinity point": - proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, randZ: bool, gen: RandomGen) = for _ in 0 ..< Iters: var r{.noInit.}: EC let P = rng.random_point(EC, randZ, gen) @@ -113,7 +113,7 @@ proc run_EC_addition_tests*( test(ec, randZ = true, gen = Long01Sequence) test "EC " & G1_or_G2 & " add is commutative": - proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, randZ: bool, gen: RandomGen) = for _ in 0 ..< Iters: var r0{.noInit.}, r1{.noInit.}: EC let P = rng.random_point(EC, randZ, gen) @@ -131,7 +131,7 @@ proc run_EC_addition_tests*( test(ec, randZ = true, gen = Long01Sequence) test "EC " & G1_or_G2 & " add is associative": - proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, randZ: bool, gen: RandomGen) = for _ in 0 ..< Iters: let a = rng.random_point(EC, randZ, gen) let b = rng.random_point(EC, randZ, gen) @@ -180,7 +180,7 @@ proc run_EC_addition_tests*( test(ec, randZ = true, gen = Long01Sequence) test "EC " & G1_or_G2 & " double and EC " & G1_or_G2 & " add are consistent": - proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, randZ: bool, gen: RandomGen) = for _ in 0 ..< Iters: let a = rng.random_point(EC, randZ, gen) @@ -220,7 +220,7 @@ proc run_EC_mul_sanity_tests*( suite testSuiteDesc & " - " & $ec & " - [" & $WordBitwidth & "-bit mode]": test "EC " & G1_or_G2 & " mul [0]P == Inf": - proc test(EC: typedesc, bits: static int, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, bits: static int, randZ: bool, gen: RandomGen) = for _ in 0 ..< ItersMul: let a = rng.random_point(EC, randZ, gen) @@ -243,7 +243,7 @@ proc run_EC_mul_sanity_tests*( test(ec, bits = ec.F.C.getCurveOrderBitwidth(), randZ = true, gen = Long01Sequence) test "EC " & G1_or_G2 & " mul [1]P == P": - proc test(EC: typedesc, bits: static int, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, bits: static int, randZ: bool, gen: RandomGen) = for _ in 0 ..< ItersMul: let a = rng.random_point(EC, randZ, gen) @@ -269,7 +269,7 @@ proc run_EC_mul_sanity_tests*( test(ec, bits = ec.F.C.getCurveOrderBitwidth(), randZ = true, gen = Long01Sequence) test "EC " & G1_or_G2 & " mul [2]P == P.double()": - proc test(EC: typedesc, bits: static int, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, bits: static int, randZ: bool, gen: RandomGen) = for _ in 0 ..< ItersMul: let a = rng.random_point(EC, randZ, gen) @@ -319,7 +319,7 @@ proc run_EC_mul_distributive_tests*( suite testSuiteDesc & " - " & $ec & " - [" & $WordBitwidth & "-bit mode]": test "EC " & G1_or_G2 & " mul is distributive over EC add": - proc test(EC: typedesc, bits: static int, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, bits: static int, randZ: bool, gen: RandomGen) = for _ in 0 ..< ItersMul: let a = rng.random_point(EC, randZ, gen) let b = rng.random_point(EC, randZ, gen) @@ -388,7 +388,7 @@ proc run_EC_mul_vs_ref_impl*( suite testSuiteDesc & " - " & $ec & " - [" & $WordBitwidth & "-bit mode]": test "EC " & G1_or_G2 & " mul constant-time is equivalent to a simple double-and-add algorithm": - proc test(EC: typedesc, bits: static int, randZ: static bool, gen: static RandomGen) = + proc test(EC: typedesc, bits: static int, randZ: bool, gen: RandomGen) = for _ in 0 ..< ItersMul: let a = rng.random_point(EC, randZ, gen) diff --git a/tests/t_ec_wstrass_prj_g1_add_double.nim b/tests/t_ec_wstrass_prj_g1_add_double.nim index eec6af8..c8d7826 100644 --- a/tests/t_ec_wstrass_prj_g1_add_double.nim +++ b/tests/t_ec_wstrass_prj_g1_add_double.nim @@ -19,7 +19,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 8 run_EC_addition_tests( ec = ECP_SWei_Proj[Fp[BN254_Snarks]], diff --git a/tests/t_ec_wstrass_prj_g1_mul_distri.nim b/tests/t_ec_wstrass_prj_g1_mul_distri.nim index 16e1bee..6e669e7 100644 --- a/tests/t_ec_wstrass_prj_g1_mul_distri.nim +++ b/tests/t_ec_wstrass_prj_g1_mul_distri.nim @@ -20,7 +20,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_distributive_tests( diff --git a/tests/t_ec_wstrass_prj_g1_mul_sanity.nim b/tests/t_ec_wstrass_prj_g1_mul_sanity.nim index b5eaff3..007d711 100644 --- a/tests/t_ec_wstrass_prj_g1_mul_sanity.nim +++ b/tests/t_ec_wstrass_prj_g1_mul_sanity.nim @@ -20,7 +20,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_sanity_tests( diff --git a/tests/t_ec_wstrass_prj_g1_mul_vs_ref.nim b/tests/t_ec_wstrass_prj_g1_mul_vs_ref.nim index 49661d1..2d16ba2 100644 --- a/tests/t_ec_wstrass_prj_g1_mul_vs_ref.nim +++ b/tests/t_ec_wstrass_prj_g1_mul_vs_ref.nim @@ -20,7 +20,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_vs_ref_impl( diff --git a/tests/t_ec_wstrass_prj_g2_add_double_bls12_381.nim b/tests/t_ec_wstrass_prj_g2_add_double_bls12_381.nim index 7d33697..bb2f130 100644 --- a/tests/t_ec_wstrass_prj_g2_add_double_bls12_381.nim +++ b/tests/t_ec_wstrass_prj_g2_add_double_bls12_381.nim @@ -20,7 +20,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 8 run_EC_addition_tests( ec = ECP_SWei_Proj[Fp2[BLS12_381]], diff --git a/tests/t_ec_wstrass_prj_g2_add_double_bn254_snarks.nim b/tests/t_ec_wstrass_prj_g2_add_double_bn254_snarks.nim index 4c024f7..7d26404 100644 --- a/tests/t_ec_wstrass_prj_g2_add_double_bn254_snarks.nim +++ b/tests/t_ec_wstrass_prj_g2_add_double_bn254_snarks.nim @@ -20,7 +20,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 8 run_EC_addition_tests( ec = ECP_SWei_Proj[Fp2[BN254_Snarks]], diff --git a/tests/t_ec_wstrass_prj_g2_mul_distri_bls12_381.nim b/tests/t_ec_wstrass_prj_g2_mul_distri_bls12_381.nim index 4c76ea7..ba10ca2 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_distri_bls12_381.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_distri_bls12_381.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_distributive_tests( diff --git a/tests/t_ec_wstrass_prj_g2_mul_distri_bn254_snarks.nim b/tests/t_ec_wstrass_prj_g2_mul_distri_bn254_snarks.nim index 2199b1c..37cc43e 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_distri_bn254_snarks.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_distri_bn254_snarks.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_distributive_tests( diff --git a/tests/t_ec_wstrass_prj_g2_mul_sanity_bls12_381.nim b/tests/t_ec_wstrass_prj_g2_mul_sanity_bls12_381.nim index 8888b2f..a15c6e2 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_sanity_bls12_381.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_sanity_bls12_381.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_sanity_tests( diff --git a/tests/t_ec_wstrass_prj_g2_mul_sanity_bn254_snarks.nim b/tests/t_ec_wstrass_prj_g2_mul_sanity_bn254_snarks.nim index a0fe1f8..a713fad 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_sanity_bn254_snarks.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_sanity_bn254_snarks.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_sanity_tests( diff --git a/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bls12_381.nim b/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bls12_381.nim index ccd5fd3..64758fb 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bls12_381.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bls12_381.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_vs_ref_impl( diff --git a/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bn254_snarks.nim b/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bn254_snarks.nim index 3bbe774..acb49a5 100644 --- a/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bn254_snarks.nim +++ b/tests/t_ec_wstrass_prj_g2_mul_vs_ref_bn254_snarks.nim @@ -21,7 +21,7 @@ import ./t_ec_template const - Iters = 24 + Iters = 12 ItersMul = Iters div 4 run_EC_mul_vs_ref_impl( diff --git a/tests/t_finite_fields_powinv.nim b/tests/t_finite_fields_powinv.nim index 83a9744..2608b91 100644 --- a/tests/t_finite_fields_powinv.nim +++ b/tests/t_finite_fields_powinv.nim @@ -20,7 +20,7 @@ import static: doAssert defined(testingCurves), "This modules requires the -d:testingCurves compile option" -const Iters = 24 +const Iters = 8 var rng: RngState let seed = uint32(getTime().toUnix() and (1'i64 shl 32 - 1)) # unixTime mod 2^32 diff --git a/tests/t_finite_fields_sqrt.nim b/tests/t_finite_fields_sqrt.nim index 10b3271..2f42f82 100644 --- a/tests/t_finite_fields_sqrt.nim +++ b/tests/t_finite_fields_sqrt.nim @@ -17,7 +17,7 @@ import ../helpers/prng_unsafe -const Iters = 24 +const Iters = 8 var rng: RngState let seed = uint32(getTime().toUnix() and (1'i64 shl 32 - 1)) # unixTime mod 2^32 diff --git a/tests/t_fp12_bls12_377.nim b/tests/t_fp12_bls12_377.nim index ae3c47f..c300d55 100644 --- a/tests/t_fp12_bls12_377.nim +++ b/tests/t_fp12_bls12_377.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 12, - Iters = 24, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp12_" & $BLS12_377, testSuiteDesc = "𝔽p12 = 𝔽p6[w] " & $BLS12_377 diff --git a/tests/t_fp12_bls12_381.nim b/tests/t_fp12_bls12_381.nim index 9f48322..eeb0a28 100644 --- a/tests/t_fp12_bls12_381.nim +++ b/tests/t_fp12_bls12_381.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 12, - Iters = 24, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp12_" & $BLS12_381, testSuiteDesc = "𝔽p12 = 𝔽p6[w] " & $BLS12_381 diff --git a/tests/t_fp12_bn254_snarks.nim b/tests/t_fp12_bn254_snarks.nim index 69957e2..096a53b 100644 --- a/tests/t_fp12_bn254_snarks.nim +++ b/tests/t_fp12_bn254_snarks.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 12, - Iters = 24, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp12_" & $BN254_Snarks, testSuiteDesc = "𝔽p12 = 𝔽p6[w] " & $BN254_Snarks diff --git a/tests/t_fp12_exponentiation.nim b/tests/t_fp12_exponentiation.nim new file mode 100644 index 0000000..3953287 --- /dev/null +++ b/tests/t_fp12_exponentiation.nim @@ -0,0 +1,219 @@ +# 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/[tables, unittest, times], + # Internals + ../constantine/config/common, + ../constantine/[arithmetic, primitives], + ../constantine/towers, + ../constantine/config/curves, + ../constantine/io/io_towers, + # Test utilities + ../helpers/[prng_unsafe, static_for] + +const + Iters = 2 + TestCurves = [ + BN254_Snarks, + BLS12_381 + ] + +type + RandomGen = enum + Uniform + HighHammingWeight + Long01Sequence + +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 "test_fp12_exponentiation xoshiro512** seed: ", seed + +func random_elem(rng: var RngState, F: typedesc, gen: RandomGen): F {.inline, noInit.} = + if gen == Uniform: + result = rng.random_unsafe(F) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(F) + else: + result = rng.random_long01Seq(F) + +proc test_sameBaseProduct(C: static Curve, gen: RandomGen) = + ## xᴬ xᴮ = xᴬ⁺ᴮ - product of power + let x = rng.random_elem(Fp12[C], gen) + + var a = rng.random_elem(BigInt[128], gen) + var b = rng.random_elem(BigInt[128], gen) + # div by 2 to ensure their sum doesn't overflow + # 128 bits + a.div2() + b.div2() + + var xa = x + xa.powUnsafeExponent(a, window = 3) + + var xb = x + xb.powUnsafeExponent(b, window = 3) + + var xapb = x + var apb: BigInt[128] + discard apb.sum(a, b) + xapb.powUnsafeExponent(apb, window = 3) + + xa *= xb + check: bool(xa == xapb) + +proc test_powpow(C: static Curve, gen: RandomGen) = + ## (xᴬ)ᴮ = xᴬᴮ - power of power + var x = rng.random_elem(Fp12[C], gen) + + var a = rng.random_elem(BigInt[128], gen) + var b = rng.random_elem(BigInt[128], gen) + + var ab: BigInt[256] + ab.prod(a, b) + + var y = x + + x.powUnsafeExponent(a, window = 3) + x.powUnsafeExponent(b, window = 3) + + y.powUnsafeExponent(ab, window = 3) + check: bool(x == y) + +proc test_powprod(C: static Curve, gen: RandomGen) = + ## (xy)ᴬ = xᴬyᴬ - power of product + var x = rng.random_elem(Fp12[C], gen) + var y = rng.random_elem(Fp12[C], gen) + + let a = rng.random_elem(BigInt[128], gen) + + var xy{.noInit.}: Fp12[C] + xy.prod(x, y) + + xy.powUnsafeExponent(a, window=3) + + x.powUnsafeExponent(a, window=3) + y.powUnsafeExponent(a, window=3) + + x *= y + + check: bool(x == xy) + +proc test_pow0(C: static Curve, gen: RandomGen) = + ## x⁰ = 1 + var x = rng.random_elem(Fp12[C], gen) + var a: BigInt[128] # 0-init + + x.powUnsafeExponent(a, window=3) + check: bool x.isOne() + +proc test_0pow0(C: static Curve, gen: RandomGen) = + ## 0⁰ = 1 + var x: Fp12[C] # 0-init + var a: BigInt[128] # 0-init + + x.powUnsafeExponent(a, window=3) + check: bool x.isOne() + +proc test_powinv(C: static Curve, gen: RandomGen) = + ## xᴬ / xᴮ = xᴬ⁻ᴮ - quotient of power + let x = rng.random_elem(Fp12[C], gen) + + var a = rng.random_elem(BigInt[128], gen) + var b = rng.random_elem(BigInt[128], gen) + # div by 2 to ensure their sum doesn't overflow + # 128 bits + a.div2() + b.div2() + # Ensure a > b + cswap(a, b, a < b) + + var xa = x + xa.powUnsafeExponent(a, window = 3) + + var xb = x + xb.powUnsafeExponent(b, window = 3) + + xb.inv() + xa *= xb + + var xamb = x + var amb: BigInt[128] + discard amb.diff(a, b) + xamb.powUnsafeExponent(amb, window = 3) + + check: bool(xa == xamb) + +proc test_invpow(C: static Curve, gen: RandomGen) = + ## (x / y)ᴬ = xᴬ / yᴬ - power of quotient + let x = rng.random_elem(Fp12[C], gen) + let y = rng.random_elem(Fp12[C], gen) + + var a = rng.random_elem(BigInt[128], gen) + + var xa = x + xa.powUnsafeExponent(a, window = 3) + + var ya = y + ya.powUnsafeExponent(a, window = 3) + ya.inv() + xa *= ya + + var xqya = x + var invy = y + invy.inv() + xqya *= invy + xqya.powUnsafeExponent(a, window = 3) + + check: bool(xa == xqya) + +suite "Exponentiation in 𝔽p12" & " [" & $WordBitwidth & "-bit mode]": + staticFor(curve, TestCurves): + test "xᴬ xᴮ = xᴬ⁺ᴮ on " & $curve: + test_sameBaseProduct(curve, gen = Uniform) + test_sameBaseProduct(curve, gen = HighHammingWeight) + test_sameBaseProduct(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "(xᴬ)ᴮ = xᴬᴮ on " & $curve: + test_powpow(curve, gen = Uniform) + test_powpow(curve, gen = HighHammingWeight) + test_powpow(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "(xy)ᴬ = xᴬyᴬ on " & $curve: + test_powprod(curve, gen = Uniform) + test_powprod(curve, gen = HighHammingWeight) + test_powprod(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "x⁰ = 1 on " & $curve: + test_pow0(curve, gen = Uniform) + test_pow0(curve, gen = HighHammingWeight) + test_pow0(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "0⁰ = 1 on " & $curve: + test_0pow0(curve, gen = Uniform) + test_0pow0(curve, gen = HighHammingWeight) + test_0pow0(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "xᴬ / xᴮ = xᴬ⁻ᴮ on " & $curve: + test_powinv(curve, gen = Uniform) + test_powinv(curve, gen = HighHammingWeight) + test_powinv(curve, gen = Long01Sequence) + + staticFor(curve, TestCurves): + test "(x / y)ᴬ = xᴬ / yᴬ on " & $curve: + test_invpow(curve, gen = Uniform) + test_invpow(curve, gen = HighHammingWeight) + test_invpow(curve, gen = Long01Sequence) diff --git a/tests/t_fp2_frobenius.nim b/tests/t_fp2_frobenius.nim new file mode 100644 index 0000000..ea9530d --- /dev/null +++ b/tests/t_fp2_frobenius.nim @@ -0,0 +1,33 @@ +# 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 + # Internals + ../constantine/towers, + ../constantine/config/curves, + # Test utilities + ./t_fp_tower_frobenius_template + +const TestCurves = [ + # BN254_Nogami + BN254_Snarks, + BLS12_377, + BLS12_381, + # BN446 + # FKM12_447 + # BLS12_461 + # BN462 + ] + +runFrobeniusTowerTests( + ExtDegree = 2, + Iters = 8, + TestCurves = TestCurves, + moduleName = "test_fp2_frobenius", + testSuiteDesc = "𝔽p2 Frobenius map: Frobenius(a, k) = a^(p^k) (mod p²)" +) diff --git a/tests/t_fp2_sqrt.nim b/tests/t_fp2_sqrt.nim index b10122a..b9f30aa 100644 --- a/tests/t_fp2_sqrt.nim +++ b/tests/t_fp2_sqrt.nim @@ -16,9 +16,20 @@ import ../constantine/config/curves, ../constantine/io/io_towers, # Test utilities - ../helpers/prng_unsafe + ../helpers/[prng_unsafe, static_for] -const Iters = 24 +const + Iters = 8 + TestCurves = [ + BN254_Snarks, + BLS12_381 + ] + +type + RandomGen = enum + Uniform + HighHammingWeight + Long01Sequence var rng: RngState let seed = uint32(getTime().toUnix() and (1'i64 shl 32 - 1)) # unixTime mod 2^32 @@ -26,33 +37,43 @@ rng.seed(seed) echo "\n------------------------------------------------------\n" echo "test_fp2_sqrt xoshiro512** seed: ", seed -proc randomSqrtCheck_p3mod4(C: static Curve) = - test "[𝔽p2] Random square root check for p ≡ 3 (mod 4) on " & $Curve(C): - for _ in 0 ..< Iters: - let a = rng.random_unsafe(Fp2[C]) - var na{.noInit.}: Fp2[C] - na.neg(a) +func random_elem(rng: var RngState, F: typedesc, gen: RandomGen): F {.inline, noInit.} = + if gen == Uniform: + result = rng.random_unsafe(F) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(F) + else: + result = rng.random_long01Seq(F) - var a2 = a - var na2 = na - a2.square() - na2.square() - check: - bool a2 == na2 - bool a2.isSquare() +proc randomSqrtCheck_p3mod4(C: static Curve, gen: RandomGen) = + for _ in 0 ..< Iters: + let a = rng.random_elem(Fp2[C], gen) + var na{.noInit.}: Fp2[C] + na.neg(a) - var r, s = a2 - # r.sqrt() - let ok = s.sqrt_if_square() - check: - bool ok - # bool(r == s) - bool(s == a or s == na) + var a2 = a + var na2 = na + a2.square() + na2.square() + check: + bool a2 == na2 + bool a2.isSquare() + + var r, s = a2 + # r.sqrt() + let ok = s.sqrt_if_square() + check: + bool ok + # bool(r == s) + bool(s == a or s == na) proc main() = suite "Modular square root" & " [" & $WordBitwidth & "-bit mode]": - randomSqrtCheck_p3mod4 BN254_Snarks - randomSqrtCheck_p3mod4 BLS12_381 + staticFor(curve, TestCurves): + test "[𝔽p2] Random square root check for p ≡ 3 (mod 4) on " & $curve: + randomSqrtCheck_p3mod4(curve, gen = Uniform) + randomSqrtCheck_p3mod4(curve, gen = HighHammingWeight) + randomSqrtCheck_p3mod4(curve, gen = Long01Sequence) suite "Modular square root - 32-bit bugs highlighted by property-based testing " & " [" & $WordBitwidth & "-bit mode]": test "sqrt_if_square invalid square BLS12_381 - #64": diff --git a/tests/t_fp4_frobenius.nim b/tests/t_fp4_frobenius.nim new file mode 100644 index 0000000..08329a3 --- /dev/null +++ b/tests/t_fp4_frobenius.nim @@ -0,0 +1,33 @@ +# 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 + # Internals + ../constantine/towers, + ../constantine/config/curves, + # Test utilities + ./t_fp_tower_frobenius_template + +const TestCurves = [ + # BN254_Nogami + # BN254_Snarks, + # BLS12_377, + BLS12_381, + # BN446 + # FKM12_447 + # BLS12_461 + # BN462 + ] + +runFrobeniusTowerTests( + ExtDegree = 4, + Iters = 8, + TestCurves = TestCurves, + moduleName = "test_fp4_frobenius", + testSuiteDesc = "𝔽p4 Frobenius map: Frobenius(a, k) = a^(p^k) (mod p⁴)" +) diff --git a/tests/t_fp6_bls12_377.nim b/tests/t_fp6_bls12_377.nim index 2f79766..bd4c804 100644 --- a/tests/t_fp6_bls12_377.nim +++ b/tests/t_fp6_bls12_377.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 6, - Iters = 24, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp6_" & $BLS12_377, testSuiteDesc = "𝔽p6 = 𝔽p2[v] " & $BLS12_377 diff --git a/tests/t_fp6_bls12_381.nim b/tests/t_fp6_bls12_381.nim index 6672874..183cecf 100644 --- a/tests/t_fp6_bls12_381.nim +++ b/tests/t_fp6_bls12_381.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 6, - Iters = 24, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp6_" & $BLS12_381, testSuiteDesc = "𝔽p6 = 𝔽p2[v] " & $BLS12_381 diff --git a/tests/t_fp6_bn254_snarks.nim b/tests/t_fp6_bn254_snarks.nim index 104ebea..a078f37 100644 --- a/tests/t_fp6_bn254_snarks.nim +++ b/tests/t_fp6_bn254_snarks.nim @@ -19,7 +19,7 @@ const TestCurves = [ runTowerTests( ExtDegree = 6, - Iters = 128, + Iters = 12, TestCurves = TestCurves, moduleName = "test_fp6_" & $BN254_Snarks, testSuiteDesc = "𝔽p6 = 𝔽p2[w] " & $BN254_Snarks diff --git a/tests/t_fp_tower_frobenius_template.nim b/tests/t_fp_tower_frobenius_template.nim new file mode 100644 index 0000000..939b5d8 --- /dev/null +++ b/tests/t_fp_tower_frobenius_template.nim @@ -0,0 +1,115 @@ +# 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. + +# ############################################################ +# +# Template tests for towered extension fields +# +# ############################################################ + + +import + # Standard library + std/[unittest, times], + # Internals + ../constantine/towers, + ../constantine/config/[common, curves], + ../constantine/arithmetic, + ../constantine/isogeny/frobenius, + # Test utilities + ../helpers/[prng_unsafe, static_for] + +echo "\n------------------------------------------------------\n" + +template ExtField(degree: static int, curve: static Curve): untyped = + when degree == 2: + Fp2[curve] + elif degree == 4: + Fp4[curve] + elif degree == 6: + Fp6[curve] + elif degree == 12: + Fp12[curve] + else: + {.error: "Unconfigured extension degree".} + +type + RandomGen = enum + Uniform + HighHammingWeight + Long01Sequence + +func random_elem(rng: var RngState, F: typedesc, gen: RandomGen): F {.inline, noInit.} = + if gen == Uniform: + result = rng.random_unsafe(F) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(F) + else: + result = rng.random_long01Seq(F) + +proc runFrobeniusTowerTests*[N]( + ExtDegree: static int, + Iters: static int, + TestCurves: static array[N, Curve], + moduleName: string, + testSuiteDesc: string + ) = + # Random seed for reproducibility + var rng: RngState + let seed = uint32(getTime().toUnix() and (1'i64 shl 32 - 1)) # unixTime mod 2^32 + rng.seed(seed) + echo moduleName, " xoshiro512** seed: ", seed + + suite testSuiteDesc & " [" & $WordBitwidth & "-bit mode]": + test "Frobenius(a) = a^p (mod p^" & $ExtDegree & ")": + proc test(Field: typedesc, Iters: static int, gen: RandomGen) = + for _ in 0 ..< Iters: + var a = rng.random_elem(Field, gen) + var fa {.noInit.}: typeof(a) + fa.frobenius_map(a, k = 1) + + a.powUnsafeExponent(Field.C.Mod, window = 3) + check: bool(a == fa) + + staticFor(curve, TestCurves): + test(ExtField(ExtDegree, curve), Iters, gen = Uniform) + test(ExtField(ExtDegree, curve), Iters, gen = HighHammingWeight) + test(ExtField(ExtDegree, curve), Iters, gen = Long01Sequence) + + # test "Frobenius(a, 2) = a^(p^2) (mod p^" & $ExtDegree & ")": + # proc test(Field: typedesc, Iters: static int, gen: RandomGen) = + # for _ in 0 ..< Iters: + # var a = rng.random_elem(Field, gen) + # var fa {.noInit.}: typeof(a) + # fa.frobenius_map(a, k = 2) + + # a.powUnsafeExponent(Field.C.Mod, window = 3) + # a.powUnsafeExponent(Field.C.Mod, window = 3) + # check: bool(a == fa) + + # staticFor(curve, TestCurves): + # test(ExtField(ExtDegree, curve), Iters, gen = Uniform) + # test(ExtField(ExtDegree, curve), Iters, gen = HighHammingWeight) + # test(ExtField(ExtDegree, curve), Iters, gen = Long01Sequence) + + # test "Frobenius(a, 3) = a^(p^3) (mod p^" & $ExtDegree & ")": + # proc test(Field: typedesc, Iters: static int, gen: RandomGen) = + # for _ in 0 ..< Iters: + # var a = rng.random_elem(Field, gen) + # var fa {.noInit.}: typeof(a) + # fa.frobenius_map(a, k = 3) + + # a.powUnsafeExponent(Field.C.Mod, window = 3) + # a.powUnsafeExponent(Field.C.Mod, window = 3) + # a.powUnsafeExponent(Field.C.Mod, window = 3) + # check: bool(a == fa) + + # staticFor(curve, TestCurves): + # test(ExtField(ExtDegree, curve), Iters, gen = Uniform) + # test(ExtField(ExtDegree, curve), Iters, gen = HighHammingWeight) + # test(ExtField(ExtDegree, curve), Iters, gen = Long01Sequence) diff --git a/tests/t_fp_tower_template.nim b/tests/t_fp_tower_template.nim index c51947f..a9fa4e4 100644 --- a/tests/t_fp_tower_template.nim +++ b/tests/t_fp_tower_template.nim @@ -41,8 +41,8 @@ type HighHammingWeight Long01Sequence -func random_elem(rng: var RngState, F: typedesc, gen: static RandomGen): F {.inline, noInit.} = - when gen == Uniform: +func random_elem(rng: var RngState, F: typedesc, gen: RandomGen): F {.inline, noInit.} = + if gen == Uniform: result = rng.random_unsafe(F) elif gen == HighHammingWeight: result = rng.random_highHammingWeight(F) @@ -77,7 +77,7 @@ proc runTowerTests*[N]( test(ExtField(ExtDegree, curve)) test "Addition, substraction negation are consistent": - proc test(Field: typedesc, Iters: static int, gen: static RandomGen) = + proc test(Field: typedesc, Iters: static int, gen: RandomGen) = # Try to exercise all code paths for in-place/out-of-place add/sum/sub/diff/double/neg # (1 - (-a) - b + (-a) - 2a) + (2a + 2b + (-b)) == 1 var accum {.noInit.}, One {.noInit.}, a{.noInit.}, na{.noInit.}, b{.noInit.}, nb{.noInit.}, a2 {.noInit.}, b2 {.noInit.}: Field diff --git a/tests/t_pairing_bls12_381_line_functions.nim b/tests/t_pairing_bls12_381_line_functions.nim new file mode 100644 index 0000000..bf4edc0 --- /dev/null +++ b/tests/t_pairing_bls12_381_line_functions.nim @@ -0,0 +1,114 @@ +# 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/[tables, unittest, times], + # Internals + ../constantine/config/common, + ../constantine/[arithmetic, primitives], + ../constantine/towers, + ../constantine/config/curves, + ../constantine/io/io_towers, + ../constantine/elliptic/[ + ec_weierstrass_affine, + ec_weierstrass_projective, + ec_scalar_mul], + ../constantine/pairing/lines_projective, + # Test utilities + ../helpers/[prng_unsafe, static_for] + +const + Iters = 4 + TestCurves = [ + BLS12_381 + ] + +type + RandomGen = enum + Uniform + HighHammingWeight + Long01Sequence + +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 "test_pairing_bls12_381_line_functions xoshiro512** seed: ", seed + +func random_point*(rng: var RngState, EC: typedesc, gen: RandomGen): EC {.noInit.} = + if gen == Uniform: + result = rng.random_unsafe(EC) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(EC) + else: + result = rng.random_long01Seq(EC) + +func random_point*(rng: var RngState, EC: typedesc, randZ: bool, gen: RandomGen): EC {.noInit.} = + if not randZ: + if gen == Uniform: + result = rng.random_unsafe(EC) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(EC) + else: + result = rng.random_long01Seq(EC) + else: + if gen == Uniform: + result = rng.random_unsafe_with_randZ(EC) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight_with_randZ(EC) + else: + result = rng.random_long01Seq_with_randZ(EC) + +suite "Pairing - Line Functions on BLS12-381": + test "Line double - lt,t(P)": + proc test_line_double(C: static Curve, randZ: bool, gen: RandomGen) = + for _ in 0 ..< Iters: + let P = rng.random_point(ECP_SWei_Aff[Fp[C]], gen) + var T = rng.random_point(ECP_SWei_Proj[Fp2[C]], randZ, gen) + let Q = rng.random_point(ECP_SWei_Proj[Fp2[C]], randZ, gen) + var l: Line[Fp2[C], C.getSexticTwist()] + + var T2: typeof(Q) + T2.double(T) + l.line_double(T, P) + + doAssert: bool(T == T2) + + staticFor(curve, TestCurves): + test_line_double(curve, randZ = false, gen = Uniform) + test_line_double(curve, randZ = true, gen = Uniform) + test_line_double(curve, randZ = false, gen = HighHammingWeight) + test_line_double(curve, randZ = true, gen = HighHammingWeight) + test_line_double(curve, randZ = false, gen = Long01Sequence) + test_line_double(curve, randZ = true, gen = Long01Sequence) + + test "Line add - lt,q(P)": + proc test_line_add(C: static Curve, randZ: bool, gen: RandomGen) = + for _ in 0 ..< Iters: + let P = rng.random_point(ECP_SWei_Aff[Fp[C]], gen) + let Q = rng.random_point(ECP_SWei_Proj[Fp2[C]], randZ, gen) + var T = rng.random_point(ECP_SWei_Proj[Fp2[C]], randZ, gen) + var l: Line[Fp2[C], C.getSexticTwist()] + + var TQ{.noInit.}: typeof(T) + TQ.sum(T, Q) + + var Qaff{.noInit.}: ECP_SWei_Aff[Fp2[C]] + Qaff.affineFromProjective(Q) + l.line_add(T, Qaff, P) + + doAssert: bool(T == TQ) + + staticFor(curve, TestCurves): + test_line_add(curve, randZ = false, gen = Uniform) + test_line_add(curve, randZ = true, gen = Uniform) + test_line_add(curve, randZ = false, gen = HighHammingWeight) + test_line_add(curve, randZ = true, gen = HighHammingWeight) + test_line_add(curve, randZ = false, gen = Long01Sequence) + test_line_add(curve, randZ = true, gen = Long01Sequence) diff --git a/tests/t_pairing_bls12_381_optate.nim b/tests/t_pairing_bls12_381_optate.nim new file mode 100644 index 0000000..6b824f8 --- /dev/null +++ b/tests/t_pairing_bls12_381_optate.nim @@ -0,0 +1,15 @@ +# 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/config/curves, + ../constantine/pairing/pairing_bls12, + # Test utilities + ./t_pairing_template + +runPairingTests(4, BLS12_381, pairing_bls12_reference) diff --git a/tests/t_pairing_bn254_nogami_optate.nim b/tests/t_pairing_bn254_nogami_optate.nim new file mode 100644 index 0000000..a32b939 --- /dev/null +++ b/tests/t_pairing_bn254_nogami_optate.nim @@ -0,0 +1,15 @@ +# 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/config/curves, + ../constantine/pairing/pairing_bn, + # Test utilities + ./t_pairing_template + +runPairingTests(4, BN254_Nogami, pairing_bn_reference) diff --git a/tests/t_pairing_bn254_snarks_optate.nim b/tests/t_pairing_bn254_snarks_optate.nim new file mode 100644 index 0000000..07946d0 --- /dev/null +++ b/tests/t_pairing_bn254_snarks_optate.nim @@ -0,0 +1,15 @@ +# 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/config/curves, + ../constantine/pairing/pairing_bn, + # Test utilities + ./t_pairing_template + +runPairingTests(4, BN254_Snarks, pairing_bn_reference) diff --git a/tests/t_pairing_fp12_sparse.nim b/tests/t_pairing_fp12_sparse.nim new file mode 100644 index 0000000..671d050 --- /dev/null +++ b/tests/t_pairing_fp12_sparse.nim @@ -0,0 +1,143 @@ +# 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/[tables, unittest, times], + # Internals + ../constantine/config/common, + ../constantine/[arithmetic, primitives], + ../constantine/towers, + ../constantine/config/curves, + ../constantine/io/io_towers, + ../constantine/pairing/[ + lines_projective, + gt_fp12 + ], + # Test utilities + ../helpers/[prng_unsafe, static_for] + +const + Iters = 8 + TestCurves = [ + # BN254_Snarks, + BLS12_381 + ] + +type + RandomGen = enum + Uniform + HighHammingWeight + Long01Sequence + +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 "test_pairing_fp12_sparse xoshiro512** seed: ", seed + +func random_elem(rng: var RngState, F: typedesc, gen: RandomGen): F {.inline, noInit.} = + if gen == Uniform: + result = rng.random_unsafe(F) + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(F) + else: + result = rng.random_long01Seq(F) + +suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent with dense 𝔽p12 mul": + test "Dense 𝔽p6 by Sparse 0y0": + proc test_fp6_0y0(C: static Curve, gen: static RandomGen) = + for _ in 0 ..< Iters: + let a = rng.random_elem(Fp6[C], gen) + let y = rng.random_elem(Fp2[C], gen) + let b = Fp6[C](c1: y) + + var r {.noInit.}, r2 {.noInit.}: Fp6[C] + + r.prod(a, b) + r2.mul_sparse_by_0y0(a, y) + + check: bool(r == r2) + + staticFor(curve, TestCurves): + test_fp6_0y0(curve, gen = Uniform) + test_fp6_0y0(curve, gen = HighHammingWeight) + test_fp6_0y0(curve, gen = Long01Sequence) + + test "Dense 𝔽p6 by Sparse xy0": + proc test_fp6_0y0(C: static Curve, gen: static RandomGen) = + for _ in 0 ..< Iters: + let a = rng.random_elem(Fp6[C], gen) + let x = rng.random_elem(Fp2[C], gen) + let y = rng.random_elem(Fp2[C], gen) + let b = Fp6[C](c0: x, c1: y) + let line = Line[Fp2[C], M_twist](x: x, y: y) + + var r {.noInit.}, r2 {.noInit.}: Fp6[C] + + r.prod(a, b) + r2.mul_by_line_xy0(a, line) + + check: bool(r == r2) + + staticFor(curve, TestCurves): + test_fp6_0y0(curve, gen = Uniform) + test_fp6_0y0(curve, gen = HighHammingWeight) + test_fp6_0y0(curve, gen = Long01Sequence) + + when Fp12[BN254_Snarks]().c0.typeof is Fp6: + test "Sparse 𝔽p12 resulting from xy00z0 line function": + proc test_fp12_xy00z0(C: static Curve, gen: static RandomGen) = + for _ in 0 ..< Iters: + var a = rng.random_elem(Fp12[C], gen) + var a2 = a + + var x = rng.random_elem(Fp2[C], gen) + var y = rng.random_elem(Fp2[C], gen) + var z = rng.random_elem(Fp2[C], gen) + + let line = Line[Fp2[C], Mtwist](x: x, y: y, z: z) + let b = Fp12[C]( + c0: Fp6[C](c0: x, c1: y), + c1: Fp6[C](c1: z) + ) + + a *= b + a2.mul_sparse_by_line_xy00z0(line) + + check: bool(a == a2) + + staticFor(curve, TestCurves): + test_fp12_xy00z0(curve, gen = Uniform) + test_fp12_xy00z0(curve, gen = HighHammingWeight) + test_fp12_xy00z0(curve, gen = Long01Sequence) + + test "Sparse 𝔽p12 resulting from xyz000 line function": + proc test_fp12_xyz000(C: static Curve, gen: static RandomGen) = + for _ in 0 ..< Iters: + var a = rng.random_elem(Fp12[C], gen) + var a2 = a + + var x = rng.random_elem(Fp2[C], gen) + var y = rng.random_elem(Fp2[C], gen) + var z = rng.random_elem(Fp2[C], gen) + + let line = Line[Fp2[C], Mtwist](x: x, y: y, z: z) + let b = Fp12[C]( + c0: Fp6[C](c0: x, c1: y, c2: z) + ) + + a *= b + a2.mul_sparse_by_line_xyz000(line) + + check: bool(a == a2) + + staticFor(curve, TestCurves): + test_fp12_xyz000(curve, gen = Uniform) + test_fp12_xyz000(curve, gen = HighHammingWeight) + test_fp12_xyz000(curve, gen = Long01Sequence) diff --git a/tests/t_pairing_template.nim b/tests/t_pairing_template.nim new file mode 100644 index 0000000..39155d0 --- /dev/null +++ b/tests/t_pairing_template.nim @@ -0,0 +1,92 @@ +# 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/unittest, times, + # Internals + ../constantine/[arithmetic, primitives], + ../constantine/towers, + ../constantine/config/curves, + ../constantine/elliptic/ec_weierstrass_projective, + ../constantine/hash_to_curve/cofactors, + # Test utilities + ../helpers/[prng_unsafe, static_for] + +export + prng_unsafe, times, unittest, + ec_weierstrass_projective, arithmetic, towers, + primitives + +type + RandomGen* = enum + Uniform + HighHammingWeight + Long01Sequence + +func random_point*(rng: var RngState, EC: typedesc, randZ: bool, gen: RandomGen): EC {.noInit.} = + if not randZ: + if gen == Uniform: + result = rng.random_unsafe(EC) + result.clearCofactorReference() + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight(EC) + result.clearCofactorReference() + else: + result = rng.random_long01Seq(EC) + result.clearCofactorReference() + else: + if gen == Uniform: + result = rng.random_unsafe_with_randZ(EC) + result.clearCofactorReference() + elif gen == HighHammingWeight: + result = rng.random_highHammingWeight_with_randZ(EC) + result.clearCofactorReference() + else: + result = rng.random_long01Seq_with_randZ(EC) + result.clearCofactorReference() + +template runPairingTests*(Iters: static int, C: static Curve, pairing_fn: untyped): untyped {.dirty.}= + var rng: RngState + let timeseed = uint32(toUnix(getTime()) and (1'i64 shl 32 - 1)) # unixTime mod 2^32 + seed(rng, timeseed) + echo "\n------------------------------------------------------\n" + echo "test_pairing_",$C,"_optate xoshiro512** seed: ", timeseed + + proc test_bilinearity_double_impl(randZ: bool, gen: RandomGen) = + for _ in 0 ..< Iters: + let P = rng.random_point(ECP_SWei_Proj[Fp[C]], randZ, gen) + let Q = rng.random_point(ECP_SWei_Proj[Fp2[C]], randZ, gen) + var P2: typeof(P) + var Q2: typeof(Q) + + var r {.noInit.}, r2 {.noInit.}, r3 {.noInit.}: Fp12[C] + + P2.double(P) + Q2.double(Q) + + r.pairing_fn(P, Q) + r.square() + r2.pairing_fn(P2, Q) + r3.pairing_fn(P, Q2) + + check: + bool(not r.isZero()) + bool(not r.isOne()) + bool(r == r2) + bool(r == r3) + bool(r2 == r3) + + suite "Pairing - Optimal Ate on " & $C: + test "Bilinearity e([2]P, Q) = e(P, [2]Q) = e(P, Q)^2": + test_bilinearity_double_impl(randZ = false, gen = Uniform) + test_bilinearity_double_impl(randZ = true, gen = Uniform) + test_bilinearity_double_impl(randZ = false, gen = HighHammingWeight) + test_bilinearity_double_impl(randZ = true, gen = HighHammingWeight) + test_bilinearity_double_impl(randZ = false, gen = Long01Sequence) + test_bilinearity_double_impl(randZ = true, gen = Long01Sequence)