Optimize Miller Loop and prepare Multi-pairing (#159)

* Pairing with affine: align API to BLST and Gurvy and common use-case.

* Implement multi-pairing / aggregate verif for BLS12-381 (+2% pairing perf)

* Generalize the optimized miller loop for single pairing

* Immplement the miller loop addchain for BLS12-377

* Miller addition chain for BN254-Nogami

* no Miller adchain for BN254-Snarks

* Update the line test with new tower https://github.com/mratsim/constantine/pull/153

* Somewhat sparse for Fp2 M-Twist

* Implement line by line multiplication for Fp12 D-Twist

* Somewhat sparse Mul for Fp12 D-Twist

* Finish the sparse and somewhat sparse multiplications
This commit is contained in:
Mamy Ratsimbazafy 2021-02-14 13:06:57 +01:00 committed by GitHub
parent 0e43c12095
commit 9ac9862401
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
23 changed files with 1448 additions and 180 deletions

View File

@ -25,7 +25,7 @@ import
# ############################################################
const Iters = 50
const Iters = 1000
const AvailableCurves = [
BLS12_377,
]
@ -37,6 +37,10 @@ proc main() =
lineDoubleBench(curve, Iters)
lineAddBench(curve, Iters)
mulFp12byLine_xyz000_Bench(curve, Iters)
mulLinebyLine_xyz000_Bench(curve, Iters)
mulFp12by_abcdefghij00_Bench(curve, Iters)
mulFp12_by_2lines_v1_xyz000_Bench(curve, Iters)
mulFp12_by_2lines_v2_xyz000_Bench(curve, Iters)
separator()
finalExpEasyBench(curve, Iters)
finalExpHardBLS12Bench(curve, Iters)

View File

@ -25,7 +25,7 @@ import
# ############################################################
const Iters = 50
const Iters = 1000
const AvailableCurves = [
BLS12_381,
]
@ -37,6 +37,10 @@ proc main() =
lineDoubleBench(curve, Iters)
lineAddBench(curve, Iters)
mulFp12byLine_xy000z_Bench(curve, Iters)
mulLinebyLine_xy000z_Bench(curve, Iters)
mulFp12by_abcd00efghij_Bench(curve, Iters)
mulFp12_by_2lines_v1_xy000z_Bench(curve, Iters)
mulFp12_by_2lines_v2_xy000z_Bench(curve, Iters)
separator()
finalExpEasyBench(curve, Iters)
finalExpHardBLS12Bench(curve, Iters)

View File

@ -25,7 +25,7 @@ import
# ############################################################
const Iters = 50
const Iters = 1000
const AvailableCurves = [
BN254_Nogami,
]
@ -37,6 +37,10 @@ proc main() =
lineDoubleBench(curve, Iters)
lineAddBench(curve, Iters)
mulFp12byLine_xyz000_Bench(curve, Iters)
mulLinebyLine_xyz000_Bench(curve, Iters)
mulFp12by_abcdefghij00_Bench(curve, Iters)
mulFp12_by_2lines_v1_xyz000_Bench(curve, Iters)
mulFp12_by_2lines_v2_xyz000_Bench(curve, Iters)
separator()
finalExpEasyBench(curve, Iters)
finalExpHardBNBench(curve, Iters)

View File

@ -25,7 +25,7 @@ import
# ############################################################
const Iters = 50
const Iters = 1000
const AvailableCurves = [
BN254_Snarks,
]
@ -37,6 +37,10 @@ proc main() =
lineDoubleBench(curve, Iters)
lineAddBench(curve, Iters)
mulFp12byLine_xyz000_Bench(curve, Iters)
mulLinebyLine_xyz000_Bench(curve, Iters)
mulFp12by_abcdefghij00_Bench(curve, Iters)
mulFp12_by_2lines_v1_xyz000_Bench(curve, Iters)
mulFp12_by_2lines_v2_xyz000_Bench(curve, Iters)
separator()
finalExpEasyBench(curve, Iters)
finalExpHardBNBench(curve, Iters)

View File

@ -27,10 +27,12 @@ import
pairing_bls12,
pairing_bn
],
../constantine/curves/zoo_pairings,
# Helpers
../helpers/prng_unsafe,
./bench_blueprint
export zoo_pairings # generic sandwich https://github.com/nim-lang/Nim/issues/11225
export notes
proc separator*() = separator(132)
@ -46,6 +48,14 @@ template bench(op: string, C: static Curve, iters: int, body: untyped): untyped
measure(iters, startTime, stopTime, startClk, stopClk, body)
report(op, $C, startTime, stopTime, startClk, stopClk, iters)
func clearCofactorReference[F; Tw: static Twisted](
ec: var ECP_ShortW_Aff[F, Tw]) =
# For now we don't have any affine operation defined
var t {.noInit.}: ECP_ShortW_Prj[F, Tw]
t.projectiveFromAffine(ec)
t.clearCofactorReference()
ec.affineFromProjective(t)
func random_point*(rng: var RngState, EC: typedesc): EC {.noInit.} =
result = rng.random_unsafe(EC)
result.clearCofactorReference()
@ -53,34 +63,24 @@ func random_point*(rng: var RngState, EC: typedesc): EC {.noInit.} =
proc lineDoubleBench*(C: static Curve, iters: int) =
var line: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
var Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Paff.affineFromProjective(P)
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
bench("Line double", C, iters):
line.line_double(T, Paff)
line.line_double(T, P)
proc lineAddBench*(C: static Curve, iters: int) =
var line: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let
P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
var
Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Qaff: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Aff[Fp2[C], OnTwist])
bench("Line add", C, iters):
line.line_add(T, Qaff, Paff)
line.line_add(T, Q, P)
proc mulFp12byLine_xyz000_Bench*(C: static Curve, iters: int) =
var line: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
var Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Paff.affineFromProjective(P)
line.line_double(T, Paff)
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
line.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("Mul 𝔽p12 by line xyz000", C, iters):
@ -89,45 +89,116 @@ proc mulFp12byLine_xyz000_Bench*(C: static Curve, iters: int) =
proc mulFp12byLine_xy000z_Bench*(C: static Curve, iters: int) =
var line: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
var Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Paff.affineFromProjective(P)
line.line_double(T, Paff)
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
line.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("Mul 𝔽p12 by line xy000z", C, iters):
f.mul_sparse_by_line_xy000z(line)
proc mulLinebyLine_xyz000_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("Mul line xyz000 by line xyz000", C, iters):
f.mul_xyz000_xyz000_into_abcdefghij00(l0, l1)
proc mulLinebyLine_xy000z_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("Mul line xy000z by line xy000z", C, iters):
f.mul_xy000z_xy000z_into_abcd00efghij(l0, l1)
proc mulFp12by_abcdefghij00_Bench*(C: static Curve, iters: int) =
var f = rng.random_unsafe(Fp12[C])
let g = rng.random_unsafe(Fp12[C])
bench("Mul 𝔽p12 by abcdefghij00", C, iters):
f.mul_sparse_by_abcdefghij00(g)
proc mulFp12by_abcd00efghij_Bench*(C: static Curve, iters: int) =
var f = rng.random_unsafe(Fp12[C])
let g = rng.random_unsafe(Fp12[C])
bench("Mul 𝔽p12 by abcd00efghij", C, iters):
f.mul_sparse_by_abcd00efghij(g)
proc mulFp12_by_2lines_v1_xyz000_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("mulFp12 by 2 lines v1", C, iters):
f.mul_sparse_by_line_xyz000(l0)
f.mul_sparse_by_line_xyz000(l1)
proc mulFp12_by_2lines_v2_xyz000_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("mulFp12 by 2 lines v2", C, iters):
var f2 {.noInit.}: Fp12[C]
f2.mul_xyz000_xyz000_into_abcdefghij00(l0, l1)
f.mul_sparse_by_abcdefghij00(f2)
proc mulFp12_by_2lines_v1_xy000z_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("mulFp12 by 2 lines v1", C, iters):
f.mul_sparse_by_line_xy000z(l0)
f.mul_sparse_by_line_xy000z(l1)
proc mulFp12_by_2lines_v2_xy000z_Bench*(C: static Curve, iters: int) =
var l0, l1: Line[Fp2[C]]
var T = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
let P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
l0.line_double(T, P)
l1.line_double(T, P)
var f = rng.random_unsafe(Fp12[C])
bench("mulFp12 by 2 lines v2", C, iters):
var f2 {.noInit.}: Fp12[C]
f2.mul_xy000z_xy000z_into_abcd00efghij(l0, l1)
f.mul_sparse_by_abcd00efghij(f2)
proc millerLoopBLS12Bench*(C: static Curve, iters: int) =
let
P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
var
Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Qaff: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Aff[Fp2[C], OnTwist])
var f: Fp12[C]
bench("Miller Loop BLS12", C, iters):
f.millerLoopGenericBLS12(Paff, Qaff)
f.millerLoopGenericBLS12(P, Q)
proc millerLoopBNBench*(C: static Curve, iters: int) =
let
P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
var
Paff: ECP_ShortW_Aff[Fp[C], NotOnTwist]
Qaff: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Aff[Fp2[C], OnTwist])
var f: Fp12[C]
bench("Miller Loop BN", C, iters):
f.millerLoopGenericBN(Paff, Qaff)
f.millerLoopGenericBN(P, Q)
proc finalExpEasyBench*(C: static Curve, iters: int) =
var r = rng.random_unsafe(Fp12[C])
@ -160,20 +231,18 @@ proc finalExpBNBench*(C: static Curve, iters: int) =
proc pairingBLS12Bench*(C: static Curve, iters: int) =
let
P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Aff[Fp2[C], OnTwist])
var f: Fp12[C]
bench("Pairing BLS12", C, iters):
f.pairing_bls12(P, Q)
proc pairingBNBench*(C: static Curve, iters: int) =
let
P = rng.random_point(ECP_ShortW_Prj[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Prj[Fp2[C], OnTwist])
P = rng.random_point(ECP_ShortW_Aff[Fp[C], NotOnTwist])
Q = rng.random_point(ECP_ShortW_Aff[Fp2[C], OnTwist])
var f: Fp12[C]
bench("Pairing BN", C, iters):
f.pairing_bn(P, Q)

View File

@ -142,7 +142,7 @@ const testDesc: seq[tuple[path: string, useGMP: bool]] = @[
# ----------------------------------------------------------
# ("tests/t_pairing_bls12_377_line_functions.nim", false),
# ("tests/t_pairing_bls12_381_line_functions.nim", false),
# ("tests/t_pairing_mul_fp12_by_lines.nim", false),
("tests/t_pairing_mul_fp12_by_lines.nim", false),
# ("tests/t_pairing_cyclotomic_fp12.nim", false),
("tests/t_pairing_bn254_nogami_optate.nim", false),
("tests/t_pairing_bn254_snarks_optate.nim", false),

View File

@ -7,10 +7,11 @@
# at your option. This file may not be copied, modified, or distributed except according to those terms.
import
../config/[curves, type_bigint],
../config/[curves, type_bigint, type_ff],
../io/io_bigints,
../towers,
../pairing/cyclotomic_fp12
../elliptic/[ec_shortweierstrass_affine, ec_shortweierstrass_projective],
../pairing/[cyclotomic_fp12, miller_loops]
# Slow generic implementation
# ------------------------------------------------------------
@ -27,13 +28,38 @@ const BLS12_377_pairing_finalexponent* = block:
# (p^12 - 1) / r * 3
BigInt[4271].fromHex"0x518fe3a450394da01ed0ec73865aed18d4251c557c299312d07b5d31105598be5439b32fda943a26e8d85c306e6c1941dd3f9d646d87211c240f5489c67b1a8663c49da97a2880dc48213527e51d370acd05663ffda035ca31c4ba994c89d66c0c97066502f8ef19bb008e047c24cf96e02493f4683ffdc39075cc1c01df9fd0ec1dc0419176c010ac1a83b777201a77f8dab474e99c59ae840de7362f7c231d500aecc1eb52616067540d419f7f9fbfd22831919b4ac04960703d9753698941c95aa2d2a04f4bf26de9d191661a013cbb09227c09424595e2639ae94d35ce708bdec2c10628eb4f981945698ef049502d2a71994fab9898c028c73dd021f13208590be27e78f0f18a88f5ffe40157a9e9fef5aa229c0aa7fdb16a887af2c4a486258bf11fb1a5d945707a89d7bf8f67e5bb28f76a460d9a1e660cbbe91bfc456b8789d5bae1dba8cbef5b03bcd0ea30f6a7b45218292b2bf3b20ed5937cb5e2250eee395821805c6383d0286c7423beb42e79f85dab2a36df8fd154f2d89e5e9aaadaaa00e0a29ecc6e329195761d6063e0a2e136a3fb7671c9134c970a8588a7f3144642a10a5af77c105f5e90987f28c6604c5dcb604c02f7d642f7f819eea6fadb8aace7c4e146a17dab2c644d4372c6979845f261b4a20cd88a20325e0c0fc806bd9f60a8502fa8f466b6919311e232e06fd6a861cb5dc24d69274c7e631cac6b93e0254460d445a0000012b53b000000000000"
# Addition chain
# Addition chains
# ------------------------------------------------------------
#
# u = 0x8508c00000000001
# Ate BLS |u|
# hex: 0x8508c00000000001
# bin: 0b1000010100001000110000000000000000000000000000000000000000000001
#
# 71 naive operations to build a naive addition chain
# or 68 with optimizations, though unsure how to use them in the Miller Loop
# and for multi-pairing it would use too much temporaries anyway.
func millerLoopAddchain*(
f: var Fp12[BLS12_377],
Q: ECP_ShortW_Aff[Fp2[BLS12_377], OnTwist],
P: ECP_ShortW_Aff[Fp[BLS12_377], NotOnTwist]
) =
## Miller Loop for BLS12-377 curve
## Computes f{u,Q}(P) with u the BLS curve parameter
var T {.noInit.}: ECP_ShortW_Prj[Fp2[BLS12_377], OnTwist]
f.miller_init_double_then_add(T, Q, P, 5) # 0b100001
f.miller_accum_double_then_add(T, Q, P, 2) # 0b10000101
f.miller_accum_double_then_add(T, Q, P, 5) # 0b1000010100001
f.miller_accum_double_then_add(T, Q, P, 4) # 0b10000101000010001
f.miller_accum_double_then_add(T, Q, P, 1) # 0b100001010000100011
f.miller_accum_double_then_add(T, Q, P, 46, add = true) # 0b1000010100001000110000000000000000000000000000000000000000000001
func pow_x*(r: var Fp12[BLS12_377], a: Fp12[BLS12_377], invert = BLS12_377_pairing_ate_param_isNeg) =
## f^x with x the curve parameter
## For BLS12_377 f^-0x8508c00000000001
## Warning: The parameter is odd and needs a correction
r.cyclotomic_square(a)
r *= a
r.cyclotomic_square()
@ -53,7 +79,7 @@ func pow_x*(r: var Fp12[BLS12_377], a: Fp12[BLS12_377], invert = BLS12_377_pairi
r.cycl_sqr_repeated(10)
r *= t100011
r.cycl_sqr_repeated(46)
r.cycl_sqr_repeated(46) # TODO: Karabina's compressed squarings
r *= a
if invert:

View File

@ -7,10 +7,11 @@
# at your option. This file may not be copied, modified, or distributed except according to those terms.
import
../config/[curves, type_bigint],
../config/[curves, type_bigint, type_ff],
../io/io_bigints,
../towers,
../pairing/cyclotomic_fp12
../elliptic/[ec_shortweierstrass_affine, ec_shortweierstrass_projective],
../pairing/[cyclotomic_fp12, miller_loops]
# Slow generic implementation
# ------------------------------------------------------------
@ -27,8 +28,58 @@ const BLS12_381_pairing_finalexponent* = block:
# (p^12 - 1) / r * 3
BigInt[4316].fromHex"0x8ca592196587127a538fd40dc3e541f9dca04bb7dc671be77cf17715a2b2fe3bea73dfb468d8f473094aecb7315a664019fbd84913caba6579c08fd42009fe1bd6fcbce15eacb2cf3218a165958cb8bfdae2d2d54207282314fc0dea9d6ff3a07dbd34efb77b732ba5f994816e296a72928cfee133bdc3ca9412b984b9783d9c6aa81297ab1cd294a502304773528bbae8706979f28efa0d355b0224e2513d6e4a5d3bb4dde0523678105d9167ff1323d6e99ac312d8a7d762336370c4347bb5a7e405d6f3496b2dd38e722d4c1f3ac25e3167ec2cb543d69430c37c2f98fcdd0dd36caa9f5aa7994cec31b24ed5e515911037b376e521070d29c9d56cfa8c3574363efb20f28c19e4105ab99edd44084bd23725017931d6740bda71e5f07600ce6b407e543c4bc40bcd4c0b600e6c98003bf8548986b14d9098746dc89d154af91ad54f337b31c79222145dd3ed254fdeda0300c49ebcd2352765f533883a3513435f3ee452496f5166c25bf503bd6ec0a0679efda3b46ebf86211d458de749460d4a2a19abe6ea2accb451ab9a096b98465d044dc2a7f86c253a4ee57b6df108eff598a8dbc483bf8b74c2789939db85ffd7e0fd55b32bc26877f5be26fa7d750500ce2fab93c0cbe7336b126a5693d0c16484f37addccc7642590dbe98538990b88637e374d545d9b34b67448d0357e60280bbd8542f1f4e813caa8e8db57364b4e0cc14f35af381dd9b71ec9292b3a3f16e42362d2019e05f30"
# Addition chain
# Addition chains
# ------------------------------------------------------------
#
# u = -0xd201000000010000
# Ate BLS |u|
# hex: 0xd201000000010000
# bin: 0b1101001000000001000000000000000000000000000000010000000000000000
#
# 68 operations to build an addition chain
func millerLoopAddchain*(
f: var Fp12[BLS12_381],
Q: ECP_ShortW_Aff[Fp2[BLS12_381], OnTwist],
P: ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]
) =
## Miller Loop for BLS12-381 curve
## Computes f{u,Q}(P) with u the BLS curve parameter
var T {.noInit.}: ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]
f.miller_init_double_then_add(T, Q, P, 1) # 0b11
f.miller_accum_double_then_add(T, Q, P, 2) # 0b1101
f.miller_accum_double_then_add(T, Q, P, 3) # 0b1101001
f.miller_accum_double_then_add(T, Q, P, 9) # 0b1101001000000001
f.miller_accum_double_then_add(T, Q, P, 32) # 0b110100100000000100000000000000000000000000000001
f.miller_accum_double_then_add(T, Q, P, 16, add = false) # 0b1101001000000001000000000000000000000000000000010000000000000000
# Negative AteParam, conjugation eliminated by final exponentiation
# f.conj()
func millerLoopAddchain*[N: static int](
f: var Fp12[BLS12_381],
Qs: array[N, ECP_ShortW_Aff[Fp2[BLS12_381], OnTwist]],
Ps: array[N, ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]]
) =
## Generic Miller Loop for BLS12 curve
## Computes f{u,Q}(P) with u the BLS curve parameter
var Ts {.noInit.}: array[N, ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]]
# Ate param addition chain
# Hex: 0xd201000000010000
# Bin: 0b1101001000000001000000000000000000000000000000010000000000000000
f.miller_init_double_then_add(Ts, Qs, Ps, 1) # 0b11
f.miller_accum_double_then_add(Ts, Qs, Ps, 2) # 0b1101
f.miller_accum_double_then_add(Ts, Qs, Ps, 3) # 0b1101001
f.miller_accum_double_then_add(Ts, Qs, Ps, 9) # 0b1101001000000001
f.miller_accum_double_then_add(Ts, Qs, Ps, 32) # 0b110100100000000100000000000000000000000000000001
f.miller_accum_double_then_add(Ts, Qs, Ps, 16, add = false) # 0b1101001000000001000000000000000000000000000000010000000000000000
# TODO: what is the threshold for Karabina's compressed squarings?
func pow_xdiv2*(r: var Fp12[BLS12_381], a: Fp12[BLS12_381], invert = BLS12_381_pairing_ate_param_isNeg) =
## f^(x/2) with x the curve parameter

