xref: /openbmc/linux/lib/crypto/mpi/ec.c (revision 83146efc)
1 /* ec.c -  Elliptic Curve functions
2  * Copyright (C) 2007 Free Software Foundation, Inc.
3  * Copyright (C) 2013 g10 Code GmbH
4  *
5  * This file is part of Libgcrypt.
6  *
7  * Libgcrypt is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU Lesser General Public License as
9  * published by the Free Software Foundation; either version 2.1 of
10  * the License, or (at your option) any later version.
11  *
12  * Libgcrypt is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  * GNU Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with this program; if not, see <http://www.gnu.org/licenses/>.
19  */
20 
21 #include "mpi-internal.h"
22 #include "longlong.h"
23 
24 #define point_init(a)  mpi_point_init((a))
25 #define point_free(a)  mpi_point_free_parts((a))
26 
27 #define log_error(fmt, ...) pr_err(fmt, ##__VA_ARGS__)
28 #define log_fatal(fmt, ...) pr_err(fmt, ##__VA_ARGS__)
29 
30 #define DIM(v) (sizeof(v)/sizeof((v)[0]))
31 
32 
33 /* Create a new point option.  NBITS gives the size in bits of one
34  * coordinate; it is only used to pre-allocate some resources and
35  * might also be passed as 0 to use a default value.
36  */
37 MPI_POINT mpi_point_new(unsigned int nbits)
38 {
39 	MPI_POINT p;
40 
41 	(void)nbits;  /* Currently not used.  */
42 
43 	p = kmalloc(sizeof(*p), GFP_KERNEL);
44 	if (p)
45 		mpi_point_init(p);
46 	return p;
47 }
48 EXPORT_SYMBOL_GPL(mpi_point_new);
49 
50 /* Release the point object P.  P may be NULL. */
51 void mpi_point_release(MPI_POINT p)
52 {
53 	if (p) {
54 		mpi_point_free_parts(p);
55 		kfree(p);
56 	}
57 }
58 EXPORT_SYMBOL_GPL(mpi_point_release);
59 
60 /* Initialize the fields of a point object.  gcry_mpi_point_free_parts
61  * may be used to release the fields.
62  */
63 void mpi_point_init(MPI_POINT p)
64 {
65 	p->x = mpi_new(0);
66 	p->y = mpi_new(0);
67 	p->z = mpi_new(0);
68 }
69 EXPORT_SYMBOL_GPL(mpi_point_init);
70 
71 /* Release the parts of a point object. */
72 void mpi_point_free_parts(MPI_POINT p)
73 {
74 	mpi_free(p->x); p->x = NULL;
75 	mpi_free(p->y); p->y = NULL;
76 	mpi_free(p->z); p->z = NULL;
77 }
78 EXPORT_SYMBOL_GPL(mpi_point_free_parts);
79 
80 /* Set the value from S into D.  */
81 static void point_set(MPI_POINT d, MPI_POINT s)
82 {
83 	mpi_set(d->x, s->x);
84 	mpi_set(d->y, s->y);
85 	mpi_set(d->z, s->z);
86 }
87 
88 static void point_resize(MPI_POINT p, struct mpi_ec_ctx *ctx)
89 {
90 	size_t nlimbs = ctx->p->nlimbs;
91 
92 	mpi_resize(p->x, nlimbs);
93 	p->x->nlimbs = nlimbs;
94 	mpi_resize(p->z, nlimbs);
95 	p->z->nlimbs = nlimbs;
96 
97 	if (ctx->model != MPI_EC_MONTGOMERY) {
98 		mpi_resize(p->y, nlimbs);
99 		p->y->nlimbs = nlimbs;
100 	}
101 }
102 
103 static void point_swap_cond(MPI_POINT d, MPI_POINT s, unsigned long swap,
104 		struct mpi_ec_ctx *ctx)
105 {
106 	mpi_swap_cond(d->x, s->x, swap);
107 	if (ctx->model != MPI_EC_MONTGOMERY)
108 		mpi_swap_cond(d->y, s->y, swap);
109 	mpi_swap_cond(d->z, s->z, swap);
110 }
111 
112 
113 /* W = W mod P.  */
114 static void ec_mod(MPI w, struct mpi_ec_ctx *ec)
115 {
116 	if (ec->t.p_barrett)
117 		mpi_mod_barrett(w, w, ec->t.p_barrett);
118 	else
119 		mpi_mod(w, w, ec->p);
120 }
121 
122 static void ec_addm(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
123 {
124 	mpi_add(w, u, v);
125 	ec_mod(w, ctx);
126 }
127 
128 static void ec_subm(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ec)
129 {
130 	mpi_sub(w, u, v);
131 	while (w->sign)
132 		mpi_add(w, w, ec->p);
133 	/*ec_mod(w, ec);*/
134 }
135 
136 static void ec_mulm(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
137 {
138 	mpi_mul(w, u, v);
139 	ec_mod(w, ctx);
140 }
141 
142 /* W = 2 * U mod P.  */
143 static void ec_mul2(MPI w, MPI u, struct mpi_ec_ctx *ctx)
144 {
145 	mpi_lshift(w, u, 1);
146 	ec_mod(w, ctx);
147 }
148 
149 static void ec_powm(MPI w, const MPI b, const MPI e,
150 		struct mpi_ec_ctx *ctx)
151 {
152 	mpi_powm(w, b, e, ctx->p);
153 	/* mpi_abs(w); */
154 }
155 
156 /* Shortcut for
157  * ec_powm(B, B, mpi_const(MPI_C_TWO), ctx);
158  * for easier optimization.
159  */
160 static void ec_pow2(MPI w, const MPI b, struct mpi_ec_ctx *ctx)
161 {
162 	/* Using mpi_mul is slightly faster (at least on amd64).  */
163 	/* mpi_powm(w, b, mpi_const(MPI_C_TWO), ctx->p); */
164 	ec_mulm(w, b, b, ctx);
165 }
166 
167 /* Shortcut for
168  * ec_powm(B, B, mpi_const(MPI_C_THREE), ctx);
169  * for easier optimization.
170  */
171 static void ec_pow3(MPI w, const MPI b, struct mpi_ec_ctx *ctx)
172 {
173 	mpi_powm(w, b, mpi_const(MPI_C_THREE), ctx->p);
174 }
175 
176 static void ec_invm(MPI x, MPI a, struct mpi_ec_ctx *ctx)
177 {
178 	if (!mpi_invm(x, a, ctx->p))
179 		log_error("ec_invm: inverse does not exist:\n");
180 }
181 
182 static void mpih_set_cond(mpi_ptr_t wp, mpi_ptr_t up,
183 		mpi_size_t usize, unsigned long set)
184 {
185 	mpi_size_t i;
186 	mpi_limb_t mask = ((mpi_limb_t)0) - set;
187 	mpi_limb_t x;
188 
189 	for (i = 0; i < usize; i++) {
190 		x = mask & (wp[i] ^ up[i]);
191 		wp[i] = wp[i] ^ x;
192 	}
193 }
194 
195 /* Routines for 2^255 - 19.  */
196 
197 #define LIMB_SIZE_25519 ((256+BITS_PER_MPI_LIMB-1)/BITS_PER_MPI_LIMB)
198 
199 static void ec_addm_25519(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
200 {
201 	mpi_ptr_t wp, up, vp;
202 	mpi_size_t wsize = LIMB_SIZE_25519;
203 	mpi_limb_t n[LIMB_SIZE_25519];
204 	mpi_limb_t borrow;
205 
206 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
207 		log_bug("addm_25519: different sizes\n");
208 
209 	memset(n, 0, sizeof(n));
210 	up = u->d;
211 	vp = v->d;
212 	wp = w->d;
213 
214 	mpihelp_add_n(wp, up, vp, wsize);
215 	borrow = mpihelp_sub_n(wp, wp, ctx->p->d, wsize);
216 	mpih_set_cond(n, ctx->p->d, wsize, (borrow != 0UL));
217 	mpihelp_add_n(wp, wp, n, wsize);
218 	wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
219 }
220 
221 static void ec_subm_25519(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
222 {
223 	mpi_ptr_t wp, up, vp;
224 	mpi_size_t wsize = LIMB_SIZE_25519;
225 	mpi_limb_t n[LIMB_SIZE_25519];
226 	mpi_limb_t borrow;
227 
228 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
229 		log_bug("subm_25519: different sizes\n");
230 
231 	memset(n, 0, sizeof(n));
232 	up = u->d;
233 	vp = v->d;
234 	wp = w->d;
235 
236 	borrow = mpihelp_sub_n(wp, up, vp, wsize);
237 	mpih_set_cond(n, ctx->p->d, wsize, (borrow != 0UL));
238 	mpihelp_add_n(wp, wp, n, wsize);
239 	wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
240 }
241 
242 static void ec_mulm_25519(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
243 {
244 	mpi_ptr_t wp, up, vp;
245 	mpi_size_t wsize = LIMB_SIZE_25519;
246 	mpi_limb_t n[LIMB_SIZE_25519*2];
247 	mpi_limb_t m[LIMB_SIZE_25519+1];
248 	mpi_limb_t cy;
249 	int msb;
250 
251 	(void)ctx;
252 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
253 		log_bug("mulm_25519: different sizes\n");
254 
255 	up = u->d;
256 	vp = v->d;
257 	wp = w->d;
258 
259 	mpihelp_mul_n(n, up, vp, wsize);
260 	memcpy(wp, n, wsize * BYTES_PER_MPI_LIMB);
261 	wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
262 
263 	memcpy(m, n+LIMB_SIZE_25519-1, (wsize+1) * BYTES_PER_MPI_LIMB);
264 	mpihelp_rshift(m, m, LIMB_SIZE_25519+1, (255 % BITS_PER_MPI_LIMB));
265 
266 	memcpy(n, m, wsize * BYTES_PER_MPI_LIMB);
267 	cy = mpihelp_lshift(m, m, LIMB_SIZE_25519, 4);
268 	m[LIMB_SIZE_25519] = cy;
269 	cy = mpihelp_add_n(m, m, n, wsize);
270 	m[LIMB_SIZE_25519] += cy;
271 	cy = mpihelp_add_n(m, m, n, wsize);
272 	m[LIMB_SIZE_25519] += cy;
273 	cy = mpihelp_add_n(m, m, n, wsize);
274 	m[LIMB_SIZE_25519] += cy;
275 
276 	cy = mpihelp_add_n(wp, wp, m, wsize);
277 	m[LIMB_SIZE_25519] += cy;
278 
279 	memset(m, 0, wsize * BYTES_PER_MPI_LIMB);
280 	msb = (wp[LIMB_SIZE_25519-1] >> (255 % BITS_PER_MPI_LIMB));
281 	m[0] = (m[LIMB_SIZE_25519] * 2 + msb) * 19;
282 	wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
283 	mpihelp_add_n(wp, wp, m, wsize);
284 
285 	m[0] = 0;
286 	cy = mpihelp_sub_n(wp, wp, ctx->p->d, wsize);
287 	mpih_set_cond(m, ctx->p->d, wsize, (cy != 0UL));
288 	mpihelp_add_n(wp, wp, m, wsize);
289 }
290 
291 static void ec_mul2_25519(MPI w, MPI u, struct mpi_ec_ctx *ctx)
292 {
293 	ec_addm_25519(w, u, u, ctx);
294 }
295 
296 static void ec_pow2_25519(MPI w, const MPI b, struct mpi_ec_ctx *ctx)
297 {
298 	ec_mulm_25519(w, b, b, ctx);
299 }
300 
301 /* Routines for 2^448 - 2^224 - 1.  */
302 
303 #define LIMB_SIZE_448 ((448+BITS_PER_MPI_LIMB-1)/BITS_PER_MPI_LIMB)
304 #define LIMB_SIZE_HALF_448 ((LIMB_SIZE_448+1)/2)
305 
306 static void ec_addm_448(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
307 {
308 	mpi_ptr_t wp, up, vp;
309 	mpi_size_t wsize = LIMB_SIZE_448;
310 	mpi_limb_t n[LIMB_SIZE_448];
311 	mpi_limb_t cy;
312 
313 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
314 		log_bug("addm_448: different sizes\n");
315 
316 	memset(n, 0, sizeof(n));
317 	up = u->d;
318 	vp = v->d;
319 	wp = w->d;
320 
321 	cy = mpihelp_add_n(wp, up, vp, wsize);
322 	mpih_set_cond(n, ctx->p->d, wsize, (cy != 0UL));
323 	mpihelp_sub_n(wp, wp, n, wsize);
324 }
325 
326 static void ec_subm_448(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
327 {
328 	mpi_ptr_t wp, up, vp;
329 	mpi_size_t wsize = LIMB_SIZE_448;
330 	mpi_limb_t n[LIMB_SIZE_448];
331 	mpi_limb_t borrow;
332 
333 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
334 		log_bug("subm_448: different sizes\n");
335 
336 	memset(n, 0, sizeof(n));
337 	up = u->d;
338 	vp = v->d;
339 	wp = w->d;
340 
341 	borrow = mpihelp_sub_n(wp, up, vp, wsize);
342 	mpih_set_cond(n, ctx->p->d, wsize, (borrow != 0UL));
343 	mpihelp_add_n(wp, wp, n, wsize);
344 }
345 
346 static void ec_mulm_448(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx)
347 {
348 	mpi_ptr_t wp, up, vp;
349 	mpi_size_t wsize = LIMB_SIZE_448;
350 	mpi_limb_t n[LIMB_SIZE_448*2];
351 	mpi_limb_t a2[LIMB_SIZE_HALF_448];
352 	mpi_limb_t a3[LIMB_SIZE_HALF_448];
353 	mpi_limb_t b0[LIMB_SIZE_HALF_448];
354 	mpi_limb_t b1[LIMB_SIZE_HALF_448];
355 	mpi_limb_t cy;
356 	int i;
357 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
358 	mpi_limb_t b1_rest, a3_rest;
359 #endif
360 
361 	if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
362 		log_bug("mulm_448: different sizes\n");
363 
364 	up = u->d;
365 	vp = v->d;
366 	wp = w->d;
367 
368 	mpihelp_mul_n(n, up, vp, wsize);
369 
370 	for (i = 0; i < (wsize + 1) / 2; i++) {
371 		b0[i] = n[i];
372 		b1[i] = n[i+wsize/2];
373 		a2[i] = n[i+wsize];
374 		a3[i] = n[i+wsize+wsize/2];
375 	}
376 
377 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
378 	b0[LIMB_SIZE_HALF_448-1] &= ((mpi_limb_t)1UL << 32)-1;
379 	a2[LIMB_SIZE_HALF_448-1] &= ((mpi_limb_t)1UL << 32)-1;
380 
381 	b1_rest = 0;
382 	a3_rest = 0;
383 
384 	for (i = (wsize + 1) / 2 - 1; i >= 0; i--) {
385 		mpi_limb_t b1v, a3v;
386 		b1v = b1[i];
387 		a3v = a3[i];
388 		b1[i] = (b1_rest << 32) | (b1v >> 32);
389 		a3[i] = (a3_rest << 32) | (a3v >> 32);
390 		b1_rest = b1v & (((mpi_limb_t)1UL << 32)-1);
391 		a3_rest = a3v & (((mpi_limb_t)1UL << 32)-1);
392 	}
393 #endif
394 
395 	cy = mpihelp_add_n(b0, b0, a2, LIMB_SIZE_HALF_448);
396 	cy += mpihelp_add_n(b0, b0, a3, LIMB_SIZE_HALF_448);
397 	for (i = 0; i < (wsize + 1) / 2; i++)
398 		wp[i] = b0[i];
399 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
400 	wp[LIMB_SIZE_HALF_448-1] &= (((mpi_limb_t)1UL << 32)-1);
401 #endif
402 
403 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
404 	cy = b0[LIMB_SIZE_HALF_448-1] >> 32;
405 #endif
406 
407 	cy = mpihelp_add_1(b1, b1, LIMB_SIZE_HALF_448, cy);
408 	cy += mpihelp_add_n(b1, b1, a2, LIMB_SIZE_HALF_448);
409 	cy += mpihelp_add_n(b1, b1, a3, LIMB_SIZE_HALF_448);
410 	cy += mpihelp_add_n(b1, b1, a3, LIMB_SIZE_HALF_448);
411 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
412 	b1_rest = 0;
413 	for (i = (wsize + 1) / 2 - 1; i >= 0; i--) {
414 		mpi_limb_t b1v = b1[i];
415 		b1[i] = (b1_rest << 32) | (b1v >> 32);
416 		b1_rest = b1v & (((mpi_limb_t)1UL << 32)-1);
417 	}
418 	wp[LIMB_SIZE_HALF_448-1] |= (b1_rest << 32);
419 #endif
420 	for (i = 0; i < wsize / 2; i++)
421 		wp[i+(wsize + 1) / 2] = b1[i];
422 
423 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
424 	cy = b1[LIMB_SIZE_HALF_448-1];
425 #endif
426 
427 	memset(n, 0, wsize * BYTES_PER_MPI_LIMB);
428 
429 #if (LIMB_SIZE_HALF_448 > LIMB_SIZE_448/2)
430 	n[LIMB_SIZE_HALF_448-1] = cy << 32;
431 #else
432 	n[LIMB_SIZE_HALF_448] = cy;
433 #endif
434 	n[0] = cy;
435 	mpihelp_add_n(wp, wp, n, wsize);
436 
437 	memset(n, 0, wsize * BYTES_PER_MPI_LIMB);
438 	cy = mpihelp_sub_n(wp, wp, ctx->p->d, wsize);
439 	mpih_set_cond(n, ctx->p->d, wsize, (cy != 0UL));
440 	mpihelp_add_n(wp, wp, n, wsize);
441 }
442 
443 static void ec_mul2_448(MPI w, MPI u, struct mpi_ec_ctx *ctx)
444 {
445 	ec_addm_448(w, u, u, ctx);
446 }
447 
448 static void ec_pow2_448(MPI w, const MPI b, struct mpi_ec_ctx *ctx)
449 {
450 	ec_mulm_448(w, b, b, ctx);
451 }
452 
453 struct field_table {
454 	const char *p;
455 
456 	/* computation routines for the field.  */
457 	void (*addm)(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx);
458 	void (*subm)(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx);
459 	void (*mulm)(MPI w, MPI u, MPI v, struct mpi_ec_ctx *ctx);
460 	void (*mul2)(MPI w, MPI u, struct mpi_ec_ctx *ctx);
461 	void (*pow2)(MPI w, const MPI b, struct mpi_ec_ctx *ctx);
462 };
463 
464 static const struct field_table field_table[] = {
465 	{
466 		"0x7FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFED",
467 		ec_addm_25519,
468 		ec_subm_25519,
469 		ec_mulm_25519,
470 		ec_mul2_25519,
471 		ec_pow2_25519
472 	},
473 	{
474 		"0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFE"
475 		"FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF",
476 		ec_addm_448,
477 		ec_subm_448,
478 		ec_mulm_448,
479 		ec_mul2_448,
480 		ec_pow2_448
481 	},
482 	{ NULL, NULL, NULL, NULL, NULL, NULL },
483 };
484 
485 /* Force recomputation of all helper variables.  */
486 static void mpi_ec_get_reset(struct mpi_ec_ctx *ec)
487 {
488 	ec->t.valid.a_is_pminus3 = 0;
489 	ec->t.valid.two_inv_p = 0;
490 }
491 
492 /* Accessor for helper variable.  */
493 static int ec_get_a_is_pminus3(struct mpi_ec_ctx *ec)
494 {
495 	MPI tmp;
496 
497 	if (!ec->t.valid.a_is_pminus3) {
498 		ec->t.valid.a_is_pminus3 = 1;
499 		tmp = mpi_alloc_like(ec->p);
500 		mpi_sub_ui(tmp, ec->p, 3);
501 		ec->t.a_is_pminus3 = !mpi_cmp(ec->a, tmp);
502 		mpi_free(tmp);
503 	}
504 
505 	return ec->t.a_is_pminus3;
506 }
507 
508 /* Accessor for helper variable.  */
509 static MPI ec_get_two_inv_p(struct mpi_ec_ctx *ec)
510 {
511 	if (!ec->t.valid.two_inv_p) {
512 		ec->t.valid.two_inv_p = 1;
513 		if (!ec->t.two_inv_p)
514 			ec->t.two_inv_p = mpi_alloc(0);
515 		ec_invm(ec->t.two_inv_p, mpi_const(MPI_C_TWO), ec);
516 	}
517 	return ec->t.two_inv_p;
518 }
519 
520 static const char *const curve25519_bad_points[] = {
521 	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffed",
522 	"0x0000000000000000000000000000000000000000000000000000000000000000",
523 	"0x0000000000000000000000000000000000000000000000000000000000000001",
524 	"0x00b8495f16056286fdb1329ceb8d09da6ac49ff1fae35616aeb8413b7c7aebe0",
525 	"0x57119fd0dd4e22d8868e1c58c45c44045bef839c55b1d0b1248c50a3bc959c5f",
526 	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffec",
527 	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffee",
528 	NULL
529 };
530 
531 static const char *const curve448_bad_points[] = {
532 	"0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffe"
533 	"ffffffffffffffffffffffffffffffffffffffffffffffffffffffff",
534 	"0x00000000000000000000000000000000000000000000000000000000"
535 	"00000000000000000000000000000000000000000000000000000000",
536 	"0x00000000000000000000000000000000000000000000000000000000"
537 	"00000000000000000000000000000000000000000000000000000001",
538 	"0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffe"
539 	"fffffffffffffffffffffffffffffffffffffffffffffffffffffffe",
540 	"0xffffffffffffffffffffffffffffffffffffffffffffffffffffffff"
541 	"00000000000000000000000000000000000000000000000000000000",
542 	NULL
543 };
544 
545 static const char *const *bad_points_table[] = {
546 	curve25519_bad_points,
547 	curve448_bad_points,
548 };
549 
550 static void mpi_ec_coefficient_normalize(MPI a, MPI p)
551 {
552 	if (a->sign) {
553 		mpi_resize(a, p->nlimbs);
554 		mpihelp_sub_n(a->d, p->d, a->d, p->nlimbs);
555 		a->nlimbs = p->nlimbs;
556 		a->sign = 0;
557 	}
558 }
559 
560 /* This function initialized a context for elliptic curve based on the
561  * field GF(p).  P is the prime specifying this field, A is the first
562  * coefficient.  CTX is expected to be zeroized.
563  */
564 void mpi_ec_init(struct mpi_ec_ctx *ctx, enum gcry_mpi_ec_models model,
565 			enum ecc_dialects dialect,
566 			int flags, MPI p, MPI a, MPI b)
567 {
568 	int i;
569 	static int use_barrett = -1 /* TODO: 1 or -1 */;
570 
571 	mpi_ec_coefficient_normalize(a, p);
572 	mpi_ec_coefficient_normalize(b, p);
573 
574 	/* Fixme: Do we want to check some constraints? e.g.  a < p  */
575 
576 	ctx->model = model;
577 	ctx->dialect = dialect;
578 	ctx->flags = flags;
579 	if (dialect == ECC_DIALECT_ED25519)
580 		ctx->nbits = 256;
581 	else
582 		ctx->nbits = mpi_get_nbits(p);
583 	ctx->p = mpi_copy(p);
584 	ctx->a = mpi_copy(a);
585 	ctx->b = mpi_copy(b);
586 
587 	ctx->d = NULL;
588 	ctx->t.two_inv_p = NULL;
589 
590 	ctx->t.p_barrett = use_barrett > 0 ? mpi_barrett_init(ctx->p, 0) : NULL;
591 
592 	mpi_ec_get_reset(ctx);
593 
594 	if (model == MPI_EC_MONTGOMERY) {
595 		for (i = 0; i < DIM(bad_points_table); i++) {
596 			MPI p_candidate = mpi_scanval(bad_points_table[i][0]);
597 			int match_p = !mpi_cmp(ctx->p, p_candidate);
598 			int j;
599 
600 			mpi_free(p_candidate);
601 			if (!match_p)
602 				continue;
603 
604 			for (j = 0; i < DIM(ctx->t.scratch) && bad_points_table[i][j]; j++)
605 				ctx->t.scratch[j] = mpi_scanval(bad_points_table[i][j]);
606 		}
607 	} else {
608 		/* Allocate scratch variables.  */
609 		for (i = 0; i < DIM(ctx->t.scratch); i++)
610 			ctx->t.scratch[i] = mpi_alloc_like(ctx->p);
611 	}
612 
613 	ctx->addm = ec_addm;
614 	ctx->subm = ec_subm;
615 	ctx->mulm = ec_mulm;
616 	ctx->mul2 = ec_mul2;
617 	ctx->pow2 = ec_pow2;
618 
619 	for (i = 0; field_table[i].p; i++) {
620 		MPI f_p;
621 
622 		f_p = mpi_scanval(field_table[i].p);
623 		if (!f_p)
624 			break;
625 
626 		if (!mpi_cmp(p, f_p)) {
627 			ctx->addm = field_table[i].addm;
628 			ctx->subm = field_table[i].subm;
629 			ctx->mulm = field_table[i].mulm;
630 			ctx->mul2 = field_table[i].mul2;
631 			ctx->pow2 = field_table[i].pow2;
632 			mpi_free(f_p);
633 
634 			mpi_resize(ctx->a, ctx->p->nlimbs);
635 			ctx->a->nlimbs = ctx->p->nlimbs;
636 
637 			mpi_resize(ctx->b, ctx->p->nlimbs);
638 			ctx->b->nlimbs = ctx->p->nlimbs;
639 
640 			for (i = 0; i < DIM(ctx->t.scratch) && ctx->t.scratch[i]; i++)
641 				ctx->t.scratch[i]->nlimbs = ctx->p->nlimbs;
642 
643 			break;
644 		}
645 
646 		mpi_free(f_p);
647 	}
648 }
649 EXPORT_SYMBOL_GPL(mpi_ec_init);
650 
651 void mpi_ec_deinit(struct mpi_ec_ctx *ctx)
652 {
653 	int i;
654 
655 	mpi_barrett_free(ctx->t.p_barrett);
656 
657 	/* Domain parameter.  */
658 	mpi_free(ctx->p);
659 	mpi_free(ctx->a);
660 	mpi_free(ctx->b);
661 	mpi_point_release(ctx->G);
662 	mpi_free(ctx->n);
663 
664 	/* The key.  */
665 	mpi_point_release(ctx->Q);
666 	mpi_free(ctx->d);
667 
668 	/* Private data of ec.c.  */
669 	mpi_free(ctx->t.two_inv_p);
670 
671 	for (i = 0; i < DIM(ctx->t.scratch); i++)
672 		mpi_free(ctx->t.scratch[i]);
673 }
674 EXPORT_SYMBOL_GPL(mpi_ec_deinit);
675 
676 /* Compute the affine coordinates from the projective coordinates in
677  * POINT.  Set them into X and Y.  If one coordinate is not required,
678  * X or Y may be passed as NULL.  CTX is the usual context. Returns: 0
679  * on success or !0 if POINT is at infinity.
680  */
681 int mpi_ec_get_affine(MPI x, MPI y, MPI_POINT point, struct mpi_ec_ctx *ctx)
682 {
683 	if (!mpi_cmp_ui(point->z, 0))
684 		return -1;
685 
686 	switch (ctx->model) {
687 	case MPI_EC_WEIERSTRASS: /* Using Jacobian coordinates.  */
688 		{
689 			MPI z1, z2, z3;
690 
691 			z1 = mpi_new(0);
692 			z2 = mpi_new(0);
693 			ec_invm(z1, point->z, ctx);  /* z1 = z^(-1) mod p  */
694 			ec_mulm(z2, z1, z1, ctx);    /* z2 = z^(-2) mod p  */
695 
696 			if (x)
697 				ec_mulm(x, point->x, z2, ctx);
698 
699 			if (y) {
700 				z3 = mpi_new(0);
701 				ec_mulm(z3, z2, z1, ctx);      /* z3 = z^(-3) mod p */
702 				ec_mulm(y, point->y, z3, ctx);
703 				mpi_free(z3);
704 			}
705 
706 			mpi_free(z2);
707 			mpi_free(z1);
708 		}
709 		return 0;
710 
711 	case MPI_EC_MONTGOMERY:
712 		{
713 			if (x)
714 				mpi_set(x, point->x);
715 
716 			if (y) {
717 				log_fatal("%s: Getting Y-coordinate on %s is not supported\n",
718 						"mpi_ec_get_affine", "Montgomery");
719 				return -1;
720 			}
721 		}
722 		return 0;
723 
724 	case MPI_EC_EDWARDS:
725 		{
726 			MPI z;
727 
728 			z = mpi_new(0);
729 			ec_invm(z, point->z, ctx);
730 
731 			mpi_resize(z, ctx->p->nlimbs);
732 			z->nlimbs = ctx->p->nlimbs;
733 
734 			if (x) {
735 				mpi_resize(x, ctx->p->nlimbs);
736 				x->nlimbs = ctx->p->nlimbs;
737 				ctx->mulm(x, point->x, z, ctx);
738 			}
739 			if (y) {
740 				mpi_resize(y, ctx->p->nlimbs);
741 				y->nlimbs = ctx->p->nlimbs;
742 				ctx->mulm(y, point->y, z, ctx);
743 			}
744 
745 			mpi_free(z);
746 		}
747 		return 0;
748 
749 	default:
750 		return -1;
751 	}
752 }
753 EXPORT_SYMBOL_GPL(mpi_ec_get_affine);
754 
755 /*  RESULT = 2 * POINT  (Weierstrass version). */
756 static void dup_point_weierstrass(MPI_POINT result,
757 		MPI_POINT point, struct mpi_ec_ctx *ctx)
758 {
759 #define x3 (result->x)
760 #define y3 (result->y)
761 #define z3 (result->z)
762 #define t1 (ctx->t.scratch[0])
763 #define t2 (ctx->t.scratch[1])
764 #define t3 (ctx->t.scratch[2])
765 #define l1 (ctx->t.scratch[3])
766 #define l2 (ctx->t.scratch[4])
767 #define l3 (ctx->t.scratch[5])
768 
769 	if (!mpi_cmp_ui(point->y, 0) || !mpi_cmp_ui(point->z, 0)) {
770 		/* P_y == 0 || P_z == 0 => [1:1:0] */
771 		mpi_set_ui(x3, 1);
772 		mpi_set_ui(y3, 1);
773 		mpi_set_ui(z3, 0);
774 	} else {
775 		if (ec_get_a_is_pminus3(ctx)) {
776 			/* Use the faster case.  */
777 			/* L1 = 3(X - Z^2)(X + Z^2) */
778 			/*                          T1: used for Z^2. */
779 			/*                          T2: used for the right term. */
780 			ec_pow2(t1, point->z, ctx);
781 			ec_subm(l1, point->x, t1, ctx);
782 			ec_mulm(l1, l1, mpi_const(MPI_C_THREE), ctx);
783 			ec_addm(t2, point->x, t1, ctx);
784 			ec_mulm(l1, l1, t2, ctx);
785 		} else {
786 			/* Standard case. */
787 			/* L1 = 3X^2 + aZ^4 */
788 			/*                          T1: used for aZ^4. */
789 			ec_pow2(l1, point->x, ctx);
790 			ec_mulm(l1, l1, mpi_const(MPI_C_THREE), ctx);
791 			ec_powm(t1, point->z, mpi_const(MPI_C_FOUR), ctx);
792 			ec_mulm(t1, t1, ctx->a, ctx);
793 			ec_addm(l1, l1, t1, ctx);
794 		}
795 		/* Z3 = 2YZ */
796 		ec_mulm(z3, point->y, point->z, ctx);
797 		ec_mul2(z3, z3, ctx);
798 
799 		/* L2 = 4XY^2 */
800 		/*                              T2: used for Y2; required later. */
801 		ec_pow2(t2, point->y, ctx);
802 		ec_mulm(l2, t2, point->x, ctx);
803 		ec_mulm(l2, l2, mpi_const(MPI_C_FOUR), ctx);
804 
805 		/* X3 = L1^2 - 2L2 */
806 		/*                              T1: used for L2^2. */
807 		ec_pow2(x3, l1, ctx);
808 		ec_mul2(t1, l2, ctx);
809 		ec_subm(x3, x3, t1, ctx);
810 
811 		/* L3 = 8Y^4 */
812 		/*                              T2: taken from above. */
813 		ec_pow2(t2, t2, ctx);
814 		ec_mulm(l3, t2, mpi_const(MPI_C_EIGHT), ctx);
815 
816 		/* Y3 = L1(L2 - X3) - L3 */
817 		ec_subm(y3, l2, x3, ctx);
818 		ec_mulm(y3, y3, l1, ctx);
819 		ec_subm(y3, y3, l3, ctx);
820 	}
821 
822 #undef x3
823 #undef y3
824 #undef z3
825 #undef t1
826 #undef t2
827 #undef t3
828 #undef l1
829 #undef l2
830 #undef l3
831 }
832 
833 /*  RESULT = 2 * POINT  (Montgomery version). */
834 static void dup_point_montgomery(MPI_POINT result,
835 				MPI_POINT point, struct mpi_ec_ctx *ctx)
836 {
837 	(void)result;
838 	(void)point;
839 	(void)ctx;
840 	log_fatal("%s: %s not yet supported\n",
841 			"mpi_ec_dup_point", "Montgomery");
842 }
843 
844 /*  RESULT = 2 * POINT  (Twisted Edwards version). */
845 static void dup_point_edwards(MPI_POINT result,
846 		MPI_POINT point, struct mpi_ec_ctx *ctx)
847 {
848 #define X1 (point->x)
849 #define Y1 (point->y)
850 #define Z1 (point->z)
851 #define X3 (result->x)
852 #define Y3 (result->y)
853 #define Z3 (result->z)
854 #define B (ctx->t.scratch[0])
855 #define C (ctx->t.scratch[1])
856 #define D (ctx->t.scratch[2])
857 #define E (ctx->t.scratch[3])
858 #define F (ctx->t.scratch[4])
859 #define H (ctx->t.scratch[5])
860 #define J (ctx->t.scratch[6])
861 
862 	/* Compute: (X_3 : Y_3 : Z_3) = 2( X_1 : Y_1 : Z_1 ) */
863 
864 	/* B = (X_1 + Y_1)^2  */
865 	ctx->addm(B, X1, Y1, ctx);
866 	ctx->pow2(B, B, ctx);
867 
868 	/* C = X_1^2 */
869 	/* D = Y_1^2 */
870 	ctx->pow2(C, X1, ctx);
871 	ctx->pow2(D, Y1, ctx);
872 
873 	/* E = aC */
874 	if (ctx->dialect == ECC_DIALECT_ED25519)
875 		ctx->subm(E, ctx->p, C, ctx);
876 	else
877 		ctx->mulm(E, ctx->a, C, ctx);
878 
879 	/* F = E + D */
880 	ctx->addm(F, E, D, ctx);
881 
882 	/* H = Z_1^2 */
883 	ctx->pow2(H, Z1, ctx);
884 
885 	/* J = F - 2H */
886 	ctx->mul2(J, H, ctx);
887 	ctx->subm(J, F, J, ctx);
888 
889 	/* X_3 = (B - C - D) · J */
890 	ctx->subm(X3, B, C, ctx);
891 	ctx->subm(X3, X3, D, ctx);
892 	ctx->mulm(X3, X3, J, ctx);
893 
894 	/* Y_3 = F · (E - D) */
895 	ctx->subm(Y3, E, D, ctx);
896 	ctx->mulm(Y3, Y3, F, ctx);
897 
898 	/* Z_3 = F · J */
899 	ctx->mulm(Z3, F, J, ctx);
900 
901 #undef X1
902 #undef Y1
903 #undef Z1
904 #undef X3
905 #undef Y3
906 #undef Z3
907 #undef B
908 #undef C
909 #undef D
910 #undef E
911 #undef F
912 #undef H
913 #undef J
914 }
915 
916 /*  RESULT = 2 * POINT  */
917 static void
918 mpi_ec_dup_point(MPI_POINT result, MPI_POINT point, struct mpi_ec_ctx *ctx)
919 {
920 	switch (ctx->model) {
921 	case MPI_EC_WEIERSTRASS:
922 		dup_point_weierstrass(result, point, ctx);
923 		break;
924 	case MPI_EC_MONTGOMERY:
925 		dup_point_montgomery(result, point, ctx);
926 		break;
927 	case MPI_EC_EDWARDS:
928 		dup_point_edwards(result, point, ctx);
929 		break;
930 	}
931 }
932 
933 /* RESULT = P1 + P2  (Weierstrass version).*/
934 static void add_points_weierstrass(MPI_POINT result,
935 		MPI_POINT p1, MPI_POINT p2,
936 		struct mpi_ec_ctx *ctx)
937 {
938 #define x1 (p1->x)
939 #define y1 (p1->y)
940 #define z1 (p1->z)
941 #define x2 (p2->x)
942 #define y2 (p2->y)
943 #define z2 (p2->z)
944 #define x3 (result->x)
945 #define y3 (result->y)
946 #define z3 (result->z)
947 #define l1 (ctx->t.scratch[0])
948 #define l2 (ctx->t.scratch[1])
949 #define l3 (ctx->t.scratch[2])
950 #define l4 (ctx->t.scratch[3])
951 #define l5 (ctx->t.scratch[4])
952 #define l6 (ctx->t.scratch[5])
953 #define l7 (ctx->t.scratch[6])
954 #define l8 (ctx->t.scratch[7])
955 #define l9 (ctx->t.scratch[8])
956 #define t1 (ctx->t.scratch[9])
957 #define t2 (ctx->t.scratch[10])
958 
959 	if ((!mpi_cmp(x1, x2)) && (!mpi_cmp(y1, y2)) && (!mpi_cmp(z1, z2))) {
960 		/* Same point; need to call the duplicate function.  */
961 		mpi_ec_dup_point(result, p1, ctx);
962 	} else if (!mpi_cmp_ui(z1, 0)) {
963 		/* P1 is at infinity.  */
964 		mpi_set(x3, p2->x);
965 		mpi_set(y3, p2->y);
966 		mpi_set(z3, p2->z);
967 	} else if (!mpi_cmp_ui(z2, 0)) {
968 		/* P2 is at infinity.  */
969 		mpi_set(x3, p1->x);
970 		mpi_set(y3, p1->y);
971 		mpi_set(z3, p1->z);
972 	} else {
973 		int z1_is_one = !mpi_cmp_ui(z1, 1);
974 		int z2_is_one = !mpi_cmp_ui(z2, 1);
975 
976 		/* l1 = x1 z2^2  */
977 		/* l2 = x2 z1^2  */
978 		if (z2_is_one)
979 			mpi_set(l1, x1);
980 		else {
981 			ec_pow2(l1, z2, ctx);
982 			ec_mulm(l1, l1, x1, ctx);
983 		}
984 		if (z1_is_one)
985 			mpi_set(l2, x2);
986 		else {
987 			ec_pow2(l2, z1, ctx);
988 			ec_mulm(l2, l2, x2, ctx);
989 		}
990 		/* l3 = l1 - l2 */
991 		ec_subm(l3, l1, l2, ctx);
992 		/* l4 = y1 z2^3  */
993 		ec_powm(l4, z2, mpi_const(MPI_C_THREE), ctx);
994 		ec_mulm(l4, l4, y1, ctx);
995 		/* l5 = y2 z1^3  */
996 		ec_powm(l5, z1, mpi_const(MPI_C_THREE), ctx);
997 		ec_mulm(l5, l5, y2, ctx);
998 		/* l6 = l4 - l5  */
999 		ec_subm(l6, l4, l5, ctx);
1000 
1001 		if (!mpi_cmp_ui(l3, 0)) {
1002 			if (!mpi_cmp_ui(l6, 0)) {
1003 				/* P1 and P2 are the same - use duplicate function. */
1004 				mpi_ec_dup_point(result, p1, ctx);
1005 			} else {
1006 				/* P1 is the inverse of P2.  */
1007 				mpi_set_ui(x3, 1);
1008 				mpi_set_ui(y3, 1);
1009 				mpi_set_ui(z3, 0);
1010 			}
1011 		} else {
1012 			/* l7 = l1 + l2  */
1013 			ec_addm(l7, l1, l2, ctx);
1014 			/* l8 = l4 + l5  */
1015 			ec_addm(l8, l4, l5, ctx);
1016 			/* z3 = z1 z2 l3  */
1017 			ec_mulm(z3, z1, z2, ctx);
1018 			ec_mulm(z3, z3, l3, ctx);
1019 			/* x3 = l6^2 - l7 l3^2  */
1020 			ec_pow2(t1, l6, ctx);
1021 			ec_pow2(t2, l3, ctx);
1022 			ec_mulm(t2, t2, l7, ctx);
1023 			ec_subm(x3, t1, t2, ctx);
1024 			/* l9 = l7 l3^2 - 2 x3  */
1025 			ec_mul2(t1, x3, ctx);
1026 			ec_subm(l9, t2, t1, ctx);
1027 			/* y3 = (l9 l6 - l8 l3^3)/2  */
1028 			ec_mulm(l9, l9, l6, ctx);
1029 			ec_powm(t1, l3, mpi_const(MPI_C_THREE), ctx); /* fixme: Use saved value*/
1030 			ec_mulm(t1, t1, l8, ctx);
1031 			ec_subm(y3, l9, t1, ctx);
1032 			ec_mulm(y3, y3, ec_get_two_inv_p(ctx), ctx);
1033 		}
1034 	}
1035 
1036 #undef x1
1037 #undef y1
1038 #undef z1
1039 #undef x2
1040 #undef y2
1041 #undef z2
1042 #undef x3
1043 #undef y3
1044 #undef z3
1045 #undef l1
1046 #undef l2
1047 #undef l3
1048 #undef l4
1049 #undef l5
1050 #undef l6
1051 #undef l7
1052 #undef l8
1053 #undef l9
1054 #undef t1
1055 #undef t2
1056 }
1057 
1058 /* RESULT = P1 + P2  (Montgomery version).*/
1059 static void add_points_montgomery(MPI_POINT result,
1060 		MPI_POINT p1, MPI_POINT p2,
1061 		struct mpi_ec_ctx *ctx)
1062 {
1063 	(void)result;
1064 	(void)p1;
1065 	(void)p2;
1066 	(void)ctx;
1067 	log_fatal("%s: %s not yet supported\n",
1068 			"mpi_ec_add_points", "Montgomery");
1069 }
1070 
1071 /* RESULT = P1 + P2  (Twisted Edwards version).*/
1072 static void add_points_edwards(MPI_POINT result,
1073 		MPI_POINT p1, MPI_POINT p2,
1074 		struct mpi_ec_ctx *ctx)
1075 {
1076 #define X1 (p1->x)
1077 #define Y1 (p1->y)
1078 #define Z1 (p1->z)
1079 #define X2 (p2->x)
1080 #define Y2 (p2->y)
1081 #define Z2 (p2->z)
1082 #define X3 (result->x)
1083 #define Y3 (result->y)
1084 #define Z3 (result->z)
1085 #define A (ctx->t.scratch[0])
1086 #define B (ctx->t.scratch[1])
1087 #define C (ctx->t.scratch[2])
1088 #define D (ctx->t.scratch[3])
1089 #define E (ctx->t.scratch[4])
1090 #define F (ctx->t.scratch[5])
1091 #define G (ctx->t.scratch[6])
1092 #define tmp (ctx->t.scratch[7])
1093 
1094 	point_resize(result, ctx);
1095 
1096 	/* Compute: (X_3 : Y_3 : Z_3) = (X_1 : Y_1 : Z_1) + (X_2 : Y_2 : Z_3) */
1097 
1098 	/* A = Z1 · Z2 */
1099 	ctx->mulm(A, Z1, Z2, ctx);
1100 
1101 	/* B = A^2 */
1102 	ctx->pow2(B, A, ctx);
1103 
1104 	/* C = X1 · X2 */
1105 	ctx->mulm(C, X1, X2, ctx);
1106 
1107 	/* D = Y1 · Y2 */
1108 	ctx->mulm(D, Y1, Y2, ctx);
1109 
1110 	/* E = d · C · D */
1111 	ctx->mulm(E, ctx->b, C, ctx);
1112 	ctx->mulm(E, E, D, ctx);
1113 
1114 	/* F = B - E */
1115 	ctx->subm(F, B, E, ctx);
1116 
1117 	/* G = B + E */
1118 	ctx->addm(G, B, E, ctx);
1119 
1120 	/* X_3 = A · F · ((X_1 + Y_1) · (X_2 + Y_2) - C - D) */
1121 	ctx->addm(tmp, X1, Y1, ctx);
1122 	ctx->addm(X3, X2, Y2, ctx);
1123 	ctx->mulm(X3, X3, tmp, ctx);
1124 	ctx->subm(X3, X3, C, ctx);
1125 	ctx->subm(X3, X3, D, ctx);
1126 	ctx->mulm(X3, X3, F, ctx);
1127 	ctx->mulm(X3, X3, A, ctx);
1128 
1129 	/* Y_3 = A · G · (D - aC) */
1130 	if (ctx->dialect == ECC_DIALECT_ED25519) {
1131 		ctx->addm(Y3, D, C, ctx);
1132 	} else {
1133 		ctx->mulm(Y3, ctx->a, C, ctx);
1134 		ctx->subm(Y3, D, Y3, ctx);
1135 	}
1136 	ctx->mulm(Y3, Y3, G, ctx);
1137 	ctx->mulm(Y3, Y3, A, ctx);
1138 
1139 	/* Z_3 = F · G */
1140 	ctx->mulm(Z3, F, G, ctx);
1141 
1142 
1143 #undef X1
1144 #undef Y1
1145 #undef Z1
1146 #undef X2
1147 #undef Y2
1148 #undef Z2
1149 #undef X3
1150 #undef Y3
1151 #undef Z3
1152 #undef A
1153 #undef B
1154 #undef C
1155 #undef D
1156 #undef E
1157 #undef F
1158 #undef G
1159 #undef tmp
1160 }
1161 
1162 /* Compute a step of Montgomery Ladder (only use X and Z in the point).
1163  * Inputs:  P1, P2, and x-coordinate of DIF = P1 - P1.
1164  * Outputs: PRD = 2 * P1 and  SUM = P1 + P2.
1165  */
1166 static void montgomery_ladder(MPI_POINT prd, MPI_POINT sum,
1167 		MPI_POINT p1, MPI_POINT p2, MPI dif_x,
1168 		struct mpi_ec_ctx *ctx)
1169 {
1170 	ctx->addm(sum->x, p2->x, p2->z, ctx);
1171 	ctx->subm(p2->z, p2->x, p2->z, ctx);
1172 	ctx->addm(prd->x, p1->x, p1->z, ctx);
1173 	ctx->subm(p1->z, p1->x, p1->z, ctx);
1174 	ctx->mulm(p2->x, p1->z, sum->x, ctx);
1175 	ctx->mulm(p2->z, prd->x, p2->z, ctx);
1176 	ctx->pow2(p1->x, prd->x, ctx);
1177 	ctx->pow2(p1->z, p1->z, ctx);
1178 	ctx->addm(sum->x, p2->x, p2->z, ctx);
1179 	ctx->subm(p2->z, p2->x, p2->z, ctx);
1180 	ctx->mulm(prd->x, p1->x, p1->z, ctx);
1181 	ctx->subm(p1->z, p1->x, p1->z, ctx);
1182 	ctx->pow2(sum->x, sum->x, ctx);
1183 	ctx->pow2(sum->z, p2->z, ctx);
1184 	ctx->mulm(prd->z, p1->z, ctx->a, ctx); /* CTX->A: (a-2)/4 */
1185 	ctx->mulm(sum->z, sum->z, dif_x, ctx);
1186 	ctx->addm(prd->z, p1->x, prd->z, ctx);
1187 	ctx->mulm(prd->z, prd->z, p1->z, ctx);
1188 }
1189 
1190 /* RESULT = P1 + P2 */
1191 void mpi_ec_add_points(MPI_POINT result,
1192 		MPI_POINT p1, MPI_POINT p2,
1193 		struct mpi_ec_ctx *ctx)
1194 {
1195 	switch (ctx->model) {
1196 	case MPI_EC_WEIERSTRASS:
1197 		add_points_weierstrass(result, p1, p2, ctx);
1198 		break;
1199 	case MPI_EC_MONTGOMERY:
1200 		add_points_montgomery(result, p1, p2, ctx);
1201 		break;
1202 	case MPI_EC_EDWARDS:
1203 		add_points_edwards(result, p1, p2, ctx);
1204 		break;
1205 	}
1206 }
1207 EXPORT_SYMBOL_GPL(mpi_ec_add_points);
1208 
1209 /* Scalar point multiplication - the main function for ECC.  If takes
1210  * an integer SCALAR and a POINT as well as the usual context CTX.
1211  * RESULT will be set to the resulting point.
1212  */
1213 void mpi_ec_mul_point(MPI_POINT result,
1214 			MPI scalar, MPI_POINT point,
1215 			struct mpi_ec_ctx *ctx)
1216 {
1217 	MPI x1, y1, z1, k, h, yy;
1218 	unsigned int i, loops;
1219 	struct gcry_mpi_point p1, p2, p1inv;
1220 
1221 	if (ctx->model == MPI_EC_EDWARDS) {
1222 		/* Simple left to right binary method.  Algorithm 3.27 from
1223 		 * {author={Hankerson, Darrel and Menezes, Alfred J. and Vanstone, Scott},
1224 		 *  title = {Guide to Elliptic Curve Cryptography},
1225 		 *  year = {2003}, isbn = {038795273X},
1226 		 *  url = {http://www.cacr.math.uwaterloo.ca/ecc/},
1227 		 *  publisher = {Springer-Verlag New York, Inc.}}
1228 		 */
1229 		unsigned int nbits;
1230 		int j;
1231 
1232 		if (mpi_cmp(scalar, ctx->p) >= 0)
1233 			nbits = mpi_get_nbits(scalar);
1234 		else
1235 			nbits = mpi_get_nbits(ctx->p);
1236 
1237 		mpi_set_ui(result->x, 0);
1238 		mpi_set_ui(result->y, 1);
1239 		mpi_set_ui(result->z, 1);
1240 		point_resize(point, ctx);
1241 
1242 		point_resize(result, ctx);
1243 		point_resize(point, ctx);
1244 
1245 		for (j = nbits-1; j >= 0; j--) {
1246 			mpi_ec_dup_point(result, result, ctx);
1247 			if (mpi_test_bit(scalar, j))
1248 				mpi_ec_add_points(result, result, point, ctx);
1249 		}
1250 		return;
1251 	} else if (ctx->model == MPI_EC_MONTGOMERY) {
1252 		unsigned int nbits;
1253 		int j;
1254 		struct gcry_mpi_point p1_, p2_;
1255 		MPI_POINT q1, q2, prd, sum;
1256 		unsigned long sw;
1257 		mpi_size_t rsize;
1258 
1259 		/* Compute scalar point multiplication with Montgomery Ladder.
1260 		 * Note that we don't use Y-coordinate in the points at all.
1261 		 * RESULT->Y will be filled by zero.
1262 		 */
1263 
1264 		nbits = mpi_get_nbits(scalar);
1265 		point_init(&p1);
1266 		point_init(&p2);
1267 		point_init(&p1_);
1268 		point_init(&p2_);
1269 		mpi_set_ui(p1.x, 1);
1270 		mpi_free(p2.x);
1271 		p2.x = mpi_copy(point->x);
1272 		mpi_set_ui(p2.z, 1);
1273 
1274 		point_resize(&p1, ctx);
1275 		point_resize(&p2, ctx);
1276 		point_resize(&p1_, ctx);
1277 		point_resize(&p2_, ctx);
1278 
1279 		mpi_resize(point->x, ctx->p->nlimbs);
1280 		point->x->nlimbs = ctx->p->nlimbs;
1281 
1282 		q1 = &p1;
1283 		q2 = &p2;
1284 		prd = &p1_;
1285 		sum = &p2_;
1286 
1287 		for (j = nbits-1; j >= 0; j--) {
1288 			MPI_POINT t;
1289 
1290 			sw = mpi_test_bit(scalar, j);
1291 			point_swap_cond(q1, q2, sw, ctx);
1292 			montgomery_ladder(prd, sum, q1, q2, point->x, ctx);
1293 			point_swap_cond(prd, sum, sw, ctx);
1294 			t = q1;  q1 = prd;  prd = t;
1295 			t = q2;  q2 = sum;  sum = t;
1296 		}
1297 
1298 		mpi_clear(result->y);
1299 		sw = (nbits & 1);
1300 		point_swap_cond(&p1, &p1_, sw, ctx);
1301 
1302 		rsize = p1.z->nlimbs;
1303 		MPN_NORMALIZE(p1.z->d, rsize);
1304 		if (rsize == 0) {
1305 			mpi_set_ui(result->x, 1);
1306 			mpi_set_ui(result->z, 0);
1307 		} else {
1308 			z1 = mpi_new(0);
1309 			ec_invm(z1, p1.z, ctx);
1310 			ec_mulm(result->x, p1.x, z1, ctx);
1311 			mpi_set_ui(result->z, 1);
1312 			mpi_free(z1);
1313 		}
1314 
1315 		point_free(&p1);
1316 		point_free(&p2);
1317 		point_free(&p1_);
1318 		point_free(&p2_);
1319 		return;
1320 	}
1321 
1322 	x1 = mpi_alloc_like(ctx->p);
1323 	y1 = mpi_alloc_like(ctx->p);
1324 	h  = mpi_alloc_like(ctx->p);
1325 	k  = mpi_copy(scalar);
1326 	yy = mpi_copy(point->y);
1327 
1328 	if (mpi_has_sign(k)) {
1329 		k->sign = 0;
1330 		ec_invm(yy, yy, ctx);
1331 	}
1332 
1333 	if (!mpi_cmp_ui(point->z, 1)) {
1334 		mpi_set(x1, point->x);
1335 		mpi_set(y1, yy);
1336 	} else {
1337 		MPI z2, z3;
1338 
1339 		z2 = mpi_alloc_like(ctx->p);
1340 		z3 = mpi_alloc_like(ctx->p);
1341 		ec_mulm(z2, point->z, point->z, ctx);
1342 		ec_mulm(z3, point->z, z2, ctx);
1343 		ec_invm(z2, z2, ctx);
1344 		ec_mulm(x1, point->x, z2, ctx);
1345 		ec_invm(z3, z3, ctx);
1346 		ec_mulm(y1, yy, z3, ctx);
1347 		mpi_free(z2);
1348 		mpi_free(z3);
1349 	}
1350 	z1 = mpi_copy(mpi_const(MPI_C_ONE));
1351 
1352 	mpi_mul(h, k, mpi_const(MPI_C_THREE)); /* h = 3k */
1353 	loops = mpi_get_nbits(h);
1354 	if (loops < 2) {
1355 		/* If SCALAR is zero, the above mpi_mul sets H to zero and thus
1356 		 * LOOPs will be zero.  To avoid an underflow of I in the main
1357 		 * loop we set LOOP to 2 and the result to (0,0,0).
1358 		 */
1359 		loops = 2;
1360 		mpi_clear(result->x);
1361 		mpi_clear(result->y);
1362 		mpi_clear(result->z);
1363 	} else {
1364 		mpi_set(result->x, point->x);
1365 		mpi_set(result->y, yy);
1366 		mpi_set(result->z, point->z);
1367 	}
1368 	mpi_free(yy); yy = NULL;
1369 
1370 	p1.x = x1; x1 = NULL;
1371 	p1.y = y1; y1 = NULL;
1372 	p1.z = z1; z1 = NULL;
1373 	point_init(&p2);
1374 	point_init(&p1inv);
1375 
1376 	/* Invert point: y = p - y mod p  */
1377 	point_set(&p1inv, &p1);
1378 	ec_subm(p1inv.y, ctx->p, p1inv.y, ctx);
1379 
1380 	for (i = loops-2; i > 0; i--) {
1381 		mpi_ec_dup_point(result, result, ctx);
1382 		if (mpi_test_bit(h, i) == 1 && mpi_test_bit(k, i) == 0) {
1383 			point_set(&p2, result);
1384 			mpi_ec_add_points(result, &p2, &p1, ctx);
1385 		}
1386 		if (mpi_test_bit(h, i) == 0 && mpi_test_bit(k, i) == 1) {
1387 			point_set(&p2, result);
1388 			mpi_ec_add_points(result, &p2, &p1inv, ctx);
1389 		}
1390 	}
1391 
1392 	point_free(&p1);
1393 	point_free(&p2);
1394 	point_free(&p1inv);
1395 	mpi_free(h);
1396 	mpi_free(k);
1397 }
1398 EXPORT_SYMBOL_GPL(mpi_ec_mul_point);
1399 
1400 /* Return true if POINT is on the curve described by CTX.  */
1401 int mpi_ec_curve_point(MPI_POINT point, struct mpi_ec_ctx *ctx)
1402 {
1403 	int res = 0;
1404 	MPI x, y, w;
1405 
1406 	x = mpi_new(0);
1407 	y = mpi_new(0);
1408 	w = mpi_new(0);
1409 
1410 	/* Check that the point is in range.  This needs to be done here and
1411 	 * not after conversion to affine coordinates.
1412 	 */
1413 	if (mpi_cmpabs(point->x, ctx->p) >= 0)
1414 		goto leave;
1415 	if (mpi_cmpabs(point->y, ctx->p) >= 0)
1416 		goto leave;
1417 	if (mpi_cmpabs(point->z, ctx->p) >= 0)
1418 		goto leave;
1419 
1420 	switch (ctx->model) {
1421 	case MPI_EC_WEIERSTRASS:
1422 		{
1423 			MPI xxx;
1424 
1425 			if (mpi_ec_get_affine(x, y, point, ctx))
1426 				goto leave;
1427 
1428 			xxx = mpi_new(0);
1429 
1430 			/* y^2 == x^3 + a·x + b */
1431 			ec_pow2(y, y, ctx);
1432 
1433 			ec_pow3(xxx, x, ctx);
1434 			ec_mulm(w, ctx->a, x, ctx);
1435 			ec_addm(w, w, ctx->b, ctx);
1436 			ec_addm(w, w, xxx, ctx);
1437 
1438 			if (!mpi_cmp(y, w))
1439 				res = 1;
1440 
1441 			mpi_free(xxx);
1442 		}
1443 		break;
1444 
1445 	case MPI_EC_MONTGOMERY:
1446 		{
1447 #define xx y
1448 			/* With Montgomery curve, only X-coordinate is valid. */
1449 			if (mpi_ec_get_affine(x, NULL, point, ctx))
1450 				goto leave;
1451 
1452 			/* The equation is: b * y^2 == x^3 + a · x^2 + x */
1453 			/* We check if right hand is quadratic residue or not by
1454 			 * Euler's criterion.
1455 			 */
1456 			/* CTX->A has (a-2)/4 and CTX->B has b^-1 */
1457 			ec_mulm(w, ctx->a, mpi_const(MPI_C_FOUR), ctx);
1458 			ec_addm(w, w, mpi_const(MPI_C_TWO), ctx);
1459 			ec_mulm(w, w, x, ctx);
1460 			ec_pow2(xx, x, ctx);
1461 			ec_addm(w, w, xx, ctx);
1462 			ec_addm(w, w, mpi_const(MPI_C_ONE), ctx);
1463 			ec_mulm(w, w, x, ctx);
1464 			ec_mulm(w, w, ctx->b, ctx);
1465 #undef xx
1466 			/* Compute Euler's criterion: w^(p-1)/2 */
1467 #define p_minus1 y
1468 			ec_subm(p_minus1, ctx->p, mpi_const(MPI_C_ONE), ctx);
1469 			mpi_rshift(p_minus1, p_minus1, 1);
1470 			ec_powm(w, w, p_minus1, ctx);
1471 
1472 			res = !mpi_cmp_ui(w, 1);
1473 #undef p_minus1
1474 		}
1475 		break;
1476 
1477 	case MPI_EC_EDWARDS:
1478 		{
1479 			if (mpi_ec_get_affine(x, y, point, ctx))
1480 				goto leave;
1481 
1482 			mpi_resize(w, ctx->p->nlimbs);
1483 			w->nlimbs = ctx->p->nlimbs;
1484 
1485 			/* a · x^2 + y^2 - 1 - b · x^2 · y^2 == 0 */
1486 			ctx->pow2(x, x, ctx);
1487 			ctx->pow2(y, y, ctx);
1488 			if (ctx->dialect == ECC_DIALECT_ED25519)
1489 				ctx->subm(w, ctx->p, x, ctx);
1490 			else
1491 				ctx->mulm(w, ctx->a, x, ctx);
1492 			ctx->addm(w, w, y, ctx);
1493 			ctx->mulm(x, x, y, ctx);
1494 			ctx->mulm(x, x, ctx->b, ctx);
1495 			ctx->subm(w, w, x, ctx);
1496 			if (!mpi_cmp_ui(w, 1))
1497 				res = 1;
1498 		}
1499 		break;
1500 	}
1501 
1502 leave:
1503 	mpi_free(w);
1504 	mpi_free(x);
1505 	mpi_free(y);
1506 
1507 	return res;
1508 }
1509 EXPORT_SYMBOL_GPL(mpi_ec_curve_point);
1510