1/* Copyright (c) 2018, Google Inc.
2 *
3 * Permission to use, copy, modify, and/or distribute this software for any
4 * purpose with or without fee is hereby granted, provided that the above
5 * copyright notice and this permission notice appear in all copies.
6 *
7 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
8 * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
9 * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY
10 * SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
11 * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION
12 * OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
13 * CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. */
14
15package main
16
17import (
18	"crypto/aes"
19	"crypto/cipher"
20	"crypto/elliptic"
21	"crypto/rand"
22	"fmt"
23	"io"
24	"math/big"
25)
26
27var (
28	p256                  elliptic.Curve
29	zero, one, p, R, Rinv *big.Int
30	deterministicRand     io.Reader
31)
32
33type coordinates int
34
35const (
36	affine coordinates = iota
37	jacobian
38)
39
40func init() {
41	p256 = elliptic.P256()
42
43	zero = new(big.Int)
44	one = new(big.Int).SetInt64(1)
45
46	p = p256.Params().P
47
48	R = new(big.Int)
49	R.SetBit(R, 256, 1)
50	R.Mod(R, p)
51
52	Rinv = new(big.Int).ModInverse(R, p)
53
54	deterministicRand = newDeterministicRand()
55}
56
57func modMul(z, x, y *big.Int) *big.Int {
58	z.Mul(x, y)
59	return z.Mod(z, p)
60}
61
62func toMontgomery(z, x *big.Int) *big.Int {
63	return modMul(z, x, R)
64}
65
66func fromMontgomery(z, x *big.Int) *big.Int {
67	return modMul(z, x, Rinv)
68}
69
70func isAffineInfinity(x, y *big.Int) bool {
71	// Infinity, in affine coordinates, is represented as (0, 0) by
72	// both Go and p256-x86_64-asm.pl.
73	return x.Sign() == 0 && y.Sign() == 0
74}
75
76func randNonZeroInt(max *big.Int) *big.Int {
77	for {
78		r, err := rand.Int(deterministicRand, max)
79		if err != nil {
80			panic(err)
81		}
82		if r.Sign() != 0 {
83			return r
84		}
85	}
86}
87
88func randPoint() (x, y *big.Int) {
89	k := randNonZeroInt(p256.Params().N)
90	return p256.ScalarBaseMult(k.Bytes())
91}
92
93func toJacobian(xIn, yIn *big.Int) (x, y, z *big.Int) {
94	if isAffineInfinity(xIn, yIn) {
95		// The Jacobian representation of infinity has Z = 0. Depending
96		// on the implementation, X and Y may be further constrained.
97		// Generalizing the curve equation to Jacobian coordinates for
98		// non-zero Z gives:
99		//
100		//   y² = x³ - 3x + b, where x = X/Z² and y = Y/Z³
101		//   Y² = X³ + aXZ⁴ + bZ⁶
102		//
103		// Taking that formula at Z = 0 gives Y² = X³. This constraint
104		// allows removing a special case in the point-on-curve check.
105		//
106		// BoringSSL, however, historically generated infinities with
107		// arbitrary X and Y and include the special case. We also have
108		// not verified that add and double preserve this
109		// property. Thus, generate test vectors with unrelated X and Y,
110		// to test that p256-x86_64-asm.pl correctly handles
111		// unconstrained representations of infinity.
112		x = randNonZeroInt(p)
113		y = randNonZeroInt(p)
114		z = zero
115		return
116	}
117
118	z = randNonZeroInt(p)
119
120	// X = xZ²
121	y = modMul(new(big.Int), z, z)
122	x = modMul(new(big.Int), xIn, y)
123
124	// Y = yZ³
125	modMul(y, y, z)
126	modMul(y, y, yIn)
127	return
128}
129
130func printMontgomery(name string, a *big.Int) {
131	a = toMontgomery(new(big.Int), a)
132	fmt.Printf("%s = %064x\n", name, a)
133}
134
135func printTestCase(ax, ay *big.Int, aCoord coordinates, bx, by *big.Int, bCoord coordinates) {
136	rx, ry := p256.Add(ax, ay, bx, by)
137
138	var az *big.Int
139	if aCoord == jacobian {
140		ax, ay, az = toJacobian(ax, ay)
141	} else if isAffineInfinity(ax, ay) {
142		az = zero
143	} else {
144		az = one
145	}
146
147	var bz *big.Int
148	if bCoord == jacobian {
149		bx, by, bz = toJacobian(bx, by)
150	} else if isAffineInfinity(bx, by) {
151		bz = zero
152	} else {
153		bz = one
154	}
155
156	fmt.Printf("Test = PointAdd\n")
157	printMontgomery("A.X", ax)
158	printMontgomery("A.Y", ay)
159	printMontgomery("A.Z", az)
160	printMontgomery("B.X", bx)
161	printMontgomery("B.Y", by)
162	printMontgomery("B.Z", bz)
163	printMontgomery("Result.X", rx)
164	printMontgomery("Result.Y", ry)
165	fmt.Printf("\n")
166}
167
168func main() {
169	fmt.Printf("# ∞ + ∞ = ∞.\n")
170	printTestCase(zero, zero, affine, zero, zero, affine)
171
172	fmt.Printf("# ∞ + ∞ = ∞, with an alternate representation of ∞.\n")
173	printTestCase(zero, zero, jacobian, zero, zero, jacobian)
174
175	gx, gy := p256.Params().Gx, p256.Params().Gy
176	fmt.Printf("# g + ∞ = g.\n")
177	printTestCase(gx, gy, affine, zero, zero, affine)
178
179	fmt.Printf("# g + ∞ = g, with an alternate representation of ∞.\n")
180	printTestCase(gx, gy, affine, zero, zero, jacobian)
181
182	fmt.Printf("# g + -g = ∞.\n")
183	minusGy := new(big.Int).Sub(p, gy)
184	printTestCase(gx, gy, affine, gx, minusGy, affine)
185
186	fmt.Printf("# Test some random Jacobian sums.\n")
187	for i := 0; i < 4; i++ {
188		ax, ay := randPoint()
189		bx, by := randPoint()
190		printTestCase(ax, ay, jacobian, bx, by, jacobian)
191	}
192
193	fmt.Printf("# Test some random Jacobian doublings.\n")
194	for i := 0; i < 4; i++ {
195		ax, ay := randPoint()
196		printTestCase(ax, ay, jacobian, ax, ay, jacobian)
197	}
198
199	fmt.Printf("# Test some random affine sums.\n")
200	for i := 0; i < 4; i++ {
201		ax, ay := randPoint()
202		bx, by := randPoint()
203		printTestCase(ax, ay, affine, bx, by, affine)
204	}
205
206	fmt.Printf("# Test some random affine doublings.\n")
207	for i := 0; i < 4; i++ {
208		ax, ay := randPoint()
209		printTestCase(ax, ay, affine, ax, ay, affine)
210	}
211}
212
213type deterministicRandom struct {
214	stream cipher.Stream
215}
216
217func newDeterministicRand() io.Reader {
218	block, err := aes.NewCipher(make([]byte, 128/8))
219	if err != nil {
220		panic(err)
221	}
222	stream := cipher.NewCTR(block, make([]byte, block.BlockSize()))
223	return &deterministicRandom{stream}
224}
225
226func (r *deterministicRandom) Read(b []byte) (n int, err error) {
227	for i := range b {
228		b[i] = 0
229	}
230	r.stream.XORKeyStream(b, b)
231	return len(b), nil
232}
233