View File

@ -7,10 +7,11 @@
# at your option. This file may not be copied, modified, or distributed except according to those terms.
import
../config/[curves, type_bigint],
../config/[curves, type_bigint, type_ff],
../io/io_bigints,
../towers,
../pairing/cyclotomic_fp12
../elliptic/[ec_shortweierstrass_affine, ec_shortweierstrass_projective],
../pairing/[cyclotomic_fp12, miller_loops]
# Slow generic implementation
# ------------------------------------------------------------
@ -27,14 +28,40 @@ const BN254_Nogami_pairing_finalexponent* = block:
# (p^12 - 1) / r
BigInt[2786].fromHex"0x2928fbb36b391596ee3fe4cbe857330da83e46fedf04d235a4a8daf5ff9f6eabcb4e3f20aa06f0a0d96b24f9af0cbbce750d61627dcbf5fec9139b8f1c46c86b49b4f8a202af26e4504f2c0f56570e9bd5b94c403f385d1908556486e24b396ddc2cdf13d06542f84fe8e82ccbad7b7423fc1ef4e8cc73d605e3e867c0a75f45ea7f6356d9846ce35d5a34f30396938818ad41914b97b99c289a7259b5d2e09477a77bd3c409b19f19e893f8ade90b0aed1b5fc8a07a3cebb41d4e9eee96b21a832ddb1e93e113edfb704fa532848c18593cd0ee90444a1b3499a800177ea38bdec62ec5191f2b6bbee449722f98d2173ad33077545c2ad10347e125a56fb40f086e9a4e62ad336a72c8b202ac3c1473d73b93d93dc0795ca0ca39226e7b4c1bb92f99248ec0806e0ad70744e9f2238736790f5185ea4c70808442a7d530c6ccd56b55a6973867ec6c73599bbd020bbe105da9c6b5c009ad8946cd6f0"
# Addition chain
# Addition chains
# ------------------------------------------------------------
#
# u = -0x4080000000000001
# Ate BN |6u+2|
# hex: 0x18300000000000004
# bin: 0x11000001100000000000000000000000000000000000000000000000000000100
func millerLoopAddchain*(
f: var Fp12[BN254_Nogami],
Q: ECP_ShortW_Aff[Fp2[BN254_Nogami], OnTwist],
P: ECP_ShortW_Aff[Fp[BN254_Nogami], NotOnTwist]
) =
## Miller Loop for BN254-Nogami curve
## Computes f{6u+2,Q}(P) with u the BLS curve parameter
var T {.noInit.}: ECP_ShortW_Prj[Fp2[BN254_Nogami], OnTwist]
f.miller_init_double_then_add(T, Q, P, 1) # 0b11
f.miller_accum_double_then_add(T, Q, P, 6) # 0b11000001
f.miller_accum_double_then_add(T, Q, P, 1) # 0b110000011
f.miller_accum_double_then_add(T, Q, P, 54) # 0b110000011000000000000000000000000000000000000000000000000000001
f.miller_accum_double_then_add(T, Q, P, 2, add = false) # 0b11000001100000000000000000000000000000000000000000000000000000100
# Negative AteParam
f.conj()
# Ate pairing for BN curves needs adjustment after basic Miller loop
f.millerCorrectionBN(T, Q, P, BN254_Nogami_pairing_ate_param_isNeg)
func pow_u*(r: var Fp12[BN254_Nogami], a: Fp12[BN254_Nogami], invert = BN254_Nogami_pairing_ate_param_isNeg) =
## f^u with u the curve parameter
## For BN254_Nogami f^-0x4080000000000001
r = a
r.cycl_sqr_repeated(7)
r.cyclotomic_square(a)
r.cycl_sqr_repeated(6)
r *= a
r.cycl_sqr_repeated(55)
r *= a

