Optimize ScalarMult using endomorphism

This implements a speedup to ScalarMult using the endomorphism available to secp256k1.

Note the constants lambda, beta, a1, b1, a2 and b2 are from here:

https://bitcointalk.org/index.php?topic=3238.0

For whatever reason, rounding to the nearest integer does not work. I had to add 3 regardless so the number of bytes for each K is 17, which does work.

Preliminary tests indicate a speedup of between 17%-20% (BenchScalarMult).

More speedup can probably be achieved once splitK uses something more like what fieldVal uses. Unfortunately, the prime for this math is the order of G (N), not P.

Note the NAF optimization was specifically not done as that's the purview of another issue.

This closes #1

btcec.go

@@ -39,6 +39,12 @@ type KoblitzCurve struct {
 	*elliptic.CurveParams
 	q          *big.Int
 	H          int // cofactor of the curve.
+	lambda     *big.Int
+	beta       *fieldVal
+	a1         *big.Int
+	b1         *big.Int
+	a2         *big.Int
+	b2         *big.Int
 	bytePoints *[32][256][3]fieldVal
 }
 
@@ -594,45 +600,90 @@ func (curve *KoblitzCurve) Double(x1, y1 *big.Int) (*big.Int, *big.Int) {
 	return curve.fieldJacobianToBigAffine(fx3, fy3, fz3)
 }
 
+// splitK returns a balanced length-two representation of k
+func (curve *KoblitzCurve) splitK(k []byte) ([]byte, []byte) {
+
+	bigIntK, c1, c2, tmp1, tmp2, k1, k2 := big.NewInt(0), big.NewInt(0), big.NewInt(0), big.NewInt(0), big.NewInt(0), big.NewInt(0), big.NewInt(0)
+	bigIntK.SetBytes(k)
+
+	c1.Mul(curve.b2, bigIntK)
+	c1.Div(c1, curve.N)
+	c1.Add(c1, big.NewInt(3))
+	c2.Mul(curve.b1, bigIntK)
+	c2.Div(c2, curve.N)
+	c2.Add(c2, big.NewInt(3))
+	tmp1.Mul(c1, curve.a1)
+	tmp2.Mul(c2, curve.a2)
+	k1.Sub(bigIntK, tmp1)
+	k1.Add(k1, tmp2)
+	tmp1.Mul(c1, curve.b1)
+	tmp2.Mul(c2, curve.b2)
+	k2.Sub(tmp2, tmp1)
+
+	return k1.Bytes(), k2.Bytes()
+}
+
 // ScalarMult returns k*(Bx, By) where k is a big endian integer.
 // Part of the elliptic.Curve interface.
 func (curve *KoblitzCurve) ScalarMult(Bx, By *big.Int, k []byte) (*big.Int, *big.Int) {
-	// This uses the left to right binary method for point multiplication:
+
+	newK := curve.moduloReduce(k)
+	k1, k2 := curve.splitK(newK)
+	m := len(k1)
+	if len(k2) > m {
+		m = len(k2)
+	}
 
 	// Point Q = ∞ (point at infinity).
 	qx, qy, qz := new(fieldVal), new(fieldVal), new(fieldVal)
 
-	// Point P = the point to multiply the scalar with.
-	px, py := curve.bigAffineToField(Bx, By)
-	pz := new(fieldVal).SetInt(1)
-
-	// Double and add as necessary depending on the bits set in the scalar.
-	for _, byteVal := range k {
+	// Point P1 = point for k1 to multiply against
+	p1x, p1y := curve.bigAffineToField(Bx, By)
+	p1z := new(fieldVal).SetInt(1)
+
+	// Point P2 = point for k2 to multiply against
+	p2x := new(fieldVal).Set(p1x).Mul(curve.beta)
+	p2y := new(fieldVal).Set(p1y)
+	p2z := new(fieldVal).SetInt(1)
+
+	// use the left to right binary method on both k1 and k2
+	var byteVal1, byteVal2 byte
+	for i := 0; i < m; i++ {
+		if i < len(k1) {
+			byteVal1 = k1[i]
+		} else {
+			byteVal1 = 0
+		}
+		if i < len(k2) {
+			byteVal2 = k2[i]
+		} else {
+			byteVal2 = 0
+		}
 		for bitNum := 0; bitNum < 8; bitNum++ {
 			// Q = 2*Q
 			curve.doubleJacobian(qx, qy, qz, qx, qy, qz)
-			if byteVal&0x80 == 0x80 {
-				// Q = Q + P
-				curve.addJacobian(qx, qy, qz, px, py, pz, qx,
-					qy, qz)
+			if byteVal1&0x80 == 0x80 {
+				// Q = Q + P1
+				curve.addJacobian(qx, qy, qz, p1x, p1y, p1z, qx, qy, qz)
+			}
+			if byteVal2&0x80 == 0x80 {
+				// Q = Q + P2
+				curve.addJacobian(qx, qy, qz, p2x, p2y, p2z, qx, qy, qz)
 			}
-			byteVal <<= 1
+			byteVal1 <<= 1
+			byteVal2 <<= 1
 		}
 	}
 
 	// Convert the Jacobian coordinate field values back to affine big.Ints.
 	return curve.fieldJacobianToBigAffine(qx, qy, qz)
 }
 
-// ScalarBaseMult returns k*G where G is the base point of the group and k is a
-// big endian integer.
-// Part of the elliptic.Curve interface.
-func (curve *KoblitzCurve) ScalarBaseMult(k []byte) (*big.Int, *big.Int) {
-
+func (curve *KoblitzCurve) moduloReduce(k []byte) []byte {
 	var newK []byte
 	// Since the order of G is curve.N, we can use a much smaller number
 	// by doing modulo curve.N
-	if len(k) > len(curve.bytePoints) {
+	if len(k) > curve.BitSize/8 {
 		// reduce k by performing modulo curve.N
 		tmpK := big.NewInt(0)
 		tmpK.SetBytes(k)
@@ -641,6 +692,15 @@ func (curve *KoblitzCurve) ScalarBaseMult(k []byte) (*big.Int, *big.Int) {
 	} else {
 		newK = k
 	}
+	return newK
+}
+
+// ScalarBaseMult returns k*G where G is the base point of the group and k is a
+// big endian integer.
+// Part of the elliptic.Curve interface.
+func (curve *KoblitzCurve) ScalarBaseMult(k []byte) (*big.Int, *big.Int) {
+
+	newK := curve.moduloReduce(k)
 
 	diff := len(curve.bytePoints) - len(newK)
 
@@ -685,6 +745,16 @@ func initS256() {
 	secp256k1.H = 1
 	secp256k1.q = new(big.Int).Div(new(big.Int).Add(secp256k1.P,
 		big.NewInt(1)), big.NewInt(4))
+
+	secp256k1.lambda, _ = new(big.Int).SetString("5363AD4CC05C30E0A5261C028812645A122E22EA20816678DF02967C1B23BD72", 16)
+
+	secp256k1.beta = new(fieldVal).SetHex("7AE96A2B657C07106E64479EAC3434E99CF0497512F58995C1396C28719501EE")
+
+	secp256k1.a1, _ = new(big.Int).SetString("3086D221A7D46BCDE86C90E49284EB15", 16)
+	secp256k1.b1, _ = new(big.Int).SetString("-E4437ED6010E88286F547FA90ABFE4C3", 16)
+	secp256k1.a2, _ = new(big.Int).SetString("114CA50F7A8E2F3F657C1108D9D44CFD8", 16)
+	secp256k1.b2, _ = new(big.Int).SetString("3086D221A7D46BCDE86C90E49284EB15", 16)
+
 	secp256k1.bytePoints = &secp256k1BytePoints
 }
 

btcec_test.go

@@ -577,6 +577,33 @@ func TestBaseMultVerify(t *testing.T) {
 	}
 }
 
+func TestScalarMult(t *testing.T) {
+	// Strategy for this test:
+	// Get a random exponent from the generator point at first
+	// This creates a new point which is used in the next iteration
+	// Use another random exponent on the new point.
+	// We use BaseMult to verify by multiplying the exponents together (mod N)
+	s256 := btcec.S256()
+	x, y := s256.Gx, s256.Gy
+	exponent := big.NewInt(1)
+	for i := 0; i < 1024; i++ {
+		data := make([]byte, 32)
+		_, err := rand.Read(data)
+		if err != nil {
+			t.Fatalf("failed to read random data at %d", i)
+			continue
+		}
+		x, y = s256.ScalarMult(x, y, data)
+		exponent.Mul(exponent, new(big.Int).SetBytes(data))
+		exponent.Mod(exponent, s256.N)
+		xWant, yWant := s256.ScalarBaseMult(exponent.Bytes())
+
+		if x.Cmp(xWant) != 0 || y.Cmp(yWant) != 0 {
+			t.Errorf("%d: bad output for %X: got (%X, %X), want (%X, %X)", i, data, x, y, xWant, yWant)
+		}
+	}
+}
+
 //TODO: test more curves?
 func BenchmarkBaseMult(b *testing.B) {
 	b.ResetTimer()