2 * Copyright 2011-2018 The OpenSSL Project Authors. All Rights Reserved.
4 * Licensed under the Apache License 2.0 (the "License"). You may not use
5 * this file except in compliance with the License. You can obtain a copy
6 * in the file LICENSE in the source distribution or at
7 * https://www.openssl.org/source/license.html
10 /* Copyright 2011 Google Inc.
12 * Licensed under the Apache License, Version 2.0 (the "License");
14 * you may not use this file except in compliance with the License.
15 * You may obtain a copy of the License at
17 * http://www.apache.org/licenses/LICENSE-2.0
19 * Unless required by applicable law or agreed to in writing, software
20 * distributed under the License is distributed on an "AS IS" BASIS,
21 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
22 * See the License for the specific language governing permissions and
23 * limitations under the License.
27 * A 64-bit implementation of the NIST P-256 elliptic curve point multiplication
29 * OpenSSL integration was taken from Emilia Kasper's work in ecp_nistp224.c.
30 * Otherwise based on Emilia's P224 work, which was inspired by my curve25519
31 * work which got its smarts from Daniel J. Bernstein's work on the same.
34 #include <openssl/opensslconf.h>
35 #ifdef OPENSSL_NO_EC_NISTP_64_GCC_128
36 NON_EMPTY_TRANSLATION_UNIT
41 # include <openssl/err.h>
44 # if defined(__SIZEOF_INT128__) && __SIZEOF_INT128__==16
45 /* even with gcc, the typedef won't work for 32-bit platforms */
46 typedef __uint128_t uint128_t; /* nonstandard; implemented by gcc on 64-bit
48 typedef __int128_t int128_t;
50 # error "Your compiler doesn't appear to support 128-bit integer types"
58 * The underlying field. P256 operates over GF(2^256-2^224+2^192+2^96-1). We
59 * can serialise an element of this field into 32 bytes. We call this an
63 typedef u8 felem_bytearray[32];
66 * These are the parameters of P256, taken from FIPS 186-3, page 86. These
67 * values are big-endian.
69 static const felem_bytearray nistp256_curve_params[5] = {
70 {0xff, 0xff, 0xff, 0xff, 0x00, 0x00, 0x00, 0x01, /* p */
71 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
72 0x00, 0x00, 0x00, 0x00, 0xff, 0xff, 0xff, 0xff,
73 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff},
74 {0xff, 0xff, 0xff, 0xff, 0x00, 0x00, 0x00, 0x01, /* a = -3 */
75 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
76 0x00, 0x00, 0x00, 0x00, 0xff, 0xff, 0xff, 0xff,
77 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xfc}, /* b */
78 {0x5a, 0xc6, 0x35, 0xd8, 0xaa, 0x3a, 0x93, 0xe7,
79 0xb3, 0xeb, 0xbd, 0x55, 0x76, 0x98, 0x86, 0xbc,
80 0x65, 0x1d, 0x06, 0xb0, 0xcc, 0x53, 0xb0, 0xf6,
81 0x3b, 0xce, 0x3c, 0x3e, 0x27, 0xd2, 0x60, 0x4b},
82 {0x6b, 0x17, 0xd1, 0xf2, 0xe1, 0x2c, 0x42, 0x47, /* x */
83 0xf8, 0xbc, 0xe6, 0xe5, 0x63, 0xa4, 0x40, 0xf2,
84 0x77, 0x03, 0x7d, 0x81, 0x2d, 0xeb, 0x33, 0xa0,
85 0xf4, 0xa1, 0x39, 0x45, 0xd8, 0x98, 0xc2, 0x96},
86 {0x4f, 0xe3, 0x42, 0xe2, 0xfe, 0x1a, 0x7f, 0x9b, /* y */
87 0x8e, 0xe7, 0xeb, 0x4a, 0x7c, 0x0f, 0x9e, 0x16,
88 0x2b, 0xce, 0x33, 0x57, 0x6b, 0x31, 0x5e, 0xce,
89 0xcb, 0xb6, 0x40, 0x68, 0x37, 0xbf, 0x51, 0xf5}
93 * The representation of field elements.
94 * ------------------------------------
96 * We represent field elements with either four 128-bit values, eight 128-bit
97 * values, or four 64-bit values. The field element represented is:
98 * v[0]*2^0 + v[1]*2^64 + v[2]*2^128 + v[3]*2^192 (mod p)
100 * v[0]*2^0 + v[1]*2^64 + v[2]*2^128 + ... + v[8]*2^512 (mod p)
102 * 128-bit values are called 'limbs'. Since the limbs are spaced only 64 bits
103 * apart, but are 128-bits wide, the most significant bits of each limb overlap
104 * with the least significant bits of the next.
106 * A field element with four limbs is an 'felem'. One with eight limbs is a
109 * A field element with four, 64-bit values is called a 'smallfelem'. Small
110 * values are used as intermediate values before multiplication.
115 typedef uint128_t limb;
116 typedef limb felem[NLIMBS];
117 typedef limb longfelem[NLIMBS * 2];
118 typedef u64 smallfelem[NLIMBS];
120 /* This is the value of the prime as four 64-bit words, little-endian. */
121 static const u64 kPrime[4] =
122 { 0xfffffffffffffffful, 0xffffffff, 0, 0xffffffff00000001ul };
123 static const u64 bottom63bits = 0x7ffffffffffffffful;
126 * bin32_to_felem takes a little-endian byte array and converts it into felem
127 * form. This assumes that the CPU is little-endian.
129 static void bin32_to_felem(felem out, const u8 in[32])
131 out[0] = *((u64 *)&in[0]);
132 out[1] = *((u64 *)&in[8]);
133 out[2] = *((u64 *)&in[16]);
134 out[3] = *((u64 *)&in[24]);
138 * smallfelem_to_bin32 takes a smallfelem and serialises into a little
139 * endian, 32 byte array. This assumes that the CPU is little-endian.
141 static void smallfelem_to_bin32(u8 out[32], const smallfelem in)
143 *((u64 *)&out[0]) = in[0];
144 *((u64 *)&out[8]) = in[1];
145 *((u64 *)&out[16]) = in[2];
146 *((u64 *)&out[24]) = in[3];
149 /* To preserve endianness when using BN_bn2bin and BN_bin2bn */
150 static void flip_endian(u8 *out, const u8 *in, unsigned len)
153 for (i = 0; i < len; ++i)
154 out[i] = in[len - 1 - i];
157 /* BN_to_felem converts an OpenSSL BIGNUM into an felem */
158 static int BN_to_felem(felem out, const BIGNUM *bn)
160 felem_bytearray b_in;
161 felem_bytearray b_out;
164 num_bytes = BN_num_bytes(bn);
165 if (num_bytes > sizeof(b_out)) {
166 ECerr(EC_F_BN_TO_FELEM, EC_R_BIGNUM_OUT_OF_RANGE);
169 if (BN_is_negative(bn)) {
170 ECerr(EC_F_BN_TO_FELEM, EC_R_BIGNUM_OUT_OF_RANGE);
173 num_bytes = BN_bn2binpad(bn, b_in, sizeof(b_in));
174 flip_endian(b_out, b_in, num_bytes);
175 bin32_to_felem(out, b_out);
179 /* felem_to_BN converts an felem into an OpenSSL BIGNUM */
180 static BIGNUM *smallfelem_to_BN(BIGNUM *out, const smallfelem in)
182 felem_bytearray b_in, b_out;
183 smallfelem_to_bin32(b_in, in);
184 flip_endian(b_out, b_in, sizeof(b_out));
185 return BN_bin2bn(b_out, sizeof(b_out), out);
193 static void smallfelem_one(smallfelem out)
201 static void smallfelem_assign(smallfelem out, const smallfelem in)
209 static void felem_assign(felem out, const felem in)
217 /* felem_sum sets out = out + in. */
218 static void felem_sum(felem out, const felem in)
226 /* felem_small_sum sets out = out + in. */
227 static void felem_small_sum(felem out, const smallfelem in)
235 /* felem_scalar sets out = out * scalar */
236 static void felem_scalar(felem out, const u64 scalar)
244 /* longfelem_scalar sets out = out * scalar */
245 static void longfelem_scalar(longfelem out, const u64 scalar)
257 # define two105m41m9 (((limb)1) << 105) - (((limb)1) << 41) - (((limb)1) << 9)
258 # define two105 (((limb)1) << 105)
259 # define two105m41p9 (((limb)1) << 105) - (((limb)1) << 41) + (((limb)1) << 9)
261 /* zero105 is 0 mod p */
262 static const felem zero105 =
263 { two105m41m9, two105, two105m41p9, two105m41p9 };
266 * smallfelem_neg sets |out| to |-small|
268 * out[i] < out[i] + 2^105
270 static void smallfelem_neg(felem out, const smallfelem small)
272 /* In order to prevent underflow, we subtract from 0 mod p. */
273 out[0] = zero105[0] - small[0];
274 out[1] = zero105[1] - small[1];
275 out[2] = zero105[2] - small[2];
276 out[3] = zero105[3] - small[3];
280 * felem_diff subtracts |in| from |out|
284 * out[i] < out[i] + 2^105
286 static void felem_diff(felem out, const felem in)
289 * In order to prevent underflow, we add 0 mod p before subtracting.
291 out[0] += zero105[0];
292 out[1] += zero105[1];
293 out[2] += zero105[2];
294 out[3] += zero105[3];
302 # define two107m43m11 (((limb)1) << 107) - (((limb)1) << 43) - (((limb)1) << 11)
303 # define two107 (((limb)1) << 107)
304 # define two107m43p11 (((limb)1) << 107) - (((limb)1) << 43) + (((limb)1) << 11)
306 /* zero107 is 0 mod p */
307 static const felem zero107 =
308 { two107m43m11, two107, two107m43p11, two107m43p11 };
311 * An alternative felem_diff for larger inputs |in|
312 * felem_diff_zero107 subtracts |in| from |out|
316 * out[i] < out[i] + 2^107
318 static void felem_diff_zero107(felem out, const felem in)
321 * In order to prevent underflow, we add 0 mod p before subtracting.
323 out[0] += zero107[0];
324 out[1] += zero107[1];
325 out[2] += zero107[2];
326 out[3] += zero107[3];
335 * longfelem_diff subtracts |in| from |out|
339 * out[i] < out[i] + 2^70 + 2^40
341 static void longfelem_diff(longfelem out, const longfelem in)
343 static const limb two70m8p6 =
344 (((limb) 1) << 70) - (((limb) 1) << 8) + (((limb) 1) << 6);
345 static const limb two70p40 = (((limb) 1) << 70) + (((limb) 1) << 40);
346 static const limb two70 = (((limb) 1) << 70);
347 static const limb two70m40m38p6 =
348 (((limb) 1) << 70) - (((limb) 1) << 40) - (((limb) 1) << 38) +
350 static const limb two70m6 = (((limb) 1) << 70) - (((limb) 1) << 6);
352 /* add 0 mod p to avoid underflow */
356 out[3] += two70m40m38p6;
362 /* in[i] < 7*2^67 < 2^70 - 2^40 - 2^38 + 2^6 */
373 # define two64m0 (((limb)1) << 64) - 1
374 # define two110p32m0 (((limb)1) << 110) + (((limb)1) << 32) - 1
375 # define two64m46 (((limb)1) << 64) - (((limb)1) << 46)
376 # define two64m32 (((limb)1) << 64) - (((limb)1) << 32)
378 /* zero110 is 0 mod p */
379 static const felem zero110 = { two64m0, two110p32m0, two64m46, two64m32 };
382 * felem_shrink converts an felem into a smallfelem. The result isn't quite
383 * minimal as the value may be greater than p.
390 static void felem_shrink(smallfelem out, const felem in)
395 static const u64 kPrime3Test = 0x7fffffff00000001ul; /* 2^63 - 2^32 + 1 */
398 tmp[3] = zero110[3] + in[3] + ((u64)(in[2] >> 64));
401 tmp[2] = zero110[2] + (u64)in[2];
402 tmp[0] = zero110[0] + in[0];
403 tmp[1] = zero110[1] + in[1];
404 /* tmp[0] < 2**110, tmp[1] < 2^111, tmp[2] < 2**65 */
407 * We perform two partial reductions where we eliminate the high-word of
408 * tmp[3]. We don't update the other words till the end.
410 a = tmp[3] >> 64; /* a < 2^46 */
411 tmp[3] = (u64)tmp[3];
413 tmp[3] += ((limb) a) << 32;
417 a = tmp[3] >> 64; /* a < 2^15 */
418 b += a; /* b < 2^46 + 2^15 < 2^47 */
419 tmp[3] = (u64)tmp[3];
421 tmp[3] += ((limb) a) << 32;
422 /* tmp[3] < 2^64 + 2^47 */
425 * This adjusts the other two words to complete the two partial
429 tmp[1] -= (((limb) b) << 32);
432 * In order to make space in tmp[3] for the carry from 2 -> 3, we
433 * conditionally subtract kPrime if tmp[3] is large enough.
435 high = (u64)(tmp[3] >> 64);
436 /* As tmp[3] < 2^65, high is either 1 or 0 */
440 * all ones if the high word of tmp[3] is 1
441 * all zeros if the high word of tmp[3] if 0
444 mask = 0 - (low >> 63);
447 * all ones if the MSB of low is 1
448 * all zeros if the MSB of low if 0
452 /* if low was greater than kPrime3Test then the MSB is zero */
454 low = 0 - (low >> 63);
457 * all ones if low was > kPrime3Test
458 * all zeros if low was <= kPrime3Test
460 mask = (mask & low) | high;
461 tmp[0] -= mask & kPrime[0];
462 tmp[1] -= mask & kPrime[1];
463 /* kPrime[2] is zero, so omitted */
464 tmp[3] -= mask & kPrime[3];
465 /* tmp[3] < 2**64 - 2**32 + 1 */
467 tmp[1] += ((u64)(tmp[0] >> 64));
468 tmp[0] = (u64)tmp[0];
469 tmp[2] += ((u64)(tmp[1] >> 64));
470 tmp[1] = (u64)tmp[1];
471 tmp[3] += ((u64)(tmp[2] >> 64));
472 tmp[2] = (u64)tmp[2];
481 /* smallfelem_expand converts a smallfelem to an felem */
482 static void smallfelem_expand(felem out, const smallfelem in)
491 * smallfelem_square sets |out| = |small|^2
495 * out[i] < 7 * 2^64 < 2^67
497 static void smallfelem_square(longfelem out, const smallfelem small)
502 a = ((uint128_t) small[0]) * small[0];
508 a = ((uint128_t) small[0]) * small[1];
515 a = ((uint128_t) small[0]) * small[2];
522 a = ((uint128_t) small[0]) * small[3];
528 a = ((uint128_t) small[1]) * small[2];
535 a = ((uint128_t) small[1]) * small[1];
541 a = ((uint128_t) small[1]) * small[3];
548 a = ((uint128_t) small[2]) * small[3];
556 a = ((uint128_t) small[2]) * small[2];
562 a = ((uint128_t) small[3]) * small[3];
570 * felem_square sets |out| = |in|^2
574 * out[i] < 7 * 2^64 < 2^67
576 static void felem_square(longfelem out, const felem in)
579 felem_shrink(small, in);
580 smallfelem_square(out, small);
584 * smallfelem_mul sets |out| = |small1| * |small2|
589 * out[i] < 7 * 2^64 < 2^67
591 static void smallfelem_mul(longfelem out, const smallfelem small1,
592 const smallfelem small2)
597 a = ((uint128_t) small1[0]) * small2[0];
603 a = ((uint128_t) small1[0]) * small2[1];
609 a = ((uint128_t) small1[1]) * small2[0];
615 a = ((uint128_t) small1[0]) * small2[2];
621 a = ((uint128_t) small1[1]) * small2[1];
627 a = ((uint128_t) small1[2]) * small2[0];
633 a = ((uint128_t) small1[0]) * small2[3];
639 a = ((uint128_t) small1[1]) * small2[2];
645 a = ((uint128_t) small1[2]) * small2[1];
651 a = ((uint128_t) small1[3]) * small2[0];
657 a = ((uint128_t) small1[1]) * small2[3];
663 a = ((uint128_t) small1[2]) * small2[2];
669 a = ((uint128_t) small1[3]) * small2[1];
675 a = ((uint128_t) small1[2]) * small2[3];
681 a = ((uint128_t) small1[3]) * small2[2];
687 a = ((uint128_t) small1[3]) * small2[3];
695 * felem_mul sets |out| = |in1| * |in2|
700 * out[i] < 7 * 2^64 < 2^67
702 static void felem_mul(longfelem out, const felem in1, const felem in2)
704 smallfelem small1, small2;
705 felem_shrink(small1, in1);
706 felem_shrink(small2, in2);
707 smallfelem_mul(out, small1, small2);
711 * felem_small_mul sets |out| = |small1| * |in2|
716 * out[i] < 7 * 2^64 < 2^67
718 static void felem_small_mul(longfelem out, const smallfelem small1,
722 felem_shrink(small2, in2);
723 smallfelem_mul(out, small1, small2);
726 # define two100m36m4 (((limb)1) << 100) - (((limb)1) << 36) - (((limb)1) << 4)
727 # define two100 (((limb)1) << 100)
728 # define two100m36p4 (((limb)1) << 100) - (((limb)1) << 36) + (((limb)1) << 4)
729 /* zero100 is 0 mod p */
730 static const felem zero100 =
731 { two100m36m4, two100, two100m36p4, two100m36p4 };
734 * Internal function for the different flavours of felem_reduce.
735 * felem_reduce_ reduces the higher coefficients in[4]-in[7].
737 * out[0] >= in[6] + 2^32*in[6] + in[7] + 2^32*in[7]
738 * out[1] >= in[7] + 2^32*in[4]
739 * out[2] >= in[5] + 2^32*in[5]
740 * out[3] >= in[4] + 2^32*in[5] + 2^32*in[6]
742 * out[0] <= out[0] + in[4] + 2^32*in[5]
743 * out[1] <= out[1] + in[5] + 2^33*in[6]
744 * out[2] <= out[2] + in[7] + 2*in[6] + 2^33*in[7]
745 * out[3] <= out[3] + 2^32*in[4] + 3*in[7]
747 static void felem_reduce_(felem out, const longfelem in)
750 /* combine common terms from below */
751 c = in[4] + (in[5] << 32);
759 /* the remaining terms */
760 /* 256: [(0,1),(96,-1),(192,-1),(224,1)] */
761 out[1] -= (in[4] << 32);
762 out[3] += (in[4] << 32);
764 /* 320: [(32,1),(64,1),(128,-1),(160,-1),(224,-1)] */
765 out[2] -= (in[5] << 32);
767 /* 384: [(0,-1),(32,-1),(96,2),(128,2),(224,-1)] */
769 out[0] -= (in[6] << 32);
770 out[1] += (in[6] << 33);
771 out[2] += (in[6] * 2);
772 out[3] -= (in[6] << 32);
774 /* 448: [(0,-1),(32,-1),(64,-1),(128,1),(160,2),(192,3)] */
776 out[0] -= (in[7] << 32);
777 out[2] += (in[7] << 33);
778 out[3] += (in[7] * 3);
782 * felem_reduce converts a longfelem into an felem.
783 * To be called directly after felem_square or felem_mul.
785 * in[0] < 2^64, in[1] < 3*2^64, in[2] < 5*2^64, in[3] < 7*2^64
786 * in[4] < 7*2^64, in[5] < 5*2^64, in[6] < 3*2^64, in[7] < 2*64
790 static void felem_reduce(felem out, const longfelem in)
792 out[0] = zero100[0] + in[0];
793 out[1] = zero100[1] + in[1];
794 out[2] = zero100[2] + in[2];
795 out[3] = zero100[3] + in[3];
797 felem_reduce_(out, in);
800 * out[0] > 2^100 - 2^36 - 2^4 - 3*2^64 - 3*2^96 - 2^64 - 2^96 > 0
801 * out[1] > 2^100 - 2^64 - 7*2^96 > 0
802 * out[2] > 2^100 - 2^36 + 2^4 - 5*2^64 - 5*2^96 > 0
803 * out[3] > 2^100 - 2^36 + 2^4 - 7*2^64 - 5*2^96 - 3*2^96 > 0
805 * out[0] < 2^100 + 2^64 + 7*2^64 + 5*2^96 < 2^101
806 * out[1] < 2^100 + 3*2^64 + 5*2^64 + 3*2^97 < 2^101
807 * out[2] < 2^100 + 5*2^64 + 2^64 + 3*2^65 + 2^97 < 2^101
808 * out[3] < 2^100 + 7*2^64 + 7*2^96 + 3*2^64 < 2^101
813 * felem_reduce_zero105 converts a larger longfelem into an felem.
819 static void felem_reduce_zero105(felem out, const longfelem in)
821 out[0] = zero105[0] + in[0];
822 out[1] = zero105[1] + in[1];
823 out[2] = zero105[2] + in[2];
824 out[3] = zero105[3] + in[3];
826 felem_reduce_(out, in);
829 * out[0] > 2^105 - 2^41 - 2^9 - 2^71 - 2^103 - 2^71 - 2^103 > 0
830 * out[1] > 2^105 - 2^71 - 2^103 > 0
831 * out[2] > 2^105 - 2^41 + 2^9 - 2^71 - 2^103 > 0
832 * out[3] > 2^105 - 2^41 + 2^9 - 2^71 - 2^103 - 2^103 > 0
834 * out[0] < 2^105 + 2^71 + 2^71 + 2^103 < 2^106
835 * out[1] < 2^105 + 2^71 + 2^71 + 2^103 < 2^106
836 * out[2] < 2^105 + 2^71 + 2^71 + 2^71 + 2^103 < 2^106
837 * out[3] < 2^105 + 2^71 + 2^103 + 2^71 < 2^106
842 * subtract_u64 sets *result = *result - v and *carry to one if the
843 * subtraction underflowed.
845 static void subtract_u64(u64 *result, u64 *carry, u64 v)
847 uint128_t r = *result;
849 *carry = (r >> 64) & 1;
854 * felem_contract converts |in| to its unique, minimal representation. On
855 * entry: in[i] < 2^109
857 static void felem_contract(smallfelem out, const felem in)
860 u64 all_equal_so_far = 0, result = 0, carry;
862 felem_shrink(out, in);
863 /* small is minimal except that the value might be > p */
867 * We are doing a constant time test if out >= kPrime. We need to compare
868 * each u64, from most-significant to least significant. For each one, if
869 * all words so far have been equal (m is all ones) then a non-equal
870 * result is the answer. Otherwise we continue.
872 for (i = 3; i < 4; i--) {
874 uint128_t a = ((uint128_t) kPrime[i]) - out[i];
876 * if out[i] > kPrime[i] then a will underflow and the high 64-bits
879 result |= all_equal_so_far & ((u64)(a >> 64));
882 * if kPrime[i] == out[i] then |equal| will be all zeros and the
883 * decrement will make it all ones.
885 equal = kPrime[i] ^ out[i];
887 equal &= equal << 32;
888 equal &= equal << 16;
893 equal = 0 - (equal >> 63);
895 all_equal_so_far &= equal;
899 * if all_equal_so_far is still all ones then the two values are equal
900 * and so out >= kPrime is true.
902 result |= all_equal_so_far;
904 /* if out >= kPrime then we subtract kPrime. */
905 subtract_u64(&out[0], &carry, result & kPrime[0]);
906 subtract_u64(&out[1], &carry, carry);
907 subtract_u64(&out[2], &carry, carry);
908 subtract_u64(&out[3], &carry, carry);
910 subtract_u64(&out[1], &carry, result & kPrime[1]);
911 subtract_u64(&out[2], &carry, carry);
912 subtract_u64(&out[3], &carry, carry);
914 subtract_u64(&out[2], &carry, result & kPrime[2]);
915 subtract_u64(&out[3], &carry, carry);
917 subtract_u64(&out[3], &carry, result & kPrime[3]);
920 static void smallfelem_square_contract(smallfelem out, const smallfelem in)
925 smallfelem_square(longtmp, in);
926 felem_reduce(tmp, longtmp);
927 felem_contract(out, tmp);
930 static void smallfelem_mul_contract(smallfelem out, const smallfelem in1,
931 const smallfelem in2)
936 smallfelem_mul(longtmp, in1, in2);
937 felem_reduce(tmp, longtmp);
938 felem_contract(out, tmp);
942 * felem_is_zero returns a limb with all bits set if |in| == 0 (mod p) and 0
947 static limb smallfelem_is_zero(const smallfelem small)
952 u64 is_zero = small[0] | small[1] | small[2] | small[3];
954 is_zero &= is_zero << 32;
955 is_zero &= is_zero << 16;
956 is_zero &= is_zero << 8;
957 is_zero &= is_zero << 4;
958 is_zero &= is_zero << 2;
959 is_zero &= is_zero << 1;
960 is_zero = 0 - (is_zero >> 63);
962 is_p = (small[0] ^ kPrime[0]) |
963 (small[1] ^ kPrime[1]) |
964 (small[2] ^ kPrime[2]) | (small[3] ^ kPrime[3]);
972 is_p = 0 - (is_p >> 63);
977 result |= ((limb) is_zero) << 64;
981 static int smallfelem_is_zero_int(const void *small)
983 return (int)(smallfelem_is_zero(small) & ((limb) 1));
987 * felem_inv calculates |out| = |in|^{-1}
989 * Based on Fermat's Little Theorem:
991 * a^{p-1} = 1 (mod p)
992 * a^{p-2} = a^{-1} (mod p)
994 static void felem_inv(felem out, const felem in)
997 /* each e_I will hold |in|^{2^I - 1} */
998 felem e2, e4, e8, e16, e32, e64;
1002 felem_square(tmp, in);
1003 felem_reduce(ftmp, tmp); /* 2^1 */
1004 felem_mul(tmp, in, ftmp);
1005 felem_reduce(ftmp, tmp); /* 2^2 - 2^0 */
1006 felem_assign(e2, ftmp);
1007 felem_square(tmp, ftmp);
1008 felem_reduce(ftmp, tmp); /* 2^3 - 2^1 */
1009 felem_square(tmp, ftmp);
1010 felem_reduce(ftmp, tmp); /* 2^4 - 2^2 */
1011 felem_mul(tmp, ftmp, e2);
1012 felem_reduce(ftmp, tmp); /* 2^4 - 2^0 */
1013 felem_assign(e4, ftmp);
1014 felem_square(tmp, ftmp);
1015 felem_reduce(ftmp, tmp); /* 2^5 - 2^1 */
1016 felem_square(tmp, ftmp);
1017 felem_reduce(ftmp, tmp); /* 2^6 - 2^2 */
1018 felem_square(tmp, ftmp);
1019 felem_reduce(ftmp, tmp); /* 2^7 - 2^3 */
1020 felem_square(tmp, ftmp);
1021 felem_reduce(ftmp, tmp); /* 2^8 - 2^4 */
1022 felem_mul(tmp, ftmp, e4);
1023 felem_reduce(ftmp, tmp); /* 2^8 - 2^0 */
1024 felem_assign(e8, ftmp);
1025 for (i = 0; i < 8; i++) {
1026 felem_square(tmp, ftmp);
1027 felem_reduce(ftmp, tmp);
1029 felem_mul(tmp, ftmp, e8);
1030 felem_reduce(ftmp, tmp); /* 2^16 - 2^0 */
1031 felem_assign(e16, ftmp);
1032 for (i = 0; i < 16; i++) {
1033 felem_square(tmp, ftmp);
1034 felem_reduce(ftmp, tmp);
1036 felem_mul(tmp, ftmp, e16);
1037 felem_reduce(ftmp, tmp); /* 2^32 - 2^0 */
1038 felem_assign(e32, ftmp);
1039 for (i = 0; i < 32; i++) {
1040 felem_square(tmp, ftmp);
1041 felem_reduce(ftmp, tmp);
1043 felem_assign(e64, ftmp);
1044 felem_mul(tmp, ftmp, in);
1045 felem_reduce(ftmp, tmp); /* 2^64 - 2^32 + 2^0 */
1046 for (i = 0; i < 192; i++) {
1047 felem_square(tmp, ftmp);
1048 felem_reduce(ftmp, tmp);
1049 } /* 2^256 - 2^224 + 2^192 */
1051 felem_mul(tmp, e64, e32);
1052 felem_reduce(ftmp2, tmp); /* 2^64 - 2^0 */
1053 for (i = 0; i < 16; i++) {
1054 felem_square(tmp, ftmp2);
1055 felem_reduce(ftmp2, tmp);
1057 felem_mul(tmp, ftmp2, e16);
1058 felem_reduce(ftmp2, tmp); /* 2^80 - 2^0 */
1059 for (i = 0; i < 8; i++) {
1060 felem_square(tmp, ftmp2);
1061 felem_reduce(ftmp2, tmp);
1063 felem_mul(tmp, ftmp2, e8);
1064 felem_reduce(ftmp2, tmp); /* 2^88 - 2^0 */
1065 for (i = 0; i < 4; i++) {
1066 felem_square(tmp, ftmp2);
1067 felem_reduce(ftmp2, tmp);
1069 felem_mul(tmp, ftmp2, e4);
1070 felem_reduce(ftmp2, tmp); /* 2^92 - 2^0 */
1071 felem_square(tmp, ftmp2);
1072 felem_reduce(ftmp2, tmp); /* 2^93 - 2^1 */
1073 felem_square(tmp, ftmp2);
1074 felem_reduce(ftmp2, tmp); /* 2^94 - 2^2 */
1075 felem_mul(tmp, ftmp2, e2);
1076 felem_reduce(ftmp2, tmp); /* 2^94 - 2^0 */
1077 felem_square(tmp, ftmp2);
1078 felem_reduce(ftmp2, tmp); /* 2^95 - 2^1 */
1079 felem_square(tmp, ftmp2);
1080 felem_reduce(ftmp2, tmp); /* 2^96 - 2^2 */
1081 felem_mul(tmp, ftmp2, in);
1082 felem_reduce(ftmp2, tmp); /* 2^96 - 3 */
1084 felem_mul(tmp, ftmp2, ftmp);
1085 felem_reduce(out, tmp); /* 2^256 - 2^224 + 2^192 + 2^96 - 3 */
1088 static void smallfelem_inv_contract(smallfelem out, const smallfelem in)
1092 smallfelem_expand(tmp, in);
1093 felem_inv(tmp, tmp);
1094 felem_contract(out, tmp);
1101 * Building on top of the field operations we have the operations on the
1102 * elliptic curve group itself. Points on the curve are represented in Jacobian
1107 * point_double calculates 2*(x_in, y_in, z_in)
1109 * The method is taken from:
1110 * http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#doubling-dbl-2001-b
1112 * Outputs can equal corresponding inputs, i.e., x_out == x_in is allowed.
1113 * while x_out == y_in is not (maybe this works, but it's not tested).
1116 point_double(felem x_out, felem y_out, felem z_out,
1117 const felem x_in, const felem y_in, const felem z_in)
1119 longfelem tmp, tmp2;
1120 felem delta, gamma, beta, alpha, ftmp, ftmp2;
1121 smallfelem small1, small2;
1123 felem_assign(ftmp, x_in);
1124 /* ftmp[i] < 2^106 */
1125 felem_assign(ftmp2, x_in);
1126 /* ftmp2[i] < 2^106 */
1129 felem_square(tmp, z_in);
1130 felem_reduce(delta, tmp);
1131 /* delta[i] < 2^101 */
1134 felem_square(tmp, y_in);
1135 felem_reduce(gamma, tmp);
1136 /* gamma[i] < 2^101 */
1137 felem_shrink(small1, gamma);
1139 /* beta = x*gamma */
1140 felem_small_mul(tmp, small1, x_in);
1141 felem_reduce(beta, tmp);
1142 /* beta[i] < 2^101 */
1144 /* alpha = 3*(x-delta)*(x+delta) */
1145 felem_diff(ftmp, delta);
1146 /* ftmp[i] < 2^105 + 2^106 < 2^107 */
1147 felem_sum(ftmp2, delta);
1148 /* ftmp2[i] < 2^105 + 2^106 < 2^107 */
1149 felem_scalar(ftmp2, 3);
1150 /* ftmp2[i] < 3 * 2^107 < 2^109 */
1151 felem_mul(tmp, ftmp, ftmp2);
1152 felem_reduce(alpha, tmp);
1153 /* alpha[i] < 2^101 */
1154 felem_shrink(small2, alpha);
1156 /* x' = alpha^2 - 8*beta */
1157 smallfelem_square(tmp, small2);
1158 felem_reduce(x_out, tmp);
1159 felem_assign(ftmp, beta);
1160 felem_scalar(ftmp, 8);
1161 /* ftmp[i] < 8 * 2^101 = 2^104 */
1162 felem_diff(x_out, ftmp);
1163 /* x_out[i] < 2^105 + 2^101 < 2^106 */
1165 /* z' = (y + z)^2 - gamma - delta */
1166 felem_sum(delta, gamma);
1167 /* delta[i] < 2^101 + 2^101 = 2^102 */
1168 felem_assign(ftmp, y_in);
1169 felem_sum(ftmp, z_in);
1170 /* ftmp[i] < 2^106 + 2^106 = 2^107 */
1171 felem_square(tmp, ftmp);
1172 felem_reduce(z_out, tmp);
1173 felem_diff(z_out, delta);
1174 /* z_out[i] < 2^105 + 2^101 < 2^106 */
1176 /* y' = alpha*(4*beta - x') - 8*gamma^2 */
1177 felem_scalar(beta, 4);
1178 /* beta[i] < 4 * 2^101 = 2^103 */
1179 felem_diff_zero107(beta, x_out);
1180 /* beta[i] < 2^107 + 2^103 < 2^108 */
1181 felem_small_mul(tmp, small2, beta);
1182 /* tmp[i] < 7 * 2^64 < 2^67 */
1183 smallfelem_square(tmp2, small1);
1184 /* tmp2[i] < 7 * 2^64 */
1185 longfelem_scalar(tmp2, 8);
1186 /* tmp2[i] < 8 * 7 * 2^64 = 7 * 2^67 */
1187 longfelem_diff(tmp, tmp2);
1188 /* tmp[i] < 2^67 + 2^70 + 2^40 < 2^71 */
1189 felem_reduce_zero105(y_out, tmp);
1190 /* y_out[i] < 2^106 */
1194 * point_double_small is the same as point_double, except that it operates on
1198 point_double_small(smallfelem x_out, smallfelem y_out, smallfelem z_out,
1199 const smallfelem x_in, const smallfelem y_in,
1200 const smallfelem z_in)
1202 felem felem_x_out, felem_y_out, felem_z_out;
1203 felem felem_x_in, felem_y_in, felem_z_in;
1205 smallfelem_expand(felem_x_in, x_in);
1206 smallfelem_expand(felem_y_in, y_in);
1207 smallfelem_expand(felem_z_in, z_in);
1208 point_double(felem_x_out, felem_y_out, felem_z_out,
1209 felem_x_in, felem_y_in, felem_z_in);
1210 felem_shrink(x_out, felem_x_out);
1211 felem_shrink(y_out, felem_y_out);
1212 felem_shrink(z_out, felem_z_out);
1215 /* copy_conditional copies in to out iff mask is all ones. */
1216 static void copy_conditional(felem out, const felem in, limb mask)
1219 for (i = 0; i < NLIMBS; ++i) {
1220 const limb tmp = mask & (in[i] ^ out[i]);
1225 /* copy_small_conditional copies in to out iff mask is all ones. */
1226 static void copy_small_conditional(felem out, const smallfelem in, limb mask)
1229 const u64 mask64 = mask;
1230 for (i = 0; i < NLIMBS; ++i) {
1231 out[i] = ((limb) (in[i] & mask64)) | (out[i] & ~mask);
1236 * point_add calculates (x1, y1, z1) + (x2, y2, z2)
1238 * The method is taken from:
1239 * http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#addition-add-2007-bl,
1240 * adapted for mixed addition (z2 = 1, or z2 = 0 for the point at infinity).
1242 * This function includes a branch for checking whether the two input points
1243 * are equal, (while not equal to the point at infinity). This case never
1244 * happens during single point multiplication, so there is no timing leak for
1245 * ECDH or ECDSA signing.
1247 static void point_add(felem x3, felem y3, felem z3,
1248 const felem x1, const felem y1, const felem z1,
1249 const int mixed, const smallfelem x2,
1250 const smallfelem y2, const smallfelem z2)
1252 felem ftmp, ftmp2, ftmp3, ftmp4, ftmp5, ftmp6, x_out, y_out, z_out;
1253 longfelem tmp, tmp2;
1254 smallfelem small1, small2, small3, small4, small5;
1255 limb x_equal, y_equal, z1_is_zero, z2_is_zero;
1257 felem_shrink(small3, z1);
1259 z1_is_zero = smallfelem_is_zero(small3);
1260 z2_is_zero = smallfelem_is_zero(z2);
1262 /* ftmp = z1z1 = z1**2 */
1263 smallfelem_square(tmp, small3);
1264 felem_reduce(ftmp, tmp);
1265 /* ftmp[i] < 2^101 */
1266 felem_shrink(small1, ftmp);
1269 /* ftmp2 = z2z2 = z2**2 */
1270 smallfelem_square(tmp, z2);
1271 felem_reduce(ftmp2, tmp);
1272 /* ftmp2[i] < 2^101 */
1273 felem_shrink(small2, ftmp2);
1275 felem_shrink(small5, x1);
1277 /* u1 = ftmp3 = x1*z2z2 */
1278 smallfelem_mul(tmp, small5, small2);
1279 felem_reduce(ftmp3, tmp);
1280 /* ftmp3[i] < 2^101 */
1282 /* ftmp5 = z1 + z2 */
1283 felem_assign(ftmp5, z1);
1284 felem_small_sum(ftmp5, z2);
1285 /* ftmp5[i] < 2^107 */
1287 /* ftmp5 = (z1 + z2)**2 - (z1z1 + z2z2) = 2z1z2 */
1288 felem_square(tmp, ftmp5);
1289 felem_reduce(ftmp5, tmp);
1290 /* ftmp2 = z2z2 + z1z1 */
1291 felem_sum(ftmp2, ftmp);
1292 /* ftmp2[i] < 2^101 + 2^101 = 2^102 */
1293 felem_diff(ftmp5, ftmp2);
1294 /* ftmp5[i] < 2^105 + 2^101 < 2^106 */
1296 /* ftmp2 = z2 * z2z2 */
1297 smallfelem_mul(tmp, small2, z2);
1298 felem_reduce(ftmp2, tmp);
1300 /* s1 = ftmp2 = y1 * z2**3 */
1301 felem_mul(tmp, y1, ftmp2);
1302 felem_reduce(ftmp6, tmp);
1303 /* ftmp6[i] < 2^101 */
1306 * We'll assume z2 = 1 (special case z2 = 0 is handled later)
1309 /* u1 = ftmp3 = x1*z2z2 */
1310 felem_assign(ftmp3, x1);
1311 /* ftmp3[i] < 2^106 */
1314 felem_assign(ftmp5, z1);
1315 felem_scalar(ftmp5, 2);
1316 /* ftmp5[i] < 2*2^106 = 2^107 */
1318 /* s1 = ftmp2 = y1 * z2**3 */
1319 felem_assign(ftmp6, y1);
1320 /* ftmp6[i] < 2^106 */
1324 smallfelem_mul(tmp, x2, small1);
1325 felem_reduce(ftmp4, tmp);
1327 /* h = ftmp4 = u2 - u1 */
1328 felem_diff_zero107(ftmp4, ftmp3);
1329 /* ftmp4[i] < 2^107 + 2^101 < 2^108 */
1330 felem_shrink(small4, ftmp4);
1332 x_equal = smallfelem_is_zero(small4);
1334 /* z_out = ftmp5 * h */
1335 felem_small_mul(tmp, small4, ftmp5);
1336 felem_reduce(z_out, tmp);
1337 /* z_out[i] < 2^101 */
1339 /* ftmp = z1 * z1z1 */
1340 smallfelem_mul(tmp, small1, small3);
1341 felem_reduce(ftmp, tmp);
1343 /* s2 = tmp = y2 * z1**3 */
1344 felem_small_mul(tmp, y2, ftmp);
1345 felem_reduce(ftmp5, tmp);
1347 /* r = ftmp5 = (s2 - s1)*2 */
1348 felem_diff_zero107(ftmp5, ftmp6);
1349 /* ftmp5[i] < 2^107 + 2^107 = 2^108 */
1350 felem_scalar(ftmp5, 2);
1351 /* ftmp5[i] < 2^109 */
1352 felem_shrink(small1, ftmp5);
1353 y_equal = smallfelem_is_zero(small1);
1355 if (x_equal && y_equal && !z1_is_zero && !z2_is_zero) {
1356 point_double(x3, y3, z3, x1, y1, z1);
1360 /* I = ftmp = (2h)**2 */
1361 felem_assign(ftmp, ftmp4);
1362 felem_scalar(ftmp, 2);
1363 /* ftmp[i] < 2*2^108 = 2^109 */
1364 felem_square(tmp, ftmp);
1365 felem_reduce(ftmp, tmp);
1367 /* J = ftmp2 = h * I */
1368 felem_mul(tmp, ftmp4, ftmp);
1369 felem_reduce(ftmp2, tmp);
1371 /* V = ftmp4 = U1 * I */
1372 felem_mul(tmp, ftmp3, ftmp);
1373 felem_reduce(ftmp4, tmp);
1375 /* x_out = r**2 - J - 2V */
1376 smallfelem_square(tmp, small1);
1377 felem_reduce(x_out, tmp);
1378 felem_assign(ftmp3, ftmp4);
1379 felem_scalar(ftmp4, 2);
1380 felem_sum(ftmp4, ftmp2);
1381 /* ftmp4[i] < 2*2^101 + 2^101 < 2^103 */
1382 felem_diff(x_out, ftmp4);
1383 /* x_out[i] < 2^105 + 2^101 */
1385 /* y_out = r(V-x_out) - 2 * s1 * J */
1386 felem_diff_zero107(ftmp3, x_out);
1387 /* ftmp3[i] < 2^107 + 2^101 < 2^108 */
1388 felem_small_mul(tmp, small1, ftmp3);
1389 felem_mul(tmp2, ftmp6, ftmp2);
1390 longfelem_scalar(tmp2, 2);
1391 /* tmp2[i] < 2*2^67 = 2^68 */
1392 longfelem_diff(tmp, tmp2);
1393 /* tmp[i] < 2^67 + 2^70 + 2^40 < 2^71 */
1394 felem_reduce_zero105(y_out, tmp);
1395 /* y_out[i] < 2^106 */
1397 copy_small_conditional(x_out, x2, z1_is_zero);
1398 copy_conditional(x_out, x1, z2_is_zero);
1399 copy_small_conditional(y_out, y2, z1_is_zero);
1400 copy_conditional(y_out, y1, z2_is_zero);
1401 copy_small_conditional(z_out, z2, z1_is_zero);
1402 copy_conditional(z_out, z1, z2_is_zero);
1403 felem_assign(x3, x_out);
1404 felem_assign(y3, y_out);
1405 felem_assign(z3, z_out);
1409 * point_add_small is the same as point_add, except that it operates on
1412 static void point_add_small(smallfelem x3, smallfelem y3, smallfelem z3,
1413 smallfelem x1, smallfelem y1, smallfelem z1,
1414 smallfelem x2, smallfelem y2, smallfelem z2)
1416 felem felem_x3, felem_y3, felem_z3;
1417 felem felem_x1, felem_y1, felem_z1;
1418 smallfelem_expand(felem_x1, x1);
1419 smallfelem_expand(felem_y1, y1);
1420 smallfelem_expand(felem_z1, z1);
1421 point_add(felem_x3, felem_y3, felem_z3, felem_x1, felem_y1, felem_z1, 0,
1423 felem_shrink(x3, felem_x3);
1424 felem_shrink(y3, felem_y3);
1425 felem_shrink(z3, felem_z3);
1429 * Base point pre computation
1430 * --------------------------
1432 * Two different sorts of precomputed tables are used in the following code.
1433 * Each contain various points on the curve, where each point is three field
1434 * elements (x, y, z).
1436 * For the base point table, z is usually 1 (0 for the point at infinity).
1437 * This table has 2 * 16 elements, starting with the following:
1438 * index | bits | point
1439 * ------+---------+------------------------------
1442 * 2 | 0 0 1 0 | 2^64G
1443 * 3 | 0 0 1 1 | (2^64 + 1)G
1444 * 4 | 0 1 0 0 | 2^128G
1445 * 5 | 0 1 0 1 | (2^128 + 1)G
1446 * 6 | 0 1 1 0 | (2^128 + 2^64)G
1447 * 7 | 0 1 1 1 | (2^128 + 2^64 + 1)G
1448 * 8 | 1 0 0 0 | 2^192G
1449 * 9 | 1 0 0 1 | (2^192 + 1)G
1450 * 10 | 1 0 1 0 | (2^192 + 2^64)G
1451 * 11 | 1 0 1 1 | (2^192 + 2^64 + 1)G
1452 * 12 | 1 1 0 0 | (2^192 + 2^128)G
1453 * 13 | 1 1 0 1 | (2^192 + 2^128 + 1)G
1454 * 14 | 1 1 1 0 | (2^192 + 2^128 + 2^64)G
1455 * 15 | 1 1 1 1 | (2^192 + 2^128 + 2^64 + 1)G
1456 * followed by a copy of this with each element multiplied by 2^32.
1458 * The reason for this is so that we can clock bits into four different
1459 * locations when doing simple scalar multiplies against the base point,
1460 * and then another four locations using the second 16 elements.
1462 * Tables for other points have table[i] = iG for i in 0 .. 16. */
1464 /* gmul is the table of precomputed base points */
1465 static const smallfelem gmul[2][16][3] = {
1469 {{0xf4a13945d898c296, 0x77037d812deb33a0, 0xf8bce6e563a440f2,
1470 0x6b17d1f2e12c4247},
1471 {0xcbb6406837bf51f5, 0x2bce33576b315ece, 0x8ee7eb4a7c0f9e16,
1472 0x4fe342e2fe1a7f9b},
1474 {{0x90e75cb48e14db63, 0x29493baaad651f7e, 0x8492592e326e25de,
1475 0x0fa822bc2811aaa5},
1476 {0xe41124545f462ee7, 0x34b1a65050fe82f5, 0x6f4ad4bcb3df188b,
1477 0xbff44ae8f5dba80d},
1479 {{0x93391ce2097992af, 0xe96c98fd0d35f1fa, 0xb257c0de95e02789,
1480 0x300a4bbc89d6726f},
1481 {0xaa54a291c08127a0, 0x5bb1eeada9d806a5, 0x7f1ddb25ff1e3c6f,
1482 0x72aac7e0d09b4644},
1484 {{0x57c84fc9d789bd85, 0xfc35ff7dc297eac3, 0xfb982fd588c6766e,
1485 0x447d739beedb5e67},
1486 {0x0c7e33c972e25b32, 0x3d349b95a7fae500, 0xe12e9d953a4aaff7,
1487 0x2d4825ab834131ee},
1489 {{0x13949c932a1d367f, 0xef7fbd2b1a0a11b7, 0xddc6068bb91dfc60,
1490 0xef9519328a9c72ff},
1491 {0x196035a77376d8a8, 0x23183b0895ca1740, 0xc1ee9807022c219c,
1492 0x611e9fc37dbb2c9b},
1494 {{0xcae2b1920b57f4bc, 0x2936df5ec6c9bc36, 0x7dea6482e11238bf,
1495 0x550663797b51f5d8},
1496 {0x44ffe216348a964c, 0x9fb3d576dbdefbe1, 0x0afa40018d9d50e5,
1497 0x157164848aecb851},
1499 {{0xe48ecafffc5cde01, 0x7ccd84e70d715f26, 0xa2e8f483f43e4391,
1500 0xeb5d7745b21141ea},
1501 {0xcac917e2731a3479, 0x85f22cfe2844b645, 0x0990e6a158006cee,
1502 0xeafd72ebdbecc17b},
1504 {{0x6cf20ffb313728be, 0x96439591a3c6b94a, 0x2736ff8344315fc5,
1505 0xa6d39677a7849276},
1506 {0xf2bab833c357f5f4, 0x824a920c2284059b, 0x66b8babd2d27ecdf,
1507 0x674f84749b0b8816},
1509 {{0x2df48c04677c8a3e, 0x74e02f080203a56b, 0x31855f7db8c7fedb,
1510 0x4e769e7672c9ddad},
1511 {0xa4c36165b824bbb0, 0xfb9ae16f3b9122a5, 0x1ec0057206947281,
1512 0x42b99082de830663},
1514 {{0x6ef95150dda868b9, 0xd1f89e799c0ce131, 0x7fdc1ca008a1c478,
1515 0x78878ef61c6ce04d},
1516 {0x9c62b9121fe0d976, 0x6ace570ebde08d4f, 0xde53142c12309def,
1517 0xb6cb3f5d7b72c321},
1519 {{0x7f991ed2c31a3573, 0x5b82dd5bd54fb496, 0x595c5220812ffcae,
1520 0x0c88bc4d716b1287},
1521 {0x3a57bf635f48aca8, 0x7c8181f4df2564f3, 0x18d1b5b39c04e6aa,
1522 0xdd5ddea3f3901dc6},
1524 {{0xe96a79fb3e72ad0c, 0x43a0a28c42ba792f, 0xefe0a423083e49f3,
1525 0x68f344af6b317466},
1526 {0xcdfe17db3fb24d4a, 0x668bfc2271f5c626, 0x604ed93c24d67ff3,
1527 0x31b9c405f8540a20},
1529 {{0xd36b4789a2582e7f, 0x0d1a10144ec39c28, 0x663c62c3edbad7a0,
1530 0x4052bf4b6f461db9},
1531 {0x235a27c3188d25eb, 0xe724f33999bfcc5b, 0x862be6bd71d70cc8,
1532 0xfecf4d5190b0fc61},
1534 {{0x74346c10a1d4cfac, 0xafdf5cc08526a7a4, 0x123202a8f62bff7a,
1535 0x1eddbae2c802e41a},
1536 {0x8fa0af2dd603f844, 0x36e06b7e4c701917, 0x0c45f45273db33a0,
1537 0x43104d86560ebcfc},
1539 {{0x9615b5110d1d78e5, 0x66b0de3225c4744b, 0x0a4a46fb6aaf363a,
1540 0xb48e26b484f7a21c},
1541 {0x06ebb0f621a01b2d, 0xc004e4048b7b0f98, 0x64131bcdfed6f668,
1542 0xfac015404d4d3dab},
1547 {{0x3a5a9e22185a5943, 0x1ab919365c65dfb6, 0x21656b32262c71da,
1548 0x7fe36b40af22af89},
1549 {0xd50d152c699ca101, 0x74b3d5867b8af212, 0x9f09f40407dca6f1,
1550 0xe697d45825b63624},
1552 {{0xa84aa9397512218e, 0xe9a521b074ca0141, 0x57880b3a18a2e902,
1553 0x4a5b506612a677a6},
1554 {0x0beada7a4c4f3840, 0x626db15419e26d9d, 0xc42604fbe1627d40,
1555 0xeb13461ceac089f1},
1557 {{0xf9faed0927a43281, 0x5e52c4144103ecbc, 0xc342967aa815c857,
1558 0x0781b8291c6a220a},
1559 {0x5a8343ceeac55f80, 0x88f80eeee54a05e3, 0x97b2a14f12916434,
1560 0x690cde8df0151593},
1562 {{0xaee9c75df7f82f2a, 0x9e4c35874afdf43a, 0xf5622df437371326,
1563 0x8a535f566ec73617},
1564 {0xc5f9a0ac223094b7, 0xcde533864c8c7669, 0x37e02819085a92bf,
1565 0x0455c08468b08bd7},
1567 {{0x0c0a6e2c9477b5d9, 0xf9a4bf62876dc444, 0x5050a949b6cdc279,
1568 0x06bada7ab77f8276},
1569 {0xc8b4aed1ea48dac9, 0xdebd8a4b7ea1070f, 0x427d49101366eb70,
1570 0x5b476dfd0e6cb18a},
1572 {{0x7c5c3e44278c340a, 0x4d54606812d66f3b, 0x29a751b1ae23c5d8,
1573 0x3e29864e8a2ec908},
1574 {0x142d2a6626dbb850, 0xad1744c4765bd780, 0x1f150e68e322d1ed,
1575 0x239b90ea3dc31e7e},
1577 {{0x78c416527a53322a, 0x305dde6709776f8e, 0xdbcab759f8862ed4,
1578 0x820f4dd949f72ff7},
1579 {0x6cc544a62b5debd4, 0x75be5d937b4e8cc4, 0x1b481b1b215c14d3,
1580 0x140406ec783a05ec},
1582 {{0x6a703f10e895df07, 0xfd75f3fa01876bd8, 0xeb5b06e70ce08ffe,
1583 0x68f6b8542783dfee},
1584 {0x90c76f8a78712655, 0xcf5293d2f310bf7f, 0xfbc8044dfda45028,
1585 0xcbe1feba92e40ce6},
1587 {{0xe998ceea4396e4c1, 0xfc82ef0b6acea274, 0x230f729f2250e927,
1588 0xd0b2f94d2f420109},
1589 {0x4305adddb38d4966, 0x10b838f8624c3b45, 0x7db2636658954e7a,
1590 0x971459828b0719e5},
1592 {{0x4bd6b72623369fc9, 0x57f2929e53d0b876, 0xc2d5cba4f2340687,
1593 0x961610004a866aba},
1594 {0x49997bcd2e407a5e, 0x69ab197d92ddcb24, 0x2cf1f2438fe5131c,
1595 0x7acb9fadcee75e44},
1597 {{0x254e839423d2d4c0, 0xf57f0c917aea685b, 0xa60d880f6f75aaea,
1598 0x24eb9acca333bf5b},
1599 {0xe3de4ccb1cda5dea, 0xfeef9341c51a6b4f, 0x743125f88bac4c4d,
1600 0x69f891c5acd079cc},
1602 {{0xeee44b35702476b5, 0x7ed031a0e45c2258, 0xb422d1e7bd6f8514,
1603 0xe51f547c5972a107},
1604 {0xa25bcd6fc9cf343d, 0x8ca922ee097c184e, 0xa62f98b3a9fe9a06,
1605 0x1c309a2b25bb1387},
1607 {{0x9295dbeb1967c459, 0xb00148833472c98e, 0xc504977708011828,
1608 0x20b87b8aa2c4e503},
1609 {0x3063175de057c277, 0x1bd539338fe582dd, 0x0d11adef5f69a044,
1610 0xf5c6fa49919776be},
1612 {{0x8c944e760fd59e11, 0x3876cba1102fad5f, 0xa454c3fad83faa56,
1613 0x1ed7d1b9332010b9},
1614 {0xa1011a270024b889, 0x05e4d0dcac0cd344, 0x52b520f0eb6a2a24,
1615 0x3a2b03f03217257a},
1617 {{0xf20fc2afdf1d043d, 0xf330240db58d5a62, 0xfc7d229ca0058c3b,
1618 0x15fee545c78dd9f6},
1619 {0x501e82885bc98cda, 0x41ef80e5d046ac04, 0x557d9f49461210fb,
1620 0x4ab5b6b2b8753f81},
1625 * select_point selects the |idx|th point from a precomputation table and
1628 static void select_point(const u64 idx, unsigned int size,
1629 const smallfelem pre_comp[16][3], smallfelem out[3])
1632 u64 *outlimbs = &out[0][0];
1634 memset(out, 0, sizeof(*out) * 3);
1636 for (i = 0; i < size; i++) {
1637 const u64 *inlimbs = (u64 *)&pre_comp[i][0][0];
1644 for (j = 0; j < NLIMBS * 3; j++)
1645 outlimbs[j] |= inlimbs[j] & mask;
1649 /* get_bit returns the |i|th bit in |in| */
1650 static char get_bit(const felem_bytearray in, int i)
1652 if ((i < 0) || (i >= 256))
1654 return (in[i >> 3] >> (i & 7)) & 1;
1658 * Interleaved point multiplication using precomputed point multiples: The
1659 * small point multiples 0*P, 1*P, ..., 17*P are in pre_comp[], the scalars
1660 * in scalars[]. If g_scalar is non-NULL, we also add this multiple of the
1661 * generator, using certain (large) precomputed multiples in g_pre_comp.
1662 * Output point (X, Y, Z) is stored in x_out, y_out, z_out
1664 static void batch_mul(felem x_out, felem y_out, felem z_out,
1665 const felem_bytearray scalars[],
1666 const unsigned num_points, const u8 *g_scalar,
1667 const int mixed, const smallfelem pre_comp[][17][3],
1668 const smallfelem g_pre_comp[2][16][3])
1671 unsigned num, gen_mul = (g_scalar != NULL);
1677 /* set nq to the point at infinity */
1678 memset(nq, 0, sizeof(nq));
1681 * Loop over all scalars msb-to-lsb, interleaving additions of multiples
1682 * of the generator (two in each of the last 32 rounds) and additions of
1683 * other points multiples (every 5th round).
1685 skip = 1; /* save two point operations in the first
1687 for (i = (num_points ? 255 : 31); i >= 0; --i) {
1690 point_double(nq[0], nq[1], nq[2], nq[0], nq[1], nq[2]);
1692 /* add multiples of the generator */
1693 if (gen_mul && (i <= 31)) {
1694 /* first, look 32 bits upwards */
1695 bits = get_bit(g_scalar, i + 224) << 3;
1696 bits |= get_bit(g_scalar, i + 160) << 2;
1697 bits |= get_bit(g_scalar, i + 96) << 1;
1698 bits |= get_bit(g_scalar, i + 32);
1699 /* select the point to add, in constant time */
1700 select_point(bits, 16, g_pre_comp[1], tmp);
1703 /* Arg 1 below is for "mixed" */
1704 point_add(nq[0], nq[1], nq[2],
1705 nq[0], nq[1], nq[2], 1, tmp[0], tmp[1], tmp[2]);
1707 smallfelem_expand(nq[0], tmp[0]);
1708 smallfelem_expand(nq[1], tmp[1]);
1709 smallfelem_expand(nq[2], tmp[2]);
1713 /* second, look at the current position */
1714 bits = get_bit(g_scalar, i + 192) << 3;
1715 bits |= get_bit(g_scalar, i + 128) << 2;
1716 bits |= get_bit(g_scalar, i + 64) << 1;
1717 bits |= get_bit(g_scalar, i);
1718 /* select the point to add, in constant time */
1719 select_point(bits, 16, g_pre_comp[0], tmp);
1720 /* Arg 1 below is for "mixed" */
1721 point_add(nq[0], nq[1], nq[2],
1722 nq[0], nq[1], nq[2], 1, tmp[0], tmp[1], tmp[2]);
1725 /* do other additions every 5 doublings */
1726 if (num_points && (i % 5 == 0)) {
1727 /* loop over all scalars */
1728 for (num = 0; num < num_points; ++num) {
1729 bits = get_bit(scalars[num], i + 4) << 5;
1730 bits |= get_bit(scalars[num], i + 3) << 4;
1731 bits |= get_bit(scalars[num], i + 2) << 3;
1732 bits |= get_bit(scalars[num], i + 1) << 2;
1733 bits |= get_bit(scalars[num], i) << 1;
1734 bits |= get_bit(scalars[num], i - 1);
1735 ec_GFp_nistp_recode_scalar_bits(&sign, &digit, bits);
1738 * select the point to add or subtract, in constant time
1740 select_point(digit, 17, pre_comp[num], tmp);
1741 smallfelem_neg(ftmp, tmp[1]); /* (X, -Y, Z) is the negative
1743 copy_small_conditional(ftmp, tmp[1], (((limb) sign) - 1));
1744 felem_contract(tmp[1], ftmp);
1747 point_add(nq[0], nq[1], nq[2],
1748 nq[0], nq[1], nq[2],
1749 mixed, tmp[0], tmp[1], tmp[2]);
1751 smallfelem_expand(nq[0], tmp[0]);
1752 smallfelem_expand(nq[1], tmp[1]);
1753 smallfelem_expand(nq[2], tmp[2]);
1759 felem_assign(x_out, nq[0]);
1760 felem_assign(y_out, nq[1]);
1761 felem_assign(z_out, nq[2]);
1764 /* Precomputation for the group generator. */
1765 struct nistp256_pre_comp_st {
1766 smallfelem g_pre_comp[2][16][3];
1767 CRYPTO_REF_COUNT references;
1768 CRYPTO_RWLOCK *lock;
1771 const EC_METHOD *EC_GFp_nistp256_method(void)
1773 static const EC_METHOD ret = {
1774 EC_FLAGS_DEFAULT_OCT,
1775 NID_X9_62_prime_field,
1776 ec_GFp_nistp256_group_init,
1777 ec_GFp_simple_group_finish,
1778 ec_GFp_simple_group_clear_finish,
1779 ec_GFp_nist_group_copy,
1780 ec_GFp_nistp256_group_set_curve,
1781 ec_GFp_simple_group_get_curve,
1782 ec_GFp_simple_group_get_degree,
1783 ec_group_simple_order_bits,
1784 ec_GFp_simple_group_check_discriminant,
1785 ec_GFp_simple_point_init,
1786 ec_GFp_simple_point_finish,
1787 ec_GFp_simple_point_clear_finish,
1788 ec_GFp_simple_point_copy,
1789 ec_GFp_simple_point_set_to_infinity,
1790 ec_GFp_simple_set_Jprojective_coordinates_GFp,
1791 ec_GFp_simple_get_Jprojective_coordinates_GFp,
1792 ec_GFp_simple_point_set_affine_coordinates,
1793 ec_GFp_nistp256_point_get_affine_coordinates,
1794 0 /* point_set_compressed_coordinates */ ,
1799 ec_GFp_simple_invert,
1800 ec_GFp_simple_is_at_infinity,
1801 ec_GFp_simple_is_on_curve,
1803 ec_GFp_simple_make_affine,
1804 ec_GFp_simple_points_make_affine,
1805 ec_GFp_nistp256_points_mul,
1806 ec_GFp_nistp256_precompute_mult,
1807 ec_GFp_nistp256_have_precompute_mult,
1808 ec_GFp_nist_field_mul,
1809 ec_GFp_nist_field_sqr,
1811 ec_GFp_simple_field_inv,
1812 0 /* field_encode */ ,
1813 0 /* field_decode */ ,
1814 0, /* field_set_to_one */
1815 ec_key_simple_priv2oct,
1816 ec_key_simple_oct2priv,
1817 0, /* set private */
1818 ec_key_simple_generate_key,
1819 ec_key_simple_check_key,
1820 ec_key_simple_generate_public_key,
1823 ecdh_simple_compute_key,
1824 ecdsa_simple_sign_setup,
1825 ecdsa_simple_sign_sig,
1826 ecdsa_simple_verify_sig,
1827 0, /* field_inverse_mod_ord */
1828 0, /* blind_coordinates */
1830 0, /* ladder_step */
1837 /******************************************************************************/
1839 * FUNCTIONS TO MANAGE PRECOMPUTATION
1842 static NISTP256_PRE_COMP *nistp256_pre_comp_new(void)
1844 NISTP256_PRE_COMP *ret = OPENSSL_zalloc(sizeof(*ret));
1847 ECerr(EC_F_NISTP256_PRE_COMP_NEW, ERR_R_MALLOC_FAILURE);
1851 ret->references = 1;
1853 ret->lock = CRYPTO_THREAD_lock_new();
1854 if (ret->lock == NULL) {
1855 ECerr(EC_F_NISTP256_PRE_COMP_NEW, ERR_R_MALLOC_FAILURE);
1862 NISTP256_PRE_COMP *EC_nistp256_pre_comp_dup(NISTP256_PRE_COMP *p)
1866 CRYPTO_UP_REF(&p->references, &i, p->lock);
1870 void EC_nistp256_pre_comp_free(NISTP256_PRE_COMP *pre)
1877 CRYPTO_DOWN_REF(&pre->references, &i, pre->lock);
1878 REF_PRINT_COUNT("EC_nistp256", x);
1881 REF_ASSERT_ISNT(i < 0);
1883 CRYPTO_THREAD_lock_free(pre->lock);
1887 /******************************************************************************/
1889 * OPENSSL EC_METHOD FUNCTIONS
1892 int ec_GFp_nistp256_group_init(EC_GROUP *group)
1895 ret = ec_GFp_simple_group_init(group);
1896 group->a_is_minus3 = 1;
1900 int ec_GFp_nistp256_group_set_curve(EC_GROUP *group, const BIGNUM *p,
1901 const BIGNUM *a, const BIGNUM *b,
1905 BIGNUM *curve_p, *curve_a, *curve_b;
1907 BN_CTX *new_ctx = NULL;
1910 ctx = new_ctx = BN_CTX_new();
1916 curve_p = BN_CTX_get(ctx);
1917 curve_a = BN_CTX_get(ctx);
1918 curve_b = BN_CTX_get(ctx);
1919 if (curve_b == NULL)
1921 BN_bin2bn(nistp256_curve_params[0], sizeof(felem_bytearray), curve_p);
1922 BN_bin2bn(nistp256_curve_params[1], sizeof(felem_bytearray), curve_a);
1923 BN_bin2bn(nistp256_curve_params[2], sizeof(felem_bytearray), curve_b);
1924 if ((BN_cmp(curve_p, p)) || (BN_cmp(curve_a, a)) || (BN_cmp(curve_b, b))) {
1925 ECerr(EC_F_EC_GFP_NISTP256_GROUP_SET_CURVE,
1926 EC_R_WRONG_CURVE_PARAMETERS);
1929 group->field_mod_func = BN_nist_mod_256;
1930 ret = ec_GFp_simple_group_set_curve(group, p, a, b, ctx);
1934 BN_CTX_free(new_ctx);
1940 * Takes the Jacobian coordinates (X, Y, Z) of a point and returns (X', Y') =
1943 int ec_GFp_nistp256_point_get_affine_coordinates(const EC_GROUP *group,
1944 const EC_POINT *point,
1945 BIGNUM *x, BIGNUM *y,
1948 felem z1, z2, x_in, y_in;
1949 smallfelem x_out, y_out;
1952 if (EC_POINT_is_at_infinity(group, point)) {
1953 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1954 EC_R_POINT_AT_INFINITY);
1957 if ((!BN_to_felem(x_in, point->X)) || (!BN_to_felem(y_in, point->Y)) ||
1958 (!BN_to_felem(z1, point->Z)))
1961 felem_square(tmp, z2);
1962 felem_reduce(z1, tmp);
1963 felem_mul(tmp, x_in, z1);
1964 felem_reduce(x_in, tmp);
1965 felem_contract(x_out, x_in);
1967 if (!smallfelem_to_BN(x, x_out)) {
1968 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1973 felem_mul(tmp, z1, z2);
1974 felem_reduce(z1, tmp);
1975 felem_mul(tmp, y_in, z1);
1976 felem_reduce(y_in, tmp);
1977 felem_contract(y_out, y_in);
1979 if (!smallfelem_to_BN(y, y_out)) {
1980 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1988 /* points below is of size |num|, and tmp_smallfelems is of size |num+1| */
1989 static void make_points_affine(size_t num, smallfelem points[][3],
1990 smallfelem tmp_smallfelems[])
1993 * Runs in constant time, unless an input is the point at infinity (which
1994 * normally shouldn't happen).
1996 ec_GFp_nistp_points_make_affine_internal(num,
2000 (void (*)(void *))smallfelem_one,
2001 smallfelem_is_zero_int,
2002 (void (*)(void *, const void *))
2004 (void (*)(void *, const void *))
2005 smallfelem_square_contract,
2007 (void *, const void *,
2009 smallfelem_mul_contract,
2010 (void (*)(void *, const void *))
2011 smallfelem_inv_contract,
2012 /* nothing to contract */
2013 (void (*)(void *, const void *))
2018 * Computes scalar*generator + \sum scalars[i]*points[i], ignoring NULL
2019 * values Result is stored in r (r can equal one of the inputs).
2021 int ec_GFp_nistp256_points_mul(const EC_GROUP *group, EC_POINT *r,
2022 const BIGNUM *scalar, size_t num,
2023 const EC_POINT *points[],
2024 const BIGNUM *scalars[], BN_CTX *ctx)
2029 BIGNUM *x, *y, *z, *tmp_scalar;
2030 felem_bytearray g_secret;
2031 felem_bytearray *secrets = NULL;
2032 smallfelem (*pre_comp)[17][3] = NULL;
2033 smallfelem *tmp_smallfelems = NULL;
2034 felem_bytearray tmp;
2035 unsigned i, num_bytes;
2036 int have_pre_comp = 0;
2037 size_t num_points = num;
2038 smallfelem x_in, y_in, z_in;
2039 felem x_out, y_out, z_out;
2040 NISTP256_PRE_COMP *pre = NULL;
2041 const smallfelem(*g_pre_comp)[16][3] = NULL;
2042 EC_POINT *generator = NULL;
2043 const EC_POINT *p = NULL;
2044 const BIGNUM *p_scalar = NULL;
2047 x = BN_CTX_get(ctx);
2048 y = BN_CTX_get(ctx);
2049 z = BN_CTX_get(ctx);
2050 tmp_scalar = BN_CTX_get(ctx);
2051 if (tmp_scalar == NULL)
2054 if (scalar != NULL) {
2055 pre = group->pre_comp.nistp256;
2057 /* we have precomputation, try to use it */
2058 g_pre_comp = (const smallfelem(*)[16][3])pre->g_pre_comp;
2060 /* try to use the standard precomputation */
2061 g_pre_comp = &gmul[0];
2062 generator = EC_POINT_new(group);
2063 if (generator == NULL)
2065 /* get the generator from precomputation */
2066 if (!smallfelem_to_BN(x, g_pre_comp[0][1][0]) ||
2067 !smallfelem_to_BN(y, g_pre_comp[0][1][1]) ||
2068 !smallfelem_to_BN(z, g_pre_comp[0][1][2])) {
2069 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2072 if (!EC_POINT_set_Jprojective_coordinates_GFp(group,
2076 if (0 == EC_POINT_cmp(group, generator, group->generator, ctx))
2077 /* precomputation matches generator */
2081 * we don't have valid precomputation: treat the generator as a
2086 if (num_points > 0) {
2087 if (num_points >= 3) {
2089 * unless we precompute multiples for just one or two points,
2090 * converting those into affine form is time well spent
2094 secrets = OPENSSL_malloc(sizeof(*secrets) * num_points);
2095 pre_comp = OPENSSL_malloc(sizeof(*pre_comp) * num_points);
2098 OPENSSL_malloc(sizeof(*tmp_smallfelems) * (num_points * 17 + 1));
2099 if ((secrets == NULL) || (pre_comp == NULL)
2100 || (mixed && (tmp_smallfelems == NULL))) {
2101 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_MALLOC_FAILURE);
2106 * we treat NULL scalars as 0, and NULL points as points at infinity,
2107 * i.e., they contribute nothing to the linear combination
2109 memset(secrets, 0, sizeof(*secrets) * num_points);
2110 memset(pre_comp, 0, sizeof(*pre_comp) * num_points);
2111 for (i = 0; i < num_points; ++i) {
2114 * we didn't have a valid precomputation, so we pick the
2118 p = EC_GROUP_get0_generator(group);
2121 /* the i^th point */
2124 p_scalar = scalars[i];
2126 if ((p_scalar != NULL) && (p != NULL)) {
2127 /* reduce scalar to 0 <= scalar < 2^256 */
2128 if ((BN_num_bits(p_scalar) > 256)
2129 || (BN_is_negative(p_scalar))) {
2131 * this is an unusual input, and we don't guarantee
2134 if (!BN_nnmod(tmp_scalar, p_scalar, group->order, ctx)) {
2135 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2138 num_bytes = BN_bn2binpad(tmp_scalar, tmp, sizeof(tmp));
2140 num_bytes = BN_bn2binpad(p_scalar, tmp, sizeof(tmp));
2141 flip_endian(secrets[i], tmp, num_bytes);
2142 /* precompute multiples */
2143 if ((!BN_to_felem(x_out, p->X)) ||
2144 (!BN_to_felem(y_out, p->Y)) ||
2145 (!BN_to_felem(z_out, p->Z)))
2147 felem_shrink(pre_comp[i][1][0], x_out);
2148 felem_shrink(pre_comp[i][1][1], y_out);
2149 felem_shrink(pre_comp[i][1][2], z_out);
2150 for (j = 2; j <= 16; ++j) {
2152 point_add_small(pre_comp[i][j][0], pre_comp[i][j][1],
2153 pre_comp[i][j][2], pre_comp[i][1][0],
2154 pre_comp[i][1][1], pre_comp[i][1][2],
2155 pre_comp[i][j - 1][0],
2156 pre_comp[i][j - 1][1],
2157 pre_comp[i][j - 1][2]);
2159 point_double_small(pre_comp[i][j][0],
2162 pre_comp[i][j / 2][0],
2163 pre_comp[i][j / 2][1],
2164 pre_comp[i][j / 2][2]);
2170 make_points_affine(num_points * 17, pre_comp[0], tmp_smallfelems);
2173 /* the scalar for the generator */
2174 if ((scalar != NULL) && (have_pre_comp)) {
2175 memset(g_secret, 0, sizeof(g_secret));
2176 /* reduce scalar to 0 <= scalar < 2^256 */
2177 if ((BN_num_bits(scalar) > 256) || (BN_is_negative(scalar))) {
2179 * this is an unusual input, and we don't guarantee
2182 if (!BN_nnmod(tmp_scalar, scalar, group->order, ctx)) {
2183 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2186 num_bytes = BN_bn2binpad(tmp_scalar, tmp, sizeof(tmp));
2188 num_bytes = BN_bn2binpad(scalar, tmp, sizeof(tmp));
2189 flip_endian(g_secret, tmp, num_bytes);
2190 /* do the multiplication with generator precomputation */
2191 batch_mul(x_out, y_out, z_out,
2192 (const felem_bytearray(*))secrets, num_points,
2194 mixed, (const smallfelem(*)[17][3])pre_comp, g_pre_comp);
2196 /* do the multiplication without generator precomputation */
2197 batch_mul(x_out, y_out, z_out,
2198 (const felem_bytearray(*))secrets, num_points,
2199 NULL, mixed, (const smallfelem(*)[17][3])pre_comp, NULL);
2200 /* reduce the output to its unique minimal representation */
2201 felem_contract(x_in, x_out);
2202 felem_contract(y_in, y_out);
2203 felem_contract(z_in, z_out);
2204 if ((!smallfelem_to_BN(x, x_in)) || (!smallfelem_to_BN(y, y_in)) ||
2205 (!smallfelem_to_BN(z, z_in))) {
2206 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2209 ret = EC_POINT_set_Jprojective_coordinates_GFp(group, r, x, y, z, ctx);
2213 EC_POINT_free(generator);
2214 OPENSSL_free(secrets);
2215 OPENSSL_free(pre_comp);
2216 OPENSSL_free(tmp_smallfelems);
2220 int ec_GFp_nistp256_precompute_mult(EC_GROUP *group, BN_CTX *ctx)
2223 NISTP256_PRE_COMP *pre = NULL;
2226 EC_POINT *generator = NULL;
2227 smallfelem tmp_smallfelems[32];
2228 felem x_tmp, y_tmp, z_tmp;
2230 BN_CTX *new_ctx = NULL;
2233 /* throw away old precomputation */
2234 EC_pre_comp_free(group);
2238 ctx = new_ctx = BN_CTX_new();
2244 x = BN_CTX_get(ctx);
2245 y = BN_CTX_get(ctx);
2248 /* get the generator */
2249 if (group->generator == NULL)
2251 generator = EC_POINT_new(group);
2252 if (generator == NULL)
2254 BN_bin2bn(nistp256_curve_params[3], sizeof(felem_bytearray), x);
2255 BN_bin2bn(nistp256_curve_params[4], sizeof(felem_bytearray), y);
2256 if (!EC_POINT_set_affine_coordinates(group, generator, x, y, ctx))
2258 if ((pre = nistp256_pre_comp_new()) == NULL)
2261 * if the generator is the standard one, use built-in precomputation
2263 if (0 == EC_POINT_cmp(group, generator, group->generator, ctx)) {
2264 memcpy(pre->g_pre_comp, gmul, sizeof(pre->g_pre_comp));
2267 if ((!BN_to_felem(x_tmp, group->generator->X)) ||
2268 (!BN_to_felem(y_tmp, group->generator->Y)) ||
2269 (!BN_to_felem(z_tmp, group->generator->Z)))
2271 felem_shrink(pre->g_pre_comp[0][1][0], x_tmp);
2272 felem_shrink(pre->g_pre_comp[0][1][1], y_tmp);
2273 felem_shrink(pre->g_pre_comp[0][1][2], z_tmp);
2275 * compute 2^64*G, 2^128*G, 2^192*G for the first table, 2^32*G, 2^96*G,
2276 * 2^160*G, 2^224*G for the second one
2278 for (i = 1; i <= 8; i <<= 1) {
2279 point_double_small(pre->g_pre_comp[1][i][0], pre->g_pre_comp[1][i][1],
2280 pre->g_pre_comp[1][i][2], pre->g_pre_comp[0][i][0],
2281 pre->g_pre_comp[0][i][1],
2282 pre->g_pre_comp[0][i][2]);
2283 for (j = 0; j < 31; ++j) {
2284 point_double_small(pre->g_pre_comp[1][i][0],
2285 pre->g_pre_comp[1][i][1],
2286 pre->g_pre_comp[1][i][2],
2287 pre->g_pre_comp[1][i][0],
2288 pre->g_pre_comp[1][i][1],
2289 pre->g_pre_comp[1][i][2]);
2293 point_double_small(pre->g_pre_comp[0][2 * i][0],
2294 pre->g_pre_comp[0][2 * i][1],
2295 pre->g_pre_comp[0][2 * i][2],
2296 pre->g_pre_comp[1][i][0], pre->g_pre_comp[1][i][1],
2297 pre->g_pre_comp[1][i][2]);
2298 for (j = 0; j < 31; ++j) {
2299 point_double_small(pre->g_pre_comp[0][2 * i][0],
2300 pre->g_pre_comp[0][2 * i][1],
2301 pre->g_pre_comp[0][2 * i][2],
2302 pre->g_pre_comp[0][2 * i][0],
2303 pre->g_pre_comp[0][2 * i][1],
2304 pre->g_pre_comp[0][2 * i][2]);
2307 for (i = 0; i < 2; i++) {
2308 /* g_pre_comp[i][0] is the point at infinity */
2309 memset(pre->g_pre_comp[i][0], 0, sizeof(pre->g_pre_comp[i][0]));
2310 /* the remaining multiples */
2311 /* 2^64*G + 2^128*G resp. 2^96*G + 2^160*G */
2312 point_add_small(pre->g_pre_comp[i][6][0], pre->g_pre_comp[i][6][1],
2313 pre->g_pre_comp[i][6][2], pre->g_pre_comp[i][4][0],
2314 pre->g_pre_comp[i][4][1], pre->g_pre_comp[i][4][2],
2315 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2316 pre->g_pre_comp[i][2][2]);
2317 /* 2^64*G + 2^192*G resp. 2^96*G + 2^224*G */
2318 point_add_small(pre->g_pre_comp[i][10][0], pre->g_pre_comp[i][10][1],
2319 pre->g_pre_comp[i][10][2], pre->g_pre_comp[i][8][0],
2320 pre->g_pre_comp[i][8][1], pre->g_pre_comp[i][8][2],
2321 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2322 pre->g_pre_comp[i][2][2]);
2323 /* 2^128*G + 2^192*G resp. 2^160*G + 2^224*G */
2324 point_add_small(pre->g_pre_comp[i][12][0], pre->g_pre_comp[i][12][1],
2325 pre->g_pre_comp[i][12][2], pre->g_pre_comp[i][8][0],
2326 pre->g_pre_comp[i][8][1], pre->g_pre_comp[i][8][2],
2327 pre->g_pre_comp[i][4][0], pre->g_pre_comp[i][4][1],
2328 pre->g_pre_comp[i][4][2]);
2330 * 2^64*G + 2^128*G + 2^192*G resp. 2^96*G + 2^160*G + 2^224*G
2332 point_add_small(pre->g_pre_comp[i][14][0], pre->g_pre_comp[i][14][1],
2333 pre->g_pre_comp[i][14][2], pre->g_pre_comp[i][12][0],
2334 pre->g_pre_comp[i][12][1], pre->g_pre_comp[i][12][2],
2335 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2336 pre->g_pre_comp[i][2][2]);
2337 for (j = 1; j < 8; ++j) {
2338 /* odd multiples: add G resp. 2^32*G */
2339 point_add_small(pre->g_pre_comp[i][2 * j + 1][0],
2340 pre->g_pre_comp[i][2 * j + 1][1],
2341 pre->g_pre_comp[i][2 * j + 1][2],
2342 pre->g_pre_comp[i][2 * j][0],
2343 pre->g_pre_comp[i][2 * j][1],
2344 pre->g_pre_comp[i][2 * j][2],
2345 pre->g_pre_comp[i][1][0],
2346 pre->g_pre_comp[i][1][1],
2347 pre->g_pre_comp[i][1][2]);
2350 make_points_affine(31, &(pre->g_pre_comp[0][1]), tmp_smallfelems);
2353 SETPRECOMP(group, nistp256, pre);
2359 EC_POINT_free(generator);
2361 BN_CTX_free(new_ctx);
2363 EC_nistp256_pre_comp_free(pre);
2367 int ec_GFp_nistp256_have_precompute_mult(const EC_GROUP *group)
2369 return HAVEPRECOMP(group, nistp256);