View File

@ -29,6 +29,15 @@ const BN254_Snarks_pairing_finalexponent* = block:
# Addition chain
# ------------------------------------------------------------
#
# u = 0x44e992b44a6909f1
# Ate BN |6u+2|
# hex: 0x19d797039be763ba8
# bin: 0x11001110101111001011100000011100110111110011101100011101110101000
#
# We don't define an addition chain for the Miller loop
# it would requires saving accumulators to actually save
# operations compared to NAF, and can we combine the saved EC[Fp2] accumulators?
func pow_u*(r: var Fp12[BN254_Snarks], a: Fp12[BN254_Snarks], invert = BN254_Snarks_pairing_ate_param_isNeg) =
## f^u with u the curve parameter

View File

@ -21,4 +21,4 @@ macro pairing*(C: static Curve, value: untyped): untyped =
## Get pairing related constants
return bindSym($C & "_pairing_" & $value)
export pow_x, pow_xdiv2, pow_u
export pow_x, pow_xdiv2, pow_u, millerLoopAddchain

View File

@ -85,6 +85,14 @@
Sylvain Duquesne and Loubna Ghammam, 2015\
https://eprint.iacr.org/2015/192
- Software Implementation,
Aranha, Dominguez Perez, A. Mrabet, Schwabe,
Guide to Pairing-Based Cryptography, 2015
- Physical Attacks,
N. El Mrabet, Goubin, Guilley, Fournier, Jauvart, Moreau, Rauzy, Rondepierre,
Guide to Pairing-Based Cryptography, 2015
- A taxonomy of pairings, their security, their complexity\
Razvan Barbulescu, Nadia El Mrabet, and Loubna Ghammam, 2019\
https://hal.archives-ouvertes.fr/hal-02129868/file/2019-485.pdf

View File

