1 /* crypto/ec/ecp_nistp256.c */
3 * Written by Adam Langley (Google) for the OpenSSL project
5 /* Copyright 2011 Google Inc.
7 * Licensed under the Apache License, Version 2.0 (the "License");
9 * you may not use this file except in compliance with the License.
10 * You may obtain a copy of the License at
12 * http://www.apache.org/licenses/LICENSE-2.0
14 * Unless required by applicable law or agreed to in writing, software
15 * distributed under the License is distributed on an "AS IS" BASIS,
16 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17 * See the License for the specific language governing permissions and
18 * limitations under the License.
22 * A 64-bit implementation of the NIST P-256 elliptic curve point multiplication
24 * OpenSSL integration was taken from Emilia Kasper's work in ecp_nistp224.c.
25 * Otherwise based on Emilia's P224 work, which was inspired by my curve25519
26 * work which got its smarts from Daniel J. Bernstein's work on the same.
29 #include <openssl/opensslconf.h>
30 #ifndef OPENSSL_NO_EC_NISTP_64_GCC_128
32 # ifndef OPENSSL_SYS_VMS
35 # include <inttypes.h>
39 # include <openssl/err.h>
42 # if defined(__GNUC__) && (__GNUC__ > 3 || (__GNUC__ == 3 && __GNUC_MINOR__ >= 1))
43 /* even with gcc, the typedef won't work for 32-bit platforms */
44 typedef __uint128_t uint128_t; /* nonstandard; implemented by gcc on 64-bit
46 typedef __int128_t int128_t;
48 # error "Need GCC 3.1 or later to define type uint128_t"
57 * The underlying field. P256 operates over GF(2^256-2^224+2^192+2^96-1). We
58 * can serialise an element of this field into 32 bytes. We call this an
62 typedef u8 felem_bytearray[32];
65 * These are the parameters of P256, taken from FIPS 186-3, page 86. These
66 * values are big-endian.
68 static const felem_bytearray nistp256_curve_params[5] = {
69 {0xff, 0xff, 0xff, 0xff, 0x00, 0x00, 0x00, 0x01, /* p */
70 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
71 0x00, 0x00, 0x00, 0x00, 0xff, 0xff, 0xff, 0xff,
72 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff},
73 {0xff, 0xff, 0xff, 0xff, 0x00, 0x00, 0x00, 0x01, /* a = -3 */
74 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
75 0x00, 0x00, 0x00, 0x00, 0xff, 0xff, 0xff, 0xff,
76 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xfc}, /* b */
77 {0x5a, 0xc6, 0x35, 0xd8, 0xaa, 0x3a, 0x93, 0xe7,
78 0xb3, 0xeb, 0xbd, 0x55, 0x76, 0x98, 0x86, 0xbc,
79 0x65, 0x1d, 0x06, 0xb0, 0xcc, 0x53, 0xb0, 0xf6,
80 0x3b, 0xce, 0x3c, 0x3e, 0x27, 0xd2, 0x60, 0x4b},
81 {0x6b, 0x17, 0xd1, 0xf2, 0xe1, 0x2c, 0x42, 0x47, /* x */
82 0xf8, 0xbc, 0xe6, 0xe5, 0x63, 0xa4, 0x40, 0xf2,
83 0x77, 0x03, 0x7d, 0x81, 0x2d, 0xeb, 0x33, 0xa0,
84 0xf4, 0xa1, 0x39, 0x45, 0xd8, 0x98, 0xc2, 0x96},
85 {0x4f, 0xe3, 0x42, 0xe2, 0xfe, 0x1a, 0x7f, 0x9b, /* y */
86 0x8e, 0xe7, 0xeb, 0x4a, 0x7c, 0x0f, 0x9e, 0x16,
87 0x2b, 0xce, 0x33, 0x57, 0x6b, 0x31, 0x5e, 0xce,
88 0xcb, 0xb6, 0x40, 0x68, 0x37, 0xbf, 0x51, 0xf5}
92 * The representation of field elements.
93 * ------------------------------------
95 * We represent field elements with either four 128-bit values, eight 128-bit
96 * values, or four 64-bit values. The field element represented is:
97 * v[0]*2^0 + v[1]*2^64 + v[2]*2^128 + v[3]*2^192 (mod p)
99 * v[0]*2^0 + v[1]*2^64 + v[2]*2^128 + ... + v[8]*2^512 (mod p)
101 * 128-bit values are called 'limbs'. Since the limbs are spaced only 64 bits
102 * apart, but are 128-bits wide, the most significant bits of each limb overlap
103 * with the least significant bits of the next.
105 * A field element with four limbs is an 'felem'. One with eight limbs is a
108 * A field element with four, 64-bit values is called a 'smallfelem'. Small
109 * values are used as intermediate values before multiplication.
114 typedef uint128_t limb;
115 typedef limb felem[NLIMBS];
116 typedef limb longfelem[NLIMBS * 2];
117 typedef u64 smallfelem[NLIMBS];
119 /* This is the value of the prime as four 64-bit words, little-endian. */
120 static const u64 kPrime[4] =
121 { 0xfffffffffffffffful, 0xffffffff, 0, 0xffffffff00000001ul };
122 static const u64 bottom63bits = 0x7ffffffffffffffful;
125 * bin32_to_felem takes a little-endian byte array and converts it into felem
126 * form. This assumes that the CPU is little-endian.
128 static void bin32_to_felem(felem out, const u8 in[32])
130 out[0] = *((u64 *)&in[0]);
131 out[1] = *((u64 *)&in[8]);
132 out[2] = *((u64 *)&in[16]);
133 out[3] = *((u64 *)&in[24]);
137 * smallfelem_to_bin32 takes a smallfelem and serialises into a little
138 * endian, 32 byte array. This assumes that the CPU is little-endian.
140 static void smallfelem_to_bin32(u8 out[32], const smallfelem in)
142 *((u64 *)&out[0]) = in[0];
143 *((u64 *)&out[8]) = in[1];
144 *((u64 *)&out[16]) = in[2];
145 *((u64 *)&out[24]) = in[3];
148 /* To preserve endianness when using BN_bn2bin and BN_bin2bn */
149 static void flip_endian(u8 *out, const u8 *in, unsigned len)
152 for (i = 0; i < len; ++i)
153 out[i] = in[len - 1 - i];
156 /* BN_to_felem converts an OpenSSL BIGNUM into an felem */
157 static int BN_to_felem(felem out, const BIGNUM *bn)
159 felem_bytearray b_in;
160 felem_bytearray b_out;
163 /* BN_bn2bin eats leading zeroes */
164 memset(b_out, 0, sizeof b_out);
165 num_bytes = BN_num_bytes(bn);
166 if (num_bytes > sizeof b_out) {
167 ECerr(EC_F_BN_TO_FELEM, EC_R_BIGNUM_OUT_OF_RANGE);
170 if (BN_is_negative(bn)) {
171 ECerr(EC_F_BN_TO_FELEM, EC_R_BIGNUM_OUT_OF_RANGE);
174 num_bytes = BN_bn2bin(bn, b_in);
175 flip_endian(b_out, b_in, num_bytes);
176 bin32_to_felem(out, b_out);
180 /* felem_to_BN converts an felem into an OpenSSL BIGNUM */
181 static BIGNUM *smallfelem_to_BN(BIGNUM *out, const smallfelem in)
183 felem_bytearray b_in, b_out;
184 smallfelem_to_bin32(b_in, in);
185 flip_endian(b_out, b_in, sizeof b_out);
186 return BN_bin2bn(b_out, sizeof b_out, out);
194 static void smallfelem_one(smallfelem out)
202 static void smallfelem_assign(smallfelem out, const smallfelem in)
210 static void felem_assign(felem out, const felem in)
218 /* felem_sum sets out = out + in. */
219 static void felem_sum(felem out, const felem in)
227 /* felem_small_sum sets out = out + in. */
228 static void felem_small_sum(felem out, const smallfelem in)
236 /* felem_scalar sets out = out * scalar */
237 static void felem_scalar(felem out, const u64 scalar)
245 /* longfelem_scalar sets out = out * scalar */
246 static void longfelem_scalar(longfelem out, const u64 scalar)
258 # define two105m41m9 (((limb)1) << 105) - (((limb)1) << 41) - (((limb)1) << 9)
259 # define two105 (((limb)1) << 105)
260 # define two105m41p9 (((limb)1) << 105) - (((limb)1) << 41) + (((limb)1) << 9)
262 /* zero105 is 0 mod p */
263 static const felem zero105 =
264 { two105m41m9, two105, two105m41p9, two105m41p9 };
267 * smallfelem_neg sets |out| to |-small|
269 * out[i] < out[i] + 2^105
271 static void smallfelem_neg(felem out, const smallfelem small)
273 /* In order to prevent underflow, we subtract from 0 mod p. */
274 out[0] = zero105[0] - small[0];
275 out[1] = zero105[1] - small[1];
276 out[2] = zero105[2] - small[2];
277 out[3] = zero105[3] - small[3];
281 * felem_diff subtracts |in| from |out|
285 * out[i] < out[i] + 2^105
287 static void felem_diff(felem out, const felem in)
290 * In order to prevent underflow, we add 0 mod p before subtracting.
292 out[0] += zero105[0];
293 out[1] += zero105[1];
294 out[2] += zero105[2];
295 out[3] += zero105[3];
303 # define two107m43m11 (((limb)1) << 107) - (((limb)1) << 43) - (((limb)1) << 11)
304 # define two107 (((limb)1) << 107)
305 # define two107m43p11 (((limb)1) << 107) - (((limb)1) << 43) + (((limb)1) << 11)
307 /* zero107 is 0 mod p */
308 static const felem zero107 =
309 { two107m43m11, two107, two107m43p11, two107m43p11 };
312 * An alternative felem_diff for larger inputs |in|
313 * felem_diff_zero107 subtracts |in| from |out|
317 * out[i] < out[i] + 2^107
319 static void felem_diff_zero107(felem out, const felem in)
322 * In order to prevent underflow, we add 0 mod p before subtracting.
324 out[0] += zero107[0];
325 out[1] += zero107[1];
326 out[2] += zero107[2];
327 out[3] += zero107[3];
336 * longfelem_diff subtracts |in| from |out|
340 * out[i] < out[i] + 2^70 + 2^40
342 static void longfelem_diff(longfelem out, const longfelem in)
344 static const limb two70m8p6 =
345 (((limb) 1) << 70) - (((limb) 1) << 8) + (((limb) 1) << 6);
346 static const limb two70p40 = (((limb) 1) << 70) + (((limb) 1) << 40);
347 static const limb two70 = (((limb) 1) << 70);
348 static const limb two70m40m38p6 =
349 (((limb) 1) << 70) - (((limb) 1) << 40) - (((limb) 1) << 38) +
351 static const limb two70m6 = (((limb) 1) << 70) - (((limb) 1) << 6);
353 /* add 0 mod p to avoid underflow */
357 out[3] += two70m40m38p6;
363 /* in[i] < 7*2^67 < 2^70 - 2^40 - 2^38 + 2^6 */
374 # define two64m0 (((limb)1) << 64) - 1
375 # define two110p32m0 (((limb)1) << 110) + (((limb)1) << 32) - 1
376 # define two64m46 (((limb)1) << 64) - (((limb)1) << 46)
377 # define two64m32 (((limb)1) << 64) - (((limb)1) << 32)
379 /* zero110 is 0 mod p */
380 static const felem zero110 = { two64m0, two110p32m0, two64m46, two64m32 };
383 * felem_shrink converts an felem into a smallfelem. The result isn't quite
384 * minimal as the value may be greater than p.
391 static void felem_shrink(smallfelem out, const felem in)
396 static const u64 kPrime3Test = 0x7fffffff00000001ul; /* 2^63 - 2^32 + 1 */
399 tmp[3] = zero110[3] + in[3] + ((u64)(in[2] >> 64));
402 tmp[2] = zero110[2] + (u64)in[2];
403 tmp[0] = zero110[0] + in[0];
404 tmp[1] = zero110[1] + in[1];
405 /* tmp[0] < 2**110, tmp[1] < 2^111, tmp[2] < 2**65 */
408 * We perform two partial reductions where we eliminate the high-word of
409 * tmp[3]. We don't update the other words till the end.
411 a = tmp[3] >> 64; /* a < 2^46 */
412 tmp[3] = (u64)tmp[3];
414 tmp[3] += ((limb) a) << 32;
418 a = tmp[3] >> 64; /* a < 2^15 */
419 b += a; /* b < 2^46 + 2^15 < 2^47 */
420 tmp[3] = (u64)tmp[3];
422 tmp[3] += ((limb) a) << 32;
423 /* tmp[3] < 2^64 + 2^47 */
426 * This adjusts the other two words to complete the two partial
430 tmp[1] -= (((limb) b) << 32);
433 * In order to make space in tmp[3] for the carry from 2 -> 3, we
434 * conditionally subtract kPrime if tmp[3] is large enough.
437 /* As tmp[3] < 2^65, high is either 1 or 0 */
442 * all ones if the high word of tmp[3] is 1
443 * all zeros if the high word of tmp[3] if 0 */
448 * all ones if the MSB of low is 1
449 * all zeros if the MSB of low if 0 */
452 /* if low was greater than kPrime3Test then the MSB is zero */
457 * all ones if low was > kPrime3Test
458 * all zeros if low was <= kPrime3Test */
459 mask = (mask & low) | high;
460 tmp[0] -= mask & kPrime[0];
461 tmp[1] -= mask & kPrime[1];
462 /* kPrime[2] is zero, so omitted */
463 tmp[3] -= mask & kPrime[3];
464 /* tmp[3] < 2**64 - 2**32 + 1 */
466 tmp[1] += ((u64)(tmp[0] >> 64));
467 tmp[0] = (u64)tmp[0];
468 tmp[2] += ((u64)(tmp[1] >> 64));
469 tmp[1] = (u64)tmp[1];
470 tmp[3] += ((u64)(tmp[2] >> 64));
471 tmp[2] = (u64)tmp[2];
480 /* smallfelem_expand converts a smallfelem to an felem */
481 static void smallfelem_expand(felem out, const smallfelem in)
490 * smallfelem_square sets |out| = |small|^2
494 * out[i] < 7 * 2^64 < 2^67
496 static void smallfelem_square(longfelem out, const smallfelem small)
501 a = ((uint128_t) small[0]) * small[0];
507 a = ((uint128_t) small[0]) * small[1];
514 a = ((uint128_t) small[0]) * small[2];
521 a = ((uint128_t) small[0]) * small[3];
527 a = ((uint128_t) small[1]) * small[2];
534 a = ((uint128_t) small[1]) * small[1];
540 a = ((uint128_t) small[1]) * small[3];
547 a = ((uint128_t) small[2]) * small[3];
555 a = ((uint128_t) small[2]) * small[2];
561 a = ((uint128_t) small[3]) * small[3];
569 * felem_square sets |out| = |in|^2
573 * out[i] < 7 * 2^64 < 2^67
575 static void felem_square(longfelem out, const felem in)
578 felem_shrink(small, in);
579 smallfelem_square(out, small);
583 * smallfelem_mul sets |out| = |small1| * |small2|
588 * out[i] < 7 * 2^64 < 2^67
590 static void smallfelem_mul(longfelem out, const smallfelem small1,
591 const smallfelem small2)
596 a = ((uint128_t) small1[0]) * small2[0];
602 a = ((uint128_t) small1[0]) * small2[1];
608 a = ((uint128_t) small1[1]) * small2[0];
614 a = ((uint128_t) small1[0]) * small2[2];
620 a = ((uint128_t) small1[1]) * small2[1];
626 a = ((uint128_t) small1[2]) * small2[0];
632 a = ((uint128_t) small1[0]) * small2[3];
638 a = ((uint128_t) small1[1]) * small2[2];
644 a = ((uint128_t) small1[2]) * small2[1];
650 a = ((uint128_t) small1[3]) * small2[0];
656 a = ((uint128_t) small1[1]) * small2[3];
662 a = ((uint128_t) small1[2]) * small2[2];
668 a = ((uint128_t) small1[3]) * small2[1];
674 a = ((uint128_t) small1[2]) * small2[3];
680 a = ((uint128_t) small1[3]) * small2[2];
686 a = ((uint128_t) small1[3]) * small2[3];
694 * felem_mul sets |out| = |in1| * |in2|
699 * out[i] < 7 * 2^64 < 2^67
701 static void felem_mul(longfelem out, const felem in1, const felem in2)
703 smallfelem small1, small2;
704 felem_shrink(small1, in1);
705 felem_shrink(small2, in2);
706 smallfelem_mul(out, small1, small2);
710 * felem_small_mul sets |out| = |small1| * |in2|
715 * out[i] < 7 * 2^64 < 2^67
717 static void felem_small_mul(longfelem out, const smallfelem small1,
721 felem_shrink(small2, in2);
722 smallfelem_mul(out, small1, small2);
725 # define two100m36m4 (((limb)1) << 100) - (((limb)1) << 36) - (((limb)1) << 4)
726 # define two100 (((limb)1) << 100)
727 # define two100m36p4 (((limb)1) << 100) - (((limb)1) << 36) + (((limb)1) << 4)
728 /* zero100 is 0 mod p */
729 static const felem zero100 =
730 { two100m36m4, two100, two100m36p4, two100m36p4 };
733 * Internal function for the different flavours of felem_reduce.
734 * felem_reduce_ reduces the higher coefficients in[4]-in[7].
736 * out[0] >= in[6] + 2^32*in[6] + in[7] + 2^32*in[7]
737 * out[1] >= in[7] + 2^32*in[4]
738 * out[2] >= in[5] + 2^32*in[5]
739 * out[3] >= in[4] + 2^32*in[5] + 2^32*in[6]
741 * out[0] <= out[0] + in[4] + 2^32*in[5]
742 * out[1] <= out[1] + in[5] + 2^33*in[6]
743 * out[2] <= out[2] + in[7] + 2*in[6] + 2^33*in[7]
744 * out[3] <= out[3] + 2^32*in[4] + 3*in[7]
746 static void felem_reduce_(felem out, const longfelem in)
749 /* combine common terms from below */
750 c = in[4] + (in[5] << 32);
758 /* the remaining terms */
759 /* 256: [(0,1),(96,-1),(192,-1),(224,1)] */
760 out[1] -= (in[4] << 32);
761 out[3] += (in[4] << 32);
763 /* 320: [(32,1),(64,1),(128,-1),(160,-1),(224,-1)] */
764 out[2] -= (in[5] << 32);
766 /* 384: [(0,-1),(32,-1),(96,2),(128,2),(224,-1)] */
768 out[0] -= (in[6] << 32);
769 out[1] += (in[6] << 33);
770 out[2] += (in[6] * 2);
771 out[3] -= (in[6] << 32);
773 /* 448: [(0,-1),(32,-1),(64,-1),(128,1),(160,2),(192,3)] */
775 out[0] -= (in[7] << 32);
776 out[2] += (in[7] << 33);
777 out[3] += (in[7] * 3);
781 * felem_reduce converts a longfelem into an felem.
782 * To be called directly after felem_square or felem_mul.
784 * in[0] < 2^64, in[1] < 3*2^64, in[2] < 5*2^64, in[3] < 7*2^64
785 * in[4] < 7*2^64, in[5] < 5*2^64, in[6] < 3*2^64, in[7] < 2*64
789 static void felem_reduce(felem out, const longfelem in)
791 out[0] = zero100[0] + in[0];
792 out[1] = zero100[1] + in[1];
793 out[2] = zero100[2] + in[2];
794 out[3] = zero100[3] + in[3];
796 felem_reduce_(out, in);
799 * out[0] > 2^100 - 2^36 - 2^4 - 3*2^64 - 3*2^96 - 2^64 - 2^96 > 0
800 * out[1] > 2^100 - 2^64 - 7*2^96 > 0
801 * out[2] > 2^100 - 2^36 + 2^4 - 5*2^64 - 5*2^96 > 0
802 * out[3] > 2^100 - 2^36 + 2^4 - 7*2^64 - 5*2^96 - 3*2^96 > 0
804 * out[0] < 2^100 + 2^64 + 7*2^64 + 5*2^96 < 2^101
805 * out[1] < 2^100 + 3*2^64 + 5*2^64 + 3*2^97 < 2^101
806 * out[2] < 2^100 + 5*2^64 + 2^64 + 3*2^65 + 2^97 < 2^101
807 * out[3] < 2^100 + 7*2^64 + 7*2^96 + 3*2^64 < 2^101
812 * felem_reduce_zero105 converts a larger longfelem into an felem.
818 static void felem_reduce_zero105(felem out, const longfelem in)
820 out[0] = zero105[0] + in[0];
821 out[1] = zero105[1] + in[1];
822 out[2] = zero105[2] + in[2];
823 out[3] = zero105[3] + in[3];
825 felem_reduce_(out, in);
828 * out[0] > 2^105 - 2^41 - 2^9 - 2^71 - 2^103 - 2^71 - 2^103 > 0
829 * out[1] > 2^105 - 2^71 - 2^103 > 0
830 * out[2] > 2^105 - 2^41 + 2^9 - 2^71 - 2^103 > 0
831 * out[3] > 2^105 - 2^41 + 2^9 - 2^71 - 2^103 - 2^103 > 0
833 * out[0] < 2^105 + 2^71 + 2^71 + 2^103 < 2^106
834 * out[1] < 2^105 + 2^71 + 2^71 + 2^103 < 2^106
835 * out[2] < 2^105 + 2^71 + 2^71 + 2^71 + 2^103 < 2^106
836 * out[3] < 2^105 + 2^71 + 2^103 + 2^71 < 2^106
841 * subtract_u64 sets *result = *result - v and *carry to one if the
842 * subtraction underflowed.
844 static void subtract_u64(u64 *result, u64 *carry, u64 v)
846 uint128_t r = *result;
848 *carry = (r >> 64) & 1;
853 * felem_contract converts |in| to its unique, minimal representation. On
854 * entry: in[i] < 2^109
856 static void felem_contract(smallfelem out, const felem in)
859 u64 all_equal_so_far = 0, result = 0, carry;
861 felem_shrink(out, in);
862 /* small is minimal except that the value might be > p */
866 * We are doing a constant time test if out >= kPrime. We need to compare
867 * each u64, from most-significant to least significant. For each one, if
868 * all words so far have been equal (m is all ones) then a non-equal
869 * result is the answer. Otherwise we continue.
871 for (i = 3; i < 4; i--) {
873 uint128_t a = ((uint128_t) kPrime[i]) - out[i];
875 * if out[i] > kPrime[i] then a will underflow and the high 64-bits
878 result |= all_equal_so_far & ((u64)(a >> 64));
881 * if kPrime[i] == out[i] then |equal| will be all zeros and the
882 * decrement will make it all ones.
884 equal = kPrime[i] ^ out[i];
886 equal &= equal << 32;
887 equal &= equal << 16;
892 equal = ((s64) equal) >> 63;
894 all_equal_so_far &= equal;
898 * if all_equal_so_far is still all ones then the two values are equal
899 * and so out >= kPrime is true.
901 result |= all_equal_so_far;
903 /* if out >= kPrime then we subtract kPrime. */
904 subtract_u64(&out[0], &carry, result & kPrime[0]);
905 subtract_u64(&out[1], &carry, carry);
906 subtract_u64(&out[2], &carry, carry);
907 subtract_u64(&out[3], &carry, carry);
909 subtract_u64(&out[1], &carry, result & kPrime[1]);
910 subtract_u64(&out[2], &carry, carry);
911 subtract_u64(&out[3], &carry, carry);
913 subtract_u64(&out[2], &carry, result & kPrime[2]);
914 subtract_u64(&out[3], &carry, carry);
916 subtract_u64(&out[3], &carry, result & kPrime[3]);
919 static void smallfelem_square_contract(smallfelem out, const smallfelem in)
924 smallfelem_square(longtmp, in);
925 felem_reduce(tmp, longtmp);
926 felem_contract(out, tmp);
929 static void smallfelem_mul_contract(smallfelem out, const smallfelem in1,
930 const smallfelem in2)
935 smallfelem_mul(longtmp, in1, in2);
936 felem_reduce(tmp, longtmp);
937 felem_contract(out, tmp);
941 * felem_is_zero returns a limb with all bits set if |in| == 0 (mod p) and 0
946 static limb smallfelem_is_zero(const smallfelem small)
951 u64 is_zero = small[0] | small[1] | small[2] | small[3];
953 is_zero &= is_zero << 32;
954 is_zero &= is_zero << 16;
955 is_zero &= is_zero << 8;
956 is_zero &= is_zero << 4;
957 is_zero &= is_zero << 2;
958 is_zero &= is_zero << 1;
959 is_zero = ((s64) is_zero) >> 63;
961 is_p = (small[0] ^ kPrime[0]) |
962 (small[1] ^ kPrime[1]) |
963 (small[2] ^ kPrime[2]) | (small[3] ^ kPrime[3]);
971 is_p = ((s64) is_p) >> 63;
976 result |= ((limb) is_zero) << 64;
980 static int smallfelem_is_zero_int(const smallfelem small)
982 return (int)(smallfelem_is_zero(small) & ((limb) 1));
986 * felem_inv calculates |out| = |in|^{-1}
988 * Based on Fermat's Little Theorem:
990 * a^{p-1} = 1 (mod p)
991 * a^{p-2} = a^{-1} (mod p)
993 static void felem_inv(felem out, const felem in)
996 /* each e_I will hold |in|^{2^I - 1} */
997 felem e2, e4, e8, e16, e32, e64;
1001 felem_square(tmp, in);
1002 felem_reduce(ftmp, tmp); /* 2^1 */
1003 felem_mul(tmp, in, ftmp);
1004 felem_reduce(ftmp, tmp); /* 2^2 - 2^0 */
1005 felem_assign(e2, ftmp);
1006 felem_square(tmp, ftmp);
1007 felem_reduce(ftmp, tmp); /* 2^3 - 2^1 */
1008 felem_square(tmp, ftmp);
1009 felem_reduce(ftmp, tmp); /* 2^4 - 2^2 */
1010 felem_mul(tmp, ftmp, e2);
1011 felem_reduce(ftmp, tmp); /* 2^4 - 2^0 */
1012 felem_assign(e4, ftmp);
1013 felem_square(tmp, ftmp);
1014 felem_reduce(ftmp, tmp); /* 2^5 - 2^1 */
1015 felem_square(tmp, ftmp);
1016 felem_reduce(ftmp, tmp); /* 2^6 - 2^2 */
1017 felem_square(tmp, ftmp);
1018 felem_reduce(ftmp, tmp); /* 2^7 - 2^3 */
1019 felem_square(tmp, ftmp);
1020 felem_reduce(ftmp, tmp); /* 2^8 - 2^4 */
1021 felem_mul(tmp, ftmp, e4);
1022 felem_reduce(ftmp, tmp); /* 2^8 - 2^0 */
1023 felem_assign(e8, ftmp);
1024 for (i = 0; i < 8; i++) {
1025 felem_square(tmp, ftmp);
1026 felem_reduce(ftmp, tmp);
1028 felem_mul(tmp, ftmp, e8);
1029 felem_reduce(ftmp, tmp); /* 2^16 - 2^0 */
1030 felem_assign(e16, ftmp);
1031 for (i = 0; i < 16; i++) {
1032 felem_square(tmp, ftmp);
1033 felem_reduce(ftmp, tmp);
1035 felem_mul(tmp, ftmp, e16);
1036 felem_reduce(ftmp, tmp); /* 2^32 - 2^0 */
1037 felem_assign(e32, ftmp);
1038 for (i = 0; i < 32; i++) {
1039 felem_square(tmp, ftmp);
1040 felem_reduce(ftmp, tmp);
1042 felem_assign(e64, ftmp);
1043 felem_mul(tmp, ftmp, in);
1044 felem_reduce(ftmp, tmp); /* 2^64 - 2^32 + 2^0 */
1045 for (i = 0; i < 192; i++) {
1046 felem_square(tmp, ftmp);
1047 felem_reduce(ftmp, tmp);
1048 } /* 2^256 - 2^224 + 2^192 */
1050 felem_mul(tmp, e64, e32);
1051 felem_reduce(ftmp2, tmp); /* 2^64 - 2^0 */
1052 for (i = 0; i < 16; i++) {
1053 felem_square(tmp, ftmp2);
1054 felem_reduce(ftmp2, tmp);
1056 felem_mul(tmp, ftmp2, e16);
1057 felem_reduce(ftmp2, tmp); /* 2^80 - 2^0 */
1058 for (i = 0; i < 8; i++) {
1059 felem_square(tmp, ftmp2);
1060 felem_reduce(ftmp2, tmp);
1062 felem_mul(tmp, ftmp2, e8);
1063 felem_reduce(ftmp2, tmp); /* 2^88 - 2^0 */
1064 for (i = 0; i < 4; i++) {
1065 felem_square(tmp, ftmp2);
1066 felem_reduce(ftmp2, tmp);
1068 felem_mul(tmp, ftmp2, e4);
1069 felem_reduce(ftmp2, tmp); /* 2^92 - 2^0 */
1070 felem_square(tmp, ftmp2);
1071 felem_reduce(ftmp2, tmp); /* 2^93 - 2^1 */
1072 felem_square(tmp, ftmp2);
1073 felem_reduce(ftmp2, tmp); /* 2^94 - 2^2 */
1074 felem_mul(tmp, ftmp2, e2);
1075 felem_reduce(ftmp2, tmp); /* 2^94 - 2^0 */
1076 felem_square(tmp, ftmp2);
1077 felem_reduce(ftmp2, tmp); /* 2^95 - 2^1 */
1078 felem_square(tmp, ftmp2);
1079 felem_reduce(ftmp2, tmp); /* 2^96 - 2^2 */
1080 felem_mul(tmp, ftmp2, in);
1081 felem_reduce(ftmp2, tmp); /* 2^96 - 3 */
1083 felem_mul(tmp, ftmp2, ftmp);
1084 felem_reduce(out, tmp); /* 2^256 - 2^224 + 2^192 + 2^96 - 3 */
1087 static void smallfelem_inv_contract(smallfelem out, const smallfelem in)
1091 smallfelem_expand(tmp, in);
1092 felem_inv(tmp, tmp);
1093 felem_contract(out, tmp);
1100 * Building on top of the field operations we have the operations on the
1101 * elliptic curve group itself. Points on the curve are represented in Jacobian
1105 * point_double calculates 2*(x_in, y_in, z_in)
1107 * The method is taken from:
1108 * http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#doubling-dbl-2001-b
1110 * Outputs can equal corresponding inputs, i.e., x_out == x_in is allowed.
1111 * while x_out == y_in is not (maybe this works, but it's not tested). */
1113 point_double(felem x_out, felem y_out, felem z_out,
1114 const felem x_in, const felem y_in, const felem z_in)
1116 longfelem tmp, tmp2;
1117 felem delta, gamma, beta, alpha, ftmp, ftmp2;
1118 smallfelem small1, small2;
1120 felem_assign(ftmp, x_in);
1121 /* ftmp[i] < 2^106 */
1122 felem_assign(ftmp2, x_in);
1123 /* ftmp2[i] < 2^106 */
1126 felem_square(tmp, z_in);
1127 felem_reduce(delta, tmp);
1128 /* delta[i] < 2^101 */
1131 felem_square(tmp, y_in);
1132 felem_reduce(gamma, tmp);
1133 /* gamma[i] < 2^101 */
1134 felem_shrink(small1, gamma);
1136 /* beta = x*gamma */
1137 felem_small_mul(tmp, small1, x_in);
1138 felem_reduce(beta, tmp);
1139 /* beta[i] < 2^101 */
1141 /* alpha = 3*(x-delta)*(x+delta) */
1142 felem_diff(ftmp, delta);
1143 /* ftmp[i] < 2^105 + 2^106 < 2^107 */
1144 felem_sum(ftmp2, delta);
1145 /* ftmp2[i] < 2^105 + 2^106 < 2^107 */
1146 felem_scalar(ftmp2, 3);
1147 /* ftmp2[i] < 3 * 2^107 < 2^109 */
1148 felem_mul(tmp, ftmp, ftmp2);
1149 felem_reduce(alpha, tmp);
1150 /* alpha[i] < 2^101 */
1151 felem_shrink(small2, alpha);
1153 /* x' = alpha^2 - 8*beta */
1154 smallfelem_square(tmp, small2);
1155 felem_reduce(x_out, tmp);
1156 felem_assign(ftmp, beta);
1157 felem_scalar(ftmp, 8);
1158 /* ftmp[i] < 8 * 2^101 = 2^104 */
1159 felem_diff(x_out, ftmp);
1160 /* x_out[i] < 2^105 + 2^101 < 2^106 */
1162 /* z' = (y + z)^2 - gamma - delta */
1163 felem_sum(delta, gamma);
1164 /* delta[i] < 2^101 + 2^101 = 2^102 */
1165 felem_assign(ftmp, y_in);
1166 felem_sum(ftmp, z_in);
1167 /* ftmp[i] < 2^106 + 2^106 = 2^107 */
1168 felem_square(tmp, ftmp);
1169 felem_reduce(z_out, tmp);
1170 felem_diff(z_out, delta);
1171 /* z_out[i] < 2^105 + 2^101 < 2^106 */
1173 /* y' = alpha*(4*beta - x') - 8*gamma^2 */
1174 felem_scalar(beta, 4);
1175 /* beta[i] < 4 * 2^101 = 2^103 */
1176 felem_diff_zero107(beta, x_out);
1177 /* beta[i] < 2^107 + 2^103 < 2^108 */
1178 felem_small_mul(tmp, small2, beta);
1179 /* tmp[i] < 7 * 2^64 < 2^67 */
1180 smallfelem_square(tmp2, small1);
1181 /* tmp2[i] < 7 * 2^64 */
1182 longfelem_scalar(tmp2, 8);
1183 /* tmp2[i] < 8 * 7 * 2^64 = 7 * 2^67 */
1184 longfelem_diff(tmp, tmp2);
1185 /* tmp[i] < 2^67 + 2^70 + 2^40 < 2^71 */
1186 felem_reduce_zero105(y_out, tmp);
1187 /* y_out[i] < 2^106 */
1191 * point_double_small is the same as point_double, except that it operates on
1195 point_double_small(smallfelem x_out, smallfelem y_out, smallfelem z_out,
1196 const smallfelem x_in, const smallfelem y_in,
1197 const smallfelem z_in)
1199 felem felem_x_out, felem_y_out, felem_z_out;
1200 felem felem_x_in, felem_y_in, felem_z_in;
1202 smallfelem_expand(felem_x_in, x_in);
1203 smallfelem_expand(felem_y_in, y_in);
1204 smallfelem_expand(felem_z_in, z_in);
1205 point_double(felem_x_out, felem_y_out, felem_z_out,
1206 felem_x_in, felem_y_in, felem_z_in);
1207 felem_shrink(x_out, felem_x_out);
1208 felem_shrink(y_out, felem_y_out);
1209 felem_shrink(z_out, felem_z_out);
1212 /* copy_conditional copies in to out iff mask is all ones. */
1213 static void copy_conditional(felem out, const felem in, limb mask)
1216 for (i = 0; i < NLIMBS; ++i) {
1217 const limb tmp = mask & (in[i] ^ out[i]);
1222 /* copy_small_conditional copies in to out iff mask is all ones. */
1223 static void copy_small_conditional(felem out, const smallfelem in, limb mask)
1226 const u64 mask64 = mask;
1227 for (i = 0; i < NLIMBS; ++i) {
1228 out[i] = ((limb) (in[i] & mask64)) | (out[i] & ~mask);
1233 * point_add calcuates (x1, y1, z1) + (x2, y2, z2)
1235 * The method is taken from:
1236 * http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#addition-add-2007-bl,
1237 * adapted for mixed addition (z2 = 1, or z2 = 0 for the point at infinity).
1239 * This function includes a branch for checking whether the two input points
1240 * are equal, (while not equal to the point at infinity). This case never
1241 * happens during single point multiplication, so there is no timing leak for
1242 * ECDH or ECDSA signing. */
1243 static void point_add(felem x3, felem y3, felem z3,
1244 const felem x1, const felem y1, const felem z1,
1245 const int mixed, const smallfelem x2,
1246 const smallfelem y2, const smallfelem z2)
1248 felem ftmp, ftmp2, ftmp3, ftmp4, ftmp5, ftmp6, x_out, y_out, z_out;
1249 longfelem tmp, tmp2;
1250 smallfelem small1, small2, small3, small4, small5;
1251 limb x_equal, y_equal, z1_is_zero, z2_is_zero;
1253 felem_shrink(small3, z1);
1255 z1_is_zero = smallfelem_is_zero(small3);
1256 z2_is_zero = smallfelem_is_zero(z2);
1258 /* ftmp = z1z1 = z1**2 */
1259 smallfelem_square(tmp, small3);
1260 felem_reduce(ftmp, tmp);
1261 /* ftmp[i] < 2^101 */
1262 felem_shrink(small1, ftmp);
1265 /* ftmp2 = z2z2 = z2**2 */
1266 smallfelem_square(tmp, z2);
1267 felem_reduce(ftmp2, tmp);
1268 /* ftmp2[i] < 2^101 */
1269 felem_shrink(small2, ftmp2);
1271 felem_shrink(small5, x1);
1273 /* u1 = ftmp3 = x1*z2z2 */
1274 smallfelem_mul(tmp, small5, small2);
1275 felem_reduce(ftmp3, tmp);
1276 /* ftmp3[i] < 2^101 */
1278 /* ftmp5 = z1 + z2 */
1279 felem_assign(ftmp5, z1);
1280 felem_small_sum(ftmp5, z2);
1281 /* ftmp5[i] < 2^107 */
1283 /* ftmp5 = (z1 + z2)**2 - (z1z1 + z2z2) = 2z1z2 */
1284 felem_square(tmp, ftmp5);
1285 felem_reduce(ftmp5, tmp);
1286 /* ftmp2 = z2z2 + z1z1 */
1287 felem_sum(ftmp2, ftmp);
1288 /* ftmp2[i] < 2^101 + 2^101 = 2^102 */
1289 felem_diff(ftmp5, ftmp2);
1290 /* ftmp5[i] < 2^105 + 2^101 < 2^106 */
1292 /* ftmp2 = z2 * z2z2 */
1293 smallfelem_mul(tmp, small2, z2);
1294 felem_reduce(ftmp2, tmp);
1296 /* s1 = ftmp2 = y1 * z2**3 */
1297 felem_mul(tmp, y1, ftmp2);
1298 felem_reduce(ftmp6, tmp);
1299 /* ftmp6[i] < 2^101 */
1302 * We'll assume z2 = 1 (special case z2 = 0 is handled later)
1305 /* u1 = ftmp3 = x1*z2z2 */
1306 felem_assign(ftmp3, x1);
1307 /* ftmp3[i] < 2^106 */
1310 felem_assign(ftmp5, z1);
1311 felem_scalar(ftmp5, 2);
1312 /* ftmp5[i] < 2*2^106 = 2^107 */
1314 /* s1 = ftmp2 = y1 * z2**3 */
1315 felem_assign(ftmp6, y1);
1316 /* ftmp6[i] < 2^106 */
1320 smallfelem_mul(tmp, x2, small1);
1321 felem_reduce(ftmp4, tmp);
1323 /* h = ftmp4 = u2 - u1 */
1324 felem_diff_zero107(ftmp4, ftmp3);
1325 /* ftmp4[i] < 2^107 + 2^101 < 2^108 */
1326 felem_shrink(small4, ftmp4);
1328 x_equal = smallfelem_is_zero(small4);
1330 /* z_out = ftmp5 * h */
1331 felem_small_mul(tmp, small4, ftmp5);
1332 felem_reduce(z_out, tmp);
1333 /* z_out[i] < 2^101 */
1335 /* ftmp = z1 * z1z1 */
1336 smallfelem_mul(tmp, small1, small3);
1337 felem_reduce(ftmp, tmp);
1339 /* s2 = tmp = y2 * z1**3 */
1340 felem_small_mul(tmp, y2, ftmp);
1341 felem_reduce(ftmp5, tmp);
1343 /* r = ftmp5 = (s2 - s1)*2 */
1344 felem_diff_zero107(ftmp5, ftmp6);
1345 /* ftmp5[i] < 2^107 + 2^107 = 2^108 */
1346 felem_scalar(ftmp5, 2);
1347 /* ftmp5[i] < 2^109 */
1348 felem_shrink(small1, ftmp5);
1349 y_equal = smallfelem_is_zero(small1);
1351 if (x_equal && y_equal && !z1_is_zero && !z2_is_zero) {
1352 point_double(x3, y3, z3, x1, y1, z1);
1356 /* I = ftmp = (2h)**2 */
1357 felem_assign(ftmp, ftmp4);
1358 felem_scalar(ftmp, 2);
1359 /* ftmp[i] < 2*2^108 = 2^109 */
1360 felem_square(tmp, ftmp);
1361 felem_reduce(ftmp, tmp);
1363 /* J = ftmp2 = h * I */
1364 felem_mul(tmp, ftmp4, ftmp);
1365 felem_reduce(ftmp2, tmp);
1367 /* V = ftmp4 = U1 * I */
1368 felem_mul(tmp, ftmp3, ftmp);
1369 felem_reduce(ftmp4, tmp);
1371 /* x_out = r**2 - J - 2V */
1372 smallfelem_square(tmp, small1);
1373 felem_reduce(x_out, tmp);
1374 felem_assign(ftmp3, ftmp4);
1375 felem_scalar(ftmp4, 2);
1376 felem_sum(ftmp4, ftmp2);
1377 /* ftmp4[i] < 2*2^101 + 2^101 < 2^103 */
1378 felem_diff(x_out, ftmp4);
1379 /* x_out[i] < 2^105 + 2^101 */
1381 /* y_out = r(V-x_out) - 2 * s1 * J */
1382 felem_diff_zero107(ftmp3, x_out);
1383 /* ftmp3[i] < 2^107 + 2^101 < 2^108 */
1384 felem_small_mul(tmp, small1, ftmp3);
1385 felem_mul(tmp2, ftmp6, ftmp2);
1386 longfelem_scalar(tmp2, 2);
1387 /* tmp2[i] < 2*2^67 = 2^68 */
1388 longfelem_diff(tmp, tmp2);
1389 /* tmp[i] < 2^67 + 2^70 + 2^40 < 2^71 */
1390 felem_reduce_zero105(y_out, tmp);
1391 /* y_out[i] < 2^106 */
1393 copy_small_conditional(x_out, x2, z1_is_zero);
1394 copy_conditional(x_out, x1, z2_is_zero);
1395 copy_small_conditional(y_out, y2, z1_is_zero);
1396 copy_conditional(y_out, y1, z2_is_zero);
1397 copy_small_conditional(z_out, z2, z1_is_zero);
1398 copy_conditional(z_out, z1, z2_is_zero);
1399 felem_assign(x3, x_out);
1400 felem_assign(y3, y_out);
1401 felem_assign(z3, z_out);
1405 * point_add_small is the same as point_add, except that it operates on
1408 static void point_add_small(smallfelem x3, smallfelem y3, smallfelem z3,
1409 smallfelem x1, smallfelem y1, smallfelem z1,
1410 smallfelem x2, smallfelem y2, smallfelem z2)
1412 felem felem_x3, felem_y3, felem_z3;
1413 felem felem_x1, felem_y1, felem_z1;
1414 smallfelem_expand(felem_x1, x1);
1415 smallfelem_expand(felem_y1, y1);
1416 smallfelem_expand(felem_z1, z1);
1417 point_add(felem_x3, felem_y3, felem_z3, felem_x1, felem_y1, felem_z1, 0,
1419 felem_shrink(x3, felem_x3);
1420 felem_shrink(y3, felem_y3);
1421 felem_shrink(z3, felem_z3);
1425 * Base point pre computation
1426 * --------------------------
1428 * Two different sorts of precomputed tables are used in the following code.
1429 * Each contain various points on the curve, where each point is three field
1430 * elements (x, y, z).
1432 * For the base point table, z is usually 1 (0 for the point at infinity).
1433 * This table has 2 * 16 elements, starting with the following:
1434 * index | bits | point
1435 * ------+---------+------------------------------
1438 * 2 | 0 0 1 0 | 2^64G
1439 * 3 | 0 0 1 1 | (2^64 + 1)G
1440 * 4 | 0 1 0 0 | 2^128G
1441 * 5 | 0 1 0 1 | (2^128 + 1)G
1442 * 6 | 0 1 1 0 | (2^128 + 2^64)G
1443 * 7 | 0 1 1 1 | (2^128 + 2^64 + 1)G
1444 * 8 | 1 0 0 0 | 2^192G
1445 * 9 | 1 0 0 1 | (2^192 + 1)G
1446 * 10 | 1 0 1 0 | (2^192 + 2^64)G
1447 * 11 | 1 0 1 1 | (2^192 + 2^64 + 1)G
1448 * 12 | 1 1 0 0 | (2^192 + 2^128)G
1449 * 13 | 1 1 0 1 | (2^192 + 2^128 + 1)G
1450 * 14 | 1 1 1 0 | (2^192 + 2^128 + 2^64)G
1451 * 15 | 1 1 1 1 | (2^192 + 2^128 + 2^64 + 1)G
1452 * followed by a copy of this with each element multiplied by 2^32.
1454 * The reason for this is so that we can clock bits into four different
1455 * locations when doing simple scalar multiplies against the base point,
1456 * and then another four locations using the second 16 elements.
1458 * Tables for other points have table[i] = iG for i in 0 .. 16. */
1460 /* gmul is the table of precomputed base points */
1461 static const smallfelem gmul[2][16][3] = {
1465 {{0xf4a13945d898c296, 0x77037d812deb33a0, 0xf8bce6e563a440f2,
1466 0x6b17d1f2e12c4247},
1467 {0xcbb6406837bf51f5, 0x2bce33576b315ece, 0x8ee7eb4a7c0f9e16,
1468 0x4fe342e2fe1a7f9b},
1470 {{0x90e75cb48e14db63, 0x29493baaad651f7e, 0x8492592e326e25de,
1471 0x0fa822bc2811aaa5},
1472 {0xe41124545f462ee7, 0x34b1a65050fe82f5, 0x6f4ad4bcb3df188b,
1473 0xbff44ae8f5dba80d},
1475 {{0x93391ce2097992af, 0xe96c98fd0d35f1fa, 0xb257c0de95e02789,
1476 0x300a4bbc89d6726f},
1477 {0xaa54a291c08127a0, 0x5bb1eeada9d806a5, 0x7f1ddb25ff1e3c6f,
1478 0x72aac7e0d09b4644},
1480 {{0x57c84fc9d789bd85, 0xfc35ff7dc297eac3, 0xfb982fd588c6766e,
1481 0x447d739beedb5e67},
1482 {0x0c7e33c972e25b32, 0x3d349b95a7fae500, 0xe12e9d953a4aaff7,
1483 0x2d4825ab834131ee},
1485 {{0x13949c932a1d367f, 0xef7fbd2b1a0a11b7, 0xddc6068bb91dfc60,
1486 0xef9519328a9c72ff},
1487 {0x196035a77376d8a8, 0x23183b0895ca1740, 0xc1ee9807022c219c,
1488 0x611e9fc37dbb2c9b},
1490 {{0xcae2b1920b57f4bc, 0x2936df5ec6c9bc36, 0x7dea6482e11238bf,
1491 0x550663797b51f5d8},
1492 {0x44ffe216348a964c, 0x9fb3d576dbdefbe1, 0x0afa40018d9d50e5,
1493 0x157164848aecb851},
1495 {{0xe48ecafffc5cde01, 0x7ccd84e70d715f26, 0xa2e8f483f43e4391,
1496 0xeb5d7745b21141ea},
1497 {0xcac917e2731a3479, 0x85f22cfe2844b645, 0x0990e6a158006cee,
1498 0xeafd72ebdbecc17b},
1500 {{0x6cf20ffb313728be, 0x96439591a3c6b94a, 0x2736ff8344315fc5,
1501 0xa6d39677a7849276},
1502 {0xf2bab833c357f5f4, 0x824a920c2284059b, 0x66b8babd2d27ecdf,
1503 0x674f84749b0b8816},
1505 {{0x2df48c04677c8a3e, 0x74e02f080203a56b, 0x31855f7db8c7fedb,
1506 0x4e769e7672c9ddad},
1507 {0xa4c36165b824bbb0, 0xfb9ae16f3b9122a5, 0x1ec0057206947281,
1508 0x42b99082de830663},
1510 {{0x6ef95150dda868b9, 0xd1f89e799c0ce131, 0x7fdc1ca008a1c478,
1511 0x78878ef61c6ce04d},
1512 {0x9c62b9121fe0d976, 0x6ace570ebde08d4f, 0xde53142c12309def,
1513 0xb6cb3f5d7b72c321},
1515 {{0x7f991ed2c31a3573, 0x5b82dd5bd54fb496, 0x595c5220812ffcae,
1516 0x0c88bc4d716b1287},
1517 {0x3a57bf635f48aca8, 0x7c8181f4df2564f3, 0x18d1b5b39c04e6aa,
1518 0xdd5ddea3f3901dc6},
1520 {{0xe96a79fb3e72ad0c, 0x43a0a28c42ba792f, 0xefe0a423083e49f3,
1521 0x68f344af6b317466},
1522 {0xcdfe17db3fb24d4a, 0x668bfc2271f5c626, 0x604ed93c24d67ff3,
1523 0x31b9c405f8540a20},
1525 {{0xd36b4789a2582e7f, 0x0d1a10144ec39c28, 0x663c62c3edbad7a0,
1526 0x4052bf4b6f461db9},
1527 {0x235a27c3188d25eb, 0xe724f33999bfcc5b, 0x862be6bd71d70cc8,
1528 0xfecf4d5190b0fc61},
1530 {{0x74346c10a1d4cfac, 0xafdf5cc08526a7a4, 0x123202a8f62bff7a,
1531 0x1eddbae2c802e41a},
1532 {0x8fa0af2dd603f844, 0x36e06b7e4c701917, 0x0c45f45273db33a0,
1533 0x43104d86560ebcfc},
1535 {{0x9615b5110d1d78e5, 0x66b0de3225c4744b, 0x0a4a46fb6aaf363a,
1536 0xb48e26b484f7a21c},
1537 {0x06ebb0f621a01b2d, 0xc004e4048b7b0f98, 0x64131bcdfed6f668,
1538 0xfac015404d4d3dab},
1543 {{0x3a5a9e22185a5943, 0x1ab919365c65dfb6, 0x21656b32262c71da,
1544 0x7fe36b40af22af89},
1545 {0xd50d152c699ca101, 0x74b3d5867b8af212, 0x9f09f40407dca6f1,
1546 0xe697d45825b63624},
1548 {{0xa84aa9397512218e, 0xe9a521b074ca0141, 0x57880b3a18a2e902,
1549 0x4a5b506612a677a6},
1550 {0x0beada7a4c4f3840, 0x626db15419e26d9d, 0xc42604fbe1627d40,
1551 0xeb13461ceac089f1},
1553 {{0xf9faed0927a43281, 0x5e52c4144103ecbc, 0xc342967aa815c857,
1554 0x0781b8291c6a220a},
1555 {0x5a8343ceeac55f80, 0x88f80eeee54a05e3, 0x97b2a14f12916434,
1556 0x690cde8df0151593},
1558 {{0xaee9c75df7f82f2a, 0x9e4c35874afdf43a, 0xf5622df437371326,
1559 0x8a535f566ec73617},
1560 {0xc5f9a0ac223094b7, 0xcde533864c8c7669, 0x37e02819085a92bf,
1561 0x0455c08468b08bd7},
1563 {{0x0c0a6e2c9477b5d9, 0xf9a4bf62876dc444, 0x5050a949b6cdc279,
1564 0x06bada7ab77f8276},
1565 {0xc8b4aed1ea48dac9, 0xdebd8a4b7ea1070f, 0x427d49101366eb70,
1566 0x5b476dfd0e6cb18a},
1568 {{0x7c5c3e44278c340a, 0x4d54606812d66f3b, 0x29a751b1ae23c5d8,
1569 0x3e29864e8a2ec908},
1570 {0x142d2a6626dbb850, 0xad1744c4765bd780, 0x1f150e68e322d1ed,
1571 0x239b90ea3dc31e7e},
1573 {{0x78c416527a53322a, 0x305dde6709776f8e, 0xdbcab759f8862ed4,
1574 0x820f4dd949f72ff7},
1575 {0x6cc544a62b5debd4, 0x75be5d937b4e8cc4, 0x1b481b1b215c14d3,
1576 0x140406ec783a05ec},
1578 {{0x6a703f10e895df07, 0xfd75f3fa01876bd8, 0xeb5b06e70ce08ffe,
1579 0x68f6b8542783dfee},
1580 {0x90c76f8a78712655, 0xcf5293d2f310bf7f, 0xfbc8044dfda45028,
1581 0xcbe1feba92e40ce6},
1583 {{0xe998ceea4396e4c1, 0xfc82ef0b6acea274, 0x230f729f2250e927,
1584 0xd0b2f94d2f420109},
1585 {0x4305adddb38d4966, 0x10b838f8624c3b45, 0x7db2636658954e7a,
1586 0x971459828b0719e5},
1588 {{0x4bd6b72623369fc9, 0x57f2929e53d0b876, 0xc2d5cba4f2340687,
1589 0x961610004a866aba},
1590 {0x49997bcd2e407a5e, 0x69ab197d92ddcb24, 0x2cf1f2438fe5131c,
1591 0x7acb9fadcee75e44},
1593 {{0x254e839423d2d4c0, 0xf57f0c917aea685b, 0xa60d880f6f75aaea,
1594 0x24eb9acca333bf5b},
1595 {0xe3de4ccb1cda5dea, 0xfeef9341c51a6b4f, 0x743125f88bac4c4d,
1596 0x69f891c5acd079cc},
1598 {{0xeee44b35702476b5, 0x7ed031a0e45c2258, 0xb422d1e7bd6f8514,
1599 0xe51f547c5972a107},
1600 {0xa25bcd6fc9cf343d, 0x8ca922ee097c184e, 0xa62f98b3a9fe9a06,
1601 0x1c309a2b25bb1387},
1603 {{0x9295dbeb1967c459, 0xb00148833472c98e, 0xc504977708011828,
1604 0x20b87b8aa2c4e503},
1605 {0x3063175de057c277, 0x1bd539338fe582dd, 0x0d11adef5f69a044,
1606 0xf5c6fa49919776be},
1608 {{0x8c944e760fd59e11, 0x3876cba1102fad5f, 0xa454c3fad83faa56,
1609 0x1ed7d1b9332010b9},
1610 {0xa1011a270024b889, 0x05e4d0dcac0cd344, 0x52b520f0eb6a2a24,
1611 0x3a2b03f03217257a},
1613 {{0xf20fc2afdf1d043d, 0xf330240db58d5a62, 0xfc7d229ca0058c3b,
1614 0x15fee545c78dd9f6},
1615 {0x501e82885bc98cda, 0x41ef80e5d046ac04, 0x557d9f49461210fb,
1616 0x4ab5b6b2b8753f81},
1621 * select_point selects the |idx|th point from a precomputation table and
1624 static void select_point(const u64 idx, unsigned int size,
1625 const smallfelem pre_comp[16][3], smallfelem out[3])
1628 u64 *outlimbs = &out[0][0];
1629 memset(outlimbs, 0, 3 * sizeof(smallfelem));
1631 for (i = 0; i < size; i++) {
1632 const u64 *inlimbs = (u64 *)&pre_comp[i][0][0];
1639 for (j = 0; j < NLIMBS * 3; j++)
1640 outlimbs[j] |= inlimbs[j] & mask;
1644 /* get_bit returns the |i|th bit in |in| */
1645 static char get_bit(const felem_bytearray in, int i)
1647 if ((i < 0) || (i >= 256))
1649 return (in[i >> 3] >> (i & 7)) & 1;
1653 * Interleaved point multiplication using precomputed point multiples: The
1654 * small point multiples 0*P, 1*P, ..., 17*P are in pre_comp[], the scalars
1655 * in scalars[]. If g_scalar is non-NULL, we also add this multiple of the
1656 * generator, using certain (large) precomputed multiples in g_pre_comp.
1657 * Output point (X, Y, Z) is stored in x_out, y_out, z_out
1659 static void batch_mul(felem x_out, felem y_out, felem z_out,
1660 const felem_bytearray scalars[],
1661 const unsigned num_points, const u8 *g_scalar,
1662 const int mixed, const smallfelem pre_comp[][17][3],
1663 const smallfelem g_pre_comp[2][16][3])
1666 unsigned num, gen_mul = (g_scalar != NULL);
1672 /* set nq to the point at infinity */
1673 memset(nq, 0, 3 * sizeof(felem));
1676 * Loop over all scalars msb-to-lsb, interleaving additions of multiples
1677 * of the generator (two in each of the last 32 rounds) and additions of
1678 * other points multiples (every 5th round).
1680 skip = 1; /* save two point operations in the first
1682 for (i = (num_points ? 255 : 31); i >= 0; --i) {
1685 point_double(nq[0], nq[1], nq[2], nq[0], nq[1], nq[2]);
1687 /* add multiples of the generator */
1688 if (gen_mul && (i <= 31)) {
1689 /* first, look 32 bits upwards */
1690 bits = get_bit(g_scalar, i + 224) << 3;
1691 bits |= get_bit(g_scalar, i + 160) << 2;
1692 bits |= get_bit(g_scalar, i + 96) << 1;
1693 bits |= get_bit(g_scalar, i + 32);
1694 /* select the point to add, in constant time */
1695 select_point(bits, 16, g_pre_comp[1], tmp);
1698 /* Arg 1 below is for "mixed" */
1699 point_add(nq[0], nq[1], nq[2],
1700 nq[0], nq[1], nq[2], 1, tmp[0], tmp[1], tmp[2]);
1702 smallfelem_expand(nq[0], tmp[0]);
1703 smallfelem_expand(nq[1], tmp[1]);
1704 smallfelem_expand(nq[2], tmp[2]);
1708 /* second, look at the current position */
1709 bits = get_bit(g_scalar, i + 192) << 3;
1710 bits |= get_bit(g_scalar, i + 128) << 2;
1711 bits |= get_bit(g_scalar, i + 64) << 1;
1712 bits |= get_bit(g_scalar, i);
1713 /* select the point to add, in constant time */
1714 select_point(bits, 16, g_pre_comp[0], tmp);
1715 /* Arg 1 below is for "mixed" */
1716 point_add(nq[0], nq[1], nq[2],
1717 nq[0], nq[1], nq[2], 1, tmp[0], tmp[1], tmp[2]);
1720 /* do other additions every 5 doublings */
1721 if (num_points && (i % 5 == 0)) {
1722 /* loop over all scalars */
1723 for (num = 0; num < num_points; ++num) {
1724 bits = get_bit(scalars[num], i + 4) << 5;
1725 bits |= get_bit(scalars[num], i + 3) << 4;
1726 bits |= get_bit(scalars[num], i + 2) << 3;
1727 bits |= get_bit(scalars[num], i + 1) << 2;
1728 bits |= get_bit(scalars[num], i) << 1;
1729 bits |= get_bit(scalars[num], i - 1);
1730 ec_GFp_nistp_recode_scalar_bits(&sign, &digit, bits);
1733 * select the point to add or subtract, in constant time
1735 select_point(digit, 17, pre_comp[num], tmp);
1736 smallfelem_neg(ftmp, tmp[1]); /* (X, -Y, Z) is the negative
1738 copy_small_conditional(ftmp, tmp[1], (((limb) sign) - 1));
1739 felem_contract(tmp[1], ftmp);
1742 point_add(nq[0], nq[1], nq[2],
1743 nq[0], nq[1], nq[2],
1744 mixed, tmp[0], tmp[1], tmp[2]);
1746 smallfelem_expand(nq[0], tmp[0]);
1747 smallfelem_expand(nq[1], tmp[1]);
1748 smallfelem_expand(nq[2], tmp[2]);
1754 felem_assign(x_out, nq[0]);
1755 felem_assign(y_out, nq[1]);
1756 felem_assign(z_out, nq[2]);
1759 /* Precomputation for the group generator. */
1761 smallfelem g_pre_comp[2][16][3];
1763 } NISTP256_PRE_COMP;
1765 const EC_METHOD *EC_GFp_nistp256_method(void)
1767 static const EC_METHOD ret = {
1768 EC_FLAGS_DEFAULT_OCT,
1769 NID_X9_62_prime_field,
1770 ec_GFp_nistp256_group_init,
1771 ec_GFp_simple_group_finish,
1772 ec_GFp_simple_group_clear_finish,
1773 ec_GFp_nist_group_copy,
1774 ec_GFp_nistp256_group_set_curve,
1775 ec_GFp_simple_group_get_curve,
1776 ec_GFp_simple_group_get_degree,
1777 ec_GFp_simple_group_check_discriminant,
1778 ec_GFp_simple_point_init,
1779 ec_GFp_simple_point_finish,
1780 ec_GFp_simple_point_clear_finish,
1781 ec_GFp_simple_point_copy,
1782 ec_GFp_simple_point_set_to_infinity,
1783 ec_GFp_simple_set_Jprojective_coordinates_GFp,
1784 ec_GFp_simple_get_Jprojective_coordinates_GFp,
1785 ec_GFp_simple_point_set_affine_coordinates,
1786 ec_GFp_nistp256_point_get_affine_coordinates,
1787 0 /* point_set_compressed_coordinates */ ,
1792 ec_GFp_simple_invert,
1793 ec_GFp_simple_is_at_infinity,
1794 ec_GFp_simple_is_on_curve,
1796 ec_GFp_simple_make_affine,
1797 ec_GFp_simple_points_make_affine,
1798 ec_GFp_nistp256_points_mul,
1799 ec_GFp_nistp256_precompute_mult,
1800 ec_GFp_nistp256_have_precompute_mult,
1801 ec_GFp_nist_field_mul,
1802 ec_GFp_nist_field_sqr,
1804 0 /* field_encode */ ,
1805 0 /* field_decode */ ,
1806 0 /* field_set_to_one */
1812 /******************************************************************************/
1814 * FUNCTIONS TO MANAGE PRECOMPUTATION
1817 static NISTP256_PRE_COMP *nistp256_pre_comp_new()
1819 NISTP256_PRE_COMP *ret = NULL;
1820 ret = (NISTP256_PRE_COMP *) OPENSSL_malloc(sizeof *ret);
1822 ECerr(EC_F_NISTP256_PRE_COMP_NEW, ERR_R_MALLOC_FAILURE);
1825 memset(ret->g_pre_comp, 0, sizeof(ret->g_pre_comp));
1826 ret->references = 1;
1830 static void *nistp256_pre_comp_dup(void *src_)
1832 NISTP256_PRE_COMP *src = src_;
1834 /* no need to actually copy, these objects never change! */
1835 CRYPTO_add(&src->references, 1, CRYPTO_LOCK_EC_PRE_COMP);
1840 static void nistp256_pre_comp_free(void *pre_)
1843 NISTP256_PRE_COMP *pre = pre_;
1848 i = CRYPTO_add(&pre->references, -1, CRYPTO_LOCK_EC_PRE_COMP);
1855 static void nistp256_pre_comp_clear_free(void *pre_)
1858 NISTP256_PRE_COMP *pre = pre_;
1863 i = CRYPTO_add(&pre->references, -1, CRYPTO_LOCK_EC_PRE_COMP);
1867 OPENSSL_cleanse(pre, sizeof *pre);
1871 /******************************************************************************/
1873 * OPENSSL EC_METHOD FUNCTIONS
1876 int ec_GFp_nistp256_group_init(EC_GROUP *group)
1879 ret = ec_GFp_simple_group_init(group);
1880 group->a_is_minus3 = 1;
1884 int ec_GFp_nistp256_group_set_curve(EC_GROUP *group, const BIGNUM *p,
1885 const BIGNUM *a, const BIGNUM *b,
1889 BN_CTX *new_ctx = NULL;
1890 BIGNUM *curve_p, *curve_a, *curve_b;
1893 if ((ctx = new_ctx = BN_CTX_new()) == NULL)
1896 if (((curve_p = BN_CTX_get(ctx)) == NULL) ||
1897 ((curve_a = BN_CTX_get(ctx)) == NULL) ||
1898 ((curve_b = BN_CTX_get(ctx)) == NULL))
1900 BN_bin2bn(nistp256_curve_params[0], sizeof(felem_bytearray), curve_p);
1901 BN_bin2bn(nistp256_curve_params[1], sizeof(felem_bytearray), curve_a);
1902 BN_bin2bn(nistp256_curve_params[2], sizeof(felem_bytearray), curve_b);
1903 if ((BN_cmp(curve_p, p)) || (BN_cmp(curve_a, a)) || (BN_cmp(curve_b, b))) {
1904 ECerr(EC_F_EC_GFP_NISTP256_GROUP_SET_CURVE,
1905 EC_R_WRONG_CURVE_PARAMETERS);
1908 group->field_mod_func = BN_nist_mod_256;
1909 ret = ec_GFp_simple_group_set_curve(group, p, a, b, ctx);
1912 if (new_ctx != NULL)
1913 BN_CTX_free(new_ctx);
1918 * Takes the Jacobian coordinates (X, Y, Z) of a point and returns (X', Y') =
1921 int ec_GFp_nistp256_point_get_affine_coordinates(const EC_GROUP *group,
1922 const EC_POINT *point,
1923 BIGNUM *x, BIGNUM *y,
1926 felem z1, z2, x_in, y_in;
1927 smallfelem x_out, y_out;
1930 if (EC_POINT_is_at_infinity(group, point)) {
1931 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1932 EC_R_POINT_AT_INFINITY);
1935 if ((!BN_to_felem(x_in, &point->X)) || (!BN_to_felem(y_in, &point->Y)) ||
1936 (!BN_to_felem(z1, &point->Z)))
1939 felem_square(tmp, z2);
1940 felem_reduce(z1, tmp);
1941 felem_mul(tmp, x_in, z1);
1942 felem_reduce(x_in, tmp);
1943 felem_contract(x_out, x_in);
1945 if (!smallfelem_to_BN(x, x_out)) {
1946 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1951 felem_mul(tmp, z1, z2);
1952 felem_reduce(z1, tmp);
1953 felem_mul(tmp, y_in, z1);
1954 felem_reduce(y_in, tmp);
1955 felem_contract(y_out, y_in);
1957 if (!smallfelem_to_BN(y, y_out)) {
1958 ECerr(EC_F_EC_GFP_NISTP256_POINT_GET_AFFINE_COORDINATES,
1966 /* points below is of size |num|, and tmp_smallfelems is of size |num+1| */
1967 static void make_points_affine(size_t num, smallfelem points[][3],
1968 smallfelem tmp_smallfelems[])
1971 * Runs in constant time, unless an input is the point at infinity (which
1972 * normally shouldn't happen).
1974 ec_GFp_nistp_points_make_affine_internal(num,
1978 (void (*)(void *))smallfelem_one,
1979 (int (*)(const void *))
1980 smallfelem_is_zero_int,
1981 (void (*)(void *, const void *))
1983 (void (*)(void *, const void *))
1984 smallfelem_square_contract,
1986 (void *, const void *,
1988 smallfelem_mul_contract,
1989 (void (*)(void *, const void *))
1990 smallfelem_inv_contract,
1991 /* nothing to contract */
1992 (void (*)(void *, const void *))
1997 * Computes scalar*generator + \sum scalars[i]*points[i], ignoring NULL
1998 * values Result is stored in r (r can equal one of the inputs).
2000 int ec_GFp_nistp256_points_mul(const EC_GROUP *group, EC_POINT *r,
2001 const BIGNUM *scalar, size_t num,
2002 const EC_POINT *points[],
2003 const BIGNUM *scalars[], BN_CTX *ctx)
2008 BN_CTX *new_ctx = NULL;
2009 BIGNUM *x, *y, *z, *tmp_scalar;
2010 felem_bytearray g_secret;
2011 felem_bytearray *secrets = NULL;
2012 smallfelem(*pre_comp)[17][3] = NULL;
2013 smallfelem *tmp_smallfelems = NULL;
2014 felem_bytearray tmp;
2015 unsigned i, num_bytes;
2016 int have_pre_comp = 0;
2017 size_t num_points = num;
2018 smallfelem x_in, y_in, z_in;
2019 felem x_out, y_out, z_out;
2020 NISTP256_PRE_COMP *pre = NULL;
2021 const smallfelem(*g_pre_comp)[16][3] = NULL;
2022 EC_POINT *generator = NULL;
2023 const EC_POINT *p = NULL;
2024 const BIGNUM *p_scalar = NULL;
2027 if ((ctx = new_ctx = BN_CTX_new()) == NULL)
2030 if (((x = BN_CTX_get(ctx)) == NULL) ||
2031 ((y = BN_CTX_get(ctx)) == NULL) ||
2032 ((z = BN_CTX_get(ctx)) == NULL) ||
2033 ((tmp_scalar = BN_CTX_get(ctx)) == NULL))
2036 if (scalar != NULL) {
2037 pre = EC_EX_DATA_get_data(group->extra_data,
2038 nistp256_pre_comp_dup,
2039 nistp256_pre_comp_free,
2040 nistp256_pre_comp_clear_free);
2042 /* we have precomputation, try to use it */
2043 g_pre_comp = (const smallfelem(*)[16][3])pre->g_pre_comp;
2045 /* try to use the standard precomputation */
2046 g_pre_comp = &gmul[0];
2047 generator = EC_POINT_new(group);
2048 if (generator == NULL)
2050 /* get the generator from precomputation */
2051 if (!smallfelem_to_BN(x, g_pre_comp[0][1][0]) ||
2052 !smallfelem_to_BN(y, g_pre_comp[0][1][1]) ||
2053 !smallfelem_to_BN(z, g_pre_comp[0][1][2])) {
2054 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2057 if (!EC_POINT_set_Jprojective_coordinates_GFp(group,
2061 if (0 == EC_POINT_cmp(group, generator, group->generator, ctx))
2062 /* precomputation matches generator */
2066 * we don't have valid precomputation: treat the generator as a
2071 if (num_points > 0) {
2072 if (num_points >= 3) {
2074 * unless we precompute multiples for just one or two points,
2075 * converting those into affine form is time well spent
2079 secrets = OPENSSL_malloc(num_points * sizeof(felem_bytearray));
2080 pre_comp = OPENSSL_malloc(num_points * 17 * 3 * sizeof(smallfelem));
2083 OPENSSL_malloc((num_points * 17 + 1) * sizeof(smallfelem));
2084 if ((secrets == NULL) || (pre_comp == NULL)
2085 || (mixed && (tmp_smallfelems == NULL))) {
2086 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_MALLOC_FAILURE);
2091 * we treat NULL scalars as 0, and NULL points as points at infinity,
2092 * i.e., they contribute nothing to the linear combination
2094 memset(secrets, 0, num_points * sizeof(felem_bytearray));
2095 memset(pre_comp, 0, num_points * 17 * 3 * sizeof(smallfelem));
2096 for (i = 0; i < num_points; ++i) {
2099 * we didn't have a valid precomputation, so we pick the
2103 p = EC_GROUP_get0_generator(group);
2106 /* the i^th point */
2109 p_scalar = scalars[i];
2111 if ((p_scalar != NULL) && (p != NULL)) {
2112 /* reduce scalar to 0 <= scalar < 2^256 */
2113 if ((BN_num_bits(p_scalar) > 256)
2114 || (BN_is_negative(p_scalar))) {
2116 * this is an unusual input, and we don't guarantee
2119 if (!BN_nnmod(tmp_scalar, p_scalar, &group->order, ctx)) {
2120 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2123 num_bytes = BN_bn2bin(tmp_scalar, tmp);
2125 num_bytes = BN_bn2bin(p_scalar, tmp);
2126 flip_endian(secrets[i], tmp, num_bytes);
2127 /* precompute multiples */
2128 if ((!BN_to_felem(x_out, &p->X)) ||
2129 (!BN_to_felem(y_out, &p->Y)) ||
2130 (!BN_to_felem(z_out, &p->Z)))
2132 felem_shrink(pre_comp[i][1][0], x_out);
2133 felem_shrink(pre_comp[i][1][1], y_out);
2134 felem_shrink(pre_comp[i][1][2], z_out);
2135 for (j = 2; j <= 16; ++j) {
2137 point_add_small(pre_comp[i][j][0], pre_comp[i][j][1],
2138 pre_comp[i][j][2], pre_comp[i][1][0],
2139 pre_comp[i][1][1], pre_comp[i][1][2],
2140 pre_comp[i][j - 1][0],
2141 pre_comp[i][j - 1][1],
2142 pre_comp[i][j - 1][2]);
2144 point_double_small(pre_comp[i][j][0],
2147 pre_comp[i][j / 2][0],
2148 pre_comp[i][j / 2][1],
2149 pre_comp[i][j / 2][2]);
2155 make_points_affine(num_points * 17, pre_comp[0], tmp_smallfelems);
2158 /* the scalar for the generator */
2159 if ((scalar != NULL) && (have_pre_comp)) {
2160 memset(g_secret, 0, sizeof(g_secret));
2161 /* reduce scalar to 0 <= scalar < 2^256 */
2162 if ((BN_num_bits(scalar) > 256) || (BN_is_negative(scalar))) {
2164 * this is an unusual input, and we don't guarantee
2167 if (!BN_nnmod(tmp_scalar, scalar, &group->order, ctx)) {
2168 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2171 num_bytes = BN_bn2bin(tmp_scalar, tmp);
2173 num_bytes = BN_bn2bin(scalar, tmp);
2174 flip_endian(g_secret, tmp, num_bytes);
2175 /* do the multiplication with generator precomputation */
2176 batch_mul(x_out, y_out, z_out,
2177 (const felem_bytearray(*))secrets, num_points,
2179 mixed, (const smallfelem(*)[17][3])pre_comp, g_pre_comp);
2181 /* do the multiplication without generator precomputation */
2182 batch_mul(x_out, y_out, z_out,
2183 (const felem_bytearray(*))secrets, num_points,
2184 NULL, mixed, (const smallfelem(*)[17][3])pre_comp, NULL);
2185 /* reduce the output to its unique minimal representation */
2186 felem_contract(x_in, x_out);
2187 felem_contract(y_in, y_out);
2188 felem_contract(z_in, z_out);
2189 if ((!smallfelem_to_BN(x, x_in)) || (!smallfelem_to_BN(y, y_in)) ||
2190 (!smallfelem_to_BN(z, z_in))) {
2191 ECerr(EC_F_EC_GFP_NISTP256_POINTS_MUL, ERR_R_BN_LIB);
2194 ret = EC_POINT_set_Jprojective_coordinates_GFp(group, r, x, y, z, ctx);
2198 if (generator != NULL)
2199 EC_POINT_free(generator);
2200 if (new_ctx != NULL)
2201 BN_CTX_free(new_ctx);
2202 if (secrets != NULL)
2203 OPENSSL_free(secrets);
2204 if (pre_comp != NULL)
2205 OPENSSL_free(pre_comp);
2206 if (tmp_smallfelems != NULL)
2207 OPENSSL_free(tmp_smallfelems);
2211 int ec_GFp_nistp256_precompute_mult(EC_GROUP *group, BN_CTX *ctx)
2214 NISTP256_PRE_COMP *pre = NULL;
2216 BN_CTX *new_ctx = NULL;
2218 EC_POINT *generator = NULL;
2219 smallfelem tmp_smallfelems[32];
2220 felem x_tmp, y_tmp, z_tmp;
2222 /* throw away old precomputation */
2223 EC_EX_DATA_free_data(&group->extra_data, nistp256_pre_comp_dup,
2224 nistp256_pre_comp_free,
2225 nistp256_pre_comp_clear_free);
2227 if ((ctx = new_ctx = BN_CTX_new()) == NULL)
2230 if (((x = BN_CTX_get(ctx)) == NULL) || ((y = BN_CTX_get(ctx)) == NULL))
2232 /* get the generator */
2233 if (group->generator == NULL)
2235 generator = EC_POINT_new(group);
2236 if (generator == NULL)
2238 BN_bin2bn(nistp256_curve_params[3], sizeof(felem_bytearray), x);
2239 BN_bin2bn(nistp256_curve_params[4], sizeof(felem_bytearray), y);
2240 if (!EC_POINT_set_affine_coordinates_GFp(group, generator, x, y, ctx))
2242 if ((pre = nistp256_pre_comp_new()) == NULL)
2245 * if the generator is the standard one, use built-in precomputation
2247 if (0 == EC_POINT_cmp(group, generator, group->generator, ctx)) {
2248 memcpy(pre->g_pre_comp, gmul, sizeof(pre->g_pre_comp));
2252 if ((!BN_to_felem(x_tmp, &group->generator->X)) ||
2253 (!BN_to_felem(y_tmp, &group->generator->Y)) ||
2254 (!BN_to_felem(z_tmp, &group->generator->Z)))
2256 felem_shrink(pre->g_pre_comp[0][1][0], x_tmp);
2257 felem_shrink(pre->g_pre_comp[0][1][1], y_tmp);
2258 felem_shrink(pre->g_pre_comp[0][1][2], z_tmp);
2260 * compute 2^64*G, 2^128*G, 2^192*G for the first table, 2^32*G, 2^96*G,
2261 * 2^160*G, 2^224*G for the second one
2263 for (i = 1; i <= 8; i <<= 1) {
2264 point_double_small(pre->g_pre_comp[1][i][0], pre->g_pre_comp[1][i][1],
2265 pre->g_pre_comp[1][i][2], pre->g_pre_comp[0][i][0],
2266 pre->g_pre_comp[0][i][1],
2267 pre->g_pre_comp[0][i][2]);
2268 for (j = 0; j < 31; ++j) {
2269 point_double_small(pre->g_pre_comp[1][i][0],
2270 pre->g_pre_comp[1][i][1],
2271 pre->g_pre_comp[1][i][2],
2272 pre->g_pre_comp[1][i][0],
2273 pre->g_pre_comp[1][i][1],
2274 pre->g_pre_comp[1][i][2]);
2278 point_double_small(pre->g_pre_comp[0][2 * i][0],
2279 pre->g_pre_comp[0][2 * i][1],
2280 pre->g_pre_comp[0][2 * i][2],
2281 pre->g_pre_comp[1][i][0], pre->g_pre_comp[1][i][1],
2282 pre->g_pre_comp[1][i][2]);
2283 for (j = 0; j < 31; ++j) {
2284 point_double_small(pre->g_pre_comp[0][2 * i][0],
2285 pre->g_pre_comp[0][2 * i][1],
2286 pre->g_pre_comp[0][2 * i][2],
2287 pre->g_pre_comp[0][2 * i][0],
2288 pre->g_pre_comp[0][2 * i][1],
2289 pre->g_pre_comp[0][2 * i][2]);
2292 for (i = 0; i < 2; i++) {
2293 /* g_pre_comp[i][0] is the point at infinity */
2294 memset(pre->g_pre_comp[i][0], 0, sizeof(pre->g_pre_comp[i][0]));
2295 /* the remaining multiples */
2296 /* 2^64*G + 2^128*G resp. 2^96*G + 2^160*G */
2297 point_add_small(pre->g_pre_comp[i][6][0], pre->g_pre_comp[i][6][1],
2298 pre->g_pre_comp[i][6][2], pre->g_pre_comp[i][4][0],
2299 pre->g_pre_comp[i][4][1], pre->g_pre_comp[i][4][2],
2300 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2301 pre->g_pre_comp[i][2][2]);
2302 /* 2^64*G + 2^192*G resp. 2^96*G + 2^224*G */
2303 point_add_small(pre->g_pre_comp[i][10][0], pre->g_pre_comp[i][10][1],
2304 pre->g_pre_comp[i][10][2], pre->g_pre_comp[i][8][0],
2305 pre->g_pre_comp[i][8][1], pre->g_pre_comp[i][8][2],
2306 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2307 pre->g_pre_comp[i][2][2]);
2308 /* 2^128*G + 2^192*G resp. 2^160*G + 2^224*G */
2309 point_add_small(pre->g_pre_comp[i][12][0], pre->g_pre_comp[i][12][1],
2310 pre->g_pre_comp[i][12][2], pre->g_pre_comp[i][8][0],
2311 pre->g_pre_comp[i][8][1], pre->g_pre_comp[i][8][2],
2312 pre->g_pre_comp[i][4][0], pre->g_pre_comp[i][4][1],
2313 pre->g_pre_comp[i][4][2]);
2315 * 2^64*G + 2^128*G + 2^192*G resp. 2^96*G + 2^160*G + 2^224*G
2317 point_add_small(pre->g_pre_comp[i][14][0], pre->g_pre_comp[i][14][1],
2318 pre->g_pre_comp[i][14][2], pre->g_pre_comp[i][12][0],
2319 pre->g_pre_comp[i][12][1], pre->g_pre_comp[i][12][2],
2320 pre->g_pre_comp[i][2][0], pre->g_pre_comp[i][2][1],
2321 pre->g_pre_comp[i][2][2]);
2322 for (j = 1; j < 8; ++j) {
2323 /* odd multiples: add G resp. 2^32*G */
2324 point_add_small(pre->g_pre_comp[i][2 * j + 1][0],
2325 pre->g_pre_comp[i][2 * j + 1][1],
2326 pre->g_pre_comp[i][2 * j + 1][2],
2327 pre->g_pre_comp[i][2 * j][0],
2328 pre->g_pre_comp[i][2 * j][1],
2329 pre->g_pre_comp[i][2 * j][2],
2330 pre->g_pre_comp[i][1][0],
2331 pre->g_pre_comp[i][1][1],
2332 pre->g_pre_comp[i][1][2]);
2335 make_points_affine(31, &(pre->g_pre_comp[0][1]), tmp_smallfelems);
2337 if (!EC_EX_DATA_set_data(&group->extra_data, pre, nistp256_pre_comp_dup,
2338 nistp256_pre_comp_free,
2339 nistp256_pre_comp_clear_free))
2345 if (generator != NULL)
2346 EC_POINT_free(generator);
2347 if (new_ctx != NULL)
2348 BN_CTX_free(new_ctx);
2350 nistp256_pre_comp_free(pre);
2354 int ec_GFp_nistp256_have_precompute_mult(const EC_GROUP *group)
2356 if (EC_EX_DATA_get_data(group->extra_data, nistp256_pre_comp_dup,
2357 nistp256_pre_comp_free,
2358 nistp256_pre_comp_clear_free)
2365 static void *dummy = &dummy;