src/ecmult_impl.h

   1 /**********************************************************************
   2  * Copyright (c) 2013, 2014 Pieter Wuille                             *
   3  * Distributed under the MIT software license, see the accompanying   *
   4  * file COPYING or http://www.opensource.org/licenses/mit-license.php.*
   5  **********************************************************************/
   6
   7 #ifndef SECP256K1_ECMULT_IMPL_H
   8 #define SECP256K1_ECMULT_IMPL_H
   9
  10 #include <string.h>
  11
  12 #include "group.h"
  13 #include "scalar.h"
  14 #include "ecmult.h"
  15
  16 #if defined(EXHAUSTIVE_TEST_ORDER)
  17 /* We need to lower these values for exhaustive tests because
  18  * the tables cannot have infinities in them (this breaks the
  19  * affine-isomorphism stuff which tracks z-ratios) */
  20 #  if EXHAUSTIVE_TEST_ORDER > 128
  21 #    define WINDOW_A 5
  22 #    define WINDOW_G 8
  23 #  elif EXHAUSTIVE_TEST_ORDER > 8
  24 #    define WINDOW_A 4
  25 #    define WINDOW_G 4
  26 #  else
  27 #    define WINDOW_A 2
  28 #    define WINDOW_G 2
  29 #  endif
  30 #else
  31 /* optimal for 128-bit and 256-bit exponents. */
  32 #define WINDOW_A 5
  33 /** larger numbers may result in slightly better performance, at the cost of
  34     exponentially larger precomputed tables. */
  35 #ifdef USE_ENDOMORPHISM
  36 /** Two tables for window size 15: 1.375 MiB. */
  37 #define WINDOW_G 15
  38 #else
  39 /** One table for window size 16: 1.375 MiB. */
  40 #define WINDOW_G 16
  41 #endif
  42 #endif
  43
  44 /** The number of entries a table with precomputed multiples needs to have. */
  45 #define ECMULT_TABLE_SIZE(w) (1 << ((w)-2))
  46
  47 /** Fill a table 'prej' with precomputed odd multiples of a. Prej will contain
  48  *  the values [1*a,3*a,...,(2*n-1)*a], so it space for n values. zr[0] will
  49  *  contain prej[0].z / a.z. The other zr[i] values = prej[i].z / prej[i-1].z.
  50  *  Prej's Z values are undefined, except for the last value.
  51  */
  52 static void secp256k1_ecmult_odd_multiples_table(int n, secp256k1_gej *prej, secp256k1_fe *zr, const secp256k1_gej *a) {
  53     secp256k1_gej d;
  54     secp256k1_ge a_ge, d_ge;
  55     int i;
  56
  57     VERIFY_CHECK(!a->infinity);
  58
  59     secp256k1_gej_double_var(&d, a, NULL);
  60
  61     /*
  62      * Perform the additions on an isomorphism where 'd' is affine: drop the z coordinate
  63      * of 'd', and scale the 1P starting value's x/y coordinates without changing its z.
  64      */
  65     d_ge.x = d.x;
  66     d_ge.y = d.y;
  67     d_ge.infinity = 0;
  68
  69     secp256k1_ge_set_gej_zinv(&a_ge, a, &d.z);
  70     prej[0].x = a_ge.x;
  71     prej[0].y = a_ge.y;
  72     prej[0].z = a->z;
  73     prej[0].infinity = 0;
  74
  75     zr[0] = d.z;
  76     for (i = 1; i < n; i++) {
  77         secp256k1_gej_add_ge_var(&prej[i], &prej[i-1], &d_ge, &zr[i]);
  78     }
  79
  80     /*
  81      * Each point in 'prej' has a z coordinate too small by a factor of 'd.z'. Only
  82      * the final point's z coordinate is actually used though, so just update that.
  83      */
  84     secp256k1_fe_mul(&prej[n-1].z, &prej[n-1].z, &d.z);
  85 }
  86
  87 /** Fill a table 'pre' with precomputed odd multiples of a.
  88  *
  89  *  There are two versions of this function:
  90  *  - secp256k1_ecmult_odd_multiples_table_globalz_windowa which brings its
  91  *    resulting point set to a single constant Z denominator, stores the X and Y
  92  *    coordinates as ge_storage points in pre, and stores the global Z in rz.
  93  *    It only operates on tables sized for WINDOW_A wnaf multiples.
  94  *  - secp256k1_ecmult_odd_multiples_table_storage_var, which converts its
  95  *    resulting point set to actually affine points, and stores those in pre.
  96  *    It operates on tables of any size, but uses heap-allocated temporaries.
  97  *
  98  *  To compute a*P + b*G, we compute a table for P using the first function,
  99  *  and for G using the second (which requires an inverse, but it only needs to
 100  *  happen once).
 101  */
 102 static void secp256k1_ecmult_odd_multiples_table_globalz_windowa(secp256k1_ge *pre, secp256k1_fe *globalz, const secp256k1_gej *a) {
 103     secp256k1_gej prej[ECMULT_TABLE_SIZE(WINDOW_A)];
 104     secp256k1_fe zr[ECMULT_TABLE_SIZE(WINDOW_A)];
 105
 106     /* Compute the odd multiples in Jacobian form. */
 107     secp256k1_ecmult_odd_multiples_table(ECMULT_TABLE_SIZE(WINDOW_A), prej, zr, a);
 108     /* Bring them to the same Z denominator. */
 109     secp256k1_ge_globalz_set_table_gej(ECMULT_TABLE_SIZE(WINDOW_A), pre, globalz, prej, zr);
 110 }
 111
 112 static void secp256k1_ecmult_odd_multiples_table_storage_var(int n, secp256k1_ge_storage *pre, const secp256k1_gej *a, const secp256k1_callback *cb) {
 113     secp256k1_gej *prej = (secp256k1_gej*)checked_malloc(cb, sizeof(secp256k1_gej) * n);
 114     secp256k1_ge *prea = (secp256k1_ge*)checked_malloc(cb, sizeof(secp256k1_ge) * n);
 115     secp256k1_fe *zr = (secp256k1_fe*)checked_malloc(cb, sizeof(secp256k1_fe) * n);
 116     int i;
 117
 118     /* Compute the odd multiples in Jacobian form. */
 119     secp256k1_ecmult_odd_multiples_table(n, prej, zr, a);
 120     /* Convert them in batch to affine coordinates. */
 121     secp256k1_ge_set_table_gej_var(prea, prej, zr, n);
 122     /* Convert them to compact storage form. */
 123     for (i = 0; i < n; i++) {
 124         secp256k1_ge_to_storage(&pre[i], &prea[i]);
 125     }
 126
 127     free(prea);
 128     free(prej);
 129     free(zr);
 130 }
 131
 132 /** The following two macro retrieves a particular odd multiple from a table
 133  *  of precomputed multiples. */
 134 #define ECMULT_TABLE_GET_GE(r,pre,n,w) do { \
 135     VERIFY_CHECK(((n) & 1) == 1); \
 136     VERIFY_CHECK((n) >= -((1 << ((w)-1)) - 1)); \
 137     VERIFY_CHECK((n) <=  ((1 << ((w)-1)) - 1)); \
 138     if ((n) > 0) { \
 139         *(r) = (pre)[((n)-1)/2]; \
 140     } else { \
 141         secp256k1_ge_neg((r), &(pre)[(-(n)-1)/2]); \
 142     } \
 143 } while(0)
 144
 145 #define ECMULT_TABLE_GET_GE_STORAGE(r,pre,n,w) do { \
 146     VERIFY_CHECK(((n) & 1) == 1); \
 147     VERIFY_CHECK((n) >= -((1 << ((w)-1)) - 1)); \
 148     VERIFY_CHECK((n) <=  ((1 << ((w)-1)) - 1)); \
 149     if ((n) > 0) { \
 150         secp256k1_ge_from_storage((r), &(pre)[((n)-1)/2]); \
 151     } else { \
 152         secp256k1_ge_from_storage((r), &(pre)[(-(n)-1)/2]); \
 153         secp256k1_ge_neg((r), (r)); \
 154     } \
 155 } while(0)
 156
 157 static void secp256k1_ecmult_context_init(secp256k1_ecmult_context *ctx) {
 158     ctx->pre_g = NULL;
 159 #ifdef USE_ENDOMORPHISM
 160     ctx->pre_g_128 = NULL;
 161 #endif
 162 }
 163
 164 static void secp256k1_ecmult_context_build(secp256k1_ecmult_context *ctx, const secp256k1_callback *cb) {
 165     secp256k1_gej gj;
 166
 167     if (ctx->pre_g != NULL) {
 168         return;
 169     }
 170
 171     /* get the generator */
 172     secp256k1_gej_set_ge(&gj, &secp256k1_ge_const_g);
 173
 174     ctx->pre_g = (secp256k1_ge_storage (*)[])checked_malloc(cb, sizeof((*ctx->pre_g)[0]) * ECMULT_TABLE_SIZE(WINDOW_G));
 175
 176     /* precompute the tables with odd multiples */
 177     secp256k1_ecmult_odd_multiples_table_storage_var(ECMULT_TABLE_SIZE(WINDOW_G), *ctx->pre_g, &gj, cb);
 178
 179 #ifdef USE_ENDOMORPHISM
 180     {
 181         secp256k1_gej g_128j;
 182         int i;
 183
 184         ctx->pre_g_128 = (secp256k1_ge_storage (*)[])checked_malloc(cb, sizeof((*ctx->pre_g_128)[0]) * ECMULT_TABLE_SIZE(WINDOW_G));
 185
 186         /* calculate 2^128*generator */
 187         g_128j = gj;
 188         for (i = 0; i < 128; i++) {
 189             secp256k1_gej_double_var(&g_128j, &g_128j, NULL);
 190         }
 191         secp256k1_ecmult_odd_multiples_table_storage_var(ECMULT_TABLE_SIZE(WINDOW_G), *ctx->pre_g_128, &g_128j, cb);
 192     }
 193 #endif
 194 }
 195
 196 static void secp256k1_ecmult_context_clone(secp256k1_ecmult_context *dst,
 197                                            const secp256k1_ecmult_context *src, const secp256k1_callback *cb) {
 198     if (src->pre_g == NULL) {
 199         dst->pre_g = NULL;
 200     } else {
 201         size_t size = sizeof((*dst->pre_g)[0]) * ECMULT_TABLE_SIZE(WINDOW_G);
 202         dst->pre_g = (secp256k1_ge_storage (*)[])checked_malloc(cb, size);
 203         memcpy(dst->pre_g, src->pre_g, size);
 204     }
 205 #ifdef USE_ENDOMORPHISM
 206     if (src->pre_g_128 == NULL) {
 207         dst->pre_g_128 = NULL;
 208     } else {
 209         size_t size = sizeof((*dst->pre_g_128)[0]) * ECMULT_TABLE_SIZE(WINDOW_G);
 210         dst->pre_g_128 = (secp256k1_ge_storage (*)[])checked_malloc(cb, size);
 211         memcpy(dst->pre_g_128, src->pre_g_128, size);
 212     }
 213 #endif
 214 }
 215
 216 static int secp256k1_ecmult_context_is_built(const secp256k1_ecmult_context *ctx) {
 217     return ctx->pre_g != NULL;
 218 }
 219
 220 static void secp256k1_ecmult_context_clear(secp256k1_ecmult_context *ctx) {
 221     free(ctx->pre_g);
 222 #ifdef USE_ENDOMORPHISM
 223     free(ctx->pre_g_128);
 224 #endif
 225     secp256k1_ecmult_context_init(ctx);
 226 }
 227
 228 /** Convert a number to WNAF notation. The number becomes represented by sum(2^i * wnaf[i], i=0..bits),
 229  *  with the following guarantees:
 230  *  - each wnaf[i] is either 0, or an odd integer between -(1<<(w-1) - 1) and (1<<(w-1) - 1)
 231  *  - two non-zero entries in wnaf are separated by at least w-1 zeroes.
 232  *  - the number of set values in wnaf is returned. This number is at most 256, and at most one more
 233  *    than the number of bits in the (absolute value) of the input.
 234  */
 235 static int secp256k1_ecmult_wnaf(int *wnaf, int len, const secp256k1_scalar *a, int w) {
 236     secp256k1_scalar s = *a;
 237     int last_set_bit = -1;
 238     int bit = 0;
 239     int sign = 1;
 240     int carry = 0;
 241
 242     VERIFY_CHECK(wnaf != NULL);
 243     VERIFY_CHECK(0 <= len && len <= 256);
 244     VERIFY_CHECK(a != NULL);
 245     VERIFY_CHECK(2 <= w && w <= 31);
 246
 247     memset(wnaf, 0, len * sizeof(wnaf[0]));
 248
 249     if (secp256k1_scalar_get_bits(&s, 255, 1)) {
 250         secp256k1_scalar_negate(&s, &s);
 251         sign = -1;
 252     }
 253
 254     while (bit < len) {
 255         int now;
 256         int word;
 257         if (secp256k1_scalar_get_bits(&s, bit, 1) == (unsigned int)carry) {
 258             bit++;
 259             continue;
 260         }
 261
 262         now = w;
 263         if (now > len - bit) {
 264             now = len - bit;
 265         }
 266
 267         word = secp256k1_scalar_get_bits_var(&s, bit, now) + carry;
 268
 269         carry = (word >> (w-1)) & 1;
 270         word -= carry << w;
 271
 272         wnaf[bit] = sign * word;
 273         last_set_bit = bit;
 274
 275         bit += now;
 276     }
 277 #ifdef VERIFY
 278     CHECK(carry == 0);
 279     while (bit < 256) {
 280         CHECK(secp256k1_scalar_get_bits(&s, bit++, 1) == 0);
 281     }
 282 #endif
 283     return last_set_bit + 1;
 284 }
 285
 286 static void secp256k1_ecmult(const secp256k1_ecmult_context *ctx, secp256k1_gej *r, const secp256k1_gej *a, const secp256k1_scalar *na, const secp256k1_scalar *ng) {
 287     secp256k1_ge pre_a[ECMULT_TABLE_SIZE(WINDOW_A)];
 288     secp256k1_ge tmpa;
 289     secp256k1_fe Z;
 290 #ifdef USE_ENDOMORPHISM
 291     secp256k1_ge pre_a_lam[ECMULT_TABLE_SIZE(WINDOW_A)];
 292     secp256k1_scalar na_1, na_lam;
 293     /* Splitted G factors. */
 294     secp256k1_scalar ng_1, ng_128;
 295     int wnaf_na_1[130];
 296     int wnaf_na_lam[130];
 297     int bits_na_1;
 298     int bits_na_lam;
 299     int wnaf_ng_1[129];
 300     int bits_ng_1;
 301     int wnaf_ng_128[129];
 302     int bits_ng_128;
 303 #else
 304     int wnaf_na[256];
 305     int bits_na;
 306     int wnaf_ng[256];
 307     int bits_ng;
 308 #endif
 309     int i;
 310     int bits;
 311
 312 #ifdef USE_ENDOMORPHISM
 313     /* split na into na_1 and na_lam (where na = na_1 + na_lam*lambda, and na_1 and na_lam are ~128 bit) */
 314     secp256k1_scalar_split_lambda(&na_1, &na_lam, na);
 315
 316     /* build wnaf representation for na_1 and na_lam. */
 317     bits_na_1   = secp256k1_ecmult_wnaf(wnaf_na_1,   130, &na_1,   WINDOW_A);
 318     bits_na_lam = secp256k1_ecmult_wnaf(wnaf_na_lam, 130, &na_lam, WINDOW_A);
 319     VERIFY_CHECK(bits_na_1 <= 130);
 320     VERIFY_CHECK(bits_na_lam <= 130);
 321     bits = bits_na_1;
 322     if (bits_na_lam > bits) {
 323         bits = bits_na_lam;
 324     }
 325 #else
 326     /* build wnaf representation for na. */
 327     bits_na     = secp256k1_ecmult_wnaf(wnaf_na,     256, na,      WINDOW_A);
 328     bits = bits_na;
 329 #endif
 330
 331     /* Calculate odd multiples of a.
 332      * All multiples are brought to the same Z 'denominator', which is stored
 333      * in Z. Due to secp256k1' isomorphism we can do all operations pretending
 334      * that the Z coordinate was 1, use affine addition formulae, and correct
 335      * the Z coordinate of the result once at the end.
 336      * The exception is the precomputed G table points, which are actually
 337      * affine. Compared to the base used for other points, they have a Z ratio
 338      * of 1/Z, so we can use secp256k1_gej_add_zinv_var, which uses the same
 339      * isomorphism to efficiently add with a known Z inverse.
 340      */
 341     secp256k1_ecmult_odd_multiples_table_globalz_windowa(pre_a, &Z, a);
 342
 343 #ifdef USE_ENDOMORPHISM
 344     for (i = 0; i < ECMULT_TABLE_SIZE(WINDOW_A); i++) {
 345         secp256k1_ge_mul_lambda(&pre_a_lam[i], &pre_a[i]);
 346     }
 347
 348     /* split ng into ng_1 and ng_128 (where gn = gn_1 + gn_128*2^128, and gn_1 and gn_128 are ~128 bit) */
 349     secp256k1_scalar_split_128(&ng_1, &ng_128, ng);
 350
 351     /* Build wnaf representation for ng_1 and ng_128 */
 352     bits_ng_1   = secp256k1_ecmult_wnaf(wnaf_ng_1,   129, &ng_1,   WINDOW_G);
 353     bits_ng_128 = secp256k1_ecmult_wnaf(wnaf_ng_128, 129, &ng_128, WINDOW_G);
 354     if (bits_ng_1 > bits) {
 355         bits = bits_ng_1;
 356     }
 357     if (bits_ng_128 > bits) {
 358         bits = bits_ng_128;
 359     }
 360 #else
 361     bits_ng     = secp256k1_ecmult_wnaf(wnaf_ng,     256, ng,      WINDOW_G);
 362     if (bits_ng > bits) {
 363         bits = bits_ng;
 364     }
 365 #endif
 366
 367     secp256k1_gej_set_infinity(r);
 368
 369     for (i = bits - 1; i >= 0; i--) {
 370         int n;
 371         secp256k1_gej_double_var(r, r, NULL);
 372 #ifdef USE_ENDOMORPHISM
 373         if (i < bits_na_1 && (n = wnaf_na_1[i])) {
 374             ECMULT_TABLE_GET_GE(&tmpa, pre_a, n, WINDOW_A);
 375             secp256k1_gej_add_ge_var(r, r, &tmpa, NULL);
 376         }
 377         if (i < bits_na_lam && (n = wnaf_na_lam[i])) {
 378             ECMULT_TABLE_GET_GE(&tmpa, pre_a_lam, n, WINDOW_A);
 379             secp256k1_gej_add_ge_var(r, r, &tmpa, NULL);
 380         }
 381         if (i < bits_ng_1 && (n = wnaf_ng_1[i])) {
 382             ECMULT_TABLE_GET_GE_STORAGE(&tmpa, *ctx->pre_g, n, WINDOW_G);
 383             secp256k1_gej_add_zinv_var(r, r, &tmpa, &Z);
 384         }
 385         if (i < bits_ng_128 && (n = wnaf_ng_128[i])) {
 386             ECMULT_TABLE_GET_GE_STORAGE(&tmpa, *ctx->pre_g_128, n, WINDOW_G);
 387             secp256k1_gej_add_zinv_var(r, r, &tmpa, &Z);
 388         }
 389 #else
 390         if (i < bits_na && (n = wnaf_na[i])) {
 391             ECMULT_TABLE_GET_GE(&tmpa, pre_a, n, WINDOW_A);
 392             secp256k1_gej_add_ge_var(r, r, &tmpa, NULL);
 393         }
 394         if (i < bits_ng && (n = wnaf_ng[i])) {
 395             ECMULT_TABLE_GET_GE_STORAGE(&tmpa, *ctx->pre_g, n, WINDOW_G);
 396             secp256k1_gej_add_zinv_var(r, r, &tmpa, &Z);
 397         }
 398 #endif
 399     }
 400
 401     if (!r->infinity) {
 402         secp256k1_fe_mul(&r->z, &r->z, &Z);
 403     }
 404 }
 405
 406 #endif /* SECP256K1_ECMULT_IMPL_H */