@ -11,14 +11,14 @@ import
ec_shortweierstrass_affine,
ec_shortweierstrass_projective
],
../isogeny/frobenius,
./lines_projective,
./mul_fp6_by_lines, ./mul_fp12_by_lines,
../curves/zoo_pairings
./mul_fp6_by_lines, ./mul_fp12_by_lines
# ############################################################
#
# Basic Miller Loop
#
# #
# Basic Miller Loop #
# #
# ############################################################
template basicMillerLoop*[FT, F1, F2](
@ -31,6 +31,8 @@ template basicMillerLoop*[FT, F1, F2](
ate_param_isNeg: untyped
) =
## Basic Miller loop iterations
mixin pairing # symbol from zoo_pairings
static:
doAssert FT.C == F1.C
doAssert FT.C == F2.C
@ -57,3 +59,268 @@ template basicMillerLoop*[FT, F1, F2](
# In GT, x^-1 == conjugate(x)
# Remark 7.1, chapter 7.1.1 of Guide to Pairing-Based Cryptography, El Mrabet, 2017
conj(f)
func millerCorrectionBN*[FT, F1, F2](
f: var FT,
T: var ECP_ShortW_Prj[F2, OnTwist],
Q: ECP_ShortW_Aff[F2, OnTwist],
P: ECP_ShortW_Aff[F1, NotOnTwist],
ate_param_isNeg: static bool
) =
## Ate pairing for BN curves need adjustment after basic Miller loop
static:
doAssert FT.C == F1.C
doAssert FT.C == F2.C
when ate_param_isNeg:
T.neg()
var V {.noInit.}: typeof(Q)
var line {.noInit.}: Line[F2]
V.frobenius_psi(Q)
line.line_add(T, V, P)
f.mul(line)
V.frobenius_psi(Q, 2)
V.neg()
line.line_add(T, V, P)
f.mul(line)
# ############################################################
# #
# Optimized Miller Loops #
# #
# ############################################################
#
# - Software Implementation, Algorithm 11.2 & 11.3
# Aranha, Dominguez Perez, A. Mrabet, Schwabe,
# Guide to Pairing-Based Cryptography, 2015
#
# - Physical Attacks,
# N. El Mrabet, Goubin, Guilley, Fournier, Jauvart, Moreau, Rauzy, Rondepierre,
# Guide to Pairing-Based Cryptography, 2015
#
# - Pairing Implementation Revisited
# Mike Scott, 2019
# https://eprint.iacr.org/2019/077.pdf
#
# Fault attacks:
# To limit exposure to some fault attacks (flipping bits with a laser on embedded):
# - changing the number of Miller loop iterations
# - flipping the bits in the Miller loop
# we hardcode unrolled addition chains.
# This should also contribute to performance.
#
# Multi-pairing discussion:
# Aranha & Scott proposes 2 different approaches for multi-pairing.
#
# -----
# Scott
#
# Algorithm 2: Calculate and store line functions for BLS12 curve
# Input: Q ∈ G2, P ∈ G1 , curve parameter u
# Output: An array g of blog2(u)c line functions ∈ Fp12
# 1 T ← Q
# 2 for i ← ceil(log2(u)) 1 to 0 do
# 3 g[i] ← lT,T(P), T ← 2T
# 4 if ui = 1 then
# 5 g[i] ← g[i].lT,Q(P), T ← T + Q
# 6 return g
#
# And to accumulate lines from a new (P, Q) tuple of points
#
# Algorithm 4: Accumulate another set of line functions into g
# Input: The array g, Qj ∈ G2 , Pj ∈ G1 , curve parameter u
# Output: Updated array g of ceil(log2(u)) line functions ∈ Fp12
# 1 T ← Qj
# 2 for i ← blog2 (u)c 1 to 0 do
# 3 t ← lT,T (Pj), T ← 2T
# 4 if ui = 1 then
# 5 t ← t.lT,Qj (Pj), T ← T + Qj
# 6 g[i] ← g[i].t
# 7 return g
#
# ------
# Aranha
#
# Algorithm 11.2 Explicit multipairing version of Algorithm 11.1.
# (we extract the Miller Loop part only)
# Input : P1 , P2 , . . . Pn ∈ G1 ,
# Q1 , Q2, . . . Qn ∈ G2
# Output: (we focus on the Miller Loop)
#
# Write l in binary form, l = sum(0 ..< m-1)
# f ← 1, l ← abs(AteParam)
# for j ← 1 to n do
# Tj ← Qj
# end
#
# for i = m-2 down to 0 do
# f ← f²
# for j ← 1 to n do
# f ← f gTj,Tj(Pj), Tj ← [2]Tj
# if li = 1 then
# f ← f gTj,Qj(Pj), Tj ← Tj + Qj
# end
# end
# end
#
# -----
# Assuming we have N tuples (Pj, Qj) of points j in 0 ..< N
# and I operations to do in our Miller loop:
# - I = HammingWeight(AteParam) + Bitwidth(AteParam)
# - HammingWeight(AteParam) corresponds to line additions
# - Bitwidth(AteParam) corresponds to line doublings
#
# Scott approach is to have:
# - I Fp12 accumulators `g`
# - 1 G2 accumulator `T`
# and then accumulating each (Pj, Qj) into their corresponding `g` accumulator.
#
# Aranha approach is to have:
# - 1 Fp12 accumulator `f`
# - N G2 accumulators `T`
# and accumulate N points per I.
#
# Scott approach is fully "online"/"streaming",
# while Aranha's saves space.
# For BLS12_381,
# I = 68 hence we would need 68*12*48 = 39168 bytes (381-bit needs 48 bytes)
# G2 has size 3*2*48 = 288 bytes (3 proj coordinates on Fp2)
# and we choose N (which can be 1 for single pairing or reverting to Scott approach).
#
# In actual use, "streaming pairings" are not used, pairings to compute are receive
# by batch, for example for blockchain you receive a batch of N blocks to verify from one peer.
# Furthermore, 39kB would be over L1 cache size and incurs cache misses.
# Additionally Aranha approach would make it easier to batch inversions
# using Montgomery's simultaneous inversion technique.
# Lastly, while a higher level API will need to store N (Pj, Qj) pairs for multi-pairings
# for Aranha approach, it can decide how big N is depending on hardware and/or protocol.
#
# Regarding optimizations, as the Fp12 accumulator is dense
# and lines are sparse (xyz000 or xy000z) Scott mentions the following costs:
# - squaring is 11m
# - Dense-sparse is 13m
# - sparse-sparse is 6m
# - Dense-(somewhat sparse) is 17m
# Hence when accumulating lines from multiple points:
# - 2x Dense-sparse is 26m
# - sparse-sparse then Dense-(somewhat sparse) is 23m
# a 11.5% speedup
#
# We can use Aranha approach but process lines function 2-by-2 merging them
# before merging them to the dense Fp12 accumulator.
#
# In benchmarks though, the speedup doesn't work for BN curves but does for BLS curves.
#
# For single pairings
# Unfortunately, it's BN254_Snarks which requires a lot of addition in the Miller loop.
# BLS12-377 and BLS12-381 require 6 and 7 line addition in their Miller loop,
# the saving is about 150 cycles per addition for about 1000 cycles saved.
# A full pairing is ~2M cycles so this is only 0.5% for significantly
# more maintenance and bounds analysis complexity.
#
# For multipairing it is interesting since for a BLS signature verification (double pairing)
# we would save 1000 cycles per Ate iteration so ~70000 cycles, while a Miller loop is ~800000 cycles.
# Miller Loop - single pairing
# ----------------------------------------------------------------------------
func miller_init_double_then_add*[FT, F1, F2](
f: var FT,
T: var ECP_ShortW_Prj[F2, OnTwist],
Q: ECP_ShortW_Aff[F2, OnTwist],
P: ECP_ShortW_Aff[F1, NotOnTwist],
numDoublings: static int
) =
## Start a Miller Loop with
## - `numDoubling` doublings
## - 1 add
##
## f is overwritten
## T is overwritten by Q
static:
doAssert f.c0 is Fp4
doAssert FT.C == F1.C
doAssert FT.C == F2.C
doAssert numDoublings >= 1
{.push checks: off.} # No OverflowError or IndexError allowed
var line {.noInit.}: Line[F2]
# First step: 0b10, T <- Q, f = 1 (mod p¹²), f *= line
# ----------------------------------------------------
T.projectiveFromAffine(Q)
# f.square() -> square(1)
line.line_double(T, P)
# Doubling steps: 0b10...00
# ----------------------------------------------------
# Process all doublings, the second is special cased
# as:
# - The first line is squared (sparse * sparse)
# - The second is (somewhat-sparse * sparse)
when numDoublings >= 2:
f.mul_sparse_sparse(line, line)
line.line_double(T, P)
f.mul(line)
for _ in 2 ..< numDoublings:
f.square()
line.line_double(T, P)
f.mul(line)
# Addition step: 0b10...01
# ------------------------------------------------
# If there was only a single doubling needed,
# we special case the addition as
# - The first line and second are sparse (sparse * sparse)
when numDoublings == 1:
# TODO: sparse * sparse
# f *= line <=> f = line for the first iteration
# With Fp2 -> Fp4 -> Fp12 towering and a M-Twist
# The line corresponds to a sparse xy000z Fp12
var line2 {.noInit.}: Line[F2]
line2.line_add(T, Q, P)
f.mul_sparse_sparse(line, line2)
else:
line.line_add(T, Q, P)
f.mul(line)
{.pop.} # No OverflowError or IndexError allowed
func miller_accum_double_then_add*[FT, F1, F2](
f: var FT,
T: var ECP_ShortW_Prj[F2, OnTwist],
Q: ECP_ShortW_Aff[F2, OnTwist],
P: ECP_ShortW_Aff[F1, NotOnTwist],
numDoublings: int,
add = true
) =
## Continue a Miller Loop with
## - `numDoubling` doublings
## - 1 add
##
## f and T are updated
#
# `numDoublings` and `add` can be hardcoded at compile-time
# to prevent fault attacks.
# But fault attacks only happen on embedded
# and embedded is likely to want to minimize codesize.
# What to do?
{.push checks: off.} # No OverflowError or IndexError allowed
var line {.noInit.}: Line[F2]
for _ in 0 ..< numDoublings:
f.square()
line.line_double(T, P)
f.mul(line)
if add:
line.line_add(T, Q, P)
f.mul(line)
# Miller Loop - multi-pairing
# ----------------------------------------------------------------------------

View File

@ -17,7 +17,8 @@ import
# ############################################################
#
# Sparse Multiplication
# by lines
# by lines for embedding degree 12
# and sextic twist
#
# ############################################################
@ -38,8 +39,14 @@ import
# Craig Costello, Tanja Lange, and Michael Naehrig, 2009
# https://eprint.iacr.org/2009/615.pdf
# 𝔽p12 by line - Sparse functions
# ----------------------------------------------------------------
# ############################################################
#
# 𝔽p12 by line - 𝔽p12 quadratic over 𝔽p6
#
# ############################################################
# D-Twist
# ------------------------------------------------------------
func mul_by_line_xy0*[C: static Curve](
r: var Fp6[C],
@ -77,7 +84,7 @@ func mul_sparse_by_line_xy00z0*[C: static Curve](
static:
doAssert C.getSexticTwist() == D_Twist
doAssert f.c0.typeof is Fp6, "This assumes 𝔽p12 as a quadratic extension of 𝔽p6"
doAssert f.c0 is Fp6, "This assumes 𝔽p12 as a quadratic extension of 𝔽p6"
var
v0 {.noInit.}: Fp6[C]
@ -101,6 +108,15 @@ func mul_sparse_by_line_xy00z0*[C: static Curve](
v3.c2.sum(v0.c2, v1.c1)
f.c0 = v3
# ############################################################
#
# 𝔽p12 by line - 𝔽p12 cubic over 𝔽p4
#
# ############################################################
# D-Twist
# ------------------------------------------------------------
func mul_sparse_by_line_xyz000*[C: static Curve](
f: var Fp12[C], l: Line[Fp2[C]]) =
## Sparse multiplication of an 𝔽p12 element
@ -110,7 +126,7 @@ func mul_sparse_by_line_xyz000*[C: static Curve](
static:
doAssert C.getSexticTwist() == D_Twist
doAssert f.c0.typeof is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
doAssert f.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
# In the following equations (taken from cubic extension implementation)
# a = f
@ -187,12 +203,132 @@ func mul_sparse_by_line_xyz000*[C: static Curve](
f2x.sum2xMod(f2x, V1)
f.c2.redc2x(f2x)
func mul_xyz000_xyz000_into_abcdefghij00*[C: static Curve](f: var Fp12[C], l0, l1: Line[Fp2[C]]) =
## Multiply 2 lines together
## The result is sparse in f.c1.c1
# In the following equations (taken from cubic extension implementation)
# a0 = (x0, y0)
# a1 = (z0, 0)
# a2 = ( 0, 0)
# b0 = (x1, y1)
# b1 = (z1, 0)
# b2 = ( 0, 0)
#
# v0 = a0 b0 = (x0, y0).(x1, y1)
# v1 = a1 b1 = (z0, 0).(z1, 0)
# v2 = a2 b2 = ( 0, 0).( 0, 0)
#
# r0 = ξ ((a1 + a2) * (b1 + b2) - v1 - v2) + v0
# = ξ (a1 b1 + a2 b1 - v1) + v0
# = v0
# r1 = (a0 + a1) * (b0 + b1) - v0 - v1 + ξ v2
# = (a0 + a1) * (b0 + b1) - v0 - v1
# r2 = (a0 + a2) * (b0 + b2) - v0 - v2 + v1
# = a0 b0 - v0 + v1
# = v1
static:
doAssert C.getSexticTwist() == D_Twist
doAssert f.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
var V0{.noInit.}, f2x{.noInit.}: doublePrec(Fp4[C])
var V1{.noInit.}: doublePrec(Fp2[C])
V0.prod2x_disjoint(l0.x, l0.y, l1.x, l1.y) # a0 b0 = (x0, y0).(x1, y1)
V1.prod2x(l0.z, l1.z) # a1 b1 = (z0, 0).(z1, 0)
# r1 = (a0 + a1) * (b0 + b1) - v0 - v1
f.c1.c0.sum(l0.x, l0.z) # x0 + z0
f.c1.c1.sum(l1.x, l1.z) # x1 + z1
f2x.prod2x_disjoint(f.c1.c0, l0.y, f.c1.c1, l1.y) # (x0 + z0, y0)(x1 + z1, y1) = (a0 + a1) * (b0 + b1)
f2x.diff2xMod(f2x, V0)
f2x.c0.diff2xMod(f2x.c0, V1)
f.c1.redc2x(f2x)
# r0 = v0
f.c0.redc2x(V0)
# r2 = v1
f.c2.c0.redc2x(V1)
f.c2.c1.setZero()
func mul_sparse_by_abcdefghij00*[C: static Curve](
a: var Fp12[C], b: Fp12[C]) =
## Sparse multiplication of an 𝔽p12 element
## by a sparse 𝔽p12 element abcdefghij00
## with each representing 𝔽p2 coordinate
static:
doAssert C.getSexticTwist() == D_Twist
doAssert a.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
# In the following equations (taken from cubic extension implementation)
# b0 = (b00, b01)
# b1 = (b10, b11)
# b2 = (b20, 0)
#
# v0 = a0 b0 = (f00, f01).(b00, b01)
# v1 = a1 b1 = (f10, f11).(b10, b11)
# v2 = a2 b2 = (f20, f21).(b20, 0)
#
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁ + β v₂
# r₂ = (a₀ + a₂) * (b₀ + b₂) - v₀ - v₂ + v₁
var V0 {.noInit.}, V1 {.noInit.}, V2 {.noinit.}: doublePrec(Fp4[C])
var t0 {.noInit.}, t1 {.noInit.}: Fp4[C]
var f2x{.noInit.}, g2x {.noinit.}: doublePrec(Fp4[C])
V0.prod2x(a.c0, b.c0)
V1.prod2x(a.c1, b.c1)
V2.mul2x_sparse_by_x0(a.c2, b.c2)
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
t0.sum(a.c1, a.c2)
t1.c0.sum(b.c1.c0, b.c2.c0) # b₂ = (b20, 0)
f2x.prod2x_disjoint(t0, t1.c0, b.c1.c1) # (a₁ + a₂).(b₁ + b₂)
f2x.diff2xMod(f2x, V1)
f2x.diff2xMod(f2x, V2)
f2x.prod2x(f2x, NonResidue)
f2x.sum2xMod(f2x, V0)
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁
t0.sum(a.c0, a.c1)
t1.sum(b.c0, b.c1)
g2x.prod2x(t0, t1)
g2x.diff2xMod(g2x, V0)
g2x.diff2xMod(g2x, V1)
# r₂ = (a₀ + a₂) and (b₀ + b₂)
t0.sum(a.c0, a.c2)
t1.c0.sum(b.c0.c0, b.c2.c0) # b₂ = (b20, 0)
# Now we are aliasing free
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
a.c0.redc2x(f2x)
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁ + β v₂
f2x.prod2x(V2, NonResidue)
g2x.sum2xMod(g2x, f2x)
a.c1.redc2x(g2x)
# r₂ = (a₀ + a₂) * (b₀ + b₂) - v₀ - v₂ + v₁
f2x.prod2x_disjoint(t0, t1.c0, b.c0.c1)
f2x.diff2xMod(f2x, V0)
f2x.diff2xMod(f2x, V2)
f2x.sum2xMod(f2x, V1)
a.c2.redc2x(f2x)
# M-Twist
# ------------------------------------------------------------
func mul_sparse_by_line_xy000z*[C: static Curve](
f: var Fp12[C], l: Line[Fp2[C]]) =
static:
doAssert C.getSexticTwist() == M_Twist
doAssert f.c0.typeof is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
doAssert f.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
# In the following equations (taken from cubic extension implementation)
# a = f
@ -271,6 +407,127 @@ func mul_sparse_by_line_xy000z*[C: static Curve](
f2x.sum2xMod(f2x, V2)
f.c1.redc2x(f2x)
func mul_xy000z_xy000z_into_abcd00efghij*[C: static Curve](f: var Fp12[C], l0, l1: Line[Fp2[C]]) =
## Multiply 2 lines together
## The result is sparse in f.c1.c0
# In the following equations (taken from cubic extension implementation)
# a0 = (x0, y0)
# a1 = ( 0, 0)
# a2 = ( 0, z0)
# b0 = (x1, y1)
# b1 = ( 0, 0)
# b2 = ( 0, z1)
#
# v0 = a0 b0 = (x0, y0).(x1, y1)
# v1 = a1 b1 = ( 0, 0).( 0, 0)
# v2 = a2 b2 = ( 0, z0).( 0, z1)
#
# r0 = ξ ((a1 + a2) * (b1 + b2) - v1 - v2) + v0
# = ξ (a1 b2 + a2 b2 - v2) + v0
# = v0
# r1 = (a0 + a1) * (b0 + b1) - v0 - v1 + ξ v2
# = a0 b0 + a1 b0 - v0 + ξ v2
# = ξ v2
# r2 = (a0 + a2) * (b0 + b2) - v0 - v2 + v1
# = (a0 + a2) * (b0 + b2) - v0 - v2
static:
doAssert C.getSexticTwist() == M_Twist
doAssert f.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
var V0{.noInit.}, f2x{.noInit.}: doublePrec(Fp4[C])
var V2{.noInit.}: doublePrec(Fp2[C])
V0.prod2x_disjoint(l0.x, l0.y, l1.x, l1.y) # a0 b0 = (x0, y0).(x1, y1)
V2.prod2x(l0.z, l1.z) # a2 b2 = ( 0, z0).( 0, z1)
V2.prod2x(V2, NonResidue)
# r2 = (a0 + a2) * (b0 + b2) - v0 - v2
f.c2.c0.sum(l0.y, l0.z) # y0 + z0
f.c2.c1.sum(l1.y, l1.z) # y1 + z1
f2x.prod2x_disjoint(l0.x, f.c2.c0, l1.x, f.c2.c1) # (x0, y0 + z0).(x1, y1 + z1) = (a0 + a2) * (b0 + b2)
f2x.diff2xMod(f2x, V0) # (a0 + a2) * (b0 + b2) - v0
f2x.c0.diff2xMod(f2x.c0, V2) # (a0 + a2) * (b0 + b2) - v0 - v2
f.c2.redc2x(f2x)
# r1 = ξ v2
f.c1.c1.redc2x(V2)
f.c1.c0.setZero()
# r0 = v0
f.c0.redc2x(V0)
func mul_sparse_by_abcd00efghij*[C: static Curve](
a: var Fp12[C], b: Fp12[C]) =
## Sparse multiplication of an 𝔽p12 element
## by a sparse 𝔽p12 element abcd00efghij
## with each representing 𝔽p2 coordinate
static:
doAssert C.getSexticTwist() == M_Twist
doAssert a.c0 is Fp4, "This assumes 𝔽p12 as a cubic extension of 𝔽p4"
# In the following equations (taken from cubic extension implementation)
# b0 = (b00, b01)
# b1 = ( 0, b11)
# b2 = (b20, b21)
#
# v0 = a0 b0 = (f00, f01).(b00, b01)
# v1 = a1 b1 = (f10, f11).( 0, b11)
# v2 = a2 b2 = (f20, f21).(b20, b21)
#
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁ + β v₂
# r₂ = (a₀ + a₂) * (b₀ + b₂) - v₀ - v₂ + v₁
var V0 {.noInit.}, V1 {.noInit.}, V2 {.noinit.}: doublePrec(Fp4[C])
var t0 {.noInit.}, t1 {.noInit.}: Fp4[C]
var f2x{.noInit.}, g2x {.noinit.}: doublePrec(Fp4[C])
V0.prod2x(a.c0, b.c0)
V1.mul2x_sparse_by_0y(a.c1, b.c1)
V2.prod2x(a.c2, b.c2)
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
t0.sum(a.c1, a.c2)
t1.c1.sum(b.c1.c1, b.c2.c1) # b₁ = ( 0, b11)
f2x.prod2x_disjoint(t0, b.c2.c0, t1.c1) # (a₁ + a₂).(b₁ + b₂)
f2x.diff2xMod(f2x, V1)
f2x.diff2xMod(f2x, V2)
f2x.prod2x(f2x, NonResidue)
f2x.sum2xMod(f2x, V0)
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁
t0.sum(a.c0, a.c1)
t1.c1.sum(b.c0.c1, b.c1.c1) # b₁ = ( 0, b11)
g2x.prod2x_disjoint(t0, b.c0.c0, t1.c1) # (a₀ + a₁).(b₀ + b₁)
g2x.diff2xMod(g2x, V0)
g2x.diff2xMod(g2x, V1)
# r₂ = (a₀ + a₂) and (b₀ + b₂)
t0.sum(a.c0, a.c2)
t1.sum(b.c0, b.c2)
# Now we are aliasing free
# r₀ = ξ ((a₁ + a₂)(b₁ + b₂) - v₁ - v₂) + v₀
a.c0.redc2x(f2x)
# r₁ = (a₀ + a₁) * (b₀ + b₁) - v₀ - v₁ + β v₂
f2x.prod2x(V2, NonResidue)
g2x.sum2xMod(g2x, f2x)
a.c1.redc2x(g2x)
# r₂ = (a₀ + a₂) * (b₀ + b₂) - v₀ - v₂ + v₁
f2x.prod2x(t0, t1)
f2x.diff2xMod(f2x, V0)
f2x.diff2xMod(f2x, V2)
f2x.sum2xMod(f2x, V1)
a.c2.redc2x(f2x)
# Dispatch
# ------------------------------------------------------------
func mul*[C](f: var Fp12[C], line: Line[Fp2[C]]) {.inline.} =
when C.getSexticTwist() == D_Twist:
f.mul_sparse_by_line_xyz000(line)
@ -278,3 +535,11 @@ func mul*[C](f: var Fp12[C], line: Line[Fp2[C]]) {.inline.} =
f.mul_sparse_by_line_xy000z(line)
else:
{.error: "A line function assumes that the curve has a twist".}
func mul_sparse_sparse*[C](f: var Fp12[C], line0, line1: Line[Fp2[C]]) {.inline.} =
when C.getSexticTwist() == D_Twist:
f.mul_xyz000_xyz000_into_abcdefghij00(line0, line1)
elif C.getSexticTwist() == M_Twist:
f.mul_xy000z_xy000z_into_abcd00efghij(line0, line1)
else:
{.error: "A line function assumes that the curve has a twist".}

View File

@ -17,7 +17,8 @@ import
# ############################################################
#
# Sparse Multiplication
# by lines
# by lines for embedding degree 6
# and sextic twist
#
# ############################################################

View File

@ -19,6 +19,8 @@ import
./lines_common,
./miller_loops
export zoo_pairings # generic sandwich https://github.com/nim-lang/Nim/issues/11225
# ############################################################
#
# Optimal ATE pairing for
@ -74,18 +76,14 @@ func finalExpGeneric[C: static Curve](f: var Fp12[C]) =
func pairing_bls12_reference*[C](
gt: var Fp12[C],
P: ECP_ShortW_Prj[Fp[C], NotOnTwist],
Q: ECP_ShortW_Prj[Fp2[C], OnTwist]) =
P: ECP_ShortW_Aff[Fp[C], NotOnTwist],
Q: ECP_ShortW_Aff[Fp2[C], OnTwist]) =
## Compute the optimal Ate Pairing for BLS12 curves
## Input: P ∈ G1, Q ∈ G2
## Output: e(P, Q) ∈ Gt
##
## Reference implementation
var Paff {.noInit.}: ECP_ShortW_Aff[Fp[C], NotOnTwist]
var Qaff {.noInit.}: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
gt.millerLoopGenericBLS12(Paff, Qaff)
gt.millerLoopGenericBLS12(P, Q)
gt.finalExpGeneric()
# Optimized pairing implementation
@ -148,15 +146,11 @@ func finalExpHard_BLS12*[C](f: var Fp12[C]) {.meter.} =
func pairing_bls12*[C](
gt: var Fp12[C],
P: ECP_ShortW_Prj[Fp[C], NotOnTwist],
Q: ECP_ShortW_Prj[Fp2[C], OnTwist]) {.meter.} =
P: ECP_ShortW_Aff[Fp[C], NotOnTwist],
Q: ECP_ShortW_Aff[Fp2[C], OnTwist]) {.meter.} =
## Compute the optimal Ate Pairing for BLS12 curves
## Input: P ∈ G1, Q ∈ G2
## Output: e(P, Q) ∈ Gt
var Paff {.noInit.}: ECP_ShortW_Aff[Fp[C], NotOnTwist]
var Qaff {.noInit.}: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
gt.millerLoopGenericBLS12(Paff, Qaff)
gt.millerLoopAddchain(Q, P)
gt.finalExpEasy()
gt.finalExpHard_BLS12()

View File

@ -0,0 +1,338 @@
# 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
../config/[common, curves, type_ff],
../towers,
../elliptic/[
ec_shortweierstrass_affine,
ec_shortweierstrass_projective
],
../curves/zoo_pairings,
./lines_projective, ./mul_fp12_by_lines,
./miller_loops
# ############################################################
#
# Optimal ATE pairing for
# BLS12-381
#
# ############################################################
#
# - Software Implementation, Algorithm 11.2 & 11.3
# Aranha, Dominguez Perez, A. Mrabet, Schwabe,
# Guide to Pairing-Based Cryptography, 2015
#
# - Physical Attacks,
# N. El Mrabet, Goubin, Guilley, Fournier, Jauvart, Moreau, Rauzy, Rondepierre,
# Guide to Pairing-Based Cryptography, 2015
#
# - Pairing Implementation Revisited
# Mike Scott, 2019
# https://eprint.iacr.org/2019/077.pdf
#
# Fault attacks:
# To limite exposure to some fault attacks (flipping bits with a laser on embedded):
# - changing the number of Miller loop iterations
# - flipping the bits in the Miller loop
# we hardcode unrolled addition chains.
# This should also contribute to performance.
#
# Multi-pairing discussion:
# Aranha & Scott proposes 2 different approaches for multi-pairing.
#
# -----
# Scott
#
# Algorithm 2: Calculate and store line functions for BLS12 curve
# Input: Q ∈ G2, P ∈ G1 , curve parameter u
# Output: An array g of blog2(u)c line functions ∈ Fp12
# 1 T ← Q
# 2 for i ← ceil(log2(u)) 1 to 0 do
# 3 g[i] ← lT,T(P), T ← 2T
# 4 if ui = 1 then
# 5 g[i] ← g[i].lT,Q(P), T ← T + Q
# 6 return g
#
# And to accumulate lines from a new (P, Q) tuple of points
#
# Algorithm 4: Accumulate another set of line functions into g
# Input: The array g, Qj ∈ G2 , Pj ∈ G1 , curve parameter u
# Output: Updated array g of ceil(log2(u)) line functions ∈ Fp12
# 1 T ← Qj
# 2 for i ← blog2 (u)c 1 to 0 do
# 3 t ← lT,T (Pj), T ← 2T
# 4 if ui = 1 then
# 5 t ← t.lT,Qj (Pj), T ← T + Qj
# 6 g[i] ← g[i].t
# 7 return g
#
# ------
# Aranha
#
# Algorithm 11.2 Explicit multipairing version of Algorithm 11.1.
# (we extract the Miller Loop part only)
# Input : P1 , P2 , . . . Pn ∈ G1 ,
# Q1 , Q2, . . . Qn ∈ G2
# Output: (we focus on the Miller Loop)
#
# Write l in binary form, l = sum(0 ..< m-1)
# f ← 1, l ← abs(AteParam)
# for j ← 1 to n do
# Tj ← Qj
# end
#
# for i = m-2 down to 0 do
# f ← f²
# for j ← 1 to n do
# f ← f gTj,Tj(Pj), Tj ← [2]Tj
# if li = 1 then
# f ← f gTj,Qj(Pj), Tj ← Tj + Qj
# end
# end
# end
#
# -----
# Assuming we have N tuples (Pj, Qj) of points j in 0 ..< N
# and I operations to do in our Miller loop:
# - I = HammingWeight(AteParam) + Bitwidth(AteParam)
# - HammingWeight(AteParam) corresponds to line additions
# - Bitwidth(AteParam) corresponds to line doublings
#
# Scott approach is to have:
# - I Fp12 accumulators `g`
# - 1 G2 accumulator `T`
# and then accumulating each (Pj, Qj) into their corresponding `g` accumulator.
#
# Aranha approach is to have:
# - 1 Fp12 accumulator `f`
# - N G2 accumulators `T`
# and accumulate N points per I.
#
# Scott approach is fully "online"/"streaming",
# while Aranha's saves space.
# For BLS12_381,
# I = 68 hence we would need 68*12*48 = 39168 bytes (381-bit needs 48 bytes)
# G2 has size 3*2*48 = 288 bytes (3 proj coordinates on Fp2)
# and we choose N (which can be 1 for single pairing or reverting to Scott approach).
#
# In actual use, "streaming pairings" are not used, pairings to compute are receive
# by batch, for example for blockchain you receive a batch of N blocks to verify from one peer.
# Furthermore, 39kB would be over L1 cache size and incurs cache misses.
# Additionally Aranha approach would make it easier to batch inversions
# using Montgomery's simultaneous inversion technique.
# Lastly, while a higher level API will need to store N (Pj, Qj) pairs for multi-pairings
# for Aranha approach, it can decide how big N is depending on hardware and/or protocol.
#
# Regarding optimizations, as the Fp12 accumulator is dense
# and lines are sparse (xyz000 or xy000z) Scott mentions the following costs:
# - Dense-sparse is 13m
# - sparse-sparse is 6m
# - Dense-(somewhat sparse) is 17m
# Hence when accumulating lines from multiple points:
# - 2x Dense-sparse is 26m
# - sparse-sparse then Dense-(somewhat sparse) is 23m
# a 11.5% speedup
#
# We can use Aranha approach but process lines function 2-by-2 merging them
# before merging them to the dense Fp12 accumulator
# Miller Loop
# -------------------------------------------------------------------------------------------------------
{.push raises: [].}
import
strutils,
../io/io_towers
func miller_first_iter[N: static int](
f: var Fp12[BLS12_381],
Ts: var array[N, ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]],
Qs: array[N, ECP_ShortW_Aff[Fp2[BLS12_381], OnTwist]],
Ps: array[N, ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]]
) =
## Start a Miller Loop
## This means
## - 1 doubling
## - 1 add
##
## f is overwritten
## Ts are overwritten by Qs
static:
doAssert N >= 1
doAssert f.c0 is Fp4
{.push checks: off.} # No OverflowError or IndexError allowed
var line {.noInit.}: Line[Fp2[BLS12_381]]
# First step: T <- Q, f = 1 (mod p¹²), f *= line
# ----------------------------------------------
for i in 0 ..< N:
Ts[i].projectiveFromAffine(Qs[i])
line.line_double(Ts[0], Ps[0])
# f *= line <=> f = line for the first iteration
# With Fp2 -> Fp4 -> Fp12 towering and a M-Twist
# The line corresponds to a sparse xy000z Fp12
f.c0.c0 = line.x
f.c0.c1 = line.y
f.c1.c0.setZero()
f.c1.c1.setZero()
f.c2.c0.setZero()
f.c2.c1 = line.z
when N >= 2:
line.line_double(Ts[1], Ps[1])
f.mul_sparse_by_line_xy000z(line) # TODO: sparse-sparse mul
# Sparse merge 2 by 2, starting from 2
for i in countup(2, N-1, 2):
# var f2 {.noInit.}: Fp12[BLS12_381] # TODO: sparse-sparse mul
var line2 {.noInit.}: Line[Fp2[BLS12_381]]
line.line_double(Ts[i], Ps[i])
line2.line_double(Ts[i+1], Ps[i+1])
# f2.mul_sparse_sparse(line, line2)
# f.mul_somewhat_sparse(f2)
f.mul_sparse_by_line_xy000z(line)
f.mul_sparse_by_line_xy000z(line2)
when N and 1 == 1: # N >= 2 and N is odd, there is a leftover
line.line_double(Ts[N-1], Ps[N-1])
f.mul_sparse_by_line_xy000z(line)
# 2nd step: Line addition as MSB is always 1
# ----------------------------------------------
when N >= 2: # f is dense, there are already many lines accumulated
# Sparse merge 2 by 2, starting from 0
for i in countup(0, N-1, 2):
# var f2 {.noInit.}: Fp12[BLS12_381] # TODO: sparse-sparse mul
var line2 {.noInit.}: Line[Fp2[BLS12_381]]
line.line_add(Ts[i], Qs[i], Ps[i])
line2.line_add(Ts[i+1], Qs[i+1], Ps[i+1])
# f2.mul_sparse_sparse(line, line2)
# f.mul_somewhat_sparse(f2)
f.mul_sparse_by_line_xy000z(line)
f.mul_sparse_by_line_xy000z(line2)
when N and 1 == 1: # N >= 2 and N is odd, there is a leftover
line.line_add(Ts[N-1], Qs[N-1], Ps[N-1])
f.mul_sparse_by_line_xy000z(line)
else: # N = 1, f is sparse
line.line_add(Ts[0], Qs[0], Ps[0])
# f.mul_sparse_sparse(line)
f.mul_sparse_by_line_xy000z(line)
{.pop.} # No OverflowError or IndexError allowed
func miller_accum_doublings[N: static int](
f: var Fp12[BLS12_381],
Ts: var array[N, ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]],
Ps: array[N, ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]],
numDoublings: int
) =
## Accumulate `numDoublings` Miller loop doubling steps into `f`
static: doAssert N >= 1
{.push checks: off.} # No OverflowError or IndexError allowed
var line {.noInit.}: Line[Fp2[BLS12_381]]
for _ in 0 ..< numDoublings:
f.square()
when N >= 2:
for i in countup(0, N-1, 2):
# var f2 {.noInit.}: Fp12[BLS12_381] # TODO: sparse-sparse mul
var line2 {.noInit.}: Line[Fp2[BLS12_381]]
line.line_double(Ts[i], Ps[i])
line2.line_double(Ts[i+1], Ps[i+1])
# f2.mul_sparse_sparse(line, line2)
# f.mul_somewhat_sparse(f2)
f.mul_sparse_by_line_xy000z(line)
f.mul_sparse_by_line_xy000z(line2)
when N and 1 == 1: # N >= 2 and N is odd, there is a leftover
line.line_double(Ts[N-1], Ps[N-1])
f.mul_sparse_by_line_xy000z(line)
else:
line.line_double(Ts[0], Ps[0])
f.mul_sparse_by_line_xy000z(line)
{.pop.} # No OverflowError or IndexError allowed
func miller_accum_addition[N: static int](
f: var Fp12[BLS12_381],
Ts: var array[N, ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]],
Qs: array[N, ECP_ShortW_Aff[Fp2[BLS12_381], OnTwist]],
Ps: array[N, ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]]
) =
## Accumulate a Miller loop addition step into `f`
static: doAssert N >= 1
{.push checks: off.} # No OverflowError or IndexError allowed
var line {.noInit.}: Line[Fp2[BLS12_381]]
when N >= 2:
# Sparse merge 2 by 2, starting from 0
for i in countup(0, N-1, 2):
# var f2 {.noInit.}: Fp12[BLS12_381] # TODO: sparse-sparse mul
var line2 {.noInit.}: Line[Fp2[BLS12_381]]
line.line_add(Ts[i], Qs[i], Ps[i])
line2.line_add(Ts[i+1], Qs[i+1], Ps[i+1])
# f2.mul_sparse_sparse(line, line2)
# f.mul_somewhat_sparse(f2)
f.mul_sparse_by_line_xy000z(line)
f.mul_sparse_by_line_xy000z(line2)
when N and 1 == 1: # N >= 2 and N is odd, there is a leftover
line.line_add(Ts[N-1], Qs[N-1], Ps[N-1])
f.mul_sparse_by_line_xy000z(line)
else:
line.line_add(Ts[0], Qs[0], Ps[0])
f.mul_sparse_by_line_xy000z(line)
{.pop.} # No OverflowError or IndexError allowed
func millerLoop_opt_BLS12_381*[N: static int](
f: var Fp12[BLS12_381],
Qs: array[N, ECP_ShortW_Aff[Fp2[BLS12_381], OnTwist]],
Ps: array[N, ECP_ShortW_Aff[Fp[BLS12_381], NotOnTwist]]
) {.meter.} =
## Generic Miller Loop for BLS12 curve
## Computes f{u,Q}(P) with u the BLS curve parameter
var Ts {.noInit.}: array[N, ECP_ShortW_Prj[Fp2[BLS12_381], OnTwist]]
# Ate param addition chain
# Hex: 0xd201000000010000
# Bin: 0b1101001000000001000000000000000000000000000000010000000000000000
var iter = 1'u64
f.miller_first_iter(Ts, Qs, Ps) # 0b11
f.miller_accum_doublings(Ts, Ps, 2) # 0b1100
f.miller_accum_addition(Ts, Qs, Ps) # 0b1101
f.miller_accum_doublings(Ts, Ps, 3) # 0b1101000
f.miller_accum_addition(Ts, Qs, Ps) # 0b1101001
f.miller_accum_doublings(Ts, Ps, 9) # 0b1101001000000000
f.miller_accum_addition(Ts, Qs, Ps) # 0b1101001000000001
f.miller_accum_doublings(Ts, Ps, 32) # 0b110100100000000100000000000000000000000000000000
f.miller_accum_addition(Ts, Qs, Ps) # 0b110100100000000100000000000000000000000000000001
f.miller_accum_doublings(Ts, Ps, 16) # 0b1101001000000001000000000000000000000000000000010000000000000000
# TODO: what is the threshold for Karabina's compressed squarings?

