Linux Audio

Check our new training course

Loading...
Note: File does not exist in v6.13.7.
   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 */
  37MPI_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}
  48EXPORT_SYMBOL_GPL(mpi_point_new);
  49
  50/* Release the point object P.  P may be NULL. */
  51void mpi_point_release(MPI_POINT p)
  52{
  53	if (p) {
  54		mpi_point_free_parts(p);
  55		kfree(p);
  56	}
  57}
  58EXPORT_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 */
  63void 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}
  69EXPORT_SYMBOL_GPL(mpi_point_init);
  70
  71/* Release the parts of a point object. */
  72void 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}
  78EXPORT_SYMBOL_GPL(mpi_point_free_parts);
  79
  80/* Set the value from S into D.  */
  81static 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
  88static 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
 103static 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.  */
 114static 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
 122static 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
 128static 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
 136static 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.  */
 143static 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
 149static 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 */
 160static 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 */
 171static 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
 176static 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
 182static 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
 199static 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
 221static 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
 242static 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
 291static void ec_mul2_25519(MPI w, MPI u, struct mpi_ec_ctx *ctx)
 292{
 293	ec_addm_25519(w, u, u, ctx);
 294}
 295
 296static 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
 306static 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
 326static 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
 346static 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
 443static void ec_mul2_448(MPI w, MPI u, struct mpi_ec_ctx *ctx)
 444{
 445	ec_addm_448(w, u, u, ctx);
 446}
 447
 448static 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
 453struct 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
 464static 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.  */
 486static 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.  */
 493static 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.  */
 509static 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
 520static const char *const curve25519_bad_points[] = {
 521	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffed",
 522	"0x0000000000000000000000000000000000000000000000000000000000000000",
 523	"0x0000000000000000000000000000000000000000000000000000000000000001",
 524	"0x00b8495f16056286fdb1329ceb8d09da6ac49ff1fae35616aeb8413b7c7aebe0",
 525	"0x57119fd0dd4e22d8868e1c58c45c44045bef839c55b1d0b1248c50a3bc959c5f",
 526	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffec",
 527	"0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffee",
 528	NULL
 529};
 530
 531static 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
 545static const char *const *bad_points_table[] = {
 546	curve25519_bad_points,
 547	curve448_bad_points,
 548};
 549
 550static 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 */
 564void 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}
 649EXPORT_SYMBOL_GPL(mpi_ec_init);
 650
 651void 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}
 674EXPORT_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 */
 681int 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}
 753EXPORT_SYMBOL_GPL(mpi_ec_get_affine);
 754
 755/*  RESULT = 2 * POINT  (Weierstrass version). */
 756static 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). */
 834static 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). */
 845static 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  */
 917static void
 918mpi_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).*/
 934static 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).*/
1059static 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).*/
1072static 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 */
1166static 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 */
1191void 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}
1207EXPORT_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 */
1213void 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}
1398EXPORT_SYMBOL_GPL(mpi_ec_mul_point);
1399
1400/* Return true if POINT is on the curve described by CTX.  */
1401int 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
1502leave:
1503	mpi_free(w);
1504	mpi_free(x);
1505	mpi_free(y);
1506
1507	return res;
1508}
1509EXPORT_SYMBOL_GPL(mpi_ec_curve_point);