1 // Copyright 2010 The Go Authors. All rights reserved.
2 // Copyright 2011 ThePiachu. All rights reserved.
3 // Copyright 2013-2014 The btcsuite developers
4 // Use of this source code is governed by an ISC
5 // license that can be found in the LICENSE file.
6 7 package ecc
8 9 // References:
10 // [SECG]: Recommended Elliptic Curve Domain Parameters
11 // http://www.secg.org/sec2-v2.pdf
12 //
13 // [GECC]: Guide to Elliptic Curve Cryptography (Hankerson, Menezes, Vanstone)
14 15 // This package operates, internally, on Jacobian coordinates. For a given
16 // (x, y) position on the curve, the Jacobian coordinates are (x1, y1, z1)
17 // where x = x1/z1² and y = y1/z1³. The greatest speedups come when the whole
18 // calculation can be performed within the transform (as in ScalarMult and
19 // ScalarBaseMult). But even for Add and Double, it's faster to apply and
20 // reverse the transform than to operate in affine coordinates.
21 22 import (
23 "crypto/elliptic"
24 "math/big"
25 "sync"
26 )
27 28 var (
29 // fieldOne is simply the integer 1 in field representation. It is
30 // used to avoid needing to create it multiple times during the internal
31 // arithmetic.
32 fieldOne = new(fieldVal).SetInt(1)
33 )
34 35 // KoblitzCurve supports a koblitz curve implementation that fits the ECC Curve
36 // interface from crypto/elliptic.
37 type KoblitzCurve struct {
38 *elliptic.CurveParams
39 40 // q is the value (P+1)/4 used to compute the square root of field
41 // elements.
42 q *big.Int
43 44 H int // cofactor of the curve.
45 halfOrder *big.Int // half the order N
46 47 // fieldB is the constant B of the curve as a fieldVal.
48 fieldB *fieldVal
49 50 // byteSize is simply the bit size / 8 and is provided for convenience
51 // since it is calculated repeatedly.
52 byteSize int
53 54 // bytePoints
55 bytePoints *[32][256][3]fieldVal
56 57 // The next 6 values are used specifically for endomorphism
58 // optimizations in ScalarMult.
59 60 // lambda must fulfill lambda^3 = 1 mod N where N is the order of G.
61 lambda *big.Int
62 63 // beta must fulfill beta^3 = 1 mod P where P is the prime field of the
64 // curve.
65 beta *fieldVal
66 67 // See the EndomorphismVectors in gensecp256k1.go to see how these are
68 // derived.
69 a1 *big.Int
70 b1 *big.Int
71 a2 *big.Int
72 b2 *big.Int
73 }
74 75 // Params returns the parameters for the curve.
76 func (curve *KoblitzCurve) Params() *elliptic.CurveParams {
77 return curve.CurveParams
78 }
79 80 // bigAffineToField takes an affine point (x, y) as big integers and converts
81 // it to an affine point as field values.
82 func (curve *KoblitzCurve) bigAffineToField(x, y *big.Int) (*fieldVal, *fieldVal) {
83 x3, y3 := new(fieldVal), new(fieldVal)
84 x3.SetByteSlice(x.Bytes())
85 y3.SetByteSlice(y.Bytes())
86 87 return x3, y3
88 }
89 90 // fieldJacobianToBigAffine takes a Jacobian point (x, y, z) as field values and
91 // converts it to an affine point as big integers.
92 func (curve *KoblitzCurve) fieldJacobianToBigAffine(x, y, z *fieldVal) (*big.Int, *big.Int) {
93 // Inversions are expensive and both point addition and point doubling
94 // are faster when working with points that have a z value of one. So,
95 // if the point needs to be converted to affine, go ahead and normalize
96 // the point itself at the same time as the calculation is the same.
97 var zInv, tempZ fieldVal
98 zInv.Set(z).Inverse() // zInv = Z^-1
99 tempZ.SquareVal(&zInv) // tempZ = Z^-2
100 x.Mul(&tempZ) // X = X/Z^2 (mag: 1)
101 y.Mul(tempZ.Mul(&zInv)) // Y = Y/Z^3 (mag: 1)
102 z.SetInt(1) // Z = 1 (mag: 1)
103 104 // Normalize the x and y values.
105 x.Normalize()
106 y.Normalize()
107 108 // Convert the field values for the now affine point to big.Ints.
109 x3, y3 := new(big.Int), new(big.Int)
110 x3.SetBytes(x.Bytes()[:])
111 y3.SetBytes(y.Bytes()[:])
112 return x3, y3
113 }
114 115 // IsOnCurve returns boolean if the point (x,y) is on the curve.
116 // Part of the elliptic.Curve interface. This function differs from the
117 // crypto/elliptic algorithm since a = 0 not -3.
118 func (curve *KoblitzCurve) IsOnCurve(x, y *big.Int) bool {
119 // Convert big ints to field values for faster arithmetic.
120 fx, fy := curve.bigAffineToField(x, y)
121 122 // Elliptic curve equation for secp256k1 is: y^2 = x^3 + 7
123 y2 := new(fieldVal).SquareVal(fy).Normalize()
124 result := new(fieldVal).SquareVal(fx).Mul(fx).AddInt(7).Normalize()
125 return y2.Equals(result)
126 }
127 128 // addZ1AndZ2EqualsOne adds two Jacobian points that are already known to have
129 // z values of 1 and stores the result in (x3, y3, z3). That is to say
130 // (x1, y1, 1) + (x2, y2, 1) = (x3, y3, z3). It performs faster addition than
131 // the generic add routine since less arithmetic is needed due to the ability to
132 // avoid the z value multiplications.
133 func (curve *KoblitzCurve) addZ1AndZ2EqualsOne(x1, y1, z1, x2, y2, x3, y3, z3 *fieldVal) {
134 // To compute the point addition efficiently, this implementation splits
135 // the equation into intermediate elements which are used to minimize
136 // the number of field multiplications using the method shown at:
137 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-mmadd-2007-bl
138 //
139 // In particular it performs the calculations using the following:
140 // H = X2-X1, HH = H^2, I = 4*HH, J = H*I, r = 2*(Y2-Y1), V = X1*I
141 // X3 = r^2-J-2*V, Y3 = r*(V-X3)-2*Y1*J, Z3 = 2*H
142 //
143 // This results in a cost of 4 field multiplications, 2 field squarings,
144 // 6 field additions, and 5 integer multiplications.
145 146 // When the x coordinates are the same for two points on the curve, the
147 // y coordinates either must be the same, in which case it is point
148 // doubling, or they are opposite and the result is the point at
149 // infinity per the group law for elliptic curve cryptography.
150 x1.Normalize()
151 y1.Normalize()
152 x2.Normalize()
153 y2.Normalize()
154 if x1.Equals(x2) {
155 if y1.Equals(y2) {
156 // Since x1 == x2 and y1 == y2, point doubling must be
157 // done, otherwise the addition would end up dividing
158 // by zero.
159 curve.doubleJacobian(x1, y1, z1, x3, y3, z3)
160 return
161 }
162 163 // Since x1 == x2 and y1 == -y2, the sum is the point at
164 // infinity per the group law.
165 x3.SetInt(0)
166 y3.SetInt(0)
167 z3.SetInt(0)
168 return
169 }
170 171 // Calculate X3, Y3, and Z3 according to the intermediate elements
172 // breakdown above.
173 var h, i, j, r, v fieldVal
174 var negJ, neg2V, negX3 fieldVal
175 h.Set(x1).Negate(1).Add(x2) // H = X2-X1 (mag: 3)
176 i.SquareVal(&h).MulInt(4) // I = 4*H^2 (mag: 4)
177 j.Mul2(&h, &i) // J = H*I (mag: 1)
178 r.Set(y1).Negate(1).Add(y2).MulInt(2) // r = 2*(Y2-Y1) (mag: 6)
179 v.Mul2(x1, &i) // V = X1*I (mag: 1)
180 negJ.Set(&j).Negate(1) // negJ = -J (mag: 2)
181 neg2V.Set(&v).MulInt(2).Negate(2) // neg2V = -(2*V) (mag: 3)
182 x3.Set(&r).Square().Add(&negJ).Add(&neg2V) // X3 = r^2-J-2*V (mag: 6)
183 negX3.Set(x3).Negate(6) // negX3 = -X3 (mag: 7)
184 j.Mul(y1).MulInt(2).Negate(2) // J = -(2*Y1*J) (mag: 3)
185 y3.Set(&v).Add(&negX3).Mul(&r).Add(&j) // Y3 = r*(V-X3)-2*Y1*J (mag: 4)
186 z3.Set(&h).MulInt(2) // Z3 = 2*H (mag: 6)
187 188 // Normalize the resulting field values to a magnitude of 1 as needed.
189 x3.Normalize()
190 y3.Normalize()
191 z3.Normalize()
192 }
193 194 // addZ1EqualsZ2 adds two Jacobian points that are already known to have the
195 // same z value and stores the result in (x3, y3, z3). That is to say
196 // (x1, y1, z1) + (x2, y2, z1) = (x3, y3, z3). It performs faster addition than
197 // the generic add routine since less arithmetic is needed due to the known
198 // equivalence.
199 func (curve *KoblitzCurve) addZ1EqualsZ2(x1, y1, z1, x2, y2, x3, y3, z3 *fieldVal) {
200 // To compute the point addition efficiently, this implementation splits
201 // the equation into intermediate elements which are used to minimize
202 // the number of field multiplications using a slightly modified version
203 // of the method shown at:
204 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-mmadd-2007-bl
205 //
206 // In particular it performs the calculations using the following:
207 // A = X2-X1, B = A^2, C=Y2-Y1, D = C^2, E = X1*B, F = X2*B
208 // X3 = D-E-F, Y3 = C*(E-X3)-Y1*(F-E), Z3 = Z1*A
209 //
210 // This results in a cost of 5 field multiplications, 2 field squarings,
211 // 9 field additions, and 0 integer multiplications.
212 213 // When the x coordinates are the same for two points on the curve, the
214 // y coordinates either must be the same, in which case it is point
215 // doubling, or they are opposite and the result is the point at
216 // infinity per the group law for elliptic curve cryptography.
217 x1.Normalize()
218 y1.Normalize()
219 x2.Normalize()
220 y2.Normalize()
221 if x1.Equals(x2) {
222 if y1.Equals(y2) {
223 // Since x1 == x2 and y1 == y2, point doubling must be
224 // done, otherwise the addition would end up dividing
225 // by zero.
226 curve.doubleJacobian(x1, y1, z1, x3, y3, z3)
227 return
228 }
229 230 // Since x1 == x2 and y1 == -y2, the sum is the point at
231 // infinity per the group law.
232 x3.SetInt(0)
233 y3.SetInt(0)
234 z3.SetInt(0)
235 return
236 }
237 238 // Calculate X3, Y3, and Z3 according to the intermediate elements
239 // breakdown above.
240 var a, b, c, d, e, f fieldVal
241 var negX1, negY1, negE, negX3 fieldVal
242 negX1.Set(x1).Negate(1) // negX1 = -X1 (mag: 2)
243 negY1.Set(y1).Negate(1) // negY1 = -Y1 (mag: 2)
244 a.Set(&negX1).Add(x2) // A = X2-X1 (mag: 3)
245 b.SquareVal(&a) // B = A^2 (mag: 1)
246 c.Set(&negY1).Add(y2) // C = Y2-Y1 (mag: 3)
247 d.SquareVal(&c) // D = C^2 (mag: 1)
248 e.Mul2(x1, &b) // E = X1*B (mag: 1)
249 negE.Set(&e).Negate(1) // negE = -E (mag: 2)
250 f.Mul2(x2, &b) // F = X2*B (mag: 1)
251 x3.Add2(&e, &f).Negate(3).Add(&d) // X3 = D-E-F (mag: 5)
252 negX3.Set(x3).Negate(5).Normalize() // negX3 = -X3 (mag: 1)
253 y3.Set(y1).Mul(f.Add(&negE)).Negate(3) // Y3 = -(Y1*(F-E)) (mag: 4)
254 y3.Add(e.Add(&negX3).Mul(&c)) // Y3 = C*(E-X3)+Y3 (mag: 5)
255 z3.Mul2(z1, &a) // Z3 = Z1*A (mag: 1)
256 257 // Normalize the resulting field values to a magnitude of 1 as needed.
258 x3.Normalize()
259 y3.Normalize()
260 }
261 262 // addZ2EqualsOne adds two Jacobian points when the second point is already
263 // known to have a z value of 1 (and the z value for the first point is not 1)
264 // and stores the result in (x3, y3, z3). That is to say (x1, y1, z1) +
265 // (x2, y2, 1) = (x3, y3, z3). It performs faster addition than the generic
266 // add routine since less arithmetic is needed due to the ability to avoid
267 // multiplications by the second point's z value.
268 func (curve *KoblitzCurve) addZ2EqualsOne(x1, y1, z1, x2, y2, x3, y3, z3 *fieldVal) {
269 // To compute the point addition efficiently, this implementation splits
270 // the equation into intermediate elements which are used to minimize
271 // the number of field multiplications using the method shown at:
272 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-madd-2007-bl
273 //
274 // In particular it performs the calculations using the following:
275 // Z1Z1 = Z1^2, U2 = X2*Z1Z1, S2 = Y2*Z1*Z1Z1, H = U2-X1, HH = H^2,
276 // I = 4*HH, J = H*I, r = 2*(S2-Y1), V = X1*I
277 // X3 = r^2-J-2*V, Y3 = r*(V-X3)-2*Y1*J, Z3 = (Z1+H)^2-Z1Z1-HH
278 //
279 // This results in a cost of 7 field multiplications, 4 field squarings,
280 // 9 field additions, and 4 integer multiplications.
281 282 // When the x coordinates are the same for two points on the curve, the
283 // y coordinates either must be the same, in which case it is point
284 // doubling, or they are opposite and the result is the point at
285 // infinity per the group law for elliptic curve cryptography. Since
286 // any number of Jacobian coordinates can represent the same affine
287 // point, the x and y values need to be converted to like terms. Due to
288 // the assumption made for this function that the second point has a z
289 // value of 1 (z2=1), the first point is already "converted".
290 var z1z1, u2, s2 fieldVal
291 x1.Normalize()
292 y1.Normalize()
293 z1z1.SquareVal(z1) // Z1Z1 = Z1^2 (mag: 1)
294 u2.Set(x2).Mul(&z1z1).Normalize() // U2 = X2*Z1Z1 (mag: 1)
295 s2.Set(y2).Mul(&z1z1).Mul(z1).Normalize() // S2 = Y2*Z1*Z1Z1 (mag: 1)
296 if x1.Equals(&u2) {
297 if y1.Equals(&s2) {
298 // Since x1 == x2 and y1 == y2, point doubling must be
299 // done, otherwise the addition would end up dividing
300 // by zero.
301 curve.doubleJacobian(x1, y1, z1, x3, y3, z3)
302 return
303 }
304 305 // Since x1 == x2 and y1 == -y2, the sum is the point at
306 // infinity per the group law.
307 x3.SetInt(0)
308 y3.SetInt(0)
309 z3.SetInt(0)
310 return
311 }
312 313 // Calculate X3, Y3, and Z3 according to the intermediate elements
314 // breakdown above.
315 var h, hh, i, j, r, rr, v fieldVal
316 var negX1, negY1, negX3 fieldVal
317 negX1.Set(x1).Negate(1) // negX1 = -X1 (mag: 2)
318 h.Add2(&u2, &negX1) // H = U2-X1 (mag: 3)
319 hh.SquareVal(&h) // HH = H^2 (mag: 1)
320 i.Set(&hh).MulInt(4) // I = 4 * HH (mag: 4)
321 j.Mul2(&h, &i) // J = H*I (mag: 1)
322 negY1.Set(y1).Negate(1) // negY1 = -Y1 (mag: 2)
323 r.Set(&s2).Add(&negY1).MulInt(2) // r = 2*(S2-Y1) (mag: 6)
324 rr.SquareVal(&r) // rr = r^2 (mag: 1)
325 v.Mul2(x1, &i) // V = X1*I (mag: 1)
326 x3.Set(&v).MulInt(2).Add(&j).Negate(3) // X3 = -(J+2*V) (mag: 4)
327 x3.Add(&rr) // X3 = r^2+X3 (mag: 5)
328 negX3.Set(x3).Negate(5) // negX3 = -X3 (mag: 6)
329 y3.Set(y1).Mul(&j).MulInt(2).Negate(2) // Y3 = -(2*Y1*J) (mag: 3)
330 y3.Add(v.Add(&negX3).Mul(&r)) // Y3 = r*(V-X3)+Y3 (mag: 4)
331 z3.Add2(z1, &h).Square() // Z3 = (Z1+H)^2 (mag: 1)
332 z3.Add(z1z1.Add(&hh).Negate(2)) // Z3 = Z3-(Z1Z1+HH) (mag: 4)
333 334 // Normalize the resulting field values to a magnitude of 1 as needed.
335 x3.Normalize()
336 y3.Normalize()
337 z3.Normalize()
338 }
339 340 // addGeneric adds two Jacobian points (x1, y1, z1) and (x2, y2, z2) without any
341 // assumptions about the z values of the two points and stores the result in
342 // (x3, y3, z3). That is to say (x1, y1, z1) + (x2, y2, z2) = (x3, y3, z3). It
343 // is the slowest of the add routines due to requiring the most arithmetic.
344 func (curve *KoblitzCurve) addGeneric(x1, y1, z1, x2, y2, z2, x3, y3, z3 *fieldVal) {
345 // To compute the point addition efficiently, this implementation splits
346 // the equation into intermediate elements which are used to minimize
347 // the number of field multiplications using the method shown at:
348 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-add-2007-bl
349 //
350 // In particular it performs the calculations using the following:
351 // Z1Z1 = Z1^2, Z2Z2 = Z2^2, U1 = X1*Z2Z2, U2 = X2*Z1Z1, S1 = Y1*Z2*Z2Z2
352 // S2 = Y2*Z1*Z1Z1, H = U2-U1, I = (2*H)^2, J = H*I, r = 2*(S2-S1)
353 // V = U1*I
354 // X3 = r^2-J-2*V, Y3 = r*(V-X3)-2*S1*J, Z3 = ((Z1+Z2)^2-Z1Z1-Z2Z2)*H
355 //
356 // This results in a cost of 11 field multiplications, 5 field squarings,
357 // 9 field additions, and 4 integer multiplications.
358 359 // When the x coordinates are the same for two points on the curve, the
360 // y coordinates either must be the same, in which case it is point
361 // doubling, or they are opposite and the result is the point at
362 // infinity. Since any number of Jacobian coordinates can represent the
363 // same affine point, the x and y values need to be converted to like
364 // terms.
365 var z1z1, z2z2, u1, u2, s1, s2 fieldVal
366 z1z1.SquareVal(z1) // Z1Z1 = Z1^2 (mag: 1)
367 z2z2.SquareVal(z2) // Z2Z2 = Z2^2 (mag: 1)
368 u1.Set(x1).Mul(&z2z2).Normalize() // U1 = X1*Z2Z2 (mag: 1)
369 u2.Set(x2).Mul(&z1z1).Normalize() // U2 = X2*Z1Z1 (mag: 1)
370 s1.Set(y1).Mul(&z2z2).Mul(z2).Normalize() // S1 = Y1*Z2*Z2Z2 (mag: 1)
371 s2.Set(y2).Mul(&z1z1).Mul(z1).Normalize() // S2 = Y2*Z1*Z1Z1 (mag: 1)
372 if u1.Equals(&u2) {
373 if s1.Equals(&s2) {
374 // Since x1 == x2 and y1 == y2, point doubling must be
375 // done, otherwise the addition would end up dividing
376 // by zero.
377 curve.doubleJacobian(x1, y1, z1, x3, y3, z3)
378 return
379 }
380 381 // Since x1 == x2 and y1 == -y2, the sum is the point at
382 // infinity per the group law.
383 x3.SetInt(0)
384 y3.SetInt(0)
385 z3.SetInt(0)
386 return
387 }
388 389 // Calculate X3, Y3, and Z3 according to the intermediate elements
390 // breakdown above.
391 var h, i, j, r, rr, v fieldVal
392 var negU1, negS1, negX3 fieldVal
393 negU1.Set(&u1).Negate(1) // negU1 = -U1 (mag: 2)
394 h.Add2(&u2, &negU1) // H = U2-U1 (mag: 3)
395 i.Set(&h).MulInt(2).Square() // I = (2*H)^2 (mag: 2)
396 j.Mul2(&h, &i) // J = H*I (mag: 1)
397 negS1.Set(&s1).Negate(1) // negS1 = -S1 (mag: 2)
398 r.Set(&s2).Add(&negS1).MulInt(2) // r = 2*(S2-S1) (mag: 6)
399 rr.SquareVal(&r) // rr = r^2 (mag: 1)
400 v.Mul2(&u1, &i) // V = U1*I (mag: 1)
401 x3.Set(&v).MulInt(2).Add(&j).Negate(3) // X3 = -(J+2*V) (mag: 4)
402 x3.Add(&rr) // X3 = r^2+X3 (mag: 5)
403 negX3.Set(x3).Negate(5) // negX3 = -X3 (mag: 6)
404 y3.Mul2(&s1, &j).MulInt(2).Negate(2) // Y3 = -(2*S1*J) (mag: 3)
405 y3.Add(v.Add(&negX3).Mul(&r)) // Y3 = r*(V-X3)+Y3 (mag: 4)
406 z3.Add2(z1, z2).Square() // Z3 = (Z1+Z2)^2 (mag: 1)
407 z3.Add(z1z1.Add(&z2z2).Negate(2)) // Z3 = Z3-(Z1Z1+Z2Z2) (mag: 4)
408 z3.Mul(&h) // Z3 = Z3*H (mag: 1)
409 410 // Normalize the resulting field values to a magnitude of 1 as needed.
411 x3.Normalize()
412 y3.Normalize()
413 }
414 415 // addJacobian adds the passed Jacobian points (x1, y1, z1) and (x2, y2, z2)
416 // together and stores the result in (x3, y3, z3).
417 func (curve *KoblitzCurve) addJacobian(x1, y1, z1, x2, y2, z2, x3, y3, z3 *fieldVal) {
418 // A point at infinity is the identity according to the group law for
419 // elliptic curve cryptography. Thus, ∞ + P = P and P + ∞ = P.
420 if (x1.IsZero() && y1.IsZero()) || z1.IsZero() {
421 x3.Set(x2)
422 y3.Set(y2)
423 z3.Set(z2)
424 return
425 }
426 if (x2.IsZero() && y2.IsZero()) || z2.IsZero() {
427 x3.Set(x1)
428 y3.Set(y1)
429 z3.Set(z1)
430 return
431 }
432 433 // Faster point addition can be achieved when certain assumptions are
434 // met. For example, when both points have the same z value, arithmetic
435 // on the z values can be avoided. This section thus checks for these
436 // conditions and calls an appropriate add function which is accelerated
437 // by using those assumptions.
438 z1.Normalize()
439 z2.Normalize()
440 isZ1One := z1.Equals(fieldOne)
441 isZ2One := z2.Equals(fieldOne)
442 switch {
443 case isZ1One && isZ2One:
444 curve.addZ1AndZ2EqualsOne(x1, y1, z1, x2, y2, x3, y3, z3)
445 return
446 case z1.Equals(z2):
447 curve.addZ1EqualsZ2(x1, y1, z1, x2, y2, x3, y3, z3)
448 return
449 case isZ2One:
450 curve.addZ2EqualsOne(x1, y1, z1, x2, y2, x3, y3, z3)
451 return
452 }
453 454 // None of the above assumptions are true, so fall back to generic
455 // point addition.
456 curve.addGeneric(x1, y1, z1, x2, y2, z2, x3, y3, z3)
457 }
458 459 // Add returns the sum of (x1,y1) and (x2,y2). Part of the elliptic.Curve
460 // interface.
461 func (curve *KoblitzCurve) Add(x1, y1, x2, y2 *big.Int) (*big.Int, *big.Int) {
462 // A point at infinity is the identity according to the group law for
463 // elliptic curve cryptography. Thus, ∞ + P = P and P + ∞ = P.
464 if x1.Sign() == 0 && y1.Sign() == 0 {
465 return x2, y2
466 }
467 if x2.Sign() == 0 && y2.Sign() == 0 {
468 return x1, y1
469 }
470 471 // Convert the affine coordinates from big integers to field values
472 // and do the point addition in Jacobian projective space.
473 fx1, fy1 := curve.bigAffineToField(x1, y1)
474 fx2, fy2 := curve.bigAffineToField(x2, y2)
475 fx3, fy3, fz3 := new(fieldVal), new(fieldVal), new(fieldVal)
476 fOne := new(fieldVal).SetInt(1)
477 curve.addJacobian(fx1, fy1, fOne, fx2, fy2, fOne, fx3, fy3, fz3)
478 479 // Convert the Jacobian coordinate field values back to affine big
480 // integers.
481 return curve.fieldJacobianToBigAffine(fx3, fy3, fz3)
482 }
483 484 // doubleZ1EqualsOne performs point doubling on the passed Jacobian point
485 // when the point is already known to have a z value of 1 and stores
486 // the result in (x3, y3, z3). That is to say (x3, y3, z3) = 2*(x1, y1, 1). It
487 // performs faster point doubling than the generic routine since less arithmetic
488 // is needed due to the ability to avoid multiplication by the z value.
489 func (curve *KoblitzCurve) doubleZ1EqualsOne(x1, y1, x3, y3, z3 *fieldVal) {
490 // This function uses the assumptions that z1 is 1, thus the point
491 // doubling formulas reduce to:
492 //
493 // X3 = (3*X1^2)^2 - 8*X1*Y1^2
494 // Y3 = (3*X1^2)*(4*X1*Y1^2 - X3) - 8*Y1^4
495 // Z3 = 2*Y1
496 //
497 // To compute the above efficiently, this implementation splits the
498 // equation into intermediate elements which are used to minimize the
499 // number of field multiplications in favor of field squarings which
500 // are roughly 35% faster than field multiplications with the current
501 // implementation at the time this was written.
502 //
503 // This uses a slightly modified version of the method shown at:
504 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-mdbl-2007-bl
505 //
506 // In particular it performs the calculations using the following:
507 // A = X1^2, B = Y1^2, C = B^2, D = 2*((X1+B)^2-A-C)
508 // E = 3*A, F = E^2, X3 = F-2*D, Y3 = E*(D-X3)-8*C
509 // Z3 = 2*Y1
510 //
511 // This results in a cost of 1 field multiplication, 5 field squarings,
512 // 6 field additions, and 5 integer multiplications.
513 var a, b, c, d, e, f fieldVal
514 z3.Set(y1).MulInt(2) // Z3 = 2*Y1 (mag: 2)
515 a.SquareVal(x1) // A = X1^2 (mag: 1)
516 b.SquareVal(y1) // B = Y1^2 (mag: 1)
517 c.SquareVal(&b) // C = B^2 (mag: 1)
518 b.Add(x1).Square() // B = (X1+B)^2 (mag: 1)
519 d.Set(&a).Add(&c).Negate(2) // D = -(A+C) (mag: 3)
520 d.Add(&b).MulInt(2) // D = 2*(B+D)(mag: 8)
521 e.Set(&a).MulInt(3) // E = 3*A (mag: 3)
522 f.SquareVal(&e) // F = E^2 (mag: 1)
523 x3.Set(&d).MulInt(2).Negate(16) // X3 = -(2*D) (mag: 17)
524 x3.Add(&f) // X3 = F+X3 (mag: 18)
525 f.Set(x3).Negate(18).Add(&d).Normalize() // F = D-X3 (mag: 1)
526 y3.Set(&c).MulInt(8).Negate(8) // Y3 = -(8*C) (mag: 9)
527 y3.Add(f.Mul(&e)) // Y3 = E*F+Y3 (mag: 10)
528 529 // Normalize the field values back to a magnitude of 1.
530 x3.Normalize()
531 y3.Normalize()
532 z3.Normalize()
533 }
534 535 // doubleGeneric performs point doubling on the passed Jacobian point without
536 // any assumptions about the z value and stores the result in (x3, y3, z3).
537 // That is to say (x3, y3, z3) = 2*(x1, y1, z1). It is the slowest of the point
538 // doubling routines due to requiring the most arithmetic.
539 func (curve *KoblitzCurve) doubleGeneric(x1, y1, z1, x3, y3, z3 *fieldVal) {
540 // Point doubling formula for Jacobian coordinates for the secp256k1
541 // curve:
542 // X3 = (3*X1^2)^2 - 8*X1*Y1^2
543 // Y3 = (3*X1^2)*(4*X1*Y1^2 - X3) - 8*Y1^4
544 // Z3 = 2*Y1*Z1
545 //
546 // To compute the above efficiently, this implementation splits the
547 // equation into intermediate elements which are used to minimize the
548 // number of field multiplications in favor of field squarings which
549 // are roughly 35% faster than field multiplications with the current
550 // implementation at the time this was written.
551 //
552 // This uses a slightly modified version of the method shown at:
553 // http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-dbl-2009-l
554 //
555 // In particular it performs the calculations using the following:
556 // A = X1^2, B = Y1^2, C = B^2, D = 2*((X1+B)^2-A-C)
557 // E = 3*A, F = E^2, X3 = F-2*D, Y3 = E*(D-X3)-8*C
558 // Z3 = 2*Y1*Z1
559 //
560 // This results in a cost of 1 field multiplication, 5 field squarings,
561 // 6 field additions, and 5 integer multiplications.
562 var a, b, c, d, e, f fieldVal
563 z3.Mul2(y1, z1).MulInt(2) // Z3 = 2*Y1*Z1 (mag: 2)
564 a.SquareVal(x1) // A = X1^2 (mag: 1)
565 b.SquareVal(y1) // B = Y1^2 (mag: 1)
566 c.SquareVal(&b) // C = B^2 (mag: 1)
567 b.Add(x1).Square() // B = (X1+B)^2 (mag: 1)
568 d.Set(&a).Add(&c).Negate(2) // D = -(A+C) (mag: 3)
569 d.Add(&b).MulInt(2) // D = 2*(B+D)(mag: 8)
570 e.Set(&a).MulInt(3) // E = 3*A (mag: 3)
571 f.SquareVal(&e) // F = E^2 (mag: 1)
572 x3.Set(&d).MulInt(2).Negate(16) // X3 = -(2*D) (mag: 17)
573 x3.Add(&f) // X3 = F+X3 (mag: 18)
574 f.Set(x3).Negate(18).Add(&d).Normalize() // F = D-X3 (mag: 1)
575 y3.Set(&c).MulInt(8).Negate(8) // Y3 = -(8*C) (mag: 9)
576 y3.Add(f.Mul(&e)) // Y3 = E*F+Y3 (mag: 10)
577 578 // Normalize the field values back to a magnitude of 1.
579 x3.Normalize()
580 y3.Normalize()
581 z3.Normalize()
582 }
583 584 // doubleJacobian doubles the passed Jacobian point (x1, y1, z1) and stores the
585 // result in (x3, y3, z3).
586 func (curve *KoblitzCurve) doubleJacobian(x1, y1, z1, x3, y3, z3 *fieldVal) {
587 // Doubling a point at infinity is still infinity.
588 if y1.IsZero() || z1.IsZero() {
589 x3.SetInt(0)
590 y3.SetInt(0)
591 z3.SetInt(0)
592 return
593 }
594 595 // Slightly faster point doubling can be achieved when the z value is 1
596 // by avoiding the multiplication on the z value. This section calls
597 // a point doubling function which is accelerated by using that
598 // assumption when possible.
599 if z1.Normalize().Equals(fieldOne) {
600 curve.doubleZ1EqualsOne(x1, y1, x3, y3, z3)
601 return
602 }
603 604 // Fall back to generic point doubling which works with arbitrary z
605 // values.
606 curve.doubleGeneric(x1, y1, z1, x3, y3, z3)
607 }
608 609 // Double returns 2*(x1,y1). Part of the elliptic.Curve interface.
610 func (curve *KoblitzCurve) Double(x1, y1 *big.Int) (*big.Int, *big.Int) {
611 if y1.Sign() == 0 {
612 return new(big.Int), new(big.Int)
613 }
614 615 // Convert the affine coordinates from big integers to field values
616 // and do the point doubling in Jacobian projective space.
617 fx1, fy1 := curve.bigAffineToField(x1, y1)
618 fx3, fy3, fz3 := new(fieldVal), new(fieldVal), new(fieldVal)
619 fOne := new(fieldVal).SetInt(1)
620 curve.doubleJacobian(fx1, fy1, fOne, fx3, fy3, fz3)
621 622 // Convert the Jacobian coordinate field values back to affine big
623 // integers.
624 return curve.fieldJacobianToBigAffine(fx3, fy3, fz3)
625 }
626 627 // splitK returns a balanced length-two representation of k and their signs.
628 // This is algorithm 3.74 from [GECC].
629 //
630 // One thing of note about this algorithm is that no matter what c1 and c2 are,
631 // the final equation of k = k1 + k2 * lambda (mod n) will hold. This is
632 // provable mathematically due to how a1/b1/a2/b2 are computed.
633 //
634 // c1 and c2 are chosen to minimize the max(k1,k2).
635 func (curve *KoblitzCurve) splitK(k []byte) ([]byte, []byte, int, int) {
636 // All math here is done with big.Int, which is slow.
637 // At some point, it might be useful to write something similar to
638 // fieldVal but for N instead of P as the prime field if this ends up
639 // being a bottleneck.
640 bigIntK := new(big.Int)
641 c1, c2 := new(big.Int), new(big.Int)
642 tmp1, tmp2 := new(big.Int), new(big.Int)
643 k1, k2 := new(big.Int), new(big.Int)
644 645 bigIntK.SetBytes(k)
646 // c1 = round(b2 * k / n) from step 4.
647 // Rounding isn't really necessary and costs too much, hence skipped
648 c1.Mul(curve.b2, bigIntK)
649 c1.Div(c1, curve.N)
650 // c2 = round(b1 * k / n) from step 4 (sign reversed to optimize one step)
651 // Rounding isn't really necessary and costs too much, hence skipped
652 c2.Mul(curve.b1, bigIntK)
653 c2.Div(c2, curve.N)
654 // k1 = k - c1 * a1 - c2 * a2 from step 5 (note c2's sign is reversed)
655 tmp1.Mul(c1, curve.a1)
656 tmp2.Mul(c2, curve.a2)
657 k1.Sub(bigIntK, tmp1)
658 k1.Add(k1, tmp2)
659 // k2 = - c1 * b1 - c2 * b2 from step 5 (note c2's sign is reversed)
660 tmp1.Mul(c1, curve.b1)
661 tmp2.Mul(c2, curve.b2)
662 k2.Sub(tmp2, tmp1)
663 664 // Note Bytes() throws out the sign of k1 and k2. This matters
665 // since k1 and/or k2 can be negative. Hence, we pass that
666 // back separately.
667 return k1.Bytes(), k2.Bytes(), k1.Sign(), k2.Sign()
668 }
669 670 // moduloReduce reduces k from more than 32 bytes to 32 bytes and under. This
671 // is done by doing a simple modulo curve.N. We can do this since G^N = 1 and
672 // thus any other valid point on the elliptic curve has the same order.
673 func (curve *KoblitzCurve) moduloReduce(k []byte) []byte {
674 // Since the order of G is curve.N, we can use a much smaller number
675 // by doing modulo curve.N
676 if len(k) > curve.byteSize {
677 // Reduce k by performing modulo curve.N.
678 tmpK := new(big.Int).SetBytes(k)
679 tmpK.Mod(tmpK, curve.N)
680 return tmpK.Bytes()
681 }
682 683 return k
684 }
685 686 // NAF takes a positive integer k and returns the Non-Adjacent Form (NAF) as two
687 // byte slices. The first is where 1s will be. The second is where -1s will
688 // be. NAF is convenient in that on average, only 1/3rd of its values are
689 // non-zero. This is algorithm 3.30 from [GECC].
690 //
691 // Essentially, this makes it possible to minimize the number of operations
692 // since the resulting ints returned will be at least 50% 0s.
693 func NAF(k []byte) ([]byte, []byte) {
694 // The essence of this algorithm is that whenever we have consecutive 1s
695 // in the binary, we want to put a -1 in the lowest bit and get a bunch
696 // of 0s up to the highest bit of consecutive 1s. This is due to this
697 // identity:
698 // 2^n + 2^(n-1) + 2^(n-2) + ... + 2^(n-k) = 2^(n+1) - 2^(n-k)
699 //
700 // The algorithm thus may need to go 1 more bit than the length of the
701 // bits we actually have, hence bits being 1 bit longer than was
702 // necessary. Since we need to know whether adding will cause a carry,
703 // we go from right-to-left in this addition.
704 var carry, curIsOne, nextIsOne bool
705 // these default to zero
706 retPos := make([]byte, len(k)+1)
707 retNeg := make([]byte, len(k)+1)
708 for i := len(k) - 1; i >= 0; i-- {
709 curByte := k[i]
710 for j := uint(0); j < 8; j++ {
711 curIsOne = curByte&1 == 1
712 if j == 7 {
713 if i == 0 {
714 nextIsOne = false
715 } else {
716 nextIsOne = k[i-1]&1 == 1
717 }
718 } else {
719 nextIsOne = curByte&2 == 2
720 }
721 if carry {
722 if curIsOne {
723 // This bit is 1, so continue to carry
724 // and don't need to do anything.
725 } else {
726 // We've hit a 0 after some number of
727 // 1s.
728 if nextIsOne {
729 // Start carrying again since
730 // a new sequence of 1s is
731 // starting.
732 retNeg[i+1] += 1 << j
733 } else {
734 // Stop carrying since 1s have
735 // stopped.
736 carry = false
737 retPos[i+1] += 1 << j
738 }
739 }
740 } else if curIsOne {
741 if nextIsOne {
742 // If this is the start of at least 2
743 // consecutive 1s, set the current one
744 // to -1 and start carrying.
745 retNeg[i+1] += 1 << j
746 carry = true
747 } else {
748 // This is a singleton, not consecutive
749 // 1s.
750 retPos[i+1] += 1 << j
751 }
752 }
753 curByte >>= 1
754 }
755 }
756 if carry {
757 retPos[0] = 1
758 return retPos, retNeg
759 }
760 return retPos[1:], retNeg[1:]
761 }
762 763 // ScalarMult returns k*(Bx, By) where k is a big endian integer.
764 // Part of the elliptic.Curve interface.
765 func (curve *KoblitzCurve) ScalarMult(Bx, By *big.Int, k []byte) (*big.Int, *big.Int) {
766 // Point Q = ∞ (point at infinity).
767 qx, qy, qz := new(fieldVal), new(fieldVal), new(fieldVal)
768 769 // Decompose K into k1 and k2 in order to halve the number of EC ops.
770 // See Algorithm 3.74 in [GECC].
771 k1, k2, signK1, signK2 := curve.splitK(curve.moduloReduce(k))
772 773 // The main equation here to remember is:
774 // k * P = k1 * P + k2 * ϕ(P)
775 //
776 // P1 below is P in the equation, P2 below is ϕ(P) in the equation
777 p1x, p1y := curve.bigAffineToField(Bx, By)
778 p1yNeg := new(fieldVal).NegateVal(p1y, 1)
779 p1z := new(fieldVal).SetInt(1)
780 781 // NOTE: ϕ(x,y) = (βx,y). The Jacobian z coordinate is 1, so this math
782 // goes through.
783 p2x := new(fieldVal).Mul2(p1x, curve.beta)
784 p2y := new(fieldVal).Set(p1y)
785 p2yNeg := new(fieldVal).NegateVal(p2y, 1)
786 p2z := new(fieldVal).SetInt(1)
787 788 // Flip the positive and negative values of the points as needed
789 // depending on the signs of k1 and k2. As mentioned in the equation
790 // above, each of k1 and k2 are multiplied by the respective point.
791 // Since -k * P is the same thing as k * -P, and the group law for
792 // elliptic curves states that P(x, y) = -P(x, -y), it's faster and
793 // simplifies the code to just make the point negative.
794 if signK1 == -1 {
795 p1y, p1yNeg = p1yNeg, p1y
796 }
797 if signK2 == -1 {
798 p2y, p2yNeg = p2yNeg, p2y
799 }
800 801 // NAF versions of k1 and k2 should have a lot more zeros.
802 //
803 // The Pos version of the bytes contain the +1s and the Neg versions
804 // contain the -1s.
805 k1PosNAF, k1NegNAF := NAF(k1)
806 k2PosNAF, k2NegNAF := NAF(k2)
807 k1Len := len(k1PosNAF)
808 k2Len := len(k2PosNAF)
809 810 m := k1Len
811 if m < k2Len {
812 m = k2Len
813 }
814 815 // Add left-to-right using the NAF optimization. See algorithm 3.77
816 // from [GECC]. This should be faster overall since there will be a lot
817 // more instances of 0, hence reducing the number of Jacobian additions
818 // at the cost of 1 possible extra doubling.
819 var k1BytePos, k1ByteNeg, k2BytePos, k2ByteNeg byte
820 for i := 0; i < m; i++ {
821 // Since we're going left-to-right, pad the front with 0s.
822 if i < m-k1Len {
823 k1BytePos = 0
824 k1ByteNeg = 0
825 } else {
826 k1BytePos = k1PosNAF[i-(m-k1Len)]
827 k1ByteNeg = k1NegNAF[i-(m-k1Len)]
828 }
829 if i < m-k2Len {
830 k2BytePos = 0
831 k2ByteNeg = 0
832 } else {
833 k2BytePos = k2PosNAF[i-(m-k2Len)]
834 k2ByteNeg = k2NegNAF[i-(m-k2Len)]
835 }
836 837 for j := 7; j >= 0; j-- {
838 // Q = 2 * Q
839 curve.doubleJacobian(qx, qy, qz, qx, qy, qz)
840 841 if k1BytePos&0x80 == 0x80 {
842 curve.addJacobian(qx, qy, qz, p1x, p1y, p1z,
843 qx, qy, qz)
844 } else if k1ByteNeg&0x80 == 0x80 {
845 curve.addJacobian(qx, qy, qz, p1x, p1yNeg, p1z,
846 qx, qy, qz)
847 }
848 849 if k2BytePos&0x80 == 0x80 {
850 curve.addJacobian(qx, qy, qz, p2x, p2y, p2z,
851 qx, qy, qz)
852 } else if k2ByteNeg&0x80 == 0x80 {
853 curve.addJacobian(qx, qy, qz, p2x, p2yNeg, p2z,
854 qx, qy, qz)
855 }
856 k1BytePos <<= 1
857 k1ByteNeg <<= 1
858 k2BytePos <<= 1
859 k2ByteNeg <<= 1
860 }
861 }
862 863 // Convert the Jacobian coordinate field values back to affine big.Ints.
864 return curve.fieldJacobianToBigAffine(qx, qy, qz)
865 }
866 867 // ScalarBaseMult returns k*G where G is the base point of the group and k is a
868 // big endian integer.
869 // Part of the elliptic.Curve interface.
870 func (curve *KoblitzCurve) ScalarBaseMult(k []byte) (*big.Int, *big.Int) {
871 newK := curve.moduloReduce(k)
872 diff := len(curve.bytePoints) - len(newK)
873 874 // Point Q = ∞ (point at infinity).
875 qx, qy, qz := new(fieldVal), new(fieldVal), new(fieldVal)
876 877 // curve.bytePoints has all 256 byte points for each 8-bit window. The
878 // strategy is to add up the byte points. This is best understood by
879 // expressing k in base-256 which it already sort of is.
880 // Each "digit" in the 8-bit window can be looked up using bytePoints
881 // and added together.
882 for i, byteVal := range newK {
883 p := curve.bytePoints[diff+i][byteVal]
884 curve.addJacobian(qx, qy, qz, &p[0], &p[1], &p[2], qx, qy, qz)
885 }
886 return curve.fieldJacobianToBigAffine(qx, qy, qz)
887 }
888 889 // QPlus1Div4 returns the (P+1)/4 constant for the curve for use in calculating
890 // square roots via exponentiation.
891 //
892 // DEPRECATED: The actual value returned is (P+1)/4, where as the original
893 // method name implies that this value is (((P+1)/4)+1)/4. This method is kept
894 // to maintain backwards compatibility of the API. Use Q() instead.
895 func (curve *KoblitzCurve) QPlus1Div4() *big.Int {
896 return curve.q
897 }
898 899 // Q returns the (P+1)/4 constant for the curve for use in calculating square
900 // roots via exponentiation.
901 func (curve *KoblitzCurve) Q() *big.Int {
902 return curve.q
903 }
904 905 var initonce sync.Once
906 var secp256k1 KoblitzCurve
907 908 func initAll() {
909 initS256()
910 }
911 912 // fromHex converts the passed hex string into a big integer pointer and will
913 // panic is there is an error. This is only provided for the hard-coded
914 // constants so errors in the source code can bet detected. It will only (and
915 // must only) be called for initialization purposes.
916 func fromHex(s string) *big.Int {
917 r, ok := new(big.Int).SetString(s, 16)
918 if !ok {
919 panic("invalid hex in source file: " + s)
920 }
921 return r
922 }
923 924 func initS256() {
925 // Curve parameters taken from [SECG] section 2.4.1.
926 secp256k1.CurveParams = new(elliptic.CurveParams)
927 secp256k1.P = fromHex("FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFEFFFFFC2F")
928 secp256k1.N = fromHex("FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFEBAAEDCE6AF48A03BBFD25E8CD0364141")
929 secp256k1.B = fromHex("0000000000000000000000000000000000000000000000000000000000000007")
930 secp256k1.Gx = fromHex("79BE667EF9DCBBAC55A06295CE870B07029BFCDB2DCE28D959F2815B16F81798")
931 secp256k1.Gy = fromHex("483ADA7726A3C4655DA4FBFC0E1108A8FD17B448A68554199C47D08FFB10D4B8")
932 secp256k1.BitSize = 256
933 // Curve name taken from https://safecurves.cr.yp.to/.
934 secp256k1.Name = "secp256k1"
935 secp256k1.q = new(big.Int).Div(new(big.Int).Add(secp256k1.P,
936 big.NewInt(1)), big.NewInt(4))
937 secp256k1.H = 1
938 secp256k1.halfOrder = new(big.Int).Rsh(secp256k1.N, 1)
939 secp256k1.fieldB = new(fieldVal).SetByteSlice(secp256k1.B.Bytes())
940 941 // Provided for convenience since this gets computed repeatedly.
942 secp256k1.byteSize = secp256k1.BitSize / 8
943 944 // Deserialize and set the pre-computed table used to accelerate scalar
945 // base multiplication. This is hard-coded data, so any errors are
946 // panics because it means something is wrong in the source code.
947 if err := loadS256BytePoints(); err != nil {
948 panic(err)
949 }
950 951 // Next 6 constants are from Hal Finney's bitcointalk.org post:
952 // https://bitcointalk.org/index.php?topic=3238.msg45565#msg45565
953 // May he rest in peace.
954 //
955 // They have also been independently derived from the code in the
956 // EndomorphismVectors function in gensecp256k1.go.
957 secp256k1.lambda = fromHex("5363AD4CC05C30E0A5261C028812645A122E22EA20816678DF02967C1B23BD72")
958 secp256k1.beta = new(fieldVal).SetHex("7AE96A2B657C07106E64479EAC3434E99CF0497512F58995C1396C28719501EE")
959 secp256k1.a1 = fromHex("3086D221A7D46BCDE86C90E49284EB15")
960 secp256k1.b1 = fromHex("-E4437ED6010E88286F547FA90ABFE4C3")
961 secp256k1.a2 = fromHex("114CA50F7A8E2F3F657C1108D9D44CFD8")
962 secp256k1.b2 = fromHex("3086D221A7D46BCDE86C90E49284EB15")
963 964 // Alternatively, we can use the parameters below, however, they seem
965 // to be about 8% slower.
966 // secp256k1.lambda = fromHex("AC9C52B33FA3CF1F5AD9E3FD77ED9BA4A880B9FC8EC739C2E0CFC810B51283CE")
967 // secp256k1.beta = new(fieldVal).SetHex("851695D49A83F8EF919BB86153CBCB16630FB68AED0A766A3EC693D68E6AFA40")
968 // secp256k1.a1 = fromHex("E4437ED6010E88286F547FA90ABFE4C3")
969 // secp256k1.b1 = fromHex("-3086D221A7D46BCDE86C90E49284EB15")
970 // secp256k1.a2 = fromHex("3086D221A7D46BCDE86C90E49284EB15")
971 // secp256k1.b2 = fromHex("114CA50F7A8E2F3F657C1108D9D44CFD8")
972 }
973 974 // S256 returns a Curve which implements secp256k1.
975 func S256() *KoblitzCurve {
976 initonce.Do(initAll)
977 return &secp256k1
978 }
979