View File

@ -20,6 +20,8 @@ import
./cyclotomic_fp12,
./miller_loops
export zoo_pairings # generic sandwich https://github.com/nim-lang/Nim/issues/11225
# ############################################################
#
# Optimal ATE pairing for
@ -65,19 +67,11 @@ func millerLoopGenericBN*[C](
ate_param, ate_param_isNeg
)
# Ate pairing for BN curves need adjustment after basic Miller loop
when C.pairing(ate_param_isNeg):
T.neg()
var V {.noInit.}: typeof(Q)
V.frobenius_psi(Q)
line.line_add(T, V, P)
f.mul(line)
V.frobenius_psi(Q, 2)
V.neg()
line.line_add(T, V, P)
f.mul(line)
# Ate pairing for BN curves needs adjustment after basic Miller loop
f.millerCorrectionBN(
T, Q, P,
pairing(C, ate_param_isNeg)
)
func finalExpGeneric[C: static Curve](f: var Fp12[C]) =
## A generic and slow implementation of final exponentiation
@ -86,18 +80,14 @@ func finalExpGeneric[C: static Curve](f: var Fp12[C]) =
func pairing_bn_reference*[C](
gt: var Fp12[C],
P: ECP_ShortW_Prj[Fp[C], NotOnTwist],
Q: ECP_ShortW_Prj[Fp2[C], OnTwist]) =
P: ECP_ShortW_Aff[Fp[C], NotOnTwist],
Q: ECP_ShortW_Aff[Fp2[C], OnTwist]) =
## Compute the optimal Ate Pairing for BN curves
## Input: P ∈ G1, Q ∈ G2
## Output: e(P, Q) ∈ Gt
##
## Reference implementation
var Paff {.noInit.}: ECP_ShortW_Aff[Fp[C], NotOnTwist]
var Qaff {.noInit.}: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
gt.millerLoopGenericBN(Paff, Qaff)
gt.millerLoopGenericBN(P, Q)
gt.finalExpGeneric()
# Optimized pairing implementation
@ -157,15 +147,14 @@ func finalExpHard_BN*[C: static Curve](f: var Fp12[C]) =
func pairing_bn*[C](
gt: var Fp12[C],
P: ECP_ShortW_Prj[Fp[C], NotOnTwist],
Q: ECP_ShortW_Prj[Fp2[C], OnTwist]) =
P: ECP_ShortW_Aff[Fp[C], NotOnTwist],
Q: ECP_ShortW_Aff[Fp2[C], OnTwist]) =
## Compute the optimal Ate Pairing for BLS12 curves
## Input: P ∈ G1, Q ∈ G2
## Output: e(P, Q) ∈ Gt
var Paff {.noInit.}: ECP_ShortW_Aff[Fp[C], NotOnTwist]
var Qaff {.noInit.}: ECP_ShortW_Aff[Fp2[C], OnTwist]
Paff.affineFromProjective(P)
Qaff.affineFromProjective(Q)
gt.millerLoopGenericBN(Paff, Qaff)
when C == BN254_Nogami:
gt.millerLoopAddChain(Q, P)
else:
gt.millerLoopGenericBN(P, Q)
gt.finalExpEasy()
gt.finalExpHard_BN()

