ecdsa/ecc_chudnovsky.c

   1 /* ecc_chudnovsky.c */
   2 /*
   3     This file is part of the ARM-Crypto-Lib.
   4     Copyright (C) 2006-2015 Daniel Otte (bg@nerilex.org)
   5
   6     This program is free software: you can redistribute it and/or modify
   7     it under the terms of the GNU General Public License as published by
   8     the Free Software Foundation, either version 3 of the License, or
   9     (at your option) any later version.
  10
  11     This program is distributed in the hope that it will be useful,
  12     but WITHOUT ANY WARRANTY; without even the implied warranty of
  13     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14     GNU General Public License for more details.
  15
  16     You should have received a copy of the GNU General Public License
  17     along with this program.  If not, see <http://www.gnu.org/licenses/>.
  18 */
  19
  20 #include <stdlib.h>
  21 #include <string.h>
  22 #include "bigint.h"
  23 #include "ecc.h"
  24 #include <stdlib.h>
  25 #include <string.h>
  26
  27 #if 1
  28 #include <stdio.h>
  29 #include "bigint_io.h"
  30 #include "uart_i.h"
  31 #include <avr/pgmspace.h>
  32 #endif
  33
  34 #if 0
  35 #define printf_P(...)
  36 #define bigint_print_hex(a)
  37 #undef putchar
  38 #define putchar(a)
  39 #endif
  40
  41 uint8_t ecc_chudnovsky_point_alloc(ecc_chudnovsky_point_t *p, size_t length_B){
  42     if(!(p->x.wordv = malloc(length_B))){
  43         printf_P(PSTR("DBG: XXX <%S %s %d>\n"), PSTR(__FILE__), __func__, __LINE__);
  44         return 1;
  45     }
  46     if(!(p->y.wordv = malloc(length_B))){
  47         printf_P(PSTR("DBG: XXX <%S %s %d>\n"), PSTR(__FILE__), __func__, __LINE__);
  48         free(p->x.wordv);
  49         return 1;
  50     }
  51     if(!(p->z1.wordv = malloc(length_B))){
  52         printf_P(PSTR("DBG: XXX <%S %s %d>\n"), PSTR(__FILE__), __func__, __LINE__);
  53         free(p->x.wordv);
  54         free(p->y.wordv);
  55         return 1;
  56     }
  57     if(!(p->z2.wordv = malloc(length_B))){
  58         printf_P(PSTR("DBG: XXX <%S %s %d>\n"), PSTR(__FILE__), __func__, __LINE__);
  59         free(p->x.wordv);
  60         free(p->y.wordv);
  61         free(p->z1.wordv);
  62         return 1;
  63     }
  64     if(!(p->z3.wordv = malloc(length_B))){
  65         printf_P(PSTR("DBG: XXX <%S %s %d>\n"), PSTR(__FILE__), __func__, __LINE__);
  66         free(p->x.wordv);
  67         free(p->y.wordv);
  68         free(p->z1.wordv);
  69         free(p->z2.wordv);
  70         return 1;
  71     }
  72     bigint_set_zero(&p->x);
  73     bigint_set_zero(&p->y);
  74     bigint_set_zero(&p->z1);
  75     bigint_set_zero(&p->z2);
  76     bigint_set_zero(&p->z3);
  77     return 0;
  78 }
  79
  80 void ecc_chudnovsky_point_free(ecc_chudnovsky_point_t *p){
  81     free(p->x.wordv);
  82     free(p->y.wordv);
  83     free(p->z1.wordv);
  84     free(p->z2.wordv);
  85     free(p->z3.wordv);
  86 }
  87
  88 /*
  89  * if (Y == 0)
  90  *   return POINT_AT_INFINITY
  91  * S = 4*X*Y^2
  92  * M = 3*(X + Z^2)*(X - Z^2)
  93  * X' = M^2 - 2*S
  94  * Y' = M*(S - X') - 8*Y^4
  95  * Z' = 2*Y*Z
  96  * Z'^2 = Z'^2
  97  * Z'^3 = Z'^2 * Z'
  98  * return (X', Y', Z', Z'^2, Z'^3)
  99  */
 100
 101 uint8_t ecc_affine_to_chudnovsky_point(ecc_chudnovsky_point_t *dest,
 102                                        const ecc_affine_point_t *src){
 103     if(src->y.length_W == 0){
 104         /* point at infinity */
 105         bigint_set_zero(&dest->y);
 106         return 0;
 107     }
 108     bigint_copy(&dest->x, &src->x);
 109     bigint_copy(&dest->y, &src->y);
 110     dest->z1.wordv[0] = 1;
 111     dest->z2.wordv[0] = 1;
 112     dest->z3.wordv[0] = 1;
 113     dest->z1.length_W = 1;
 114     dest->z2.length_W = 1;
 115     dest->z3.length_W = 1;
 116     dest->z1.info = 0;
 117     dest->z2.info = 0;
 118     dest->z3.info = 0;
 119     bigint_adjust(&dest->z1);
 120     bigint_adjust(&dest->z2);
 121     bigint_adjust(&dest->z3);
 122     return 0;
 123 }
 124
 125 uint8_t ecc_chudnovsky_to_affine_point(ecc_affine_point_t *dest,
 126                                        const ecc_chudnovsky_point_t *src,
 127                                        const ecc_curve_sp_t *curve){
 128     if(src->y.length_W == 0){
 129         /* point at infinity */
 130         bigint_set_zero(&dest->y);
 131         return 0;
 132     }
 133     bigint_word_t t_w[curve->p->length_W * 2];
 134     bigint_word_t z1_w[curve->p->length_W * 2];
 135     bigint_word_t z2_w[curve->p->length_W];
 136     bigint_t t, z1, z2;
 137     t.wordv =  t_w;
 138     z1.wordv = z1_w;
 139     z2.wordv = z2_w;
 140
 141     bigint_inverse(&z1, &src->z1, curve->p);
 142     bigint_square(&t, &z1);
 143     curve->reduce_p(&t);
 144     bigint_copy(&z2, &t);
 145     bigint_mul_u(&t, &src->x, &z2);
 146     curve->reduce_p(&t);
 147     bigint_copy(&dest->x, &t);
 148     bigint_mul_u(&t, &z1, &z2);
 149     curve->reduce_p(&t);
 150     bigint_mul_u(&t, &t, &src->y);
 151     curve->reduce_p(&t);
 152     bigint_copy(&dest->y, &t);
 153
 154     return 0;
 155 }
 156
 157 void ecc_chudnovsky_point_print(const ecc_chudnovsky_point_t *p){
 158     if(p->y.length_W == 0){
 159         printf_P(PSTR(" ECC point = point-at-infinity\n"));
 160         return;
 161     }
 162     printf_P(PSTR(" ECC point x  = "));
 163     bigint_print_hex(&p->x);
 164     printf_P(PSTR("\n ECC point y  = "));
 165     bigint_print_hex(&p->y);
 166     printf_P(PSTR("\n ECC point z1 = "));
 167     bigint_print_hex(&p->z1);
 168     printf_P(PSTR("\n ECC point z2 = "));
 169     bigint_print_hex(&p->z2);
 170     printf_P(PSTR("\n ECC point z3 = "));
 171     bigint_print_hex(&p->z3);
 172     putchar('\n');
 173 }
 174
 175 uint8_t ecc_chudnovsky_point_double_sp(ecc_chudnovsky_point_t *dest,
 176                                        const ecc_chudnovsky_point_t *a,
 177                                        const ecc_curve_sp_t *curve){
 178     if(a->y.length_W == 0){
 179         /* point at infinity */
 180         bigint_set_zero(&dest->y);
 181         return 0;
 182     }
 183     bigint_word_t s_w[curve->p->length_W * 2];
 184     bigint_word_t m_w[curve->p->length_W * 2];
 185     bigint_word_t t_w[curve->p->length_W * 2];
 186     bigint_t s, m, t;
 187
 188     s.wordv = s_w;
 189     m.wordv = m_w;
 190     t.wordv = t_w;
 191
 192     /* compute s*/
 193     bigint_square(&t, &a->y);
 194     curve->reduce_p(&t);
 195     bigint_mul_u(&s, &t, &a->x);
 196     curve->reduce_p(&s);
 197     bigint_shiftleft(&s, 2);
 198     curve->reduce_p(&s);
 199
 200     /* compute m */
 201     bigint_sub_u(&t, &a->x, &a->z2);
 202 //    /**/curve->reduce_p(&t);
 203     bigint_add_u(&m, &a->x, &a->z2);
 204 //    /**/curve->reduce_p(&m);
 205     bigint_mul_s(&m, &m, &t);
 206     curve->reduce_p(&m);
 207     bigint_copy(&t, &m);
 208     bigint_shiftleft(&t, 1);
 209     bigint_add_s(&m, &m, &t);
 210     curve->reduce_p(&m);
 211
 212     /* compute new z1 */
 213     bigint_mul_u(&t, &a->z1, &a->y);
 214     curve->reduce_p(&t);
 215     bigint_shiftleft(&t, 1);
 216     curve->reduce_p(&t);
 217     bigint_copy(&dest->z1, &t);
 218
 219     /* compute new x */
 220     bigint_square(&t, &m);
 221     curve->reduce_p(&t);
 222     bigint_sub_s(&t, &t, &s);
 223     bigint_sub_s(&t, &t, &s);
 224     curve->reduce_p(&t);
 225     bigint_copy(&dest->x, &t);
 226
 227     /* compute new y */
 228     bigint_sub_s(&s, &s, &t);
 229     curve->reduce_p(&s);
 230     bigint_mul_s(&s, &s, &m);
 231     curve->reduce_p(&s);
 232     bigint_square(&t, &a->y);
 233     curve->reduce_p(&t);
 234     bigint_square(&t, &t);
 235     curve->reduce_p(&t);
 236     bigint_shiftleft(&t, 3);
 237     curve->reduce_p(&t);
 238     bigint_sub_s(&s, &s, &t);
 239     curve->reduce_p(&s);
 240     bigint_copy(&dest->y, &s);
 241
 242     /* compute new z2 */
 243     bigint_square(&t, &dest->z1);
 244     curve->reduce_p(&t);
 245     bigint_copy(&dest->z2, &t);
 246
 247     /* compute new z3 */
 248     bigint_mul_u(&t, &t, &dest->z1);
 249     curve->reduce_p(&t);
 250     bigint_copy(&dest->z3, &t);
 251
 252     return 0;
 253 }
 254
 255 void ecc_chudnovsky_point_copy(ecc_chudnovsky_point_t *dest,
 256                                const ecc_chudnovsky_point_t *src){
 257     bigint_copy(&dest->x, &src->x);
 258     bigint_copy(&dest->y, &src->y);
 259     bigint_copy(&dest->z1, &src->z1);
 260     bigint_copy(&dest->z2, &src->z2);
 261     bigint_copy(&dest->z3, &src->z3);
 262 }
 263
 264 uint8_t ecc_chudnovsky_point_add_sp(ecc_chudnovsky_point_t *dest,
 265                                     const ecc_chudnovsky_point_t *a,
 266                                     const ecc_chudnovsky_point_t *b,
 267                                     const ecc_curve_sp_t *curve){
 268     if(a->y.length_W == 0){
 269         ecc_chudnovsky_point_copy(dest, b);
 270         return 0;
 271     }
 272     if(b->y.length_W == 0){
 273         ecc_chudnovsky_point_copy(dest, a);
 274         return 0;
 275     }
 276     bigint_word_t u1_w[curve->p->length_W * 2];
 277     bigint_word_t u2_w[curve->p->length_W * 2];
 278     bigint_word_t s1_w[curve->p->length_W * 2];
 279     bigint_word_t s2_w[curve->p->length_W * 2];
 280     bigint_t u1, u2, s1, s2;
 281
 282     u1.wordv = u1_w;
 283     u2.wordv = u2_w;
 284     s1.wordv = s1_w;
 285     s2.wordv = s2_w;
 286
 287     /* compute u1 */
 288     bigint_mul_u(&u1, &a->x, &b->z2);
 289     curve->reduce_p(&u1);
 290
 291     /* compute u2 */
 292     bigint_mul_u(&u2, &b->x, &a->z2);
 293     curve->reduce_p(&u2);
 294
 295     /* compute s1 */
 296     bigint_mul_u(&s1, &a->y, &b->z3);
 297     curve->reduce_p(&s1);
 298
 299     /* compute s2 */
 300     bigint_mul_u(&s2, &b->y, &a->z3);
 301     curve->reduce_p(&s2);
 302
 303     if(bigint_cmp_u(&u1, &u2) == 0){
 304         if(bigint_cmp_u(&s1, &s2)){
 305             /* point at infinity */
 306             bigint_set_zero(&dest->y);
 307             return 0;
 308         }else{
 309             /* a == b --> dest = 2*a */
 310             ecc_chudnovsky_point_double_sp(dest, a, curve);
 311         }
 312     }
 313
 314     bigint_word_t h_w[curve->p->length_W * 2];
 315     bigint_word_t r_w[curve->p->length_W * 2];
 316     bigint_t h, r;
 317
 318     h.wordv = h_w;
 319     r.wordv = r_w;
 320     /* compute h */
 321     bigint_sub_s(&h, &u2, &u1);
 322     /**/curve->reduce_p(&h);
 323
 324     /* compute r */
 325     bigint_sub_s(&r, &s2, &s1);
 326 //    /**/curve->reduce_p(&r);
 327
 328     /* compute new z */
 329     bigint_mul_u(&s2, &a->z1, &b->z1);
 330     curve->reduce_p(&s2);
 331     bigint_mul_s(&s2, &s2, &h);
 332     curve->reduce_p(&s2);
 333     bigint_copy(&dest->z1, &s2);
 334
 335     /* compute u1*h^2 and h^3 */
 336     bigint_square(&s2, &h);
 337     curve->reduce_p(&s2);
 338     bigint_mul_s(&h, &s2, &h);
 339     curve->reduce_p(&h);
 340     bigint_mul_s(&u1, &s2, &u1);
 341     curve->reduce_p(&u1);
 342
 343     /* compute new x */
 344     bigint_square(&u2, &r);
 345     curve->reduce_p(&u2);
 346     bigint_sub_s(&u2, &u2, &h);
 347     curve->reduce_p(&u2);
 348     bigint_sub_s(&u2, &u2, &u1);
 349     bigint_sub_s(&u2, &u2, &u1);
 350     curve->reduce_p(&u2);
 351     bigint_copy(&dest->x, &u2);
 352
 353     /* compute new y */
 354     bigint_sub_s(&u1, &u1, &u2);
 355     curve->reduce_p(&u1);
 356     bigint_mul_s(&s2, &u1, &r);
 357     curve->reduce_p(&s2);
 358     bigint_mul_s(&s1, &s1, &h);
 359     curve->reduce_p(&s1);
 360     bigint_sub_s(&s2, &s2, &s1);
 361     curve->reduce_p(&s2);
 362     bigint_copy(&dest->y, &s2);
 363
 364
 365     /* compute new z2 */
 366     bigint_square(&s1, &dest->z1);
 367     curve->reduce_p(&s1);
 368     bigint_copy(&dest->z2, &s1);
 369
 370     /* compute new z2 */
 371     bigint_mul_u(&s1, &s1, &dest->z1);
 372     curve->reduce_p(&s1);
 373     bigint_copy(&dest->z3, &s1);
 374
 375     return 0;
 376 }
 377
 378 uint8_t ecc_chudnovsky_double_and_add(ecc_chudnovsky_point_t *dest,
 379                                       const bigint_t *k,
 380                                       const ecc_chudnovsky_point_t *p,
 381                                       const ecc_curve_sp_t *curve){
 382     uint16_t i;
 383     uint8_t s = 0;
 384     bigint_word_t v, t;
 385
 386     for(i = k->length_W; i > 0; --i){
 387         v = 1 << (BIGINT_WORD_SIZE - 1);
 388         t = k->wordv[i - 1];
 389         do{
 390             if(s){
 391                 ecc_chudnovsky_point_double_sp(dest, dest, curve);
 392                 if(v & t){
 393                     ecc_chudnovsky_point_add_sp(dest, dest, p, curve);
 394                 }
 395             }else{
 396                 if(v & t){
 397                     s = 1;
 398                     ecc_chudnovsky_point_copy(dest, p);
 399                 }
 400             }
 401             v >>= 1;
 402         }while(v);
 403     }
 404     return 0;
 405 }
 406
 407 uint8_t bigint_to_naf(uint8_t *dest, uint16_t *length, const bigint_t *src){
 408     if(src->length_W == 0){
 409         *dest = 0;
 410         *length = 2;
 411         return 0;
 412     }
 413
 414     memset(dest, 0, src->length_W * sizeof(bigint_word_t) * 2 +1);
 415
 416     uint16_t i = 0;
 417     uint8_t t; /* 3 -> -1 ; 1 -> 1; 0 -> 0 (2 should not happen) */
 418     bigint_t k, p;
 419     bigint_word_t k_w[src->length_W + 1];
 420     bigint_word_t p_w = 1;
 421     p.wordv = &p_w;
 422     p.info = 0;
 423     p.length_W = 1;
 424     k.wordv = k_w;
 425     bigint_copy(&k, src);
 426
 427     while(k.length_W >= 1){
 428         if(k.wordv[0] & 1){
 429             t = k.wordv[0] & 3;
 430             if(t == 1){
 431                 bigint_sub_u(&k, &k, &p);
 432             }else{
 433                 bigint_add_u(&k, &k, &p);
 434             }
 435         }else{
 436             t = 0;
 437         }
 438         dest[(i * 2) / 8] |= t << ((2 * i) & 7);
 439         bigint_shiftright(&k, 1);
 440         i += 1;
 441     }
 442     *length = i;
 443     return 0;
 444 }
 445
 446 void print_naf(uint8_t *naf, uint16_t length){
 447     if(!length){
 448         return;
 449     }
 450     --length;
 451     int8_t t;
 452     do{
 453         t = (naf[(length * 2) / 8] >> ((length * 2) & 7)) & 3;
 454         switch(t & 3){
 455         case 0: putchar('0'); break;
 456         case 1: putchar('1'); break;
 457         case 3: putchar('-'); putchar('1'); break;
 458         default: putchar('E');
 459         }
 460         if(length){
 461             putchar(' ');
 462         }
 463     }while(length--);
 464 }
 465
 466 uint8_t ecc_chudnovsky_naf_multiplication(ecc_chudnovsky_point_t *dest,
 467                                           const bigint_t *k,
 468                                           const ecc_chudnovsky_point_t *p,
 469                                           const ecc_curve_sp_t *curve){
 470     if(k->length_W == 0 || p->y.length_W == 0){
 471         bigint_set_zero(&dest->y);
 472         return 0;
 473     }
 474
 475     uint8_t *t, q;
 476     uint16_t i;
 477     ecc_chudnovsky_point_t p_;
 478     bigint_word_t y_[curve->p->length_W];
 479
 480     /* p_ = -p*/
 481
 482     memcpy(&p_, p, sizeof(p_));
 483     p_.y.wordv = y_;
 484     bigint_copy(&p_.y, &p->y);
 485     p_.y.info |= BIGINT_NEG_MASK;
 486     bigint_add_s(&p_.y, &p_.y, curve->p);
 487
 488
 489     if(!(t = calloc(k->length_W * sizeof(bigint_word_t) * 2 + 1, 1))){
 490         return 1;
 491     }
 492     bigint_to_naf(t, &i, k);
 493
 494  //   printf(" naf: ");
 495  //   print_naf(t, i);
 496
 497     --i;
 498     dest->y.length_W = 0;
 499     do{
 500         q = (t[(i * 2) / 8] >> ((i * 2) & 7)) & 3;
 501         ecc_chudnovsky_point_double_sp(dest, dest, curve);
 502         if(q == 1){
 503             ecc_chudnovsky_point_add_sp(dest, dest, p, curve);
 504         }
 505         if(q == 3){
 506             ecc_chudnovsky_point_add_sp(dest, dest, &p_, curve);
 507         }
 508     }while(i--);
 509
 510     free(t);
 511
 512     return 0;
 513 }
 514
 515 uint8_t ecc_chudnovsky_multiplication(ecc_chudnovsky_point_t *dest,
 516                                       const bigint_t *k,
 517                                       const ecc_chudnovsky_point_t *p,
 518                                       const ecc_curve_sp_t *curve){
 519     return ecc_chudnovsky_naf_multiplication(dest, k, p, curve);
 520 }
 521
 522
 523
 524 uint8_t ecc_chudnovsky_multipy_and_sum(ecc_chudnovsky_point_t *dest,
 525                                       const bigint_t *k,
 526                                       const ecc_chudnovsky_point_t *p,
 527                                       const bigint_t *l,
 528                                       const ecc_chudnovsky_point_t *q,
 529                                       const ecc_curve_sp_t *curve){
 530     return ecc_chudnovsky_naf_multiplication(dest, k, p, curve);
 531 }