genprecomps.go raw
1 // Copyright 2015 The btcsuite developers
2 // Copyright (c) 2015-2022 The Decred developers
3 // Use of this source code is governed by an ISC
4 // license that can be found in the LICENSE file.
5
6 package main
7
8 // References:
9 // [GECC]: Guide to Elliptic Curve Cryptography (Hankerson, Menezes, Vanstone)
10
11 import (
12 "fmt"
13 "math/big"
14 "os"
15
16 "next.orly.dev/pkg/nostr/crypto/ec/secp256k1"
17 "next.orly.dev/pkg/lol/chk"
18 "next.orly.dev/pkg/lol/log"
19 )
20
21 // curveParams houses the secp256k1 curve parameters for convenient access.
22 var curveParams = secp256k1.Params()
23
24 // bigAffineToJacobian takes an affine point (x, y) as big integers and converts
25 // it to Jacobian point with Z=1.
26 func bigAffineToJacobian(x, y *big.Int, result *secp256k1.JacobianPoint) {
27 result.X.SetByteSlice(x.Bytes())
28 result.Y.SetByteSlice(y.Bytes())
29 result.Z.SetInt(1)
30 }
31
32 // serializedBytePoints returns a serialized byte slice which contains all possible points per
33 // 8-bit window. This is used to when generating compressedbytepoints.go.
34 func serializedBytePoints() []byte {
35 // Calculate G^(2^i) for i in 0..255. These are used to avoid recomputing
36 // them for each digit of the 8-bit windows.
37 doublingPoints := make([]secp256k1.JacobianPoint, curveParams.BitSize)
38 var q secp256k1.JacobianPoint
39 bigAffineToJacobian(curveParams.Gx, curveParams.Gy, &q)
40 for i := 0; i < curveParams.BitSize; i++ {
41 // Q = 2*Q.
42 doublingPoints[i] = q
43 secp256k1.DoubleNonConst(&q, &q)
44 }
45
46 // Separate the bits into byte-sized windows.
47 curveByteSize := curveParams.BitSize / 8
48 serialized := make([]byte, curveByteSize*256*2*32)
49 offset := 0
50 for byteNum := 0; byteNum < curveByteSize; byteNum++ {
51 // Grab the 8 bits that make up this byte from doubling points.
52 startingBit := 8 * (curveByteSize - byteNum - 1)
53 windowPoints := doublingPoints[startingBit : startingBit+8]
54
55 // Compute all points in this window, convert them to affine, and
56 // serialize them.
57 for i := 0; i < 256; i++ {
58 var point secp256k1.JacobianPoint
59 for bit := 0; bit < 8; bit++ {
60 if i>>uint(bit)&1 == 1 {
61 secp256k1.AddNonConst(&point, &windowPoints[bit], &point)
62 }
63 }
64 point.ToAffine()
65 point.X.PutBytesUnchecked(serialized[offset:])
66 offset += 32
67 point.Y.PutBytesUnchecked(serialized[offset:])
68 offset += 32
69 }
70 }
71 return serialized
72 }
73
74 // sqrt returns the square root of the provided big integer using Newton's
75 // method. It's only compiled and used during generation of pre-computed
76 // values, so speed is not a huge concern.
77 func sqrt(n *big.Int) *big.Int {
78 // Initial guess = 2^(log_2(n)/2)
79 guess := big.NewInt(2)
80 guess.Exp(guess, big.NewInt(int64(n.BitLen()/2)), nil)
81 // Now refine using Newton's method.
82 big2 := big.NewInt(2)
83 prevGuess := big.NewInt(0)
84 for {
85 prevGuess.Set(guess)
86 guess.Add(guess, new(big.Int).Div(n, guess))
87 guess.Div(guess, big2)
88 if guess.Cmp(prevGuess) == 0 {
89 break
90 }
91 }
92 return guess
93 }
94
95 // endomorphismParams houses the parameters needed to make use of the secp256k1
96 // endomorphism.
97 type endomorphismParams struct {
98 lambda *big.Int
99 beta *big.Int
100 a1, b1 *big.Int
101 a2, b2 *big.Int
102 z1, z2 *big.Int
103 }
104
105 // endomorphismVectors runs the first 3 steps of algorithm 3.74 from [GECC] to
106 // generate the linearly independent vectors needed to generate a balanced
107 // length-two representation of a multiplier such that k = k1 + k2λ (mod N) and
108 // returns them. Since the values will always be the same given the fact that N
109 // and λ are fixed, the final results can be accelerated by storing the
110 // precomputed values.
111 func endomorphismVectors(lambda *big.Int) (a1, b1, a2, b2 *big.Int) {
112 // This section uses an extended Euclidean algorithm to generate a
113 // sequence of equations:
114 // s[i] * N + t[i] * λ = r[i]
115 nSqrt := sqrt(curveParams.N)
116 u, v := new(big.Int).Set(curveParams.N), new(big.Int).Set(lambda)
117 x1, y1 := big.NewInt(1), big.NewInt(0)
118 x2, y2 := big.NewInt(0), big.NewInt(1)
119 q, r := new(big.Int), new(big.Int)
120 qu, qx1, qy1 := new(big.Int), new(big.Int), new(big.Int)
121 s, t := new(big.Int), new(big.Int)
122 ri, ti := new(big.Int), new(big.Int)
123 a1, b1, a2, b2 = new(big.Int), new(big.Int), new(big.Int), new(big.Int)
124 found, oneMore := false, false
125 for u.Sign() != 0 {
126 // q = v/u
127 q.Div(v, u)
128 // r = v - q*u
129 qu.Mul(q, u)
130 r.Sub(v, qu)
131 // s = x2 - q*x1
132 qx1.Mul(q, x1)
133 s.Sub(x2, qx1)
134 // t = y2 - q*y1
135 qy1.Mul(q, y1)
136 t.Sub(y2, qy1)
137 // v = u, u = r, x2 = x1, x1 = s, y2 = y1, y1 = t
138 v.Set(u)
139 u.Set(r)
140 x2.Set(x1)
141 x1.Set(s)
142 y2.Set(y1)
143 y1.Set(t)
144 // As soon as the remainder is less than the sqrt of n, the
145 // values of a1 and b1 are known.
146 if !found && r.Cmp(nSqrt) < 0 {
147 // When this condition executes ri and ti represent the
148 // r[i] and t[i] values such that i is the greatest
149 // index for which r >= sqrt(n). Meanwhile, the current
150 // r and t values are r[i+1] and t[i+1], respectively.
151 //
152 // a1 = r[i+1], b1 = -t[i+1]
153 a1.Set(r)
154 b1.Neg(t)
155 found = true
156 oneMore = true
157 // Skip to the next iteration so ri and ti are not
158 // modified.
159 continue
160
161 } else if oneMore {
162 // When this condition executes ri and ti still
163 // represent the r[i] and t[i] values while the current
164 // r and t are r[i+2] and t[i+2], respectively.
165 //
166 // sum1 = r[i]^2 + t[i]^2
167 rSquared := new(big.Int).Mul(ri, ri)
168 tSquared := new(big.Int).Mul(ti, ti)
169 sum1 := new(big.Int).Add(rSquared, tSquared)
170 // sum2 = r[i+2]^2 + t[i+2]^2
171 r2Squared := new(big.Int).Mul(r, r)
172 t2Squared := new(big.Int).Mul(t, t)
173 sum2 := new(big.Int).Add(r2Squared, t2Squared)
174 // if (r[i]^2 + t[i]^2) <= (r[i+2]^2 + t[i+2]^2)
175 if sum1.Cmp(sum2) <= 0 {
176 // a2 = r[i], b2 = -t[i]
177 a2.Set(ri)
178 b2.Neg(ti)
179 } else {
180 // a2 = r[i+2], b2 = -t[i+2]
181 a2.Set(r)
182 b2.Neg(t)
183 }
184 // All done.
185 break
186 }
187 ri.Set(r)
188 ti.Set(t)
189 }
190
191 return a1, b1, a2, b2
192 }
193
194 // deriveEndomorphismParams calculates and returns parameters needed to make use
195 // of the secp256k1 endomorphism.
196 func deriveEndomorphismParams() [2]endomorphismParams {
197 // roots returns the solutions of the characteristic polynomial of the
198 // secp256k1 endomorphism.
199 //
200 // The characteristic polynomial for the endomorphism is:
201 //
202 // X^2 + X + 1 ≡ 0 (mod n).
203 //
204 // Solving for X:
205 //
206 // 4X^2 + 4X + 4 ≡ 0 (mod n) | (*4, possible because gcd(4, n) = 1)
207 // (2X + 1)^2 + 3 ≡ 0 (mod n) | (factor by completing the square)
208 // (2X + 1)^2 ≡ -3 (mod n) | (-3)
209 // (2X + 1) ≡ ±sqrt(-3) (mod n) | (sqrt)
210 // 2X ≡ ±sqrt(-3) - 1 (mod n) | (-1)
211 // X ≡ (±sqrt(-3)-1)*2^-1 (mod n) | (*2^-1)
212 //
213 // So, the roots are:
214 // X1 ≡ (-(sqrt(-3)+1)*2^-1 (mod n)
215 // X2 ≡ (sqrt(-3)-1)*2^-1 (mod n)
216 //
217 // It is also worth noting that X is a cube root of unity, meaning
218 // X^3 - 1 ≡ 0 (mod n), hence it can be factored as (X - 1)(X^2 + X + 1) ≡ 0
219 // and (X1)^2 ≡ X2 (mod n) and (X2)^2 ≡ X1 (mod n).
220 roots := func(prime *big.Int) [2]big.Int {
221 var result [2]big.Int
222 one := big.NewInt(1)
223 twoInverse := new(big.Int).ModInverse(big.NewInt(2), prime)
224 negThree := new(big.Int).Neg(big.NewInt(3))
225 sqrtNegThree := new(big.Int).ModSqrt(negThree, prime)
226 sqrtNegThreePlusOne := new(big.Int).Add(sqrtNegThree, one)
227 negSqrtNegThreePlusOne := new(big.Int).Neg(sqrtNegThreePlusOne)
228 result[0].Mul(negSqrtNegThreePlusOne, twoInverse)
229 result[0].Mod(&result[0], prime)
230 sqrtNegThreeMinusOne := new(big.Int).Sub(sqrtNegThree, one)
231 result[1].Mul(sqrtNegThreeMinusOne, twoInverse)
232 result[1].Mod(&result[1], prime)
233 return result
234 }
235 // Find the λ's and β's which are the solutions for the characteristic
236 // polynomial of the secp256k1 endomorphism modulo the curve order and
237 // modulo the field order, respectively.
238 lambdas := roots(curveParams.N)
239 betas := roots(curveParams.P)
240 // Ensure the calculated roots are actually the roots of the characteristic
241 // polynomial.
242 checkRoots := func(foundRoots [2]big.Int, prime *big.Int) {
243 // X^2 + X + 1 ≡ 0 (mod p)
244 one := big.NewInt(1)
245 for i := 0; i < len(foundRoots); i++ {
246 root := &foundRoots[i]
247 result := new(big.Int).Mul(root, root)
248 result.Add(result, root)
249 result.Add(result, one)
250 result.Mod(result, prime)
251 if result.Sign() != 0 {
252 panic(
253 fmt.Sprintf(
254 "%[1]x^2 + %[1]x + 1 != 0 (mod %x)", root,
255 prime,
256 ),
257 )
258 }
259 }
260 }
261 checkRoots(lambdas, curveParams.N)
262 checkRoots(betas, curveParams.P)
263 // checkVectors ensures the passed vectors satisfy the equation:
264 // a + b*λ ≡ 0 (mod n)
265 checkVectors := func(a, b *big.Int, lambda *big.Int) {
266 result := new(big.Int).Mul(b, lambda)
267 result.Add(result, a)
268 result.Mod(result, curveParams.N)
269 if result.Sign() != 0 {
270 panic(
271 fmt.Sprintf(
272 "%x + %x*lambda != 0 (mod %x)", a, b,
273 curveParams.N,
274 ),
275 )
276 }
277 }
278 var endoParams [2]endomorphismParams
279 for i := 0; i < 2; i++ {
280 // Calculate the linearly independent vectors needed to generate a
281 // balanced length-two representation of a scalar such that
282 // k = k1 + k2*λ (mod n) for each of the solutions.
283 lambda := &lambdas[i]
284 a1, b1, a2, b2 := endomorphismVectors(lambda)
285 // Ensure the derived vectors satisfy the required equation.
286 checkVectors(a1, b1, lambda)
287 checkVectors(a2, b2, lambda)
288 // Calculate the precomputed estimates also used when generating the
289 // aforementioned decomposition.
290 //
291 // z1 = floor(b2<<320 / n)
292 // z2 = floor(((-b1)%n)<<320) / n)
293 const shift = 320
294 z1 := new(big.Int).Lsh(b2, shift)
295 z1.Div(z1, curveParams.N)
296 z2 := new(big.Int).Neg(b1)
297 z2.Lsh(z2, shift)
298 z2.Div(z2, curveParams.N)
299 params := &endoParams[i]
300 params.lambda = lambda
301 params.beta = &betas[i]
302 params.a1 = a1
303 params.b1 = b1
304 params.a2 = a2
305 params.b2 = b2
306 params.z1 = z1
307 params.z2 = z2
308 }
309 return endoParams
310 }
311
312 func main() {
313 if _, err := os.Stat(".git"); chk.T(err) {
314 fmt.Printf("File exists\n")
315 _, _ = fmt.Fprintln(
316 os.Stderr,
317 "This generator must be run with working directory at the root of"+
318 " the repository",
319 )
320 os.Exit(1)
321 }
322 serialized := serializedBytePoints()
323 embedded, err := os.Create("secp256k1/rawbytepoints.bin")
324 if err != nil {
325 log.F.Ln(err)
326 os.Exit(1)
327 }
328 n, err := embedded.Write(serialized)
329 if err != nil {
330 panic(err)
331 }
332 if n != len(serialized) {
333 fmt.Printf(
334 "failed to write all of byte points, wrote %d expected %d",
335 n, len(serialized),
336 )
337 panic("fail")
338 }
339 _ = embedded.Close()
340 }
341