View File

@ -83,3 +83,20 @@ func log2*[T: SomeUnsignedInt](n: T): T =
else:
static: doAssert sizeof(T) <= sizeof(uint32)
T(log2Impl(uint32(n)))
func hammingWeight*(x: uint32): int {.inline.} =
## Counts the set bits in integer.
# https://graphics.stanford.edu/~seander/bithacks.html#CountBitsSetParallel
var v = x
v = v - ((v shr 1) and 0x55555555)
v = (v and 0x33333333) + ((v shr 2) and 0x33333333)
cast[int](((v + (v shr 4) and 0xF0F0F0F) * 0x1010101) shr 24)
func hammingWeight*(x: uint64): int {.inline.} =
## Counts the set bits in integer.
# https://graphics.stanford.edu/~seander/bithacks.html#CountBitsSetParallel
var v = x
v = v - ((v shr 1'u64) and 0x5555555555555555'u64)
v = (v and 0x3333333333333333'u64) + ((v shr 2'u64) and 0x3333333333333333'u64)
v = (v + (v shr 4'u64) and 0x0F0F0F0F0F0F0F0F'u64)
cast[int]((v * 0x0101010101010101'u64) shr 56'u64)

View File

@ -513,22 +513,22 @@ func square2x_complex(r: var QuadraticExt2x, a: QuadraticExt) =
func prod2x_disjoint*[Fdbl, F](
r: var QuadraticExt2x[FDbl],
a: QuadraticExt[F],
a0, a1: F,
b0, b1: F) =
## Return a * (b0, b1) in r
static: doAssert Fdbl is doublePrec(F)
var V0 {.noInit.}, V1 {.noInit.}: typeof(r.c0) # Double-precision
var t0 {.noInit.}, t1 {.noInit.}: typeof(a.c0) # Single-width
var t0 {.noInit.}, t1 {.noInit.}: typeof(a0) # Single-width
# Require 2 extra bits
V0.prod2x(a.c0, b0) # v0 = a0b0
V1.prod2x(a.c1, b1) # v1 = a1b1
V0.prod2x(a0, b0) # v0 = a0b0
V1.prod2x(a1, b1) # v1 = a1b1
when F.has1extraBit():
t0.sumUnr(a.c0, a.c1)
t0.sumUnr(a0, a1)
t1.sumUnr(b0, b1)
else:
t0.sum(a.c0, a.c1)
t0.sum(a0, a1)
t1.sum(b0, b1)
r.c1.prod2x(t0, t1) # r1 = (a0 + a1)(b0 + b1)
@ -655,12 +655,19 @@ func inv2xImpl(r: var QuadraticExt, a: QuadraticExt) =
# Dispatch
# ----------------------------------------------------------------------
func prod2x_disjoint*[Fdbl, F](
r: var QuadraticExt2x[FDbl],
a: QuadraticExt[F],
b0, b1: F) =
## Return (a0, a1) * (b0, b1) in r
r.prod2x_disjoint(a.c0, a.c1, b0, b1)
func prod2x*(r: var QuadraticExt2x, a, b: QuadraticExt) =
mixin fromComplexExtension
when a.fromComplexExtension():
r.prod2x_complex(a, b)
else:
r.prod2x_disjoint(a, b.c0, b.c1)
r.prod2x_disjoint(a.c0, a.c1, b.c0, b.c1)
func square2x*(r: var QuadraticExt2x, a: QuadraticExt) =
mixin fromComplexExtension

View File

@ -166,6 +166,7 @@ The optimizations can be of algebraic, algorithmic or "implementation details" n
- [x] M-Twist
- [x] Fused line add + elliptic curve add
- [x] Fused line double + elliptic curve double
- [x] 6-way sparse multiplication line * Gₜ element
- [ ] Jacobian projective coordinates
- [ ] D-Twist
- [ ] Fused line add + elliptic curve add
@ -173,18 +174,15 @@ The optimizations can be of algebraic, algorithmic or "implementation details" n
- [ ] M-Twist
- [ ] Fused line add + elliptic curve add
- [ ] Fused line double + elliptic curve double
- [x] Sparse multiplication line * Gₜ element
- [x] D-Twist
- [x] 6-way sparse
- [ ] Pseudo 8-sparse
- [x] M-Twist
- [x] 6-way sparse
- [ ] Pseudo 8-sparse
- [x] 6-way sparse multiplication line * Gₜ element
- [ ] Affine coordinates
- [ ] 7-way sparse multiplication line * Gₜ element
- [ ] Pseudo-8 sparse multiplication line * Gₜ element
- Miller Loop
- [x] NAF recoding
- [ ] Quadruple-and-add and Octuple-and-add
- [ ] addition chain
- [x] addition chains
- Final exponentiation
- [x] Cyclotomic squaring

View File

@ -57,7 +57,7 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
for _ in 0 ..< Iters:
let a = rng.random_elem(Fp4[C], gen)
let y = rng.random_elem(Fp2[C], gen)
let b = Fp4[C](c1: y)
let b = Fp4[C](coords: [Fp2[C](), y])
var r {.noInit.}, r2 {.noInit.}: Fp4[C]
@ -76,7 +76,7 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
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)
let b = Fp6[C](coords: [Fp2[C](), y, Fp2[C]()])
var r {.noInit.}, r2 {.noInit.}: Fp6[C]
@ -91,12 +91,12 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
test_fp6_0y0(curve, gen = Long01Sequence)
test "Dense 𝔽p6 by Sparse xy0":
proc test_fp6_0y0(C: static Curve, gen: static RandomGen) =
proc test_fp6_xy0(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 b = Fp6[C](coords: [x, y, Fp2[C]()])
let line = Line[Fp2[C]](x: x, y: y)
var r {.noInit.}, r2 {.noInit.}: Fp6[C]
@ -107,9 +107,9 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
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_fp6_xy0(curve, gen = Uniform)
test_fp6_xy0(curve, gen = HighHammingWeight)
test_fp6_xy0(curve, gen = Long01Sequence)
when Fp12[BN254_Snarks]().c0.typeof is Fp6:
test "Sparse 𝔽p12/𝔽p6 resulting from xy00z0 line function":
@ -124,8 +124,8 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
let line = Line[Fp2[C]](x: x, y: y, z: z)
let b = Fp12[C](
c0: Fp6[C](c0: x, c1: y),
c1: Fp6[C](c1: z)
c0: Fp6[C](coords: [ x, y, Fp2[C]()]),
c1: Fp6[C](coords: [Fp2[C](), z, Fp2[C]()])
)
a *= b
@ -150,7 +150,7 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
let line = Line[Fp2[C]](x: x, y: y, z: z)
let b = Fp12[C](
c0: Fp6[C](c0: x, c1: y, c2: z)
c0: Fp6[C](coords: [x, y, z])
)
a *= b
@ -178,9 +178,11 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
let line = Line[Fp2[C]](x: x, y: y, z: z)
let b = Fp12[C](
c0: Fp4[C](c0: x, c1: y),
# c1
c2: Fp4[C]( c1: z),
coords: [
Fp4[C](coords: [x, y]),
Fp4[C](),
Fp4[C](coords: [Fp2[C](), z])
]
)
a *= b
@ -194,7 +196,7 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
test_fp12_xy000z(curve, gen = Long01Sequence)
test "Sparse 𝔽p12/𝔽p4 resulting from xyz000 line function (D-twist only)":
proc test_fp12_xy000z(C: static Curve, gen: static RandomGen) =
proc test_fp12_xyz000(C: static Curve, gen: static RandomGen) =
when C.getSexticTwist() == D_Twist:
for _ in 0 ..< Iters:
var a = rng.random_elem(Fp12[C], gen)
@ -206,9 +208,11 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
let line = Line[Fp2[C]](x: x, y: y, z: z)
let b = Fp12[C](
c0: Fp4[C](c0: x, c1: y),
c1: Fp4[C](c0: z ),
# c2:
coords: [
Fp4[C](coords: [x, y]),
Fp4[C](coords: [z, Fp2[C]()]),
Fp4[C]()
]
)
a *= b
@ -217,6 +221,179 @@ suite "Pairing - Sparse 𝔽p12 multiplication by line function is consistent wi
check: bool(a == a2)
staticFor(curve, TestCurves):
test_fp12_xy000z(curve, gen = Uniform)
test_fp12_xy000z(curve, gen = HighHammingWeight)
test_fp12_xy000z(curve, gen = Long01Sequence)
test_fp12_xyz000(curve, gen = Uniform)
test_fp12_xyz000(curve, gen = HighHammingWeight)
test_fp12_xyz000(curve, gen = Long01Sequence)
test "Somewhat-sparse 𝔽p12/𝔽p4 resulting from xy000z*xy000z line functions (M-twist only)":
proc test_fp12_xy000z_xy000z(C: static Curve, gen: static RandomGen) =
when C.getSexticTwist() == M_Twist:
for _ in 0 ..< Iters:
var x0 = rng.random_elem(Fp2[C], gen)
var y0 = rng.random_elem(Fp2[C], gen)
var z0 = rng.random_elem(Fp2[C], gen)
let line0 = Line[Fp2[C]](x: x0, y: y0, z: z0)
let f0 = Fp12[C](
coords: [
Fp4[C](coords: [x0, y0]),
Fp4[C](),
Fp4[C](coords: [Fp2[C](), z0])
]
)
var x1 = rng.random_elem(Fp2[C], gen)
var y1 = rng.random_elem(Fp2[C], gen)
var z1 = rng.random_elem(Fp2[C], gen)
let line1 = Line[Fp2[C]](x: x1, y: y1, z: z1)
let f1 = Fp12[C](
coords: [
Fp4[C](coords: [x1, y1]),
Fp4[C](),
Fp4[C](coords: [Fp2[C](), z1])
]
)
var r: Fp12[C]
r.prod(f0, f1)
var rl: Fp12[C]
rl.mul_xy000z_xy000z_into_abcd00efghij(line0, line1)
check: bool(r == rl)
test "Somewhat-sparse 𝔽p12/𝔽p4 resulting from xyz000*xyz000 line functions (D-twist only)":
proc test_fp12_xyz000_xyz000(C: static Curve, gen: static RandomGen) =
when C.getSexticTwist() == D_Twist:
for _ in 0 ..< Iters:
var x0 = rng.random_elem(Fp2[C], gen)
var y0 = rng.random_elem(Fp2[C], gen)
var z0 = rng.random_elem(Fp2[C], gen)
let line0 = Line[Fp2[C]](x: x0, y: y0, z: z0)
let f0 = Fp12[C](
coords: [
Fp4[C](coords: [x0, y0]),
Fp4[C](coords: [z0, Fp2[C]()]),
Fp4[C]()
]
)
var x1 = rng.random_elem(Fp2[C], gen)
var y1 = rng.random_elem(Fp2[C], gen)
var z1 = rng.random_elem(Fp2[C], gen)
let line1 = Line[Fp2[C]](x: x1, y: y1, z: z1)
let f1 = Fp12[C](
coords: [
Fp4[C](coords: [x1, y1]),
Fp4[C](coords: [z1, Fp2[C]()]),
Fp4[C]()
]
)
var r: Fp12[C]
r.prod(f0, f1)
var rl: Fp12[C]
rl.mul_xyz000_xyz000_into_abcdefghij00(line0, line1)
check: bool(r == rl)
staticFor(curve, TestCurves):
test_fp12_xyz000_xyz000(curve, gen = Uniform)
test_fp12_xyz000_xyz000(curve, gen = HighHammingWeight)
test_fp12_xyz000_xyz000(curve, gen = Long01Sequence)
test "Somewhat-sparse 𝔽p12/𝔽p4 mul by the product (xyz000*xyz000) of line functions (D-twist only)":
proc test_fp12_abcdefghij00(C: static Curve, gen: static RandomGen) =
when C.getSexticTwist() == D_Twist:
for _ in 0 ..< Iters:
var x0 = rng.random_elem(Fp2[C], gen)
var y0 = rng.random_elem(Fp2[C], gen)
var z0 = rng.random_elem(Fp2[C], gen)
let line0 = Line[Fp2[C]](x: x0, y: y0, z: z0)
let f0 = Fp12[C](
coords: [
Fp4[C](coords: [x0, y0]),
Fp4[C](coords: [z0, Fp2[C]()]),
Fp4[C]()
]
)
var x1 = rng.random_elem(Fp2[C], gen)
var y1 = rng.random_elem(Fp2[C], gen)
var z1 = rng.random_elem(Fp2[C], gen)
let line1 = Line[Fp2[C]](x: x1, y: y1, z: z1)
let f1 = Fp12[C](
coords: [
Fp4[C](coords: [x1, y1]),
Fp4[C](coords: [z1, Fp2[C]()]),
Fp4[C]()
]
)
var rl: Fp12[C]
rl.mul_xyz000_xyz000_into_abcdefghij00(line0, line1)
var f = rng.random_elem(Fp12[C], gen)
var f2 = f
f *= rl
f2.mul_sparse_by_abcdefghij00(rl)
check: bool(f == f2)
staticFor(curve, TestCurves):
test_fp12_abcdefghij00(curve, gen = Uniform)
test_fp12_abcdefghij00(curve, gen = HighHammingWeight)
test_fp12_abcdefghij00(curve, gen = Long01Sequence)
test "Somewhat-sparse 𝔽p12/𝔽p4 mul by the product (xy000z*xy000z) of line functions (M-twist only)":
proc test_fp12_abcd00efghij(C: static Curve, gen: static RandomGen) =
when C.getSexticTwist() == M_Twist:
for _ in 0 ..< Iters:
var x0 = rng.random_elem(Fp2[C], gen)
var y0 = rng.random_elem(Fp2[C], gen)
var z0 = rng.random_elem(Fp2[C], gen)
let line0 = Line[Fp2[C]](x: x0, y: y0, z: z0)
let f0 = Fp12[C](
coords: [
Fp4[C](coords: [x0, y0]),
Fp4[C](),
Fp4[C](coords: [Fp2[C](), z0])
]
)
var x1 = rng.random_elem(Fp2[C], gen)
var y1 = rng.random_elem(Fp2[C], gen)
var z1 = rng.random_elem(Fp2[C], gen)
let line1 = Line[Fp2[C]](x: x1, y: y1, z: z1)
let f1 = Fp12[C](
coords: [
Fp4[C](coords: [x1, y1]),
Fp4[C](),
Fp4[C](coords: [Fp2[C](), z1])
]
)
var rl: Fp12[C]
rl.mul_xy000z_xy000z_into_abcd00efghij(line0, line1)
var f = rng.random_elem(Fp12[C], gen)
var f2 = f
f *= rl
f2.mul_sparse_by_abcd00efghij(rl)
check: bool(f == f2)
staticFor(curve, TestCurves):
test_fp12_abcd00efghij(curve, gen = Uniform)
test_fp12_abcd00efghij(curve, gen = HighHammingWeight)
test_fp12_abcd00efghij(curve, gen = Long01Sequence)

View File

@ -14,14 +14,15 @@ import
../constantine/[arithmetic, primitives],
../constantine/towers,
../constantine/config/curves,
../constantine/elliptic/ec_shortweierstrass_projective,
../constantine/elliptic/[ec_shortweierstrass_affine, ec_shortweierstrass_projective],
../constantine/hash_to_curve/cofactors,
# Test utilities
../helpers/prng_unsafe
export
prng_unsafe, times, unittest,
ec_shortweierstrass_projective, arithmetic, towers,
ec_shortweierstrass_affine, ec_shortweierstrass_projective,
arithmetic, towers,
primitives
type
@ -30,29 +31,32 @@ type
HighHammingWeight
Long01Sequence
template affineType[F; Tw: static Twisted](
ec: ECP_ShortW_Prj[F, Tw]): type =
ECP_ShortW_Aff[F, Tw]
func clearCofactorReference[F; Tw: static Twisted](
ec: var ECP_ShortW_Aff[F, Tw]) =
# For now we don't have any affine operation defined
var t {.noInit.}: ECP_ShortW_Prj[F, Tw]
t.projectiveFromAffine(ec)
t.clearCofactorReference()
ec.affineFromProjective(t)
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()
if gen == Uniform:
result = rng.random_unsafe(EC)
result.clearCofactorReference()
elif gen == HighHammingWeight:
result = rng.random_highHammingWeight(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()
result = rng.random_long01Seq(EC)
result.clearCofactorReference()
template runPairingTests*(Iters: static int, C: static Curve, G1, G2, GT: typedesc, pairing_fn: untyped): untyped {.dirty.}=
bind affineType
var rng: RngState
let timeseed = uint32(toUnix(getTime()) and (1'i64 shl 32 - 1)) # unixTime mod 2^32
seed(rng, timeseed)
@ -71,10 +75,18 @@ template runPairingTests*(Iters: static int, C: static Curve, G1, G2, GT: typede
P2.double(P)
Q2.double(Q)
r.pairing_fn(P, Q)
var Pa {.noInit.}, Pa2 {.noInit.}: affineType(P)
var Qa {.noInit.}, Qa2 {.noInit.}: affineType(Q)
Pa.affineFromProjective(P)
Pa2.affineFromProjective(P2)
Qa.affineFromProjective(Q)
Qa2.affineFromProjective(Q2)
r.pairing_fn(Pa, Qa)
r.square()
r2.pairing_fn(P2, Q)
r3.pairing_fn(P, Q2)
r2.pairing_fn(Pa2, Qa)
r3.pairing_fn(Pa, Qa2)
doAssert bool(not r.isZero())
doAssert bool(not r.isOne())
@ -85,8 +97,5 @@ template runPairingTests*(Iters: static int, C: static Curve, G1, G2, GT: typede
suite "Pairing - Optimal Ate on " & $C & " [" & $WordBitwidth & "-bit mode]":
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)