cln16sidh: add three-point ladder

curve.go

@@ -176,20 +176,29 @@ func (xQ *ProjectivePoint) Pow3k(curve *ProjectiveCurveParameters, xP *Projectiv
 	return xQ
 }
 
+// Given x(P) and a scalar m in little-endian bytes, compute x([m]P) using the
+// Montgomery ladder.  This is described in Algorithm 8 of Costello-Smith.
+//
+// This function's execution time is dependent only on the byte-length of the
+// input scalar.  All scalars of the same input length execute in uniform time.
+// The scalar can be padded with zero bytes to ensure a uniform length.
+//
+// Safe to overlap the source with the destination.
 func (xQ *ProjectivePoint) ScalarMult(curve *ProjectiveCurveParameters, xP *ProjectivePoint, scalar []uint8) *ProjectivePoint {
 	cachedParams := curve.cachedParams()
 	var x0, x1, tmp ProjectivePoint
+
 	x0.x.One()
 	x0.z.Zero()
 	x1 = *xP
-	prevBit := uint8(0)
+
 	// Iterate over the bits of the scalar, top to bottom
+	prevBit := uint8(0)
 	for i := len(scalar) - 1; i >= 0; i-- {
 		scalarByte := scalar[i]
 		for j := 7; j >= 0; j-- {
 			bit := (scalarByte >> uint(j)) & 0x1
 			ProjectivePointConditionalSwap(&x0, &x1, (bit ^ prevBit))
-			// could avoid use of tmp by having unified double/add
 			tmp.Double(&x0, &cachedParams)
 			x1.Add(&x0, &x1, xP)
 			x0 = tmp
@@ -201,3 +210,78 @@ func (xQ *ProjectivePoint) ScalarMult(curve *ProjectiveCurveParameters, xP *Proj
 	*xQ = x0
 	return xQ
 }
+
+// Given x(P), x(Q), x(P-Q), as well as a scalar m in little-endian bytes,
+// compute x(P + [m]Q) using the "three-point ladder" of de Feo, Jao, and Plut.
+//
+// Safe to overlap the source with the destination.
+//
+// This function's execution time is dependent only on the byte-length of the
+// input scalar.  All scalars of the same input length execute in uniform time.
+// The scalar can be padded with zero bytes to ensure a uniform length.
+//
+// The algorithm, as described in de Feo-Jao-Plut, is as follows:
+//
+// (x0, x1, x2) <--- (x(O), x(Q), x(P))
+//
+// for i = |m| down to 0, indexing the bits of m:
+//     Invariant: (x0, x1, x2) == (x( [t]Q ), x( [t+1]Q ), x( P + [t]Q ))
+//          where t = m//2^i is the high bits of m, starting at i
+//     if m_i == 0:
+//         (x0, x1, x2) <--- (xDBL(x0), xADD(x1, x0, x(Q)), xADD(x2, x0, x(P)))
+//         Invariant: (x0, x1, x2) == (x( [2t]Q ), x( [2t+1]Q ), x( P + [2t]Q ))
+//                                 == (x( [t']Q ), x( [t'+1]Q ), x( P + [t']Q ))
+//              where t' = m//2^{i-1} is the high bits of m, starting at i-1
+//     if m_i == 1:
+//         (x0, x1, x2) <--- (xADD(x1, x0, x(Q)), xDBL(x1), xADD(x2, x1, x(P-Q)))
+//         Invariant: (x0, x1, x2) == (x( [2t+1]Q ), x( [2t+2]Q ), x( P + [2t+1]Q ))
+//                                 == (x( [t']Q ),   x( [t'+1]Q ), x( P + [t']Q ))
+//              where t' = m//2^{i-1} is the high bits of m, starting at i-1
+// return x2
+//
+// Notice that the roles of (x0,x1) and (x(P), x(P-Q)) swap depending on the
+// current bit of the scalar.  Instead of swapping which operations we do, we
+// can swap variable names, producing the following uniform algorithm:
+//
+// (x0, x1, x2) <--- (x(O), x(Q), x(P))
+// (y0, y1) <--- (x(P), x(P-Q))
+//
+// for i = |m| down to 0, indexing the bits of m:
+//      (x0, x1) <--- SWAP( m_{i+1} xor m_i, (x0,x1) )
+//      (y0, y1) <--- SWAP( m_{i+1} xor m_i, (y0,y1) )
+//      (x0, x1, x2) <--- ( xDBL(x0), xADD(x1,x0,x(Q)), xADD(x2, x0, y0) )
+//
+// return x2
+//
+func (xR *ProjectivePoint) ThreePointLadder(curve *ProjectiveCurveParameters, xP, xQ, xPmQ *ProjectivePoint, scalar []uint8) *ProjectivePoint {
+	cachedParams := curve.cachedParams()
+	var x0, x1, x2, y0, y1, tmp ProjectivePoint
+
+	// (x0, x1, x2) <--- (x(O), x(Q), x(P))
+	x0.x.One()
+	x0.z.Zero()
+	x1 = *xQ
+	x2 = *xP
+	// (y0, y1) <--- (x(P), x(P-Q))
+	y0 = *xP
+	y1 = *xPmQ
+
+	// Iterate over the bits of the scalar, top to bottom
+	prevBit := uint8(0)
+	for i := len(scalar) - 1; i >= 0; i-- {
+		scalarByte := scalar[i]
+		for j := 7; j >= 0; j-- {
+			bit := (scalarByte >> uint(j)) & 0x1
+			ProjectivePointConditionalSwap(&x0, &x1, (bit ^ prevBit))
+			ProjectivePointConditionalSwap(&y0, &y1, (bit ^ prevBit))
+			x2.Add(&x2, &x0, &y0) // = xADD(x2, x0, y0)
+			tmp.Double(&x0, &cachedParams)
+			x1.Add(&x1, &x0, xQ) // = xADD(x1, x0, x(Q))
+			x0 = tmp             // = xDBL(x0)
+			prevBit = bit
+		}
+	}
+
+	*xR = x2
+	return xR
+}

curve_test.go

@@ -142,6 +142,29 @@ func TestScalarMultVersusSage(t *testing.T) {
 	}
 }
 
+func TestThreePointLadderVersusSage(t *testing.T) {
+	var xP, xQ, xPmQ, xR ProjectivePoint
+
+	xP.x = ExtensionFieldElement{a: fp751Element{0xe8d05f30aac47247, 0x576ec00c55441de7, 0xbf1a8ec5fe558518, 0xd77cb17f77515881, 0x8e9852837ee73ec4, 0x8159634ad4f44a6b, 0x2e4eb5533a798c5, 0x9be8c4354d5bc849, 0xf47dc61806496b84, 0x25d0e130295120e0, 0xdbef54095f8139e3, 0x5a724f20862c}, b: fp751Element{0x3ca30d7623602e30, 0xfb281eddf45f07b7, 0xd2bf62d5901a45bc, 0xc67c9baf86306dd2, 0x4e2bd93093f538ca, 0xcfd92075c25b9cbe, 0xceafe9a3095bcbab, 0x7d928ad380c85414, 0x37c5f38b2afdc095, 0x75325899a7b779f4, 0xf130568249f20fdd, 0x178f264767d1}}
+	xP.z.One()
+
+	xQ.x = ExtensionFieldElement{a: fp751Element{0x2b71a2a93ad1e10e, 0xf0b9842a92cfb333, 0xae17373615a27f5c, 0x3039239f428330c4, 0xa0c4b735ed7dcf98, 0x6e359771ddf6af6a, 0xe986e4cac4584651, 0x8233a2b622d5518, 0xbfd67bf5f06b818b, 0xdffe38d0f5b966a6, 0xa86b36a3272ee00a, 0x193e2ea4f68f}, b: fp751Element{0x5a0f396459d9d998, 0x479f42250b1b7dda, 0x4016b57e2a15bf75, 0xc59f915203fa3749, 0xd5f90257399cf8da, 0x1fb2dadfd86dcef4, 0x600f20e6429021dc, 0x17e347d380c57581, 0xc1b0d5fa8fe3e440, 0xbcf035330ac20e8, 0x50c2eb5f6a4f03e6, 0x86b7c4571}}
+	xQ.z.One()
+
+	xPmQ.x = ExtensionFieldElement{a: fp751Element{0x4aafa9f378f7b5ff, 0x1172a683aa8eee0, 0xea518d8cbec2c1de, 0xe191bcbb63674557, 0x97bc19637b259011, 0xdbeae5c9f4a2e454, 0x78f64d1b72a42f95, 0xe71cb4ea7e181e54, 0xe4169d4c48543994, 0x6198c2286a98730f, 0xd21d675bbab1afa5, 0x2e7269fce391}, b: fp751Element{0x23355783ce1d0450, 0x683164cf4ce3d93f, 0xae6d1c4d25970fd8, 0x7807007fb80b48cf, 0xa005a62ec2bbb8a2, 0x6b5649bd016004cb, 0xbb1a13fa1330176b, 0xbf38e51087660461, 0xe577fddc5dd7b930, 0x5f38116f56947cd3, 0x3124f30b98c36fde, 0x4ca9b6e6db37}}
+	xPmQ.z.One()
+
+	xR.ThreePointLadder(&curve, &xP, &xQ, &xPmQ, mScalarBytes[:])
+
+	affine_xR := xR.toAffine()
+
+	sageAffine_xR := ExtensionFieldElement{a: fp751Element{0x729465ba800d4fd5, 0x9398015b59e514a1, 0x1a59dd6be76c748e, 0x1a7db94eb28dd55c, 0x444686e680b1b8ec, 0xcc3d4ace2a2454ff, 0x51d3dab4ec95a419, 0xc3b0f33594acac6a, 0x9598a74e7fd44f8a, 0x4fbf8c638f1c2e37, 0x844e347033052f51, 0x6cd6de3eafcf}, b: fp751Element{0x85da145412d73430, 0xd83c0e3b66eb3232, 0xd08ff2d453ec1369, 0xa64aaacfdb395b13, 0xe9cba211a20e806e, 0xa4f80b175d937cfc, 0x556ce5c64b1f7937, 0xb59b39ea2b3fdf7a, 0xc2526b869a4196b3, 0x8dad90bca9371750, 0xdfb4a30c9d9147a2, 0x346d2130629b}}
+
+	if !affine_xR.VartimeEq(&sageAffine_xR) {
+		t.Error("\nExpected\n", sageAffine_xR, "\nfound\n", affine_xR)
+	}
+}
+
 func TestPointTripleVersusAddDouble(t *testing.T) {
 	tripleEqualsAddDouble := func(curve ProjectiveCurveParameters, P ProjectivePoint) bool {