prj_pt.c
1 /* 2 * Copyright (C) 2017 - This file is part of libecc project 3 * 4 * Authors: 5 * Ryad BENADJILA <ryadbenadjila@gmail.com> 6 * Arnaud EBALARD <arnaud.ebalard@ssi.gouv.fr> 7 * Jean-Pierre FLORI <jean-pierre.flori@ssi.gouv.fr> 8 * 9 * Contributors: 10 * Nicolas VIVET <nicolas.vivet@ssi.gouv.fr> 11 * Karim KHALFALLAH <karim.khalfallah@ssi.gouv.fr> 12 * 13 * This software is licensed under a dual BSD and GPL v2 license. 14 * See LICENSE file at the root folder of the project. 15 */ 16 #include <libecc/curves/ec_shortw.h> 17 #include <libecc/curves/prj_pt.h> 18 #include <libecc/nn/nn_logical.h> 19 #include <libecc/nn/nn_add.h> 20 #include <libecc/nn/nn_rand.h> 21 #include <libecc/fp/fp_add.h> 22 #include <libecc/fp/fp_mul.h> 23 #include <libecc/fp/fp_montgomery.h> 24 #include <libecc/fp/fp_rand.h> 25 26 #define PRJ_PT_MAGIC ((word_t)(0xe1cd70babb1d5afeULL)) 27 28 /* 29 * Check given projective point has been correctly initialized (using 30 * prj_pt_init()). Returns 0 on success, -1 on error. 31 */ 32 int prj_pt_check_initialized(prj_pt_src_t in) 33 { 34 int ret; 35 36 MUST_HAVE(((in != NULL) && (in->magic == PRJ_PT_MAGIC)), ret, err); 37 ret = ec_shortw_crv_check_initialized(in->crv); 38 39 err: 40 return ret; 41 } 42 43 /* 44 * Initialize the projective point structure on given curve as the point at 45 * infinity. The function returns 0 on success, -1 on error. 46 */ 47 int prj_pt_init(prj_pt_t in, ec_shortw_crv_src_t curve) 48 { 49 int ret; 50 51 ret = ec_shortw_crv_check_initialized(curve); EG(ret, err); 52 53 MUST_HAVE((in != NULL), ret, err); 54 55 ret = fp_init(&(in->X), curve->a.ctx); EG(ret, err); 56 ret = fp_init(&(in->Y), curve->a.ctx); EG(ret, err); 57 ret = fp_init(&(in->Z), curve->a.ctx); EG(ret, err); 58 in->crv = curve; 59 in->magic = PRJ_PT_MAGIC; 60 61 err: 62 return ret; 63 } 64 65 /* 66 * Initialize the projective point structure on given curve using given 67 * coordinates. The function returns 0 on success, -1 on error. 68 */ 69 int prj_pt_init_from_coords(prj_pt_t in, 70 ec_shortw_crv_src_t curve, 71 fp_src_t xcoord, fp_src_t ycoord, fp_src_t zcoord) 72 { 73 int ret; 74 75 ret = prj_pt_init(in, curve); EG(ret, err); 76 ret = fp_copy(&(in->X), xcoord); EG(ret, err); 77 ret = fp_copy(&(in->Y), ycoord); EG(ret, err); 78 ret = fp_copy(&(in->Z), zcoord); 79 80 err: 81 return ret; 82 } 83 84 /* 85 * Uninit given projective point structure. The function returns 0 on success, 86 * -1 on error. This is an error if passed point has not already been 87 * initialized first. 88 */ 89 void prj_pt_uninit(prj_pt_t in) 90 { 91 if((in != NULL) && (in->magic == PRJ_PT_MAGIC) && (in->crv != NULL)){ 92 in->crv = NULL; 93 in->magic = WORD(0); 94 95 fp_uninit(&(in->X)); 96 fp_uninit(&(in->Y)); 97 fp_uninit(&(in->Z)); 98 } 99 100 return; 101 } 102 103 /* 104 * Checks if projective point is the point at infinity (last coordinate is 0). 105 * In that case, 'iszero' out parameter is set to 1. It is set to 0 if the 106 * point is not the point at infinity. The function returns 0 on success, -1 on 107 * error. On error, 'iszero' is not meaningful. 108 */ 109 int prj_pt_iszero(prj_pt_src_t in, int *iszero) 110 { 111 int ret; 112 113 ret = prj_pt_check_initialized(in); EG(ret, err); 114 ret = fp_iszero(&(in->Z), iszero); 115 116 err: 117 return ret; 118 } 119 120 /* 121 * Set given projective point 'out' to the point at infinity. The functions 122 * returns 0 on success, -1 on error. 123 */ 124 int prj_pt_zero(prj_pt_t out) 125 { 126 int ret; 127 128 ret = prj_pt_check_initialized(out); EG(ret, err); 129 130 ret = fp_zero(&(out->X)); EG(ret, err); 131 ret = fp_one(&(out->Y)); EG(ret, err); 132 ret = fp_zero(&(out->Z)); 133 134 err: 135 return ret; 136 } 137 138 /* 139 * Check if a projective point is indeed on its curve. The function sets 140 * 'on_curve' out parameter to 1 if the point is on the curve, 0 if not. 141 * The function returns 0 on success, -1 on error. 'on_curve' is not 142 * meaningful on error. 143 */ 144 int prj_pt_is_on_curve(prj_pt_src_t in, int *on_curve) 145 { 146 int ret, cmp; 147 148 /* In order to check that we are on the curve, we 149 * use the projective formula of the curve: 150 * 151 * Y**2 * Z = X**3 + a * X * Z**2 + b * Z**3 152 * 153 */ 154 fp X, Y, Z; 155 X.magic = Y.magic = Z.magic = WORD(0); 156 157 ret = prj_pt_check_initialized(in); EG(ret, err); 158 ret = ec_shortw_crv_check_initialized(in->crv); EG(ret, err); 159 MUST_HAVE((on_curve != NULL), ret, err); 160 161 ret = fp_init(&X, in->X.ctx); EG(ret, err); 162 ret = fp_init(&Y, in->X.ctx); EG(ret, err); 163 ret = fp_init(&Z, in->X.ctx); EG(ret, err); 164 165 /* Compute X**3 + a * X * Z**2 + b * Z**3 on one side */ 166 ret = fp_sqr(&X, &(in->X)); EG(ret, err); 167 ret = fp_mul(&X, &X, &(in->X)); EG(ret, err); 168 ret = fp_mul(&Z, &(in->X), &(in->crv->a)); EG(ret, err); 169 ret = fp_mul(&Y, &(in->crv->b), &(in->Z)); EG(ret, err); 170 ret = fp_add(&Z, &Z, &Y); EG(ret, err); 171 ret = fp_mul(&Z, &Z, &(in->Z)); EG(ret, err); 172 ret = fp_mul(&Z, &Z, &(in->Z)); EG(ret, err); 173 ret = fp_add(&X, &X, &Z); EG(ret, err); 174 175 /* Compute Y**2 * Z on the other side */ 176 ret = fp_sqr(&Y, &(in->Y)); EG(ret, err); 177 ret = fp_mul(&Y, &Y, &(in->Z)); EG(ret, err); 178 179 /* Compare the two values */ 180 ret = fp_cmp(&X, &Y, &cmp); EG(ret, err); 181 182 (*on_curve) = (!cmp); 183 184 err: 185 fp_uninit(&X); 186 fp_uninit(&Y); 187 fp_uninit(&Z); 188 189 return ret; 190 } 191 192 /* 193 * The function copies 'in' projective point to 'out'. 'out' is initialized by 194 * the function. The function returns 0 on sucess, -1 on error. 195 */ 196 int prj_pt_copy(prj_pt_t out, prj_pt_src_t in) 197 { 198 int ret; 199 200 ret = prj_pt_check_initialized(in); EG(ret, err); 201 202 ret = prj_pt_init(out, in->crv); EG(ret, err); 203 204 ret = fp_copy(&(out->X), &(in->X)); EG(ret, err); 205 ret = fp_copy(&(out->Y), &(in->Y)); EG(ret, err); 206 ret = fp_copy(&(out->Z), &(in->Z)); 207 208 err: 209 return ret; 210 } 211 212 /* 213 * Convert given projective point 'in' to affine representation in 'out'. 'out' 214 * is initialized by the function. The function returns 0 on success, -1 on 215 * error. Passing the point at infinty to the function is considered as an 216 * error. 217 */ 218 int prj_pt_to_aff(aff_pt_t out, prj_pt_src_t in) 219 { 220 int ret, iszero; 221 222 ret = prj_pt_check_initialized(in); EG(ret, err); 223 224 ret = prj_pt_iszero(in, &iszero); EG(ret, err); 225 MUST_HAVE((!iszero), ret, err); 226 227 ret = aff_pt_init(out, in->crv); EG(ret, err); 228 229 ret = fp_inv(&(out->x), &(in->Z)); EG(ret, err); 230 ret = fp_mul(&(out->y), &(in->Y), &(out->x)); EG(ret, err); 231 ret = fp_mul(&(out->x), &(in->X), &(out->x)); 232 233 err: 234 return ret; 235 } 236 237 /* 238 * Get the unique Z = 1 projective point representation ("equivalent" to affine 239 * point). The function returns 0 on success, -1 on error. 240 */ 241 int prj_pt_unique(prj_pt_t out, prj_pt_src_t in) 242 { 243 int ret, iszero; 244 245 ret = prj_pt_check_initialized(in); EG(ret, err); 246 ret = prj_pt_iszero(in, &iszero); EG(ret, err); 247 MUST_HAVE((!iszero), ret, err); 248 249 if(out == in){ 250 /* Aliasing case */ 251 fp tmp; 252 tmp.magic = WORD(0); 253 254 ret = fp_init(&tmp, (in->Z).ctx); EG(ret, err); 255 ret = fp_inv(&tmp, &(in->Z)); EG(ret, err1); 256 ret = fp_mul(&(out->Y), &(in->Y), &tmp); EG(ret, err1); 257 ret = fp_mul(&(out->X), &(in->X), &tmp); EG(ret, err1); 258 ret = fp_one(&(out->Z)); EG(ret, err1); 259 err1: 260 fp_uninit(&tmp); EG(ret, err); 261 } 262 else{ 263 ret = prj_pt_init(out, in->crv); EG(ret, err); 264 ret = fp_inv(&(out->X), &(in->Z)); EG(ret, err); 265 ret = fp_mul(&(out->Y), &(in->Y), &(out->X)); EG(ret, err); 266 ret = fp_mul(&(out->X), &(in->X), &(out->X)); EG(ret, err); 267 ret = fp_one(&(out->Z)); EG(ret, err); 268 } 269 270 271 err: 272 return ret; 273 } 274 275 /* 276 * Converts affine point 'in' to projective representation in 'out'. 'out' is 277 * initialized by the function. The function returns 0 on success, -1 on error. 278 */ 279 int ec_shortw_aff_to_prj(prj_pt_t out, aff_pt_src_t in) 280 { 281 int ret, on_curve; 282 283 ret = aff_pt_check_initialized(in); EG(ret, err); 284 285 /* The input affine point must be on the curve */ 286 ret = aff_pt_is_on_curve(in, &on_curve); EG(ret, err); 287 MUST_HAVE(on_curve, ret, err); 288 289 ret = prj_pt_init(out, in->crv); EG(ret, err); 290 ret = fp_copy(&(out->X), &(in->x)); EG(ret, err); 291 ret = fp_copy(&(out->Y), &(in->y)); EG(ret, err); 292 ret = nn_one(&(out->Z).fp_val); /* Z = 1 */ 293 294 err: 295 return ret; 296 } 297 298 /* 299 * Compare projective points 'in1' and 'in2'. On success, 'cmp' is set to 300 * the result of the comparison (0 if in1 == in2, !0 if in1 != in2). The 301 * function returns 0 on success, -1 on error. 302 */ 303 int prj_pt_cmp(prj_pt_src_t in1, prj_pt_src_t in2, int *cmp) 304 { 305 fp X1, X2, Y1, Y2; 306 int ret, x_cmp, y_cmp; 307 X1.magic = X2.magic = Y1.magic = Y2.magic = WORD(0); 308 309 MUST_HAVE((cmp != NULL), ret, err); 310 ret = prj_pt_check_initialized(in1); EG(ret, err); 311 ret = prj_pt_check_initialized(in2); EG(ret, err); 312 313 MUST_HAVE((in1->crv == in2->crv), ret, err); 314 315 ret = fp_init(&X1, (in1->X).ctx); EG(ret, err); 316 ret = fp_init(&X2, (in2->X).ctx); EG(ret, err); 317 ret = fp_init(&Y1, (in1->Y).ctx); EG(ret, err); 318 ret = fp_init(&Y2, (in2->Y).ctx); EG(ret, err); 319 320 /* 321 * Montgomery multiplication is used as it is faster than 322 * usual multiplication and the spurious multiplicative 323 * factor does not matter. 324 */ 325 ret = fp_mul_monty(&X1, &(in1->X), &(in2->Z)); EG(ret, err); 326 ret = fp_mul_monty(&X2, &(in2->X), &(in1->Z)); EG(ret, err); 327 ret = fp_mul_monty(&Y1, &(in1->Y), &(in2->Z)); EG(ret, err); 328 ret = fp_mul_monty(&Y2, &(in2->Y), &(in1->Z)); EG(ret, err); 329 330 ret = fp_mul_monty(&X1, &(in1->X), &(in2->Z)); EG(ret, err); 331 ret = fp_mul_monty(&X2, &(in2->X), &(in1->Z)); EG(ret, err); 332 ret = fp_mul_monty(&Y1, &(in1->Y), &(in2->Z)); EG(ret, err); 333 ret = fp_mul_monty(&Y2, &(in2->Y), &(in1->Z)); EG(ret, err); 334 ret = fp_cmp(&X1, &X2, &x_cmp); EG(ret, err); 335 ret = fp_cmp(&Y1, &Y2, &y_cmp); 336 337 if (!ret) { 338 (*cmp) = (x_cmp | y_cmp); 339 } 340 341 err: 342 fp_uninit(&Y2); 343 fp_uninit(&Y1); 344 fp_uninit(&X2); 345 fp_uninit(&X1); 346 347 return ret; 348 } 349 350 /* 351 * NOTE: this internal functions assumes that upper layer have checked that in1 and in2 352 * are initialized, and that cmp is not NULL. 353 */ 354 ATTRIBUTE_WARN_UNUSED_RET static inline int _prj_pt_eq_or_opp_X(prj_pt_src_t in1, prj_pt_src_t in2, int *cmp) 355 { 356 int ret; 357 fp X1, X2; 358 X1.magic = X2.magic = WORD(0); 359 360 /* 361 * Montgomery multiplication is used as it is faster than 362 * usual multiplication and the spurious multiplicative 363 * factor does not matter. 364 */ 365 ret = fp_init(&X1, (in1->X).ctx); EG(ret, err); 366 ret = fp_init(&X2, (in2->X).ctx); EG(ret, err); 367 ret = fp_mul_monty(&X1, &(in1->X), &(in2->Z)); EG(ret, err); 368 ret = fp_mul_monty(&X2, &(in2->X), &(in1->Z)); EG(ret, err); 369 ret = fp_cmp(&X1, &X2, cmp); 370 371 err: 372 fp_uninit(&X1); 373 fp_uninit(&X2); 374 375 return ret; 376 } 377 378 /* 379 * NOTE: this internal functions assumes that upper layer have checked that in1 and in2 380 * are initialized, and that eq_or_opp is not NULL. 381 */ 382 ATTRIBUTE_WARN_UNUSED_RET static inline int _prj_pt_eq_or_opp_Y(prj_pt_src_t in1, prj_pt_src_t in2, int *eq_or_opp) 383 { 384 int ret; 385 fp Y1, Y2; 386 Y1.magic = Y2.magic = WORD(0); 387 388 /* 389 * Montgomery multiplication is used as it is faster than 390 * usual multiplication and the spurious multiplicative 391 * factor does not matter. 392 */ 393 ret = fp_init(&Y1, (in1->Y).ctx); EG(ret, err); 394 ret = fp_init(&Y2, (in2->Y).ctx); EG(ret, err); 395 ret = fp_mul_monty(&Y1, &(in1->Y), &(in2->Z)); EG(ret, err); 396 ret = fp_mul_monty(&Y2, &(in2->Y), &(in1->Z)); EG(ret, err); 397 ret = fp_eq_or_opp(&Y1, &Y2, eq_or_opp); 398 399 err: 400 fp_uninit(&Y1); 401 fp_uninit(&Y2); 402 403 return ret; 404 } 405 406 /* 407 * The functions tests if given projective points 'in1' and 'in2' are equal or 408 * opposite. On success, the result of the comparison is given via 'eq_or_opp' 409 * out parameter (1 if equal or opposite, 0 otherwise). The function returns 410 * 0 on succes, -1 on error. 411 */ 412 int prj_pt_eq_or_opp(prj_pt_src_t in1, prj_pt_src_t in2, int *eq_or_opp) 413 { 414 int ret, cmp, _eq_or_opp; 415 416 ret = prj_pt_check_initialized(in1); EG(ret, err); 417 ret = prj_pt_check_initialized(in2); EG(ret, err); 418 MUST_HAVE((in1->crv == in2->crv), ret, err); 419 MUST_HAVE((eq_or_opp != NULL), ret, err); 420 421 ret = _prj_pt_eq_or_opp_X(in1, in2, &cmp); EG(ret, err); 422 ret = _prj_pt_eq_or_opp_Y(in1, in2, &_eq_or_opp); 423 424 if (!ret) { 425 (*eq_or_opp) = ((cmp == 0) & _eq_or_opp); 426 } 427 428 err: 429 return ret; 430 } 431 432 /* Compute the opposite of a projective point. Supports aliasing. 433 * Returns 0 on success, -1 on failure. 434 */ 435 int prj_pt_neg(prj_pt_t out, prj_pt_src_t in) 436 { 437 int ret; 438 439 ret = prj_pt_check_initialized(in); EG(ret, err); 440 441 if (out != in) { /* Copy point if not aliased */ 442 ret = prj_pt_init(out, in->crv); EG(ret, err); 443 ret = prj_pt_copy(out, in); EG(ret, err); 444 } 445 446 /* Then, negate Y */ 447 ret = fp_neg(&(out->Y), &(out->Y)); 448 449 err: 450 return ret; 451 } 452 453 /* 454 * Import a projective point from a buffer with the following layout; the 3 455 * coordinates (elements of Fp) are each encoded on p_len bytes, where p_len 456 * is the size of p in bytes (e.g. 66 for a prime p of 521 bits). Each 457 * coordinate is encoded in big endian. Size of buffer must exactly match 458 * 3 * p_len. The projective point is initialized by the function. 459 * 460 * The function returns 0 on success, -1 on error. 461 */ 462 int prj_pt_import_from_buf(prj_pt_t pt, 463 const u8 *pt_buf, 464 u16 pt_buf_len, ec_shortw_crv_src_t crv) 465 { 466 int on_curve, ret; 467 fp_ctx_src_t ctx; 468 u16 coord_len; 469 470 ret = ec_shortw_crv_check_initialized(crv); EG(ret, err); 471 MUST_HAVE((pt_buf != NULL) && (pt != NULL), ret, err); 472 473 ctx = crv->a.ctx; 474 coord_len = (u16)BYTECEIL(ctx->p_bitlen); 475 MUST_HAVE((pt_buf_len == (3 * coord_len)), ret, err); 476 477 ret = fp_init_from_buf(&(pt->X), ctx, pt_buf, coord_len); EG(ret, err); 478 ret = fp_init_from_buf(&(pt->Y), ctx, pt_buf + coord_len, coord_len); EG(ret, err); 479 ret = fp_init_from_buf(&(pt->Z), ctx, pt_buf + (2 * coord_len), coord_len); EG(ret, err); 480 481 /* Set the curve */ 482 pt->crv = crv; 483 484 /* Mark the point as initialized */ 485 pt->magic = PRJ_PT_MAGIC; 486 487 /* Check that the point is indeed on the provided curve, uninitialize it 488 * if this is not the case. 489 */ 490 ret = prj_pt_is_on_curve(pt, &on_curve); EG(ret, err); 491 if (!on_curve){ 492 prj_pt_uninit(pt); 493 ret = -1; 494 } 495 496 err: 497 PTR_NULLIFY(ctx); 498 499 return ret; 500 } 501 502 /* 503 * Import a projective point from an affine point buffer with the following layout; the 2 504 * coordinates (elements of Fp) are each encoded on p_len bytes, where p_len 505 * is the size of p in bytes (e.g. 66 for a prime p of 521 bits). Each 506 * coordinate is encoded in big endian. Size of buffer must exactly match 507 * 2 * p_len. The projective point is initialized by the function. 508 * 509 * The function returns 0 on success, -1 on error. 510 */ 511 int prj_pt_import_from_aff_buf(prj_pt_t pt, 512 const u8 *pt_buf, 513 u16 pt_buf_len, ec_shortw_crv_src_t crv) 514 { 515 int ret, on_curve; 516 fp_ctx_src_t ctx; 517 u16 coord_len; 518 519 ret = ec_shortw_crv_check_initialized(crv); EG(ret, err); 520 MUST_HAVE((pt_buf != NULL) && (pt != NULL), ret, err); 521 522 ctx = crv->a.ctx; 523 coord_len = (u16)BYTECEIL(ctx->p_bitlen); 524 MUST_HAVE((pt_buf_len == (2 * coord_len)), ret, err); 525 526 ret = fp_init_from_buf(&(pt->X), ctx, pt_buf, coord_len); EG(ret, err); 527 ret = fp_init_from_buf(&(pt->Y), ctx, pt_buf + coord_len, coord_len); EG(ret, err); 528 /* Z coordinate is set to 1 */ 529 ret = fp_init(&(pt->Z), ctx); EG(ret, err); 530 ret = fp_one(&(pt->Z)); EG(ret, err); 531 532 /* Set the curve */ 533 pt->crv = crv; 534 535 /* Mark the point as initialized */ 536 pt->magic = PRJ_PT_MAGIC; 537 538 /* Check that the point is indeed on the provided curve, uninitialize it 539 * if this is not the case. 540 */ 541 ret = prj_pt_is_on_curve(pt, &on_curve); EG(ret, err); 542 if (!on_curve){ 543 prj_pt_uninit(pt); 544 ret = -1; 545 } 546 547 err: 548 PTR_NULLIFY(ctx); 549 550 return ret; 551 } 552 553 554 /* Export a projective point to a buffer with the following layout; the 3 555 * coordinates (elements of Fp) are each encoded on p_len bytes, where p_len 556 * is the size of p in bytes (e.g. 66 for a prime p of 521 bits). Each 557 * coordinate is encoded in big endian. Size of buffer must exactly match 558 * 3 * p_len. 559 * 560 * The function returns 0 on success, -1 on error. 561 */ 562 int prj_pt_export_to_buf(prj_pt_src_t pt, u8 *pt_buf, u32 pt_buf_len) 563 { 564 fp_ctx_src_t ctx; 565 u16 coord_len; 566 int ret, on_curve; 567 568 ret = prj_pt_check_initialized(pt); EG(ret, err); 569 570 MUST_HAVE((pt_buf != NULL), ret, err); 571 572 /* The point to be exported must be on the curve */ 573 ret = prj_pt_is_on_curve(pt, &on_curve); EG(ret, err); 574 MUST_HAVE((on_curve), ret, err); 575 576 ctx = pt->crv->a.ctx; 577 coord_len = (u16)BYTECEIL(ctx->p_bitlen); 578 MUST_HAVE((pt_buf_len == (3 * coord_len)), ret, err); 579 580 /* Export the three coordinates */ 581 ret = fp_export_to_buf(pt_buf, coord_len, &(pt->X)); EG(ret, err); 582 ret = fp_export_to_buf(pt_buf + coord_len, coord_len, &(pt->Y)); EG(ret, err); 583 ret = fp_export_to_buf(pt_buf + (2 * coord_len), coord_len, &(pt->Z)); 584 585 err: 586 PTR_NULLIFY(ctx); 587 588 return ret; 589 } 590 591 /* 592 * Export a projective point to an affine point buffer with the following 593 * layout; the 2 coordinates (elements of Fp) are each encoded on p_len bytes, 594 * where p_len is the size of p in bytes (e.g. 66 for a prime p of 521 bits). 595 * Each coordinate is encoded in big endian. Size of buffer must exactly match 596 * 2 * p_len. 597 * 598 * The function returns 0 on success, -1 on error. 599 */ 600 int prj_pt_export_to_aff_buf(prj_pt_src_t pt, u8 *pt_buf, u32 pt_buf_len) 601 { 602 int ret, on_curve; 603 aff_pt tmp_aff; 604 tmp_aff.magic = WORD(0); 605 606 ret = prj_pt_check_initialized(pt); EG(ret, err); 607 608 MUST_HAVE((pt_buf != NULL), ret, err); 609 610 /* The point to be exported must be on the curve */ 611 ret = prj_pt_is_on_curve(pt, &on_curve); EG(ret, err); 612 MUST_HAVE((on_curve), ret, err); 613 614 /* Move to the affine unique representation */ 615 ret = prj_pt_to_aff(&tmp_aff, pt); EG(ret, err); 616 617 /* Export the affine point to the buffer */ 618 ret = aff_pt_export_to_buf(&tmp_aff, pt_buf, pt_buf_len); 619 620 err: 621 aff_pt_uninit(&tmp_aff); 622 623 return ret; 624 } 625 626 627 #ifdef NO_USE_COMPLETE_FORMULAS 628 629 /* 630 * The function is an internal one: no check is performed on parameters, 631 * this MUST be done by the caller: 632 * 633 * - in is initialized 634 * - in and out must not be aliased 635 * 636 * The function will initialize 'out'. The function returns 0 on success, -1 637 * on error. 638 */ 639 ATTRIBUTE_WARN_UNUSED_RET static int __prj_pt_dbl_monty_no_cf(prj_pt_t out, prj_pt_src_t in) 640 { 641 fp XX, ZZ, w, s, ss, sss, R, RR, B, h; 642 int ret; 643 XX.magic = ZZ.magic = w.magic = s.magic = WORD(0); 644 ss.magic = sss.magic = R.magic = WORD(0); 645 RR.magic = B.magic = h.magic = WORD(0); 646 647 ret = prj_pt_init(out, in->crv); EG(ret, err); 648 649 ret = fp_init(&XX, out->crv->a.ctx); EG(ret, err); 650 ret = fp_init(&ZZ, out->crv->a.ctx); EG(ret, err); 651 ret = fp_init(&w, out->crv->a.ctx); EG(ret, err); 652 ret = fp_init(&s, out->crv->a.ctx); EG(ret, err); 653 ret = fp_init(&ss, out->crv->a.ctx); EG(ret, err); 654 ret = fp_init(&sss, out->crv->a.ctx); EG(ret, err); 655 ret = fp_init(&R, out->crv->a.ctx); EG(ret, err); 656 ret = fp_init(&RR, out->crv->a.ctx); EG(ret, err); 657 ret = fp_init(&B, out->crv->a.ctx); EG(ret, err); 658 ret = fp_init(&h, out->crv->a.ctx); EG(ret, err); 659 660 /* XX = X1² */ 661 ret = fp_sqr_monty(&XX, &(in->X)); EG(ret, err); 662 663 /* ZZ = Z1² */ 664 ret = fp_sqr_monty(&ZZ, &(in->Z)); EG(ret, err); 665 666 /* w = a*ZZ+3*XX */ 667 ret = fp_mul_monty(&w, &(in->crv->a_monty), &ZZ); EG(ret, err); 668 ret = fp_add_monty(&w, &w, &XX); EG(ret, err); 669 ret = fp_add_monty(&w, &w, &XX); EG(ret, err); 670 ret = fp_add_monty(&w, &w, &XX); EG(ret, err); 671 672 /* s = 2*Y1*Z1 */ 673 ret = fp_mul_monty(&s, &(in->Y), &(in->Z)); EG(ret, err); 674 ret = fp_add_monty(&s, &s, &s); EG(ret, err); 675 676 /* ss = s² */ 677 ret = fp_sqr_monty(&ss, &s); EG(ret, err); 678 679 /* sss = s*ss */ 680 ret = fp_mul_monty(&sss, &s, &ss); EG(ret, err); 681 682 /* R = Y1*s */ 683 ret = fp_mul_monty(&R, &(in->Y), &s); EG(ret, err); 684 685 /* RR = R² */ 686 ret = fp_sqr_monty(&RR, &R); EG(ret, err); 687 688 /* B = (X1+R)²-XX-RR */ 689 ret = fp_add_monty(&R, &R, &(in->X)); EG(ret, err); 690 ret = fp_sqr_monty(&B, &R); EG(ret, err); 691 ret = fp_sub_monty(&B, &B, &XX); EG(ret, err); 692 ret = fp_sub_monty(&B, &B, &RR); EG(ret, err); 693 694 /* h = w²-2*B */ 695 ret = fp_sqr_monty(&h, &w); EG(ret, err); 696 ret = fp_sub_monty(&h, &h, &B); EG(ret, err); 697 ret = fp_sub_monty(&h, &h, &B); EG(ret, err); 698 699 /* X3 = h*s */ 700 ret = fp_mul_monty(&(out->X), &h, &s); EG(ret, err); 701 702 /* Y3 = w*(B-h)-2*RR */ 703 ret = fp_sub_monty(&B, &B, &h); EG(ret, err); 704 ret = fp_mul_monty(&(out->Y), &w, &B); EG(ret, err); 705 ret = fp_sub_monty(&(out->Y), &(out->Y), &RR); EG(ret, err); 706 ret = fp_sub_monty(&(out->Y), &(out->Y), &RR); EG(ret, err); 707 708 /* Z3 = sss */ 709 ret = fp_copy(&(out->Z), &sss); 710 711 err: 712 fp_uninit(&XX); 713 fp_uninit(&ZZ); 714 fp_uninit(&w); 715 fp_uninit(&s); 716 fp_uninit(&ss); 717 fp_uninit(&sss); 718 fp_uninit(&R); 719 fp_uninit(&RR); 720 fp_uninit(&B); 721 fp_uninit(&h); 722 723 return ret; 724 } 725 726 /* 727 * The function is an internal one: no check is performed on parameters, 728 * this MUST be done by the caller: 729 * 730 * - in1 and in2 are initialized 731 * - in1 and in2 are on the same curve 732 * - in1/in2 and out must not be aliased 733 * - in1 and in2 must not be equal, opposite or have identical value 734 * 735 * The function will initialize 'out'. The function returns 0 on success, -1 736 * on error. 737 */ 738 ATTRIBUTE_WARN_UNUSED_RET static int ___prj_pt_add_monty_no_cf(prj_pt_t out, 739 prj_pt_src_t in1, 740 prj_pt_src_t in2) 741 { 742 fp Y1Z2, X1Z2, Z1Z2, u, uu, v, vv, vvv, R, A; 743 int ret; 744 Y1Z2.magic = X1Z2.magic = Z1Z2.magic = u.magic = uu.magic = v.magic = WORD(0); 745 vv.magic = vvv.magic = R.magic = A.magic = WORD(0); 746 747 ret = prj_pt_init(out, in1->crv); EG(ret, err); 748 749 ret = fp_init(&Y1Z2, out->crv->a.ctx); EG(ret, err); 750 ret = fp_init(&X1Z2, out->crv->a.ctx); EG(ret, err); 751 ret = fp_init(&Z1Z2, out->crv->a.ctx); EG(ret, err); 752 ret = fp_init(&u, out->crv->a.ctx); EG(ret, err); 753 ret = fp_init(&uu, out->crv->a.ctx); EG(ret, err); 754 ret = fp_init(&v, out->crv->a.ctx); EG(ret, err); 755 ret = fp_init(&vv, out->crv->a.ctx); EG(ret, err); 756 ret = fp_init(&vvv, out->crv->a.ctx); EG(ret, err); 757 ret = fp_init(&R, out->crv->a.ctx); EG(ret, err); 758 ret = fp_init(&A, out->crv->a.ctx); EG(ret, err); 759 760 /* Y1Z2 = Y1*Z2 */ 761 ret = fp_mul_monty(&Y1Z2, &(in1->Y), &(in2->Z)); EG(ret, err); 762 763 /* X1Z2 = X1*Z2 */ 764 ret = fp_mul_monty(&X1Z2, &(in1->X), &(in2->Z)); EG(ret, err); 765 766 /* Z1Z2 = Z1*Z2 */ 767 ret = fp_mul_monty(&Z1Z2, &(in1->Z), &(in2->Z)); EG(ret, err); 768 769 /* u = Y2*Z1-Y1Z2 */ 770 ret = fp_mul_monty(&u, &(in2->Y), &(in1->Z)); EG(ret, err); 771 ret = fp_sub_monty(&u, &u, &Y1Z2); EG(ret, err); 772 773 /* uu = u² */ 774 ret = fp_sqr_monty(&uu, &u); EG(ret, err); 775 776 /* v = X2*Z1-X1Z2 */ 777 ret = fp_mul_monty(&v, &(in2->X), &(in1->Z)); EG(ret, err); 778 ret = fp_sub_monty(&v, &v, &X1Z2); EG(ret, err); 779 780 /* vv = v² */ 781 ret = fp_sqr_monty(&vv, &v); EG(ret, err); 782 783 /* vvv = v*vv */ 784 ret = fp_mul_monty(&vvv, &v, &vv); EG(ret, err); 785 786 /* R = vv*X1Z2 */ 787 ret = fp_mul_monty(&R, &vv, &X1Z2); EG(ret, err); 788 789 /* A = uu*Z1Z2-vvv-2*R */ 790 ret = fp_mul_monty(&A, &uu, &Z1Z2); EG(ret, err); 791 ret = fp_sub_monty(&A, &A, &vvv); EG(ret, err); 792 ret = fp_sub_monty(&A, &A, &R); EG(ret, err); 793 ret = fp_sub_monty(&A, &A, &R); EG(ret, err); 794 795 /* X3 = v*A */ 796 ret = fp_mul_monty(&(out->X), &v, &A); EG(ret, err); 797 798 /* Y3 = u*(R-A)-vvv*Y1Z2 */ 799 ret = fp_sub_monty(&R, &R, &A); EG(ret, err); 800 ret = fp_mul_monty(&(out->Y), &u, &R); EG(ret, err); 801 ret = fp_mul_monty(&R, &vvv, &Y1Z2); EG(ret, err); 802 ret = fp_sub_monty(&(out->Y), &(out->Y), &R); EG(ret, err); 803 804 /* Z3 = vvv*Z1Z2 */ 805 ret = fp_mul_monty(&(out->Z), &vvv, &Z1Z2); 806 807 err: 808 fp_uninit(&Y1Z2); 809 fp_uninit(&X1Z2); 810 fp_uninit(&Z1Z2); 811 fp_uninit(&u); 812 fp_uninit(&uu); 813 fp_uninit(&v); 814 fp_uninit(&vv); 815 fp_uninit(&vvv); 816 fp_uninit(&R); 817 fp_uninit(&A); 818 819 return ret; 820 } 821 822 /* 823 * Public version of the addition w/o complete formulas to handle the case 824 * where the inputs are zero or opposite. Returns 0 on success, -1 on error. 825 */ 826 ATTRIBUTE_WARN_UNUSED_RET static int __prj_pt_add_monty_no_cf(prj_pt_t out, prj_pt_src_t in1, prj_pt_src_t in2) 827 { 828 int ret, iszero, eq_or_opp, cmp; 829 830 ret = prj_pt_check_initialized(in1); EG(ret, err); 831 ret = prj_pt_check_initialized(in2); EG(ret, err); 832 MUST_HAVE((in1->crv == in2->crv), ret, err); 833 834 ret = prj_pt_iszero(in1, &iszero); EG(ret, err); 835 if (iszero) { 836 /* in1 at infinity, output in2 in all cases */ 837 ret = prj_pt_init(out, in2->crv); EG(ret, err); 838 ret = prj_pt_copy(out, in2); EG(ret, err); 839 } else { 840 /* in1 not at infinity, output in2 */ 841 ret = prj_pt_iszero(in2, &iszero); EG(ret, err); 842 if (iszero) { 843 /* in2 at infinity, output in1 */ 844 ret = prj_pt_init(out, in1->crv); EG(ret, err); 845 ret = prj_pt_copy(out, in1); EG(ret, err); 846 } else { 847 /* enither in1, nor in2 at infinity */ 848 849 /* 850 * The following test which guarantees in1 and in2 are not 851 * equal or opposite needs to be rewritten because it 852 * has a *HUGE* impact on perf (ec_self_tests run on 853 * all test vectors takes 24 times as long with this 854 * enabled). The same exists in non monty version. 855 */ 856 ret = prj_pt_eq_or_opp(in1, in2, &eq_or_opp); EG(ret, err); 857 if (eq_or_opp) { 858 /* in1 and in2 are either equal or opposite */ 859 ret = prj_pt_cmp(in1, in2, &cmp); EG(ret, err); 860 if (cmp == 0) { 861 /* in1 == in2 => doubling w/o cf */ 862 ret = __prj_pt_dbl_monty_no_cf(out, in1); EG(ret, err); 863 } else { 864 /* in1 == -in2 => output zero (point at infinity) */ 865 ret = prj_pt_init(out, in1->crv); EG(ret, err); 866 ret = prj_pt_zero(out); EG(ret, err); 867 } 868 } else { 869 /* 870 * in1 and in2 are neither 0, nor equal or 871 * opposite. Use the basic monty addition 872 * implementation w/o complete formulas. 873 */ 874 ret = ___prj_pt_add_monty_no_cf(out, in1, in2); EG(ret, err); 875 } 876 } 877 } 878 879 err: 880 return ret; 881 } 882 883 884 #else /* NO_USE_COMPLETE_FORMULAS */ 885 886 887 /* 888 * If NO_USE_COMPLETE_FORMULAS flag is not defined addition formulas from Algorithm 3 889 * of https://joostrenes.nl/publications/complete.pdf are used, otherwise 890 * http://www.hyperelliptic.org/EFD/g1p/auto-shortw-projective.html#doubling-dbl-2007-bl 891 */ 892 ATTRIBUTE_WARN_UNUSED_RET static int __prj_pt_dbl_monty_cf(prj_pt_t out, prj_pt_src_t in) 893 { 894 fp t0, t1, t2, t3; 895 int ret; 896 t0.magic = t1.magic = t2.magic = t3.magic = WORD(0); 897 898 ret = prj_pt_init(out, in->crv); EG(ret, err); 899 900 ret = fp_init(&t0, out->crv->a.ctx); EG(ret, err); 901 ret = fp_init(&t1, out->crv->a.ctx); EG(ret, err); 902 ret = fp_init(&t2, out->crv->a.ctx); EG(ret, err); 903 ret = fp_init(&t3, out->crv->a.ctx); EG(ret, err); 904 905 ret = fp_mul_monty(&t0, &in->X, &in->X); EG(ret, err); 906 ret = fp_mul_monty(&t1, &in->Y, &in->Y); EG(ret, err); 907 ret = fp_mul_monty(&t2, &in->Z, &in->Z); EG(ret, err); 908 ret = fp_mul_monty(&t3, &in->X, &in->Y); EG(ret, err); 909 ret = fp_add_monty(&t3, &t3, &t3); EG(ret, err); 910 911 ret = fp_mul_monty(&out->Z, &in->X, &in->Z); EG(ret, err); 912 ret = fp_add_monty(&out->Z, &out->Z, &out->Z); EG(ret, err); 913 ret = fp_mul_monty(&out->X, &in->crv->a_monty, &out->Z); EG(ret, err); 914 ret = fp_mul_monty(&out->Y, &in->crv->b3_monty, &t2); EG(ret, err); 915 ret = fp_add_monty(&out->Y, &out->X, &out->Y); EG(ret, err); 916 917 ret = fp_sub_monty(&out->X, &t1, &out->Y); EG(ret, err); 918 ret = fp_add_monty(&out->Y, &t1, &out->Y); EG(ret, err); 919 ret = fp_mul_monty(&out->Y, &out->X, &out->Y); EG(ret, err); 920 ret = fp_mul_monty(&out->X, &t3, &out->X); EG(ret, err); 921 ret = fp_mul_monty(&out->Z, &in->crv->b3_monty, &out->Z); EG(ret, err); 922 923 ret = fp_mul_monty(&t2, &in->crv->a_monty, &t2); EG(ret, err); 924 ret = fp_sub_monty(&t3, &t0, &t2); EG(ret, err); 925 ret = fp_mul_monty(&t3, &in->crv->a_monty, &t3); EG(ret, err); 926 ret = fp_add_monty(&t3, &t3, &out->Z); EG(ret, err); 927 ret = fp_add_monty(&out->Z, &t0, &t0); EG(ret, err); 928 929 ret = fp_add_monty(&t0, &out->Z, &t0); EG(ret, err); 930 ret = fp_add_monty(&t0, &t0, &t2); EG(ret, err); 931 ret = fp_mul_monty(&t0, &t0, &t3); EG(ret, err); 932 ret = fp_add_monty(&out->Y, &out->Y, &t0); EG(ret, err); 933 ret = fp_mul_monty(&t2, &in->Y, &in->Z); EG(ret, err); 934 935 ret = fp_add_monty(&t2, &t2, &t2); EG(ret, err); 936 ret = fp_mul_monty(&t0, &t2, &t3); EG(ret, err); 937 ret = fp_sub_monty(&out->X, &out->X, &t0); EG(ret, err); 938 ret = fp_mul_monty(&out->Z, &t2, &t1); EG(ret, err); 939 ret = fp_add_monty(&out->Z, &out->Z, &out->Z); EG(ret, err); 940 941 ret = fp_add_monty(&out->Z, &out->Z, &out->Z); 942 943 err: 944 fp_uninit(&t0); 945 fp_uninit(&t1); 946 fp_uninit(&t2); 947 fp_uninit(&t3); 948 949 return ret; 950 } 951 952 /* 953 * If NO_USE_COMPLETE_FORMULAS flag is not defined addition formulas from Algorithm 1 954 * of https://joostrenes.nl/publications/complete.pdf are used, otherwise 955 * http://www.hyperelliptic.org/EFD/g1p/auto-shortw-projective.html#addition-add-1998-cmo-2 956 */ 957 958 /* 959 * The function is an internal one: no check is performed on parameters, 960 * this MUST be done by the caller: 961 * 962 * - in1 and in2 are initialized 963 * - in1 and in2 are on the same curve 964 * - in1/in2 and out must not be aliased 965 * - in1 and in2 must not be an "exceptional" pair, i.e. (in1-in2) is not a point 966 * of order exactly 2 967 * 968 * The function will initialize 'out'. The function returns 0 on success, -1 969 * on error. 970 */ 971 ATTRIBUTE_WARN_UNUSED_RET static int __prj_pt_add_monty_cf(prj_pt_t out, 972 prj_pt_src_t in1, 973 prj_pt_src_t in2) 974 { 975 int cmp1, cmp2; 976 fp t0, t1, t2, t3, t4, t5; 977 int ret; 978 t0.magic = t1.magic = t2.magic = WORD(0); 979 t3.magic = t4.magic = t5.magic = WORD(0); 980 981 ret = prj_pt_init(out, in1->crv); EG(ret, err); 982 983 ret = fp_init(&t0, out->crv->a.ctx); EG(ret, err); 984 ret = fp_init(&t1, out->crv->a.ctx); EG(ret, err); 985 ret = fp_init(&t2, out->crv->a.ctx); EG(ret, err); 986 ret = fp_init(&t3, out->crv->a.ctx); EG(ret, err); 987 ret = fp_init(&t4, out->crv->a.ctx); EG(ret, err); 988 ret = fp_init(&t5, out->crv->a.ctx); EG(ret, err); 989 990 ret = fp_mul_monty(&t0, &in1->X, &in2->X); EG(ret, err); 991 ret = fp_mul_monty(&t1, &in1->Y, &in2->Y); EG(ret, err); 992 ret = fp_mul_monty(&t2, &in1->Z, &in2->Z); EG(ret, err); 993 ret = fp_add_monty(&t3, &in1->X, &in1->Y); EG(ret, err); 994 ret = fp_add_monty(&t4, &in2->X, &in2->Y); EG(ret, err); 995 996 ret = fp_mul_monty(&t3, &t3, &t4); EG(ret, err); 997 ret = fp_add_monty(&t4, &t0, &t1); EG(ret, err); 998 ret = fp_sub_monty(&t3, &t3, &t4); EG(ret, err); 999 ret = fp_add_monty(&t4, &in1->X, &in1->Z); EG(ret, err); 1000 ret = fp_add_monty(&t5, &in2->X, &in2->Z); EG(ret, err); 1001 1002 ret = fp_mul_monty(&t4, &t4, &t5); EG(ret, err); 1003 ret = fp_add_monty(&t5, &t0, &t2); EG(ret, err); 1004 ret = fp_sub_monty(&t4, &t4, &t5); EG(ret, err); 1005 ret = fp_add_monty(&t5, &in1->Y, &in1->Z); EG(ret, err); 1006 ret = fp_add_monty(&out->X, &in2->Y, &in2->Z); EG(ret, err); 1007 1008 ret = fp_mul_monty(&t5, &t5, &out->X); EG(ret, err); 1009 ret = fp_add_monty(&out->X, &t1, &t2); EG(ret, err); 1010 ret = fp_sub_monty(&t5, &t5, &out->X); EG(ret, err); 1011 ret = fp_mul_monty(&out->Z, &in1->crv->a_monty, &t4); EG(ret, err); 1012 ret = fp_mul_monty(&out->X, &in1->crv->b3_monty, &t2); EG(ret, err); 1013 1014 ret = fp_add_monty(&out->Z, &out->X, &out->Z); EG(ret, err); 1015 ret = fp_sub_monty(&out->X, &t1, &out->Z); EG(ret, err); 1016 ret = fp_add_monty(&out->Z, &t1, &out->Z); EG(ret, err); 1017 ret = fp_mul_monty(&out->Y, &out->X, &out->Z); EG(ret, err); 1018 ret = fp_add_monty(&t1, &t0, &t0); EG(ret, err); 1019 1020 ret = fp_add_monty(&t1, &t1, &t0); EG(ret, err); 1021 ret = fp_mul_monty(&t2, &in1->crv->a_monty, &t2); EG(ret, err); 1022 ret = fp_mul_monty(&t4, &in1->crv->b3_monty, &t4); EG(ret, err); 1023 ret = fp_add_monty(&t1, &t1, &t2); EG(ret, err); 1024 ret = fp_sub_monty(&t2, &t0, &t2); EG(ret, err); 1025 1026 ret = fp_mul_monty(&t2, &in1->crv->a_monty, &t2); EG(ret, err); 1027 ret = fp_add_monty(&t4, &t4, &t2); EG(ret, err); 1028 ret = fp_mul_monty(&t0, &t1, &t4); EG(ret, err); 1029 ret = fp_add_monty(&out->Y, &out->Y, &t0); EG(ret, err); 1030 ret = fp_mul_monty(&t0, &t5, &t4); EG(ret, err); 1031 1032 ret = fp_mul_monty(&out->X, &t3, &out->X); EG(ret, err); 1033 ret = fp_sub_monty(&out->X, &out->X, &t0); EG(ret, err); 1034 ret = fp_mul_monty(&t0, &t3, &t1); EG(ret, err); 1035 ret = fp_mul_monty(&out->Z, &t5, &out->Z); EG(ret, err); 1036 ret = fp_add_monty(&out->Z, &out->Z, &t0); 1037 1038 /* Check for "exceptional" pairs of input points with 1039 * checking if Y = Z = 0 as output (see the Bosma-Lenstra 1040 * article "Complete Systems of Two Addition Laws for 1041 * Elliptic Curves"). This should only happen on composite 1042 * order curves (i.e. not on prime order curves). 1043 * 1044 * In this case, we raise an error as the result is 1045 * not sound. This should not happen in our nominal 1046 * cases with regular signature and protocols, and if 1047 * it happens this usually means that bad points have 1048 * been injected. 1049 * 1050 * NOTE: if for some reasons you need to deal with 1051 * all the possible pairs of points including these 1052 * exceptional pairs of inputs with an order 2 difference, 1053 * you should fallback to the incomplete formulas using the 1054 * COMPLETE=0 compilation toggle. Beware that in this 1055 * case, the library will be more sensitive to 1056 * side-channel attacks. 1057 */ 1058 ret = fp_iszero(&(out->Z), &cmp1); EG(ret, err); 1059 ret = fp_iszero(&(out->Y), &cmp2); EG(ret, err); 1060 MUST_HAVE(!((cmp1) && (cmp2)), ret, err); 1061 1062 err: 1063 fp_uninit(&t0); 1064 fp_uninit(&t1); 1065 fp_uninit(&t2); 1066 fp_uninit(&t3); 1067 fp_uninit(&t4); 1068 fp_uninit(&t5); 1069 1070 return ret; 1071 } 1072 #endif /* NO_USE_COMPLETE_FORMULAS */ 1073 1074 /* 1075 * Internal function: 1076 * 1077 * - not supporting aliasing, 1078 * - requiring caller to check in parameter is initialized 1079 * 1080 * Based on library configuration, the function either use complete formulas 1081 * or not. 1082 */ 1083 static int _prj_pt_dbl_monty(prj_pt_t out, prj_pt_src_t in) 1084 { 1085 int ret; 1086 1087 #ifdef NO_USE_COMPLETE_FORMULAS 1088 int iszero; 1089 ret = prj_pt_iszero(in, &iszero); EG(ret, err); 1090 if (iszero) { 1091 ret = prj_pt_init(out, in->crv); EG(ret, err); 1092 ret = prj_pt_zero(out); 1093 } else { 1094 ret = __prj_pt_dbl_monty_no_cf(out, in); 1095 } 1096 #else 1097 ret = __prj_pt_dbl_monty_cf(out, in); EG(ret, err); 1098 #endif 1099 1100 err: 1101 return ret; 1102 } 1103 1104 /* 1105 * Internal version that peform in place doubling of given val, 1106 * by using a temporary copy. Sanity checks on parameters must 1107 * be done by caller. 1108 */ 1109 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_dbl_monty_aliased(prj_pt_t val) 1110 { 1111 prj_pt out_cpy; 1112 int ret; 1113 out_cpy.magic = WORD(0); 1114 1115 ret = _prj_pt_dbl_monty(&out_cpy, val); EG(ret, err); 1116 ret = prj_pt_copy(val, &out_cpy); 1117 1118 err: 1119 prj_pt_uninit(&out_cpy); 1120 1121 return ret; 1122 } 1123 1124 /* 1125 * Public function for projective point doubling. The function handles the init 1126 * check of 'in' parameter which must be guaranteed for internal functions. 1127 * 'out' parameter need not be initialized and can be aliased with 'in' 1128 * parameter. 1129 * 1130 * The function returns 0 on success, -1 on error. 1131 */ 1132 ATTRIBUTE_WARN_UNUSED_RET int prj_pt_dbl(prj_pt_t out, prj_pt_src_t in) 1133 { 1134 int ret; 1135 1136 ret = prj_pt_check_initialized(in); EG(ret, err); 1137 1138 if (out == in) { 1139 ret = _prj_pt_dbl_monty_aliased(out); 1140 } else { 1141 ret = _prj_pt_dbl_monty(out, in); 1142 } 1143 1144 err: 1145 return ret; 1146 } 1147 1148 /* 1149 * Internal function: 1150 * 1151 * - not supporting aliasing, 1152 * - requiring caller to check in1 and in2 parameter 1153 * 1154 * Based on library configuration, the function either use complete formulas 1155 * or not. 1156 */ 1157 ATTRIBUTE_WARN_UNUSED_RET static inline int _prj_pt_add_monty(prj_pt_t out, 1158 prj_pt_src_t in1, 1159 prj_pt_src_t in2) 1160 { 1161 #ifndef NO_USE_COMPLETE_FORMULAS 1162 return __prj_pt_add_monty_cf(out, in1, in2); 1163 #else 1164 return __prj_pt_add_monty_no_cf(out, in1, in2); 1165 #endif 1166 } 1167 1168 /* 1169 * The function is an internal one that specifically handles aliasing. No check 1170 * is performed on parameters, this MUST be done by the caller: 1171 * 1172 * - in1 and in2 are initialized 1173 * - in1 and in2 are on the same curve 1174 * 1175 * The function will initialize 'out'. The function returns 0 on success, -1 1176 * on error. 1177 */ 1178 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_add_monty_aliased(prj_pt_t out, 1179 prj_pt_src_t in1, 1180 prj_pt_src_t in2) 1181 { 1182 int ret; 1183 prj_pt out_cpy; 1184 out_cpy.magic = WORD(0); 1185 1186 ret = _prj_pt_add_monty(&out_cpy, in1, in2); EG(ret, err); 1187 ret = prj_pt_copy(out, &out_cpy); EG(ret, err); 1188 1189 err: 1190 prj_pt_uninit(&out_cpy); 1191 1192 return ret; 1193 } 1194 1195 /* 1196 * Public function for projective point addition. The function handles the 1197 * init checks of 'in1' and 'in2' parameters, along with the check they 1198 * use the same curve. This must be guaranteed for internal functions. 1199 * 'out' parameter need not be initialized and can be aliased with either 1200 * 'in1' or 'in2' parameter. 1201 * 1202 * The function returns 0 on success, -1 on error. 1203 */ 1204 int prj_pt_add(prj_pt_t out, prj_pt_src_t in1, prj_pt_src_t in2) 1205 { 1206 int ret; 1207 1208 ret = prj_pt_check_initialized(in1); EG(ret, err); 1209 ret = prj_pt_check_initialized(in2); EG(ret, err); 1210 MUST_HAVE((in1->crv == in2->crv), ret, err); 1211 1212 if ((out == in1) || (out == in2)) { 1213 ret = _prj_pt_add_monty_aliased(out, in1, in2); 1214 } else { 1215 ret = _prj_pt_add_monty(out, in1, in2); 1216 } 1217 1218 err: 1219 return ret; 1220 } 1221 1222 /*******************************************************************************/ 1223 /****** Scalar multiplication algorithms ***************************************/ 1224 /***********/ 1225 /* 1226 * The description below summarizes the following algorithms. 1227 * 1228 * Double-and-Add-Always and Montgomery Ladder masked using Itoh et al. anti-ADPA 1229 * (Address-bit DPA) countermeasure. 1230 * See "A Practical Countermeasure against Address-Bit Differential Power Analysis" 1231 * by Itoh, Izu and Takenaka for more information. 1232 * 1233 * NOTE: these masked variants of the Double-and-Add-Always and Montgomery Ladder algorithms 1234 * are used by default as Itoh et al. countermeasure has a very small impact on performance 1235 * and is inherently more robust againt DPA. The only case where we use another variant is 1236 * for devices with low memory as Itoh requires many temporary variables that consume many 1237 * temporary stack space. 1238 * 1239 * NOTE: the algorithms inherently depend on the MSB of the 1240 * scalar. In order to avoid leaking this MSB and fall into HNP (Hidden Number 1241 * Problem) issues, we use the trick described in https://eprint.iacr.org/2011/232.pdf 1242 * to have the MSB always set. However, since the scalar m might be less or bigger than 1243 * the order q of the curve, we distinguish three situations: 1244 * - The scalar m is < q (the order), in this case we compute: 1245 * - 1246 * | m' = m + (2 * q) if [log(k + q)] == [log(q)], 1247 * | m' = m + q otherwise. 1248 * - 1249 * - The scalar m is >= q and < q**2, in this case we compute: 1250 * - 1251 * | m' = m + (2 * (q**2)) if [log(k + (q**2))] == [log(q**2)], 1252 * | m' = m + (q**2) otherwise. 1253 * - 1254 * - The scalar m is >= (q**2), in this case m == m' 1255 * 1256 * => We only deal with 0 <= m < (q**2) using the countermeasure. When m >= (q**2), 1257 * we stick with m' = m, accepting MSB issues (not much can be done in this case 1258 * anyways). In the two first cases, Double-and-Add-Always and Montgomery Ladder are 1259 * performed in constant time wrt the size of the scalar m. 1260 */ 1261 /***********/ 1262 /* 1263 * Internal point blinding function: as it is internal, in is supposed to be initialized and 1264 * aliasing is NOT supported. 1265 */ 1266 ATTRIBUTE_WARN_UNUSED_RET static int _blind_projective_point(prj_pt_t out, prj_pt_src_t in) 1267 { 1268 int ret; 1269 1270 /* Random for projective coordinates masking */ 1271 /* NOTE: to limit stack usage, we reuse out->Z as a temporary 1272 * variable. This does not work if in == out, hence the check. 1273 */ 1274 MUST_HAVE((in != out), ret, err); 1275 1276 ret = prj_pt_init(out, in->crv); EG(ret, err); 1277 1278 /* Get a random value l in Fp */ 1279 ret = fp_get_random(&(out->Z), in->X.ctx); EG(ret, err); 1280 1281 /* 1282 * Blind the point with projective coordinates 1283 * (X, Y, Z) => (l*X, l*Y, l*Z) 1284 */ 1285 ret = fp_mul_monty(&(out->X), &(in->X), &(out->Z)); EG(ret, err); 1286 ret = fp_mul_monty(&(out->Y), &(in->Y), &(out->Z)); EG(ret, err); 1287 ret = fp_mul_monty(&(out->Z), &(in->Z), &(out->Z)); 1288 1289 err: 1290 return ret; 1291 } 1292 1293 /* If nothing is specified regarding the scalar multiplication algorithm, we use 1294 * the Montgomery Ladder. For the specific case of small stack devices, we release 1295 * some pressure on the stack by explicitly using double and always WITHOUT the Itoh 1296 * et al. countermeasure against A-DPA as it is quite consuming. 1297 */ 1298 #if defined(USE_SMALL_STACK) && defined(USE_MONTY_LADDER) 1299 #error "Small stack is only compatible with USE_DOUBLE_ADD_ALWAYS while USE_MONTY_LADDER has been explicitly asked!" 1300 #endif 1301 1302 #if defined(USE_SMALL_STACK) 1303 #ifndef USE_DOUBLE_ADD_ALWAYS 1304 #define USE_DOUBLE_ADD_ALWAYS 1305 #endif 1306 #endif 1307 1308 #if !defined(USE_DOUBLE_ADD_ALWAYS) && !defined(USE_MONTY_LADDER) 1309 #define USE_MONTY_LADDER 1310 #endif 1311 1312 #if defined(USE_DOUBLE_ADD_ALWAYS) && defined(USE_MONTY_LADDER) 1313 #error "You can either choose USE_DOUBLE_ADD_ALWAYS or USE_MONTY_LADDER, not both!" 1314 #endif 1315 1316 #if defined(USE_DOUBLE_ADD_ALWAYS) && !defined(USE_SMALL_STACK) 1317 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_mul_ltr_monty_dbl_add_always(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1318 { 1319 /* We use Itoh et al. notations here for T and the random r */ 1320 prj_pt T[3]; 1321 bitcnt_t mlen; 1322 u8 mbit, rbit; 1323 /* Random for masking the Double and Add Always algorithm */ 1324 nn r; 1325 /* The new scalar we will use with MSB fixed to 1 (noted m' above). 1326 * This helps dealing with constant time. 1327 */ 1328 nn m_msb_fixed; 1329 nn_src_t curve_order; 1330 nn curve_order_square; 1331 int ret, ret_ops, cmp; 1332 r.magic = m_msb_fixed.magic = curve_order_square.magic = WORD(0); 1333 T[0].magic = T[1].magic = T[2].magic = WORD(0); 1334 1335 /* Compute m' from m depending on the rule described above */ 1336 curve_order = &(in->crv->order); 1337 /* First compute q**2 */ 1338 ret = nn_sqr(&curve_order_square, curve_order); EG(ret, err); 1339 /* Then compute m' depending on m size */ 1340 ret = nn_cmp(m, curve_order, &cmp); EG(ret, err); 1341 if (cmp < 0){ 1342 bitcnt_t msb_bit_len, order_bitlen; 1343 1344 /* Case where m < q */ 1345 ret = nn_add(&m_msb_fixed, m, curve_order); EG(ret, err); 1346 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err); 1347 ret = nn_bitlen(curve_order, &order_bitlen); EG(ret, err); 1348 ret = nn_cnd_add((msb_bit_len == order_bitlen), &m_msb_fixed, 1349 &m_msb_fixed, curve_order); EG(ret, err); 1350 } else { 1351 ret = nn_cmp(m, &curve_order_square, &cmp); EG(ret, err); 1352 if (cmp < 0) { 1353 bitcnt_t msb_bit_len, curve_order_square_bitlen; 1354 1355 /* Case where m >= q and m < (q**2) */ 1356 ret = nn_add(&m_msb_fixed, m, &curve_order_square); EG(ret, err); 1357 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err); 1358 ret = nn_bitlen(&curve_order_square, &curve_order_square_bitlen); EG(ret, err); 1359 ret = nn_cnd_add((msb_bit_len == curve_order_square_bitlen), 1360 &m_msb_fixed, &m_msb_fixed, &curve_order_square); EG(ret, err); 1361 } else { 1362 /* Case where m >= (q**2) */ 1363 ret = nn_copy(&m_msb_fixed, m); EG(ret, err); 1364 } 1365 } 1366 ret = nn_bitlen(&m_msb_fixed, &mlen); EG(ret, err); 1367 MUST_HAVE(mlen != 0, ret, err); 1368 mlen--; 1369 1370 /* Hide possible internal failures for double and add 1371 * operations and perform the operation in constant time. 1372 */ 1373 ret_ops = 0; 1374 1375 /* Get a random r with the same size of m_msb_fixed */ 1376 ret = nn_get_random_len(&r, m_msb_fixed.wlen * WORD_BYTES); EG(ret, err); 1377 1378 ret = nn_getbit(&r, mlen, &rbit); EG(ret, err); 1379 1380 /* Initialize points */ 1381 ret = prj_pt_init(&T[0], in->crv); EG(ret, err); 1382 ret = prj_pt_init(&T[1], in->crv); EG(ret, err); 1383 1384 /* 1385 * T[2] = R(P) 1386 * Blind the point with projective coordinates 1387 * (X, Y, Z) => (l*X, l*Y, l*Z) 1388 */ 1389 ret = _blind_projective_point(&T[2], in); EG(ret, err); 1390 1391 /* T[r[n-1]] = T[2] */ 1392 ret = prj_pt_copy(&T[rbit], &T[2]); EG(ret, err); 1393 1394 /* Main loop of Double and Add Always */ 1395 while (mlen > 0) { 1396 u8 rbit_next; 1397 --mlen; 1398 /* rbit is r[i+1], and rbit_next is r[i] */ 1399 ret = nn_getbit(&r, mlen, &rbit_next); EG(ret, err); 1400 1401 /* mbit is m[i] */ 1402 ret = nn_getbit(&m_msb_fixed, mlen, &mbit); EG(ret, err); 1403 1404 /* Double: T[r[i+1]] = ECDBL(T[r[i+1]]) */ 1405 #ifndef NO_USE_COMPLETE_FORMULAS 1406 /* 1407 * NOTE: in case of complete formulas, we use the 1408 * addition for doubling, incurring a small performance hit 1409 * for better SCA resistance. 1410 */ 1411 ret_ops |= prj_pt_add(&T[rbit], &T[rbit], &T[rbit]); 1412 #else 1413 ret_ops |= prj_pt_dbl(&T[rbit], &T[rbit]); 1414 #endif 1415 /* Add: T[1-r[i+1]] = ECADD(T[r[i+1]],T[2]) */ 1416 ret_ops |= prj_pt_add(&T[1-rbit], &T[rbit], &T[2]); 1417 1418 /* 1419 * T[r[i]] = T[d[i] ^ r[i+1]] 1420 * NOTE: we use the low level nn_copy function here to avoid 1421 * any possible leakage on operands with prj_pt_copy 1422 */ 1423 ret = nn_copy(&(T[rbit_next].X.fp_val), &(T[mbit ^ rbit].X.fp_val)); EG(ret, err); 1424 ret = nn_copy(&(T[rbit_next].Y.fp_val), &(T[mbit ^ rbit].Y.fp_val)); EG(ret, err); 1425 ret = nn_copy(&(T[rbit_next].Z.fp_val), &(T[mbit ^ rbit].Z.fp_val)); EG(ret, err); 1426 1427 /* Update rbit */ 1428 rbit = rbit_next; 1429 } 1430 /* Output: T[r[0]] */ 1431 ret = prj_pt_copy(out, &T[rbit]); EG(ret, err); 1432 1433 /* Take into consideration our double and add errors */ 1434 ret |= ret_ops; 1435 1436 err: 1437 prj_pt_uninit(&T[0]); 1438 prj_pt_uninit(&T[1]); 1439 prj_pt_uninit(&T[2]); 1440 nn_uninit(&r); 1441 nn_uninit(&m_msb_fixed); 1442 nn_uninit(&curve_order_square); 1443 1444 PTR_NULLIFY(curve_order); 1445 1446 return ret; 1447 } 1448 #endif 1449 1450 #if defined(USE_DOUBLE_ADD_ALWAYS) && defined(USE_SMALL_STACK) 1451 /* NOTE: in small stack case where we compile for low memory devices, we do not use Itoh et al. countermeasure 1452 * as it requires too much temporary space on the stack. 1453 */ 1454 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_mul_ltr_monty_dbl_add_always(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1455 { 1456 int ret, ret_ops; 1457 1458 /* Hide possible internal failures for double and add 1459 * operations and perform the operation in constant time. 1460 */ 1461 ret_ops = 0; 1462 1463 /* Blind the input point projective coordinates */ 1464 ret = _blind_projective_point(out, in); EG(ret, err); 1465 1466 /*******************/ 1467 { 1468 bitcnt_t mlen; 1469 u8 mbit; 1470 /* The new scalar we will use with MSB fixed to 1 (noted m' above). 1471 * This helps dealing with constant time. 1472 */ 1473 nn m_msb_fixed; 1474 nn_src_t curve_order; 1475 int cmp; 1476 m_msb_fixed.magic = WORD(0); 1477 1478 { 1479 nn curve_order_square; 1480 curve_order_square.magic = WORD(0); 1481 1482 /* Compute m' from m depending on the rule described above */ 1483 curve_order = &(in->crv->order); 1484 /* First compute q**2 */ 1485 ret = nn_sqr(&curve_order_square, curve_order); EG(ret, err1); 1486 /* Then compute m' depending on m size */ 1487 ret = nn_cmp(m, curve_order, &cmp); EG(ret, err1); 1488 if (cmp < 0){ 1489 bitcnt_t msb_bit_len, order_bitlen; 1490 1491 /* Case where m < q */ 1492 ret = nn_add(&m_msb_fixed, m, curve_order); EG(ret, err1); 1493 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err1); 1494 ret = nn_bitlen(curve_order, &order_bitlen); EG(ret, err1); 1495 ret = nn_cnd_add((msb_bit_len == order_bitlen), &m_msb_fixed, 1496 &m_msb_fixed, curve_order); EG(ret, err1); 1497 } else { 1498 ret = nn_cmp(m, &curve_order_square, &cmp); EG(ret, err1); 1499 if (cmp < 0) { 1500 bitcnt_t msb_bit_len, curve_order_square_bitlen; 1501 1502 /* Case where m >= q and m < (q**2) */ 1503 ret = nn_add(&m_msb_fixed, m, &curve_order_square); EG(ret, err1); 1504 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err1); 1505 ret = nn_bitlen(&curve_order_square, &curve_order_square_bitlen); EG(ret, err1); 1506 ret = nn_cnd_add((msb_bit_len == curve_order_square_bitlen), 1507 &m_msb_fixed, &m_msb_fixed, &curve_order_square); EG(ret, err1); 1508 } else { 1509 /* Case where m >= (q**2) */ 1510 ret = nn_copy(&m_msb_fixed, m); EG(ret, err1); 1511 } 1512 } 1513 err1: 1514 nn_uninit(&curve_order_square); EG(ret, err); 1515 } 1516 1517 ret = nn_bitlen(&m_msb_fixed, &mlen); EG(ret, err); 1518 MUST_HAVE((mlen != 0), ret, err); 1519 mlen--; 1520 1521 { 1522 prj_pt dbl; 1523 dbl.magic = WORD(0); 1524 1525 /* Initialize temporary point */ 1526 ret = prj_pt_init(&dbl, in->crv); EG(ret, err2); 1527 1528 /* Main loop of Double and Add Always */ 1529 while (mlen > 0) { 1530 --mlen; 1531 /* mbit is m[i] */ 1532 ret = nn_getbit(&m_msb_fixed, mlen, &mbit); EG(ret, err2); 1533 1534 #ifndef NO_USE_COMPLETE_FORMULAS 1535 /* 1536 * NOTE: in case of complete formulas, we use the 1537 * addition for doubling, incurring a small performance hit 1538 * for better SCA resistance. 1539 */ 1540 ret_ops |= prj_pt_add(&dbl, out, out); 1541 #else 1542 ret_ops |= prj_pt_dbl(&dbl, out); 1543 #endif 1544 ret_ops |= prj_pt_add(out, &dbl, in); 1545 /* Swap */ 1546 ret = nn_cnd_swap(!mbit, &(out->X.fp_val), &(dbl.X.fp_val)); EG(ret, err2); 1547 ret = nn_cnd_swap(!mbit, &(out->Y.fp_val), &(dbl.Y.fp_val)); EG(ret, err2); 1548 ret = nn_cnd_swap(!mbit, &(out->Z.fp_val), &(dbl.Z.fp_val)); EG(ret, err2); 1549 } 1550 err2: 1551 prj_pt_uninit(&dbl); EG(ret, err); 1552 } 1553 1554 err: 1555 nn_uninit(&m_msb_fixed); 1556 1557 PTR_NULLIFY(curve_order); 1558 } 1559 1560 /* Take into consideration our double and add errors */ 1561 ret |= ret_ops; 1562 1563 return ret; 1564 } 1565 #endif 1566 1567 1568 #ifdef USE_MONTY_LADDER 1569 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_mul_ltr_monty_ladder(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1570 { 1571 /* We use Itoh et al. notations here for T and the random r */ 1572 prj_pt T[3]; 1573 bitcnt_t mlen; 1574 u8 mbit, rbit; 1575 /* Random for masking the Montgomery Ladder algorithm */ 1576 nn r; 1577 /* The new scalar we will use with MSB fixed to 1 (noted m' above). 1578 * This helps dealing with constant time. 1579 */ 1580 nn m_msb_fixed; 1581 nn_src_t curve_order; 1582 nn curve_order_square; 1583 int ret, ret_ops, cmp; 1584 r.magic = m_msb_fixed.magic = curve_order_square.magic = WORD(0); 1585 T[0].magic = T[1].magic = T[2].magic = WORD(0); 1586 1587 /* Compute m' from m depending on the rule described above */ 1588 curve_order = &(in->crv->order); 1589 1590 /* First compute q**2 */ 1591 ret = nn_sqr(&curve_order_square, curve_order); EG(ret, err); 1592 1593 /* Then compute m' depending on m size */ 1594 ret = nn_cmp(m, curve_order, &cmp); EG(ret, err); 1595 if (cmp < 0) { 1596 bitcnt_t msb_bit_len, order_bitlen; 1597 1598 /* Case where m < q */ 1599 ret = nn_add(&m_msb_fixed, m, curve_order); EG(ret, err); 1600 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err); 1601 ret = nn_bitlen(curve_order, &order_bitlen); EG(ret, err); 1602 ret = nn_cnd_add((msb_bit_len == order_bitlen), &m_msb_fixed, 1603 &m_msb_fixed, curve_order); EG(ret, err); 1604 } else { 1605 ret = nn_cmp(m, &curve_order_square, &cmp); EG(ret, err); 1606 if (cmp < 0) { 1607 bitcnt_t msb_bit_len, curve_order_square_bitlen; 1608 1609 /* Case where m >= q and m < (q**2) */ 1610 ret = nn_add(&m_msb_fixed, m, &curve_order_square); EG(ret, err); 1611 ret = nn_bitlen(&m_msb_fixed, &msb_bit_len); EG(ret, err); 1612 ret = nn_bitlen(&curve_order_square, &curve_order_square_bitlen); EG(ret, err); 1613 ret = nn_cnd_add((msb_bit_len == curve_order_square_bitlen), 1614 &m_msb_fixed, &m_msb_fixed, &curve_order_square); EG(ret, err); 1615 } else { 1616 /* Case where m >= (q**2) */ 1617 ret = nn_copy(&m_msb_fixed, m); EG(ret, err); 1618 } 1619 } 1620 1621 ret = nn_bitlen(&m_msb_fixed, &mlen); EG(ret, err); 1622 MUST_HAVE((mlen != 0), ret, err); 1623 mlen--; 1624 1625 /* Hide possible internal failures for double and add 1626 * operations and perform the operation in constant time. 1627 */ 1628 ret_ops = 0; 1629 1630 /* Get a random r with the same size of m_msb_fixed */ 1631 ret = nn_get_random_len(&r, (u16)(m_msb_fixed.wlen * WORD_BYTES)); EG(ret, err); 1632 1633 ret = nn_getbit(&r, mlen, &rbit); EG(ret, err); 1634 1635 /* Initialize points */ 1636 ret = prj_pt_init(&T[0], in->crv); EG(ret, err); 1637 ret = prj_pt_init(&T[1], in->crv); EG(ret, err); 1638 ret = prj_pt_init(&T[2], in->crv); EG(ret, err); 1639 1640 /* Initialize T[r[n-1]] to input point */ 1641 /* 1642 * Blind the point with projective coordinates 1643 * (X, Y, Z) => (l*X, l*Y, l*Z) 1644 */ 1645 ret = _blind_projective_point(&T[rbit], in); EG(ret, err); 1646 1647 /* Initialize T[1-r[n-1]] with ECDBL(T[r[n-1]])) */ 1648 #ifndef NO_USE_COMPLETE_FORMULAS 1649 /* 1650 * NOTE: in case of complete formulas, we use the 1651 * addition for doubling, incurring a small performance hit 1652 * for better SCA resistance. 1653 */ 1654 ret_ops |= prj_pt_add(&T[1-rbit], &T[rbit], &T[rbit]); 1655 #else 1656 ret_ops |= prj_pt_dbl(&T[1-rbit], &T[rbit]); 1657 #endif 1658 1659 /* Main loop of the Montgomery Ladder */ 1660 while (mlen > 0) { 1661 u8 rbit_next; 1662 --mlen; 1663 /* rbit is r[i+1], and rbit_next is r[i] */ 1664 ret = nn_getbit(&r, mlen, &rbit_next); EG(ret, err); 1665 1666 /* mbit is m[i] */ 1667 ret = nn_getbit(&m_msb_fixed, mlen, &mbit); EG(ret, err); 1668 /* Double: T[2] = ECDBL(T[d[i] ^ r[i+1]]) */ 1669 1670 #ifndef NO_USE_COMPLETE_FORMULAS 1671 /* NOTE: in case of complete formulas, we use the 1672 * addition for doubling, incurring a small performance hit 1673 * for better SCA resistance. 1674 */ 1675 ret_ops |= prj_pt_add(&T[2], &T[mbit ^ rbit], &T[mbit ^ rbit]); 1676 #else 1677 ret_ops |= prj_pt_dbl(&T[2], &T[mbit ^ rbit]); 1678 #endif 1679 1680 /* Add: T[1] = ECADD(T[0],T[1]) */ 1681 ret_ops |= prj_pt_add(&T[1], &T[0], &T[1]); 1682 1683 /* T[0] = T[2-(d[i] ^ r[i])] */ 1684 /* 1685 * NOTE: we use the low level nn_copy function here to avoid 1686 * any possible leakage on operands with prj_pt_copy 1687 */ 1688 ret = nn_copy(&(T[0].X.fp_val), &(T[2-(mbit ^ rbit_next)].X.fp_val)); EG(ret, err); 1689 ret = nn_copy(&(T[0].Y.fp_val), &(T[2-(mbit ^ rbit_next)].Y.fp_val)); EG(ret, err); 1690 ret = nn_copy(&(T[0].Z.fp_val), &(T[2-(mbit ^ rbit_next)].Z.fp_val)); EG(ret, err); 1691 1692 /* T[1] = T[1+(d[i] ^ r[i])] */ 1693 /* NOTE: we use the low level nn_copy function here to avoid 1694 * any possible leakage on operands with prj_pt_copy 1695 */ 1696 ret = nn_copy(&(T[1].X.fp_val), &(T[1+(mbit ^ rbit_next)].X.fp_val)); EG(ret, err); 1697 ret = nn_copy(&(T[1].Y.fp_val), &(T[1+(mbit ^ rbit_next)].Y.fp_val)); EG(ret, err); 1698 ret = nn_copy(&(T[1].Z.fp_val), &(T[1+(mbit ^ rbit_next)].Z.fp_val)); EG(ret, err); 1699 1700 /* Update rbit */ 1701 rbit = rbit_next; 1702 } 1703 /* Output: T[r[0]] */ 1704 ret = prj_pt_copy(out, &T[rbit]); EG(ret, err); 1705 1706 /* Take into consideration our double and add errors */ 1707 ret |= ret_ops; 1708 1709 err: 1710 prj_pt_uninit(&T[0]); 1711 prj_pt_uninit(&T[1]); 1712 prj_pt_uninit(&T[2]); 1713 nn_uninit(&r); 1714 nn_uninit(&m_msb_fixed); 1715 nn_uninit(&curve_order_square); 1716 1717 PTR_NULLIFY(curve_order); 1718 1719 return ret; 1720 } 1721 #endif 1722 1723 /* Main projective scalar multiplication function. 1724 * Depending on the preprocessing options, we use either the 1725 * Double and Add Always algorithm, or the Montgomery Ladder one. 1726 */ 1727 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_mul_ltr_monty(prj_pt_t out, nn_src_t m, prj_pt_src_t in){ 1728 #if defined(USE_DOUBLE_ADD_ALWAYS) 1729 return _prj_pt_mul_ltr_monty_dbl_add_always(out, m, in); 1730 #elif defined(USE_MONTY_LADDER) 1731 return _prj_pt_mul_ltr_monty_ladder(out, m, in); 1732 #else 1733 #error "Error: neither Double and Add Always nor Montgomery Ladder has been selected!" 1734 #endif 1735 } 1736 1737 /* version with 'm' passed via 'out'. */ 1738 ATTRIBUTE_WARN_UNUSED_RET static int _prj_pt_mul_ltr_monty_aliased(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1739 { 1740 prj_pt out_cpy; 1741 int ret; 1742 out_cpy.magic = WORD(0); 1743 1744 ret = prj_pt_init(&out_cpy, in->crv); EG(ret, err); 1745 ret = _prj_pt_mul_ltr_monty(&out_cpy, m, in); EG(ret, err); 1746 ret = prj_pt_copy(out, &out_cpy); 1747 1748 err: 1749 prj_pt_uninit(&out_cpy); 1750 return ret; 1751 } 1752 1753 /* Aliased version. This is the public main interface of our 1754 * scalar multiplication algorithm. Checks that the input point 1755 * and that the output point are on the curve are performed here 1756 * (before and after calling the core algorithm, albeit Double and 1757 * Add Always or Montgomery Ladder). 1758 */ 1759 int prj_pt_mul(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1760 { 1761 int ret, on_curve; 1762 1763 ret = prj_pt_check_initialized(in); EG(ret, err); 1764 ret = nn_check_initialized(m); EG(ret, err); 1765 1766 /* Check that the input is on the curve */ 1767 MUST_HAVE((!prj_pt_is_on_curve(in, &on_curve)) && on_curve, ret, err); 1768 1769 if (out == in) { 1770 ret = _prj_pt_mul_ltr_monty_aliased(out, m, in); EG(ret, err); 1771 } else { 1772 ret = _prj_pt_mul_ltr_monty(out, m, in); EG(ret, err); 1773 } 1774 1775 /* Check that the output is on the curve */ 1776 MUST_HAVE((!prj_pt_is_on_curve(out, &on_curve)) && on_curve, ret, err); 1777 1778 err: 1779 return ret; 1780 } 1781 1782 int prj_pt_mul_blind(prj_pt_t out, nn_src_t m, prj_pt_src_t in) 1783 { 1784 /* Blind the scalar m with (b*q) where q is the curve order. 1785 * NOTE: the curve order and the "generator" order are 1786 * usually the same (i.e. cofactor = 1) for the classical 1787 * prime fields curves. However some exceptions exist 1788 * (e.g. Wei25519 and Wei448), and in this case it is 1789 * curcial to use the curve order for a generic blinding 1790 * working on any point on the curve. 1791 */ 1792 nn b; 1793 nn_src_t q; 1794 int ret; 1795 b.magic = WORD(0); 1796 1797 ret = prj_pt_check_initialized(in); EG(ret, err); 1798 1799 q = &(in->crv->order); 1800 1801 ret = nn_init(&b, 0); EG(ret, err); 1802 1803 ret = nn_get_random_mod(&b, q); EG(ret, err); 1804 1805 ret = nn_mul(&b, &b, q); EG(ret, err); 1806 ret = nn_add(&b, &b, m); EG(ret, err); 1807 1808 /* NOTE: point blinding is performed in the lower functions */ 1809 /* NOTE: check that input and output points are on the curve are 1810 * performed in the lower functions. 1811 */ 1812 1813 /* Perform the scalar multiplication */ 1814 ret = prj_pt_mul(out, &b, in); 1815 1816 err: 1817 nn_uninit(&b); 1818 1819 PTR_NULLIFY(q); 1820 1821 return ret; 1822 } 1823 1824 /* Naive double and add scalar multiplication. 1825 * 1826 * This scalar multiplication is used on public values and is optimized with no 1827 * countermeasures, and it is usually faster as scalar can be small with few bits 1828 * to process (e.g. cofactors, etc.). 1829 * 1830 * out is initialized by the function. 1831 * 1832 * XXX: WARNING: this function must only be used on public points! 1833 * 1834 */ 1835 static int __prj_pt_unprotected_mult(prj_pt_t out, nn_src_t scalar, prj_pt_src_t public_in) 1836 { 1837 u8 expbit; 1838 bitcnt_t explen; 1839 int ret, iszero, on_curve; 1840 1841 ret = prj_pt_check_initialized(public_in); EG(ret, err); 1842 ret = nn_check_initialized(scalar); EG(ret, err); 1843 1844 /* This function does not support aliasing */ 1845 MUST_HAVE((out != public_in), ret, err); 1846 1847 /* Check that the input is on the curve */ 1848 MUST_HAVE((!prj_pt_is_on_curve(public_in, &on_curve)) && on_curve, ret, err); 1849 1850 ret = nn_iszero(scalar, &iszero); EG(ret, err); 1851 /* Multiplication by zero is the point at infinity */ 1852 if(iszero){ 1853 ret = prj_pt_zero(out); EG(ret, err); 1854 goto err; 1855 } 1856 1857 ret = nn_bitlen(scalar, &explen); EG(ret, err); 1858 /* Sanity check */ 1859 MUST_HAVE((explen > 0), ret, err); 1860 explen = (bitcnt_t)(explen - 1); 1861 ret = prj_pt_copy(out, public_in); EG(ret, err); 1862 1863 while (explen > 0) { 1864 explen = (bitcnt_t)(explen - 1); 1865 ret = nn_getbit(scalar, explen, &expbit); EG(ret, err); 1866 ret = prj_pt_dbl(out, out); EG(ret, err); 1867 if(expbit){ 1868 ret = prj_pt_add(out, out, public_in); EG(ret, err); 1869 } 1870 } 1871 1872 /* Check that the output is on the curve */ 1873 MUST_HAVE((!prj_pt_is_on_curve(out, &on_curve)) && on_curve, ret, err); 1874 1875 err: 1876 VAR_ZEROIFY(expbit); 1877 VAR_ZEROIFY(explen); 1878 1879 return ret; 1880 } 1881 1882 /* Aliased version of __prj_pt_unprotected_mult */ 1883 int _prj_pt_unprotected_mult(prj_pt_t out, nn_src_t scalar, prj_pt_src_t public_in) 1884 { 1885 int ret; 1886 1887 if(out == public_in){ 1888 prj_pt A; 1889 A.magic = WORD(0); 1890 1891 ret = prj_pt_copy(&A, public_in); EG(ret, err1); 1892 ret = __prj_pt_unprotected_mult(out, scalar, &A); 1893 err1: 1894 prj_pt_uninit(&A); 1895 goto err; 1896 } 1897 else{ 1898 ret = __prj_pt_unprotected_mult(out, scalar, public_in); 1899 } 1900 err: 1901 return ret; 1902 } 1903 /* 1904 * Check if an integer is (a multiple of) a projective point order. 1905 * 1906 * The function returns 0 on success, -1 on error. The value check is set to 1 if the projective 1907 * point has order in_isorder, 0 otherwise. The value is meaningless on error. 1908 */ 1909 int check_prj_pt_order(prj_pt_src_t in_shortw, nn_src_t in_isorder, prj_pt_sensitivity s, int *check) 1910 { 1911 int ret, iszero; 1912 prj_pt res; 1913 res.magic = WORD(0); 1914 1915 /* First sanity checks */ 1916 ret = prj_pt_check_initialized(in_shortw); EG(ret, err); 1917 ret = nn_check_initialized(in_isorder); EG(ret, err); 1918 MUST_HAVE((check != NULL), ret, err); 1919 1920 /* Then, perform the scalar multiplication */ 1921 if(s == PUBLIC_PT){ 1922 /* If this is a public point, we can use the naive scalar multiplication */ 1923 ret = _prj_pt_unprotected_mult(&res, in_isorder, in_shortw); EG(ret, err); 1924 } 1925 else{ 1926 /* If the point is private, it is sensitive and we proceed with the secure 1927 * scalar blind multiplication. 1928 */ 1929 ret = prj_pt_mul_blind(&res, in_isorder, in_shortw); EG(ret, err); 1930 } 1931 1932 /* Check if we have the point at infinity */ 1933 ret = prj_pt_iszero(&res, &iszero); EG(ret, err); 1934 (*check) = iszero; 1935 1936 err: 1937 prj_pt_uninit(&res); 1938 1939 return ret; 1940 } 1941 1942 /*****************************************************************************/ 1943 1944 /* 1945 * Map points from Edwards to short Weierstrass projective points through Montgomery (composition mapping). 1946 * Point at infinity (0, 1) -> (0, 1, 0) is treated as an exception, which is trivially not constant time. 1947 * This is OK since our mapping functions should be used at the non sensitive input and output 1948 * interfaces. 1949 * 1950 * The function returns 0 on success, -1 on error. 1951 */ 1952 int aff_pt_edwards_to_prj_pt_shortw(aff_pt_edwards_src_t in_edwards, 1953 ec_shortw_crv_src_t shortw_crv, 1954 prj_pt_t out_shortw, 1955 fp_src_t alpha_edwards) 1956 { 1957 int ret, iszero, cmp; 1958 aff_pt out_shortw_aff; 1959 fp one; 1960 out_shortw_aff.magic = one.magic = WORD(0); 1961 1962 /* Check the curves compatibility */ 1963 ret = aff_pt_edwards_check_initialized(in_edwards); EG(ret, err); 1964 ret = curve_edwards_shortw_check(in_edwards->crv, shortw_crv, alpha_edwards); EG(ret, err); 1965 1966 /* Initialize output point with curve */ 1967 ret = prj_pt_init(out_shortw, shortw_crv); EG(ret, err); 1968 1969 ret = fp_init(&one, in_edwards->x.ctx); EG(ret, err); 1970 ret = fp_one(&one); EG(ret, err); 1971 1972 /* Check if we are the point at infinity 1973 * This check induces a non consant time exception, but the current function must be called on 1974 * public data anyways. 1975 */ 1976 ret = fp_iszero(&(in_edwards->x), &iszero); EG(ret, err); 1977 ret = fp_cmp(&(in_edwards->y), &one, &cmp); EG(ret, err); 1978 if(iszero && (cmp == 0)){ 1979 ret = prj_pt_zero(out_shortw); EG(ret, err); 1980 ret = 0; 1981 goto err; 1982 } 1983 1984 /* Use the affine mapping */ 1985 ret = aff_pt_edwards_to_shortw(in_edwards, shortw_crv, &out_shortw_aff, alpha_edwards); EG(ret, err); 1986 /* And then map the short Weierstrass affine to projective coordinates */ 1987 ret = ec_shortw_aff_to_prj(out_shortw, &out_shortw_aff); 1988 1989 err: 1990 fp_uninit(&one); 1991 aff_pt_uninit(&out_shortw_aff); 1992 1993 return ret; 1994 } 1995 1996 /* 1997 * Map points from short Weierstrass projective points to Edwards through Montgomery (composition mapping). 1998 * Point at infinity with Z=0 (in projective coordinates) -> (0, 1) is treated as an exception, which is trivially not constant time. 1999 * This is OK since our mapping functions should be used at the non sensitive input and output 2000 * interfaces. 2001 * 2002 * The function returns 0 on success, -1 on error. 2003 */ 2004 int prj_pt_shortw_to_aff_pt_edwards(prj_pt_src_t in_shortw, 2005 ec_edwards_crv_src_t edwards_crv, 2006 aff_pt_edwards_t out_edwards, 2007 fp_src_t alpha_edwards) 2008 { 2009 int ret, iszero; 2010 aff_pt in_shortw_aff; 2011 in_shortw_aff.magic = WORD(0); 2012 2013 /* Check the curves compatibility */ 2014 ret = prj_pt_check_initialized(in_shortw); EG(ret, err); 2015 ret = curve_edwards_shortw_check(edwards_crv, in_shortw->crv, alpha_edwards); EG(ret, err); 2016 2017 /* Initialize output point with curve */ 2018 ret = aff_pt_init(&in_shortw_aff, in_shortw->crv); EG(ret, err); 2019 2020 /* Check if we are the point at infinity 2021 * This check induces a non consant time exception, but the current function must be called on 2022 * public data anyways. 2023 */ 2024 ret = prj_pt_iszero(in_shortw, &iszero); EG(ret, err); 2025 if(iszero){ 2026 fp zero, one; 2027 zero.magic = one.magic = WORD(0); 2028 2029 ret = fp_init(&zero, in_shortw->X.ctx); EG(ret, err1); 2030 ret = fp_init(&one, in_shortw->X.ctx); EG(ret, err1); 2031 2032 ret = fp_zero(&zero); EG(ret, err1); 2033 ret = fp_one(&one); EG(ret, err1); 2034 2035 ret = aff_pt_edwards_init_from_coords(out_edwards, edwards_crv, &zero, &one); 2036 2037 err1: 2038 fp_uninit(&zero); 2039 fp_uninit(&one); 2040 2041 goto err; 2042 } 2043 2044 /* Map projective to affine on the short Weierstrass */ 2045 ret = prj_pt_to_aff(&in_shortw_aff, in_shortw); EG(ret, err); 2046 /* Use the affine mapping */ 2047 ret = aff_pt_shortw_to_edwards(&in_shortw_aff, edwards_crv, out_edwards, alpha_edwards); 2048 2049 err: 2050 aff_pt_uninit(&in_shortw_aff); 2051 2052 return ret; 2053 } 2054 2055 /* 2056 * Map points from Montgomery to short Weierstrass projective points. 2057 * 2058 * The function returns 0 on success, -1 on error. 2059 */ 2060 int aff_pt_montgomery_to_prj_pt_shortw(aff_pt_montgomery_src_t in_montgomery, 2061 ec_shortw_crv_src_t shortw_crv, 2062 prj_pt_t out_shortw) 2063 { 2064 int ret; 2065 aff_pt out_shortw_aff; 2066 out_shortw_aff.magic = WORD(0); 2067 2068 /* Check the curves compatibility */ 2069 ret = aff_pt_montgomery_check_initialized(in_montgomery); EG(ret, err); 2070 ret = curve_montgomery_shortw_check(in_montgomery->crv, shortw_crv); EG(ret, err); 2071 2072 /* Initialize output point with curve */ 2073 ret = prj_pt_init(out_shortw, shortw_crv); EG(ret, err); 2074 2075 /* Use the affine mapping */ 2076 ret = aff_pt_montgomery_to_shortw(in_montgomery, shortw_crv, &out_shortw_aff); EG(ret, err); 2077 /* And then map the short Weierstrass affine to projective coordinates */ 2078 ret = ec_shortw_aff_to_prj(out_shortw, &out_shortw_aff); 2079 2080 err: 2081 aff_pt_uninit(&out_shortw_aff); 2082 2083 return ret; 2084 } 2085 2086 /* 2087 * Map points from short Weierstrass projective points to Montgomery. 2088 * 2089 * The function returns 0 on success, -1 on error. 2090 */ 2091 int prj_pt_shortw_to_aff_pt_montgomery(prj_pt_src_t in_shortw, ec_montgomery_crv_src_t montgomery_crv, aff_pt_montgomery_t out_montgomery) 2092 { 2093 int ret; 2094 aff_pt in_shortw_aff; 2095 in_shortw_aff.magic = WORD(0); 2096 2097 /* Check the curves compatibility */ 2098 ret = prj_pt_check_initialized(in_shortw); EG(ret, err); 2099 ret = curve_montgomery_shortw_check(montgomery_crv, in_shortw->crv); EG(ret, err); 2100 2101 /* Initialize output point with curve */ 2102 ret = aff_pt_init(&in_shortw_aff, in_shortw->crv); EG(ret, err); 2103 2104 /* Map projective to affine on the short Weierstrass */ 2105 ret = prj_pt_to_aff(&in_shortw_aff, in_shortw); EG(ret, err); 2106 /* Use the affine mapping */ 2107 ret = aff_pt_shortw_to_montgomery(&in_shortw_aff, montgomery_crv, out_montgomery); 2108 2109 err: 2110 aff_pt_uninit(&in_shortw_aff); 2111 2112 return ret; 2113 }