drandcurve.gno
11.76 Kb · 313 lines
1package bls254
2
3// RFC 9380 hash_to_field and the Shallue-van de Woestijne map to BN254 G1, as
4// used by drand's bls-bn254-unchained-on-g1 scheme.
5//
6// This is the part of the drand path that differs most from the existing SPLASH
7// beacon. SPLASH maps to G1 by try-and-increment: hash, treat the digest as x,
8// bump a counter until x^3+3 is a square. drand uses RFC 9380 properly --
9// expand_message_xmd to 96 bytes, split into two 48-byte field elements, SVDW
10// each onto the curve, add the two points. The two constructions produce
11// completely different points for the same message and are not interchangeable.
12//
13// BN254 G1 is y^2 = x^3 + 3 over Fp, so A = 0 and B = 3. The cofactor is 1, so
14// there is no cofactor clearing: every point on the curve is in G1.
15
16import (
17 "crypto/bn254"
18 "crypto/modexp"
19)
20
21// ---- SVDW constants ----
22//
23// RFC 9380 section 6.6.1 defines, for the curve y^2 = x^3 + A*x + B and a
24// non-square-ish parameter Z:
25//
26// c1 = g(Z)
27// c2 = -Z / 2
28// c3 = sqrt(-g(Z) * (3*Z^2 + 4*A)) with sgn0(c3) == 0
29// c4 = -4 * g(Z) / (3*Z^2 + 4*A)
30//
31// With A = 0, B = 3, Z = 1 that collapses to
32//
33// g(Z) = Z^3 + A*Z + B = 1 + 3 = 4
34// c1 = 4
35// c2 = -1/2 mod p, which is (p-1)/2 because p is odd
36// c3 = sqrt(-4 * 3) = sqrt(-12) mod p, negated if odd so that sgn0(c3) == 0
37// c4 = -16/3 mod p
38//
39// I derived the four values from p with python3 rather than copying them:
40//
41// p = 0x30644e72e131a029b85045b68181585d97816a916871ca8d3c208c16d87cfd47
42// Z, A, B = 1, 0, 3
43// gZ = (Z**3 + A*Z + B) % p
44// c1 = gZ
45// c2 = (-Z * pow(2, -1, p)) % p
46// d = (-gZ * (3*Z*Z + 4*A)) % p # = -12 mod p
47// c3 = pow(d, (p+1)//4, p) # p % 4 == 3, so this is a sqrt
48// assert c3*c3 % p == d
49// if c3 % 2 == 1: c3 = p - c3 # enforce sgn0(c3) == 0
50// c4 = (-4*gZ * pow(3*Z*Z + 4*A, -1, p)) % p
51//
52// and then checked them two ways before writing them down:
53//
54// - against drand/kyber v1.3.2 pairing/bn254/constants.go, whose c1..c4 are
55// stored in Montgomery form; decoding them by multiplying by R^-1 mod p
56// (R = 2^256) reproduced all four exactly.
57// - against the Solidity reference randa-mu/bls-solidity src/libraries/BLS.sol,
58// whose C1/C2/C3/C4 literals are textually identical to the hex below.
59//
60// TestSVDWConstants re-derives all four inside the VM from pBytes alone, using
61// only modexp and the limb helpers, so none of the above has to be trusted at
62// read time.
63const (
64 // c1 = g(Z) = 4
65 c1Hex = "0000000000000000000000000000000000000000000000000000000000000004"
66 // c2 = -Z/2 = (p-1)/2
67 c2Hex = "183227397098d014dc2822db40c0ac2ecbc0b548b438e5469e10460b6c3e7ea3"
68 // c3 = sqrt(-12) mod p, the even root
69 c3Hex = "00000000000000016789af3a83522eb353c98fc6b36d713d5d8d1cc5dffffffa"
70 // c4 = -16/3 mod p
71 c4Hex = "10216f7ba065e00de81ac1e7808072c9dd2b2385cd7b438469602eb24829a9bd"
72
73 // invExpHex is p-2, the exponent for a modular inverse by Fermat. modexp
74 // with this exponent is also RFC 9380's inv0: it maps 0 to 0, because
75 // 0^(p-2) mod p is 0.
76 invExpHex = "30644e72e131a029b85045b68181585d97816a916871ca8d3c208c16d87cfd45"
77 // legExpHex is (p-1)/2, the Legendre symbol exponent. Numerically equal to
78 // c2Hex; kept separate because the two mean different things and a future
79 // edit to one must not silently change the other.
80 legExpHex = "183227397098d014dc2822db40c0ac2ecbc0b548b438e5469e10460b6c3e7ea3"
81)
82
83var (
84 c1Bytes = mustDecodeHex(c1Hex)
85 c2Bytes = mustDecodeHex(c2Hex)
86 c3Bytes = mustDecodeHex(c3Hex)
87 c4Bytes = mustDecodeHex(c4Hex)
88
89 invExpBytes = mustDecodeHex(invExpHex)
90 legExpBytes = mustDecodeHex(legExpHex)
91
92 fpOne = mustDecodeHex("0000000000000000000000000000000000000000000000000000000000000001")
93 fpThree = mustDecodeHex("0000000000000000000000000000000000000000000000000000000000000003")
94)
95
96// ---- Fp helpers, all on 32-byte big-endian values reduced mod p ----
97
98func fpAdd(a, b []byte) []byte { return fromLimbs(addModP(toLimbs(a), toLimbs(b), pLimbs)) }
99func fpSub(a, b []byte) []byte { return fromLimbs(subModP(toLimbs(a), toLimbs(b))) }
100func fpNeg(a []byte) []byte { return fromLimbs(negModP(toLimbs(a))) }
101
102// fpMul and fpSqr are THE field multiply and square for the drand path. Every
103// product in the SVDW map below goes through these two names and nowhere else, so
104// the Montgomery-versus-modexp choice is made once, here, with the measurements
105// that decided it recorded at the bottom of fieldmul.gno.
106//
107// Measured on this toolchain, per operation, harness baseline subtracted:
108//
109// mulModP (4-limb Montgomery) 2,200,369 gas
110// mulModPIdentity(squaring identity) 519,791 gas <- 4.23x cheaper
111// sqrModP (Montgomery a*a) 2,200,184 gas
112// sqrModPNative (one native modexp) 69,602 gas <- 31.6x cheaper
113//
114// (mulModPIdentity was 511,291 until 2026-09-29; the 8,500 is finding 21's
115// unreduced-operand guard inside mulModPLimbs. The ranking is untouched.)
116//
117// so both names bind to the identity/native pair. mulModP and sqrModP are kept in
118// fieldmul.gno as the independent oracle the equivalence test compares against,
119// not as dead code.
120func fpMul(a, b []byte) []byte { return mulModPIdentity(a, b) }
121func fpSqr(a []byte) []byte { return sqrModPNative(a) }
122
123// fpInv0 is RFC 9380's inv0: the multiplicative inverse, with inv0(0) == 0.
124// Fermat's little theorem gives a^(p-2) == a^-1 for a != 0, and modexp maps a
125// zero base to zero, which is exactly inv0's extra clause.
126func fpInv0(a []byte) []byte { return modexp.ModExp(a, invExpBytes, pBytes) }
127
128// fpSqrt returns a square root of a when one exists. p % 4 == 3, so a^((p+1)/4)
129// is a root; callers must have established that a is a square, or check the
130// result by squaring it.
131func fpSqrt(a []byte) []byte { return modexp.ModExp(a, sqrtExpBytes, pBytes) }
132
133// fpIsSquareNonZero reports whether the Legendre symbol of a is exactly +1, i.e.
134// a is a non-zero quadratic residue.
135//
136// RFC 9380's step 15 uses is_square, which is true for 0 as well. This uses
137// "Legendre == 1", which is false for 0, because that is what both reference
138// implementations of this exact ciphersuite do: drand/kyber v1.3.2's mapToPoint
139// branches on `legendre(g(x1)) == 1`, and randa-mu/bls-solidity's BLS.sol on
140// `legendre(g(x1)) == 1`. The two differ only when g(x) == 0, which happens for
141// no u any adversary can reach in practice -- but the signers run kyber, so
142// matching kyber is what makes signatures verify, and matching the RFC letter
143// here would be a consensus bug rather than a fix.
144func fpIsSquareNonZero(a []byte) bool {
145 return bytesEqual(modexp.ModExp(a, legExpBytes, pBytes), fpOne)
146}
147
148// fpCurveRHS evaluates g(x) = x^3 + 3, the right-hand side of the BN254 G1
149// equation. Returns nil if x is not a reduced 32-byte field element.
150func fpCurveRHS(x []byte) []byte {
151 x2 := fpSqr(x)
152 if x2 == nil {
153 return nil
154 }
155 x3 := fpMul(x2, x)
156 if x3 == nil {
157 return nil
158 }
159 return fpAdd(x3, fpThree)
160}
161
162// fpSgn0 is RFC 9380's sgn0 for m == 1: the low bit of the integer
163// representation.
164func fpSgn0(a []byte) byte { return a[31] & 1 }
165
166// ---- hash_to_field ----
167
168// hashToFieldDST is RFC 9380 hash_to_field with m = 1, count = 2 and L = 48: it
169// expands msg to 96 bytes and reduces each 48-byte half mod p.
170//
171// The mod-p reduction needs no new code: modexp.ModExp(chunk, 1, p) is
172// OS2IP(chunk) mod p and returns exactly len(p) == 32 bytes, left-padded.
173//
174// ok is false only if expandMessageXMD rejects its input, i.e. a DST longer than
175// 255 bytes. It cannot be false for DrandDST.
176func hashToFieldDST(msg, dst []byte) (u0, u1 [32]byte, ok bool) {
177 b := expandMessageXMD(msg, dst, 96)
178 if len(b) != 96 {
179 return u0, u1, false
180 }
181 copy(u0[:], modexp.ModExp(b[0:48], []byte{1}, pBytes))
182 copy(u1[:], modexp.ModExp(b[48:96], []byte{1}, pBytes))
183 return u0, u1, true
184}
185
186// hashToFieldDrand is hash_to_field for drand's bn254 G1 ciphersuite: two field
187// elements, expand_message_xmd with H = legacy Keccak-256, DST = DrandDST.
188//
189// Panics only on an unreachable branch: DrandDST is a fixed 43-byte constant and
190// the output length is a fixed 96, so expandMessageXMD cannot reject. The panic
191// is there so that a future edit which breaks that invariant fails loudly
192// instead of returning two zero field elements, which would map to a valid-looking
193// point and silently accept the wrong signatures. VerifyDrand does not go through
194// this function; it uses hashToFieldDST and stays total.
195func hashToFieldDrand(msg []byte) (u0, u1 [32]byte) {
196 u0, u1, ok := hashToFieldDST(msg, []byte(DrandDST))
197 if !ok {
198 panic("bls254: hashToFieldDrand: expand_message_xmd rejected a fixed-size input")
199 }
200 return u0, u1
201}
202
203// ---- map_to_curve ----
204
205// mapToPointSVDW maps one field element to a point on BN254 G1 by the
206// Shallue-van de Woestijne method, RFC 9380 section 6.6.1 with Z = 1.
207//
208// Returns 64 bytes, x|y, each 32-byte big-endian. Returns nil if u is not a
209// reduced field element (u >= p), or -- and this should be unreachable -- if the
210// result fails the curve equation. The cofactor of BN254 G1 is 1, so the output
211// needs no cofactor clearing and is already in the prime-order subgroup.
212//
213// The step numbering in the comments is RFC 9380's.
214func mapToPointSVDW(u [32]byte) []byte {
215 ub := u[:]
216 if cmp(toLimbs(ub), pLimbs) >= 0 {
217 return nil // field element >= p
218 }
219
220 // 1-2. tv1 = u^2 * c1
221 tv1 := fpMul(fpSqr(ub), c1Bytes)
222 if tv1 == nil {
223 return nil
224 }
225 // 3. tv2 = 1 + tv1
226 tv2 := fpAdd(fpOne, tv1)
227 // 4. tv1 = 1 - tv1
228 tv1 = fpSub(fpOne, tv1)
229 // 5-6. tv3 = inv0(tv1 * tv2)
230 tv3 := fpInv0(fpMul(tv1, tv2))
231 // 7-9. tv4 = u * tv1 * tv3 * c3
232 tv4 := fpMul(fpMul(fpMul(ub, tv1), tv3), c3Bytes)
233 if tv4 == nil {
234 return nil
235 }
236 // 10. x1 = c2 - tv4
237 x1 := fpSub(c2Bytes, tv4)
238 // 16. x2 = c2 + tv4
239 x2 := fpAdd(c2Bytes, tv4)
240 // 22-26. x3 = c4 * (tv2^2 * tv3)^2 + Z, with Z = 1
241 tv8 := fpMul(fpSqr(tv2), tv3)
242 x3 := fpAdd(fpMul(c4Bytes, fpSqr(tv8)), fpOne)
243
244 // 11-21, 27-28. The first of x1, x2, x3 whose g(x) is a non-zero square
245 // wins. See fpIsSquareNonZero for why this is Legendre == 1 rather than the
246 // RFC's is_square.
247 var x, gx []byte
248 if g := fpCurveRHS(x1); fpIsSquareNonZero(g) {
249 x, gx = x1, g
250 } else if g := fpCurveRHS(x2); fpIsSquareNonZero(g) {
251 x, gx = x2, g
252 } else {
253 x, gx = x3, fpCurveRHS(x3)
254 }
255 if gx == nil {
256 return nil
257 }
258
259 // 29-33. y = sqrt(g(x))
260 y := fpSqrt(gx)
261 // 34-35. y gets u's sign.
262 if fpSgn0(ub) != fpSgn0(y) {
263 y = fpNeg(y)
264 }
265
266 // Postcondition. SVDW is supposed to land on the curve unconditionally; if
267 // it has not, something above is wrong and returning the point would hand a
268 // caller an off-curve "G1 element". Two extra multiplications to never do
269 // that. This also catches the g(x3) branch being reached with a non-square,
270 // which would make fpSqrt return a non-root.
271 if !bytesEqual(fpSqr(y), gx) {
272 return nil
273 }
274
275 out := make([]byte, 0, 64)
276 out = append(out, x...)
277 out = append(out, y...)
278 return out
279}
280
281// ---- hash_to_curve ----
282
283// HashToG1Drand maps msg to a G1 point the way drand's
284// bls-bn254-unchained-on-g1 scheme does: RFC 9380 hash_to_curve, SVDW map,
285// expand_message_xmd with H = legacy Keccak-256 and
286// DST = BLS_SIG_BN254G1_XMD:KECCAK-256_SVDW_RO_NUL_.
287//
288// Returns 64 bytes, x|y. Returns nil if anything in the construction fails,
289// which for a fixed DST means only that the native G1 addition rejected its
290// input -- it never should.
291//
292// This is NOT HashToG1. HashToG1 is the existing SPLASH beacon's
293// try-and-increment construction over a different domain tag, and the two return
294// different points for the same message. Callers must not mix them.
295func HashToG1Drand(msg []byte) []byte {
296 u0, u1, ok := hashToFieldDST(msg, []byte(DrandDST))
297 if !ok {
298 return nil
299 }
300 p0 := mapToPointSVDW(u0)
301 p1 := mapToPointSVDW(u1)
302 if p0 == nil || p1 == nil {
303 return nil
304 }
305 in := make([]byte, 0, bn254.G1AddInputSize)
306 in = append(in, p0...)
307 in = append(in, p1...)
308 sum := bn254.G1Add(in)
309 if len(sum) != bn254.G1OutputSize {
310 return nil
311 }
312 return sum
313}