bn_gf2m.c 29 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167
  1. /*
  2. * Copyright 2002-2018 The OpenSSL Project Authors. All Rights Reserved.
  3. * Copyright (c) 2002, Oracle and/or its affiliates. All rights reserved
  4. *
  5. * Licensed under the Apache License 2.0 (the "License"). You may not use
  6. * this file except in compliance with the License. You can obtain a copy
  7. * in the file LICENSE in the source distribution or at
  8. * https://www.openssl.org/source/license.html
  9. */
  10. #include <assert.h>
  11. #include <limits.h>
  12. #include <stdio.h>
  13. #include "internal/cryptlib.h"
  14. #include "bn_local.h"
  15. #ifndef OPENSSL_NO_EC2M
  16. /*
  17. * Maximum number of iterations before BN_GF2m_mod_solve_quad_arr should
  18. * fail.
  19. */
  20. # define MAX_ITERATIONS 50
  21. # define SQR_nibble(w) ((((w) & 8) << 3) \
  22. | (((w) & 4) << 2) \
  23. | (((w) & 2) << 1) \
  24. | ((w) & 1))
  25. /* Platform-specific macros to accelerate squaring. */
  26. # if defined(SIXTY_FOUR_BIT) || defined(SIXTY_FOUR_BIT_LONG)
  27. # define SQR1(w) \
  28. SQR_nibble((w) >> 60) << 56 | SQR_nibble((w) >> 56) << 48 | \
  29. SQR_nibble((w) >> 52) << 40 | SQR_nibble((w) >> 48) << 32 | \
  30. SQR_nibble((w) >> 44) << 24 | SQR_nibble((w) >> 40) << 16 | \
  31. SQR_nibble((w) >> 36) << 8 | SQR_nibble((w) >> 32)
  32. # define SQR0(w) \
  33. SQR_nibble((w) >> 28) << 56 | SQR_nibble((w) >> 24) << 48 | \
  34. SQR_nibble((w) >> 20) << 40 | SQR_nibble((w) >> 16) << 32 | \
  35. SQR_nibble((w) >> 12) << 24 | SQR_nibble((w) >> 8) << 16 | \
  36. SQR_nibble((w) >> 4) << 8 | SQR_nibble((w) )
  37. # endif
  38. # ifdef THIRTY_TWO_BIT
  39. # define SQR1(w) \
  40. SQR_nibble((w) >> 28) << 24 | SQR_nibble((w) >> 24) << 16 | \
  41. SQR_nibble((w) >> 20) << 8 | SQR_nibble((w) >> 16)
  42. # define SQR0(w) \
  43. SQR_nibble((w) >> 12) << 24 | SQR_nibble((w) >> 8) << 16 | \
  44. SQR_nibble((w) >> 4) << 8 | SQR_nibble((w) )
  45. # endif
  46. # if !defined(OPENSSL_BN_ASM_GF2m)
  47. /*
  48. * Product of two polynomials a, b each with degree < BN_BITS2 - 1, result is
  49. * a polynomial r with degree < 2 * BN_BITS - 1 The caller MUST ensure that
  50. * the variables have the right amount of space allocated.
  51. */
  52. # ifdef THIRTY_TWO_BIT
  53. static void bn_GF2m_mul_1x1(BN_ULONG *r1, BN_ULONG *r0, const BN_ULONG a,
  54. const BN_ULONG b)
  55. {
  56. register BN_ULONG h, l, s;
  57. BN_ULONG tab[8], top2b = a >> 30;
  58. register BN_ULONG a1, a2, a4;
  59. a1 = a & (0x3FFFFFFF);
  60. a2 = a1 << 1;
  61. a4 = a2 << 1;
  62. tab[0] = 0;
  63. tab[1] = a1;
  64. tab[2] = a2;
  65. tab[3] = a1 ^ a2;
  66. tab[4] = a4;
  67. tab[5] = a1 ^ a4;
  68. tab[6] = a2 ^ a4;
  69. tab[7] = a1 ^ a2 ^ a4;
  70. s = tab[b & 0x7];
  71. l = s;
  72. s = tab[b >> 3 & 0x7];
  73. l ^= s << 3;
  74. h = s >> 29;
  75. s = tab[b >> 6 & 0x7];
  76. l ^= s << 6;
  77. h ^= s >> 26;
  78. s = tab[b >> 9 & 0x7];
  79. l ^= s << 9;
  80. h ^= s >> 23;
  81. s = tab[b >> 12 & 0x7];
  82. l ^= s << 12;
  83. h ^= s >> 20;
  84. s = tab[b >> 15 & 0x7];
  85. l ^= s << 15;
  86. h ^= s >> 17;
  87. s = tab[b >> 18 & 0x7];
  88. l ^= s << 18;
  89. h ^= s >> 14;
  90. s = tab[b >> 21 & 0x7];
  91. l ^= s << 21;
  92. h ^= s >> 11;
  93. s = tab[b >> 24 & 0x7];
  94. l ^= s << 24;
  95. h ^= s >> 8;
  96. s = tab[b >> 27 & 0x7];
  97. l ^= s << 27;
  98. h ^= s >> 5;
  99. s = tab[b >> 30];
  100. l ^= s << 30;
  101. h ^= s >> 2;
  102. /* compensate for the top two bits of a */
  103. if (top2b & 01) {
  104. l ^= b << 30;
  105. h ^= b >> 2;
  106. }
  107. if (top2b & 02) {
  108. l ^= b << 31;
  109. h ^= b >> 1;
  110. }
  111. *r1 = h;
  112. *r0 = l;
  113. }
  114. # endif
  115. # if defined(SIXTY_FOUR_BIT) || defined(SIXTY_FOUR_BIT_LONG)
  116. static void bn_GF2m_mul_1x1(BN_ULONG *r1, BN_ULONG *r0, const BN_ULONG a,
  117. const BN_ULONG b)
  118. {
  119. register BN_ULONG h, l, s;
  120. BN_ULONG tab[16], top3b = a >> 61;
  121. register BN_ULONG a1, a2, a4, a8;
  122. a1 = a & (0x1FFFFFFFFFFFFFFFULL);
  123. a2 = a1 << 1;
  124. a4 = a2 << 1;
  125. a8 = a4 << 1;
  126. tab[0] = 0;
  127. tab[1] = a1;
  128. tab[2] = a2;
  129. tab[3] = a1 ^ a2;
  130. tab[4] = a4;
  131. tab[5] = a1 ^ a4;
  132. tab[6] = a2 ^ a4;
  133. tab[7] = a1 ^ a2 ^ a4;
  134. tab[8] = a8;
  135. tab[9] = a1 ^ a8;
  136. tab[10] = a2 ^ a8;
  137. tab[11] = a1 ^ a2 ^ a8;
  138. tab[12] = a4 ^ a8;
  139. tab[13] = a1 ^ a4 ^ a8;
  140. tab[14] = a2 ^ a4 ^ a8;
  141. tab[15] = a1 ^ a2 ^ a4 ^ a8;
  142. s = tab[b & 0xF];
  143. l = s;
  144. s = tab[b >> 4 & 0xF];
  145. l ^= s << 4;
  146. h = s >> 60;
  147. s = tab[b >> 8 & 0xF];
  148. l ^= s << 8;
  149. h ^= s >> 56;
  150. s = tab[b >> 12 & 0xF];
  151. l ^= s << 12;
  152. h ^= s >> 52;
  153. s = tab[b >> 16 & 0xF];
  154. l ^= s << 16;
  155. h ^= s >> 48;
  156. s = tab[b >> 20 & 0xF];
  157. l ^= s << 20;
  158. h ^= s >> 44;
  159. s = tab[b >> 24 & 0xF];
  160. l ^= s << 24;
  161. h ^= s >> 40;
  162. s = tab[b >> 28 & 0xF];
  163. l ^= s << 28;
  164. h ^= s >> 36;
  165. s = tab[b >> 32 & 0xF];
  166. l ^= s << 32;
  167. h ^= s >> 32;
  168. s = tab[b >> 36 & 0xF];
  169. l ^= s << 36;
  170. h ^= s >> 28;
  171. s = tab[b >> 40 & 0xF];
  172. l ^= s << 40;
  173. h ^= s >> 24;
  174. s = tab[b >> 44 & 0xF];
  175. l ^= s << 44;
  176. h ^= s >> 20;
  177. s = tab[b >> 48 & 0xF];
  178. l ^= s << 48;
  179. h ^= s >> 16;
  180. s = tab[b >> 52 & 0xF];
  181. l ^= s << 52;
  182. h ^= s >> 12;
  183. s = tab[b >> 56 & 0xF];
  184. l ^= s << 56;
  185. h ^= s >> 8;
  186. s = tab[b >> 60];
  187. l ^= s << 60;
  188. h ^= s >> 4;
  189. /* compensate for the top three bits of a */
  190. if (top3b & 01) {
  191. l ^= b << 61;
  192. h ^= b >> 3;
  193. }
  194. if (top3b & 02) {
  195. l ^= b << 62;
  196. h ^= b >> 2;
  197. }
  198. if (top3b & 04) {
  199. l ^= b << 63;
  200. h ^= b >> 1;
  201. }
  202. *r1 = h;
  203. *r0 = l;
  204. }
  205. # endif
  206. /*
  207. * Product of two polynomials a, b each with degree < 2 * BN_BITS2 - 1,
  208. * result is a polynomial r with degree < 4 * BN_BITS2 - 1 The caller MUST
  209. * ensure that the variables have the right amount of space allocated.
  210. */
  211. static void bn_GF2m_mul_2x2(BN_ULONG *r, const BN_ULONG a1, const BN_ULONG a0,
  212. const BN_ULONG b1, const BN_ULONG b0)
  213. {
  214. BN_ULONG m1, m0;
  215. /* r[3] = h1, r[2] = h0; r[1] = l1; r[0] = l0 */
  216. bn_GF2m_mul_1x1(r + 3, r + 2, a1, b1);
  217. bn_GF2m_mul_1x1(r + 1, r, a0, b0);
  218. bn_GF2m_mul_1x1(&m1, &m0, a0 ^ a1, b0 ^ b1);
  219. /* Correction on m1 ^= l1 ^ h1; m0 ^= l0 ^ h0; */
  220. r[2] ^= m1 ^ r[1] ^ r[3]; /* h0 ^= m1 ^ l1 ^ h1; */
  221. r[1] = r[3] ^ r[2] ^ r[0] ^ m1 ^ m0; /* l1 ^= l0 ^ h0 ^ m0; */
  222. }
  223. # else
  224. void bn_GF2m_mul_2x2(BN_ULONG *r, BN_ULONG a1, BN_ULONG a0, BN_ULONG b1,
  225. BN_ULONG b0);
  226. # endif
  227. /*
  228. * Add polynomials a and b and store result in r; r could be a or b, a and b
  229. * could be equal; r is the bitwise XOR of a and b.
  230. */
  231. int BN_GF2m_add(BIGNUM *r, const BIGNUM *a, const BIGNUM *b)
  232. {
  233. int i;
  234. const BIGNUM *at, *bt;
  235. bn_check_top(a);
  236. bn_check_top(b);
  237. if (a->top < b->top) {
  238. at = b;
  239. bt = a;
  240. } else {
  241. at = a;
  242. bt = b;
  243. }
  244. if (bn_wexpand(r, at->top) == NULL)
  245. return 0;
  246. for (i = 0; i < bt->top; i++) {
  247. r->d[i] = at->d[i] ^ bt->d[i];
  248. }
  249. for (; i < at->top; i++) {
  250. r->d[i] = at->d[i];
  251. }
  252. r->top = at->top;
  253. bn_correct_top(r);
  254. return 1;
  255. }
  256. /*-
  257. * Some functions allow for representation of the irreducible polynomials
  258. * as an int[], say p. The irreducible f(t) is then of the form:
  259. * t^p[0] + t^p[1] + ... + t^p[k]
  260. * where m = p[0] > p[1] > ... > p[k] = 0.
  261. */
  262. /* Performs modular reduction of a and store result in r. r could be a. */
  263. int BN_GF2m_mod_arr(BIGNUM *r, const BIGNUM *a, const int p[])
  264. {
  265. int j, k;
  266. int n, dN, d0, d1;
  267. BN_ULONG zz, *z;
  268. bn_check_top(a);
  269. if (p[0] == 0) {
  270. /* reduction mod 1 => return 0 */
  271. BN_zero(r);
  272. return 1;
  273. }
  274. /*
  275. * Since the algorithm does reduction in the r value, if a != r, copy the
  276. * contents of a into r so we can do reduction in r.
  277. */
  278. if (a != r) {
  279. if (!bn_wexpand(r, a->top))
  280. return 0;
  281. for (j = 0; j < a->top; j++) {
  282. r->d[j] = a->d[j];
  283. }
  284. r->top = a->top;
  285. }
  286. z = r->d;
  287. /* start reduction */
  288. dN = p[0] / BN_BITS2;
  289. for (j = r->top - 1; j > dN;) {
  290. zz = z[j];
  291. if (z[j] == 0) {
  292. j--;
  293. continue;
  294. }
  295. z[j] = 0;
  296. for (k = 1; p[k] != 0; k++) {
  297. /* reducing component t^p[k] */
  298. n = p[0] - p[k];
  299. d0 = n % BN_BITS2;
  300. d1 = BN_BITS2 - d0;
  301. n /= BN_BITS2;
  302. z[j - n] ^= (zz >> d0);
  303. if (d0)
  304. z[j - n - 1] ^= (zz << d1);
  305. }
  306. /* reducing component t^0 */
  307. n = dN;
  308. d0 = p[0] % BN_BITS2;
  309. d1 = BN_BITS2 - d0;
  310. z[j - n] ^= (zz >> d0);
  311. if (d0)
  312. z[j - n - 1] ^= (zz << d1);
  313. }
  314. /* final round of reduction */
  315. while (j == dN) {
  316. d0 = p[0] % BN_BITS2;
  317. zz = z[dN] >> d0;
  318. if (zz == 0)
  319. break;
  320. d1 = BN_BITS2 - d0;
  321. /* clear up the top d1 bits */
  322. if (d0)
  323. z[dN] = (z[dN] << d1) >> d1;
  324. else
  325. z[dN] = 0;
  326. z[0] ^= zz; /* reduction t^0 component */
  327. for (k = 1; p[k] != 0; k++) {
  328. BN_ULONG tmp_ulong;
  329. /* reducing component t^p[k] */
  330. n = p[k] / BN_BITS2;
  331. d0 = p[k] % BN_BITS2;
  332. d1 = BN_BITS2 - d0;
  333. z[n] ^= (zz << d0);
  334. if (d0 && (tmp_ulong = zz >> d1))
  335. z[n + 1] ^= tmp_ulong;
  336. }
  337. }
  338. bn_correct_top(r);
  339. return 1;
  340. }
  341. /*
  342. * Performs modular reduction of a by p and store result in r. r could be a.
  343. * This function calls down to the BN_GF2m_mod_arr implementation; this wrapper
  344. * function is only provided for convenience; for best performance, use the
  345. * BN_GF2m_mod_arr function.
  346. */
  347. int BN_GF2m_mod(BIGNUM *r, const BIGNUM *a, const BIGNUM *p)
  348. {
  349. int ret = 0;
  350. int arr[6];
  351. bn_check_top(a);
  352. bn_check_top(p);
  353. ret = BN_GF2m_poly2arr(p, arr, OSSL_NELEM(arr));
  354. if (!ret || ret > (int)OSSL_NELEM(arr)) {
  355. BNerr(BN_F_BN_GF2M_MOD, BN_R_INVALID_LENGTH);
  356. return 0;
  357. }
  358. ret = BN_GF2m_mod_arr(r, a, arr);
  359. bn_check_top(r);
  360. return ret;
  361. }
  362. /*
  363. * Compute the product of two polynomials a and b, reduce modulo p, and store
  364. * the result in r. r could be a or b; a could be b.
  365. */
  366. int BN_GF2m_mod_mul_arr(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  367. const int p[], BN_CTX *ctx)
  368. {
  369. int zlen, i, j, k, ret = 0;
  370. BIGNUM *s;
  371. BN_ULONG x1, x0, y1, y0, zz[4];
  372. bn_check_top(a);
  373. bn_check_top(b);
  374. if (a == b) {
  375. return BN_GF2m_mod_sqr_arr(r, a, p, ctx);
  376. }
  377. BN_CTX_start(ctx);
  378. if ((s = BN_CTX_get(ctx)) == NULL)
  379. goto err;
  380. zlen = a->top + b->top + 4;
  381. if (!bn_wexpand(s, zlen))
  382. goto err;
  383. s->top = zlen;
  384. for (i = 0; i < zlen; i++)
  385. s->d[i] = 0;
  386. for (j = 0; j < b->top; j += 2) {
  387. y0 = b->d[j];
  388. y1 = ((j + 1) == b->top) ? 0 : b->d[j + 1];
  389. for (i = 0; i < a->top; i += 2) {
  390. x0 = a->d[i];
  391. x1 = ((i + 1) == a->top) ? 0 : a->d[i + 1];
  392. bn_GF2m_mul_2x2(zz, x1, x0, y1, y0);
  393. for (k = 0; k < 4; k++)
  394. s->d[i + j + k] ^= zz[k];
  395. }
  396. }
  397. bn_correct_top(s);
  398. if (BN_GF2m_mod_arr(r, s, p))
  399. ret = 1;
  400. bn_check_top(r);
  401. err:
  402. BN_CTX_end(ctx);
  403. return ret;
  404. }
  405. /*
  406. * Compute the product of two polynomials a and b, reduce modulo p, and store
  407. * the result in r. r could be a or b; a could equal b. This function calls
  408. * down to the BN_GF2m_mod_mul_arr implementation; this wrapper function is
  409. * only provided for convenience; for best performance, use the
  410. * BN_GF2m_mod_mul_arr function.
  411. */
  412. int BN_GF2m_mod_mul(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  413. const BIGNUM *p, BN_CTX *ctx)
  414. {
  415. int ret = 0;
  416. const int max = BN_num_bits(p) + 1;
  417. int *arr = NULL;
  418. bn_check_top(a);
  419. bn_check_top(b);
  420. bn_check_top(p);
  421. if ((arr = OPENSSL_malloc(sizeof(*arr) * max)) == NULL)
  422. goto err;
  423. ret = BN_GF2m_poly2arr(p, arr, max);
  424. if (!ret || ret > max) {
  425. BNerr(BN_F_BN_GF2M_MOD_MUL, BN_R_INVALID_LENGTH);
  426. goto err;
  427. }
  428. ret = BN_GF2m_mod_mul_arr(r, a, b, arr, ctx);
  429. bn_check_top(r);
  430. err:
  431. OPENSSL_free(arr);
  432. return ret;
  433. }
  434. /* Square a, reduce the result mod p, and store it in a. r could be a. */
  435. int BN_GF2m_mod_sqr_arr(BIGNUM *r, const BIGNUM *a, const int p[],
  436. BN_CTX *ctx)
  437. {
  438. int i, ret = 0;
  439. BIGNUM *s;
  440. bn_check_top(a);
  441. BN_CTX_start(ctx);
  442. if ((s = BN_CTX_get(ctx)) == NULL)
  443. goto err;
  444. if (!bn_wexpand(s, 2 * a->top))
  445. goto err;
  446. for (i = a->top - 1; i >= 0; i--) {
  447. s->d[2 * i + 1] = SQR1(a->d[i]);
  448. s->d[2 * i] = SQR0(a->d[i]);
  449. }
  450. s->top = 2 * a->top;
  451. bn_correct_top(s);
  452. if (!BN_GF2m_mod_arr(r, s, p))
  453. goto err;
  454. bn_check_top(r);
  455. ret = 1;
  456. err:
  457. BN_CTX_end(ctx);
  458. return ret;
  459. }
  460. /*
  461. * Square a, reduce the result mod p, and store it in a. r could be a. This
  462. * function calls down to the BN_GF2m_mod_sqr_arr implementation; this
  463. * wrapper function is only provided for convenience; for best performance,
  464. * use the BN_GF2m_mod_sqr_arr function.
  465. */
  466. int BN_GF2m_mod_sqr(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  467. {
  468. int ret = 0;
  469. const int max = BN_num_bits(p) + 1;
  470. int *arr = NULL;
  471. bn_check_top(a);
  472. bn_check_top(p);
  473. if ((arr = OPENSSL_malloc(sizeof(*arr) * max)) == NULL)
  474. goto err;
  475. ret = BN_GF2m_poly2arr(p, arr, max);
  476. if (!ret || ret > max) {
  477. BNerr(BN_F_BN_GF2M_MOD_SQR, BN_R_INVALID_LENGTH);
  478. goto err;
  479. }
  480. ret = BN_GF2m_mod_sqr_arr(r, a, arr, ctx);
  481. bn_check_top(r);
  482. err:
  483. OPENSSL_free(arr);
  484. return ret;
  485. }
  486. /*
  487. * Invert a, reduce modulo p, and store the result in r. r could be a. Uses
  488. * Modified Almost Inverse Algorithm (Algorithm 10) from Hankerson, D.,
  489. * Hernandez, J.L., and Menezes, A. "Software Implementation of Elliptic
  490. * Curve Cryptography Over Binary Fields".
  491. */
  492. static int BN_GF2m_mod_inv_vartime(BIGNUM *r, const BIGNUM *a,
  493. const BIGNUM *p, BN_CTX *ctx)
  494. {
  495. BIGNUM *b, *c = NULL, *u = NULL, *v = NULL, *tmp;
  496. int ret = 0;
  497. bn_check_top(a);
  498. bn_check_top(p);
  499. BN_CTX_start(ctx);
  500. b = BN_CTX_get(ctx);
  501. c = BN_CTX_get(ctx);
  502. u = BN_CTX_get(ctx);
  503. v = BN_CTX_get(ctx);
  504. if (v == NULL)
  505. goto err;
  506. if (!BN_GF2m_mod(u, a, p))
  507. goto err;
  508. if (BN_is_zero(u))
  509. goto err;
  510. if (!BN_copy(v, p))
  511. goto err;
  512. # if 0
  513. if (!BN_one(b))
  514. goto err;
  515. while (1) {
  516. while (!BN_is_odd(u)) {
  517. if (BN_is_zero(u))
  518. goto err;
  519. if (!BN_rshift1(u, u))
  520. goto err;
  521. if (BN_is_odd(b)) {
  522. if (!BN_GF2m_add(b, b, p))
  523. goto err;
  524. }
  525. if (!BN_rshift1(b, b))
  526. goto err;
  527. }
  528. if (BN_abs_is_word(u, 1))
  529. break;
  530. if (BN_num_bits(u) < BN_num_bits(v)) {
  531. tmp = u;
  532. u = v;
  533. v = tmp;
  534. tmp = b;
  535. b = c;
  536. c = tmp;
  537. }
  538. if (!BN_GF2m_add(u, u, v))
  539. goto err;
  540. if (!BN_GF2m_add(b, b, c))
  541. goto err;
  542. }
  543. # else
  544. {
  545. int i;
  546. int ubits = BN_num_bits(u);
  547. int vbits = BN_num_bits(v); /* v is copy of p */
  548. int top = p->top;
  549. BN_ULONG *udp, *bdp, *vdp, *cdp;
  550. if (!bn_wexpand(u, top))
  551. goto err;
  552. udp = u->d;
  553. for (i = u->top; i < top; i++)
  554. udp[i] = 0;
  555. u->top = top;
  556. if (!bn_wexpand(b, top))
  557. goto err;
  558. bdp = b->d;
  559. bdp[0] = 1;
  560. for (i = 1; i < top; i++)
  561. bdp[i] = 0;
  562. b->top = top;
  563. if (!bn_wexpand(c, top))
  564. goto err;
  565. cdp = c->d;
  566. for (i = 0; i < top; i++)
  567. cdp[i] = 0;
  568. c->top = top;
  569. vdp = v->d; /* It pays off to "cache" *->d pointers,
  570. * because it allows optimizer to be more
  571. * aggressive. But we don't have to "cache"
  572. * p->d, because *p is declared 'const'... */
  573. while (1) {
  574. while (ubits && !(udp[0] & 1)) {
  575. BN_ULONG u0, u1, b0, b1, mask;
  576. u0 = udp[0];
  577. b0 = bdp[0];
  578. mask = (BN_ULONG)0 - (b0 & 1);
  579. b0 ^= p->d[0] & mask;
  580. for (i = 0; i < top - 1; i++) {
  581. u1 = udp[i + 1];
  582. udp[i] = ((u0 >> 1) | (u1 << (BN_BITS2 - 1))) & BN_MASK2;
  583. u0 = u1;
  584. b1 = bdp[i + 1] ^ (p->d[i + 1] & mask);
  585. bdp[i] = ((b0 >> 1) | (b1 << (BN_BITS2 - 1))) & BN_MASK2;
  586. b0 = b1;
  587. }
  588. udp[i] = u0 >> 1;
  589. bdp[i] = b0 >> 1;
  590. ubits--;
  591. }
  592. if (ubits <= BN_BITS2) {
  593. if (udp[0] == 0) /* poly was reducible */
  594. goto err;
  595. if (udp[0] == 1)
  596. break;
  597. }
  598. if (ubits < vbits) {
  599. i = ubits;
  600. ubits = vbits;
  601. vbits = i;
  602. tmp = u;
  603. u = v;
  604. v = tmp;
  605. tmp = b;
  606. b = c;
  607. c = tmp;
  608. udp = vdp;
  609. vdp = v->d;
  610. bdp = cdp;
  611. cdp = c->d;
  612. }
  613. for (i = 0; i < top; i++) {
  614. udp[i] ^= vdp[i];
  615. bdp[i] ^= cdp[i];
  616. }
  617. if (ubits == vbits) {
  618. BN_ULONG ul;
  619. int utop = (ubits - 1) / BN_BITS2;
  620. while ((ul = udp[utop]) == 0 && utop)
  621. utop--;
  622. ubits = utop * BN_BITS2 + BN_num_bits_word(ul);
  623. }
  624. }
  625. bn_correct_top(b);
  626. }
  627. # endif
  628. if (!BN_copy(r, b))
  629. goto err;
  630. bn_check_top(r);
  631. ret = 1;
  632. err:
  633. # ifdef BN_DEBUG /* BN_CTX_end would complain about the
  634. * expanded form */
  635. bn_correct_top(c);
  636. bn_correct_top(u);
  637. bn_correct_top(v);
  638. # endif
  639. BN_CTX_end(ctx);
  640. return ret;
  641. }
  642. /*-
  643. * Wrapper for BN_GF2m_mod_inv_vartime that blinds the input before calling.
  644. * This is not constant time.
  645. * But it does eliminate first order deduction on the input.
  646. */
  647. int BN_GF2m_mod_inv(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  648. {
  649. BIGNUM *b = NULL;
  650. int ret = 0;
  651. BN_CTX_start(ctx);
  652. if ((b = BN_CTX_get(ctx)) == NULL)
  653. goto err;
  654. /* generate blinding value */
  655. do {
  656. if (!BN_priv_rand_ex(b, BN_num_bits(p) - 1,
  657. BN_RAND_TOP_ANY, BN_RAND_BOTTOM_ANY, ctx))
  658. goto err;
  659. } while (BN_is_zero(b));
  660. /* r := a * b */
  661. if (!BN_GF2m_mod_mul(r, a, b, p, ctx))
  662. goto err;
  663. /* r := 1/(a * b) */
  664. if (!BN_GF2m_mod_inv_vartime(r, r, p, ctx))
  665. goto err;
  666. /* r := b/(a * b) = 1/a */
  667. if (!BN_GF2m_mod_mul(r, r, b, p, ctx))
  668. goto err;
  669. ret = 1;
  670. err:
  671. BN_CTX_end(ctx);
  672. return ret;
  673. }
  674. /*
  675. * Invert xx, reduce modulo p, and store the result in r. r could be xx.
  676. * This function calls down to the BN_GF2m_mod_inv implementation; this
  677. * wrapper function is only provided for convenience; for best performance,
  678. * use the BN_GF2m_mod_inv function.
  679. */
  680. int BN_GF2m_mod_inv_arr(BIGNUM *r, const BIGNUM *xx, const int p[],
  681. BN_CTX *ctx)
  682. {
  683. BIGNUM *field;
  684. int ret = 0;
  685. bn_check_top(xx);
  686. BN_CTX_start(ctx);
  687. if ((field = BN_CTX_get(ctx)) == NULL)
  688. goto err;
  689. if (!BN_GF2m_arr2poly(p, field))
  690. goto err;
  691. ret = BN_GF2m_mod_inv(r, xx, field, ctx);
  692. bn_check_top(r);
  693. err:
  694. BN_CTX_end(ctx);
  695. return ret;
  696. }
  697. /*
  698. * Divide y by x, reduce modulo p, and store the result in r. r could be x
  699. * or y, x could equal y.
  700. */
  701. int BN_GF2m_mod_div(BIGNUM *r, const BIGNUM *y, const BIGNUM *x,
  702. const BIGNUM *p, BN_CTX *ctx)
  703. {
  704. BIGNUM *xinv = NULL;
  705. int ret = 0;
  706. bn_check_top(y);
  707. bn_check_top(x);
  708. bn_check_top(p);
  709. BN_CTX_start(ctx);
  710. xinv = BN_CTX_get(ctx);
  711. if (xinv == NULL)
  712. goto err;
  713. if (!BN_GF2m_mod_inv(xinv, x, p, ctx))
  714. goto err;
  715. if (!BN_GF2m_mod_mul(r, y, xinv, p, ctx))
  716. goto err;
  717. bn_check_top(r);
  718. ret = 1;
  719. err:
  720. BN_CTX_end(ctx);
  721. return ret;
  722. }
  723. /*
  724. * Divide yy by xx, reduce modulo p, and store the result in r. r could be xx
  725. * * or yy, xx could equal yy. This function calls down to the
  726. * BN_GF2m_mod_div implementation; this wrapper function is only provided for
  727. * convenience; for best performance, use the BN_GF2m_mod_div function.
  728. */
  729. int BN_GF2m_mod_div_arr(BIGNUM *r, const BIGNUM *yy, const BIGNUM *xx,
  730. const int p[], BN_CTX *ctx)
  731. {
  732. BIGNUM *field;
  733. int ret = 0;
  734. bn_check_top(yy);
  735. bn_check_top(xx);
  736. BN_CTX_start(ctx);
  737. if ((field = BN_CTX_get(ctx)) == NULL)
  738. goto err;
  739. if (!BN_GF2m_arr2poly(p, field))
  740. goto err;
  741. ret = BN_GF2m_mod_div(r, yy, xx, field, ctx);
  742. bn_check_top(r);
  743. err:
  744. BN_CTX_end(ctx);
  745. return ret;
  746. }
  747. /*
  748. * Compute the bth power of a, reduce modulo p, and store the result in r. r
  749. * could be a. Uses simple square-and-multiply algorithm A.5.1 from IEEE
  750. * P1363.
  751. */
  752. int BN_GF2m_mod_exp_arr(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  753. const int p[], BN_CTX *ctx)
  754. {
  755. int ret = 0, i, n;
  756. BIGNUM *u;
  757. bn_check_top(a);
  758. bn_check_top(b);
  759. if (BN_is_zero(b))
  760. return BN_one(r);
  761. if (BN_abs_is_word(b, 1))
  762. return (BN_copy(r, a) != NULL);
  763. BN_CTX_start(ctx);
  764. if ((u = BN_CTX_get(ctx)) == NULL)
  765. goto err;
  766. if (!BN_GF2m_mod_arr(u, a, p))
  767. goto err;
  768. n = BN_num_bits(b) - 1;
  769. for (i = n - 1; i >= 0; i--) {
  770. if (!BN_GF2m_mod_sqr_arr(u, u, p, ctx))
  771. goto err;
  772. if (BN_is_bit_set(b, i)) {
  773. if (!BN_GF2m_mod_mul_arr(u, u, a, p, ctx))
  774. goto err;
  775. }
  776. }
  777. if (!BN_copy(r, u))
  778. goto err;
  779. bn_check_top(r);
  780. ret = 1;
  781. err:
  782. BN_CTX_end(ctx);
  783. return ret;
  784. }
  785. /*
  786. * Compute the bth power of a, reduce modulo p, and store the result in r. r
  787. * could be a. This function calls down to the BN_GF2m_mod_exp_arr
  788. * implementation; this wrapper function is only provided for convenience;
  789. * for best performance, use the BN_GF2m_mod_exp_arr function.
  790. */
  791. int BN_GF2m_mod_exp(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  792. const BIGNUM *p, BN_CTX *ctx)
  793. {
  794. int ret = 0;
  795. const int max = BN_num_bits(p) + 1;
  796. int *arr = NULL;
  797. bn_check_top(a);
  798. bn_check_top(b);
  799. bn_check_top(p);
  800. if ((arr = OPENSSL_malloc(sizeof(*arr) * max)) == NULL)
  801. goto err;
  802. ret = BN_GF2m_poly2arr(p, arr, max);
  803. if (!ret || ret > max) {
  804. BNerr(BN_F_BN_GF2M_MOD_EXP, BN_R_INVALID_LENGTH);
  805. goto err;
  806. }
  807. ret = BN_GF2m_mod_exp_arr(r, a, b, arr, ctx);
  808. bn_check_top(r);
  809. err:
  810. OPENSSL_free(arr);
  811. return ret;
  812. }
  813. /*
  814. * Compute the square root of a, reduce modulo p, and store the result in r.
  815. * r could be a. Uses exponentiation as in algorithm A.4.1 from IEEE P1363.
  816. */
  817. int BN_GF2m_mod_sqrt_arr(BIGNUM *r, const BIGNUM *a, const int p[],
  818. BN_CTX *ctx)
  819. {
  820. int ret = 0;
  821. BIGNUM *u;
  822. bn_check_top(a);
  823. if (p[0] == 0) {
  824. /* reduction mod 1 => return 0 */
  825. BN_zero(r);
  826. return 1;
  827. }
  828. BN_CTX_start(ctx);
  829. if ((u = BN_CTX_get(ctx)) == NULL)
  830. goto err;
  831. if (!BN_set_bit(u, p[0] - 1))
  832. goto err;
  833. ret = BN_GF2m_mod_exp_arr(r, a, u, p, ctx);
  834. bn_check_top(r);
  835. err:
  836. BN_CTX_end(ctx);
  837. return ret;
  838. }
  839. /*
  840. * Compute the square root of a, reduce modulo p, and store the result in r.
  841. * r could be a. This function calls down to the BN_GF2m_mod_sqrt_arr
  842. * implementation; this wrapper function is only provided for convenience;
  843. * for best performance, use the BN_GF2m_mod_sqrt_arr function.
  844. */
  845. int BN_GF2m_mod_sqrt(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  846. {
  847. int ret = 0;
  848. const int max = BN_num_bits(p) + 1;
  849. int *arr = NULL;
  850. bn_check_top(a);
  851. bn_check_top(p);
  852. if ((arr = OPENSSL_malloc(sizeof(*arr) * max)) == NULL)
  853. goto err;
  854. ret = BN_GF2m_poly2arr(p, arr, max);
  855. if (!ret || ret > max) {
  856. BNerr(BN_F_BN_GF2M_MOD_SQRT, BN_R_INVALID_LENGTH);
  857. goto err;
  858. }
  859. ret = BN_GF2m_mod_sqrt_arr(r, a, arr, ctx);
  860. bn_check_top(r);
  861. err:
  862. OPENSSL_free(arr);
  863. return ret;
  864. }
  865. /*
  866. * Find r such that r^2 + r = a mod p. r could be a. If no r exists returns
  867. * 0. Uses algorithms A.4.7 and A.4.6 from IEEE P1363.
  868. */
  869. int BN_GF2m_mod_solve_quad_arr(BIGNUM *r, const BIGNUM *a_, const int p[],
  870. BN_CTX *ctx)
  871. {
  872. int ret = 0, count = 0, j;
  873. BIGNUM *a, *z, *rho, *w, *w2, *tmp;
  874. bn_check_top(a_);
  875. if (p[0] == 0) {
  876. /* reduction mod 1 => return 0 */
  877. BN_zero(r);
  878. return 1;
  879. }
  880. BN_CTX_start(ctx);
  881. a = BN_CTX_get(ctx);
  882. z = BN_CTX_get(ctx);
  883. w = BN_CTX_get(ctx);
  884. if (w == NULL)
  885. goto err;
  886. if (!BN_GF2m_mod_arr(a, a_, p))
  887. goto err;
  888. if (BN_is_zero(a)) {
  889. BN_zero(r);
  890. ret = 1;
  891. goto err;
  892. }
  893. if (p[0] & 0x1) { /* m is odd */
  894. /* compute half-trace of a */
  895. if (!BN_copy(z, a))
  896. goto err;
  897. for (j = 1; j <= (p[0] - 1) / 2; j++) {
  898. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  899. goto err;
  900. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  901. goto err;
  902. if (!BN_GF2m_add(z, z, a))
  903. goto err;
  904. }
  905. } else { /* m is even */
  906. rho = BN_CTX_get(ctx);
  907. w2 = BN_CTX_get(ctx);
  908. tmp = BN_CTX_get(ctx);
  909. if (tmp == NULL)
  910. goto err;
  911. do {
  912. if (!BN_priv_rand_ex(rho, p[0], BN_RAND_TOP_ONE, BN_RAND_BOTTOM_ANY,
  913. ctx))
  914. goto err;
  915. if (!BN_GF2m_mod_arr(rho, rho, p))
  916. goto err;
  917. BN_zero(z);
  918. if (!BN_copy(w, rho))
  919. goto err;
  920. for (j = 1; j <= p[0] - 1; j++) {
  921. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  922. goto err;
  923. if (!BN_GF2m_mod_sqr_arr(w2, w, p, ctx))
  924. goto err;
  925. if (!BN_GF2m_mod_mul_arr(tmp, w2, a, p, ctx))
  926. goto err;
  927. if (!BN_GF2m_add(z, z, tmp))
  928. goto err;
  929. if (!BN_GF2m_add(w, w2, rho))
  930. goto err;
  931. }
  932. count++;
  933. } while (BN_is_zero(w) && (count < MAX_ITERATIONS));
  934. if (BN_is_zero(w)) {
  935. BNerr(BN_F_BN_GF2M_MOD_SOLVE_QUAD_ARR, BN_R_TOO_MANY_ITERATIONS);
  936. goto err;
  937. }
  938. }
  939. if (!BN_GF2m_mod_sqr_arr(w, z, p, ctx))
  940. goto err;
  941. if (!BN_GF2m_add(w, z, w))
  942. goto err;
  943. if (BN_GF2m_cmp(w, a)) {
  944. BNerr(BN_F_BN_GF2M_MOD_SOLVE_QUAD_ARR, BN_R_NO_SOLUTION);
  945. goto err;
  946. }
  947. if (!BN_copy(r, z))
  948. goto err;
  949. bn_check_top(r);
  950. ret = 1;
  951. err:
  952. BN_CTX_end(ctx);
  953. return ret;
  954. }
  955. /*
  956. * Find r such that r^2 + r = a mod p. r could be a. If no r exists returns
  957. * 0. This function calls down to the BN_GF2m_mod_solve_quad_arr
  958. * implementation; this wrapper function is only provided for convenience;
  959. * for best performance, use the BN_GF2m_mod_solve_quad_arr function.
  960. */
  961. int BN_GF2m_mod_solve_quad(BIGNUM *r, const BIGNUM *a, const BIGNUM *p,
  962. BN_CTX *ctx)
  963. {
  964. int ret = 0;
  965. const int max = BN_num_bits(p) + 1;
  966. int *arr = NULL;
  967. bn_check_top(a);
  968. bn_check_top(p);
  969. if ((arr = OPENSSL_malloc(sizeof(*arr) * max)) == NULL)
  970. goto err;
  971. ret = BN_GF2m_poly2arr(p, arr, max);
  972. if (!ret || ret > max) {
  973. BNerr(BN_F_BN_GF2M_MOD_SOLVE_QUAD, BN_R_INVALID_LENGTH);
  974. goto err;
  975. }
  976. ret = BN_GF2m_mod_solve_quad_arr(r, a, arr, ctx);
  977. bn_check_top(r);
  978. err:
  979. OPENSSL_free(arr);
  980. return ret;
  981. }
  982. /*
  983. * Convert the bit-string representation of a polynomial ( \sum_{i=0}^n a_i *
  984. * x^i) into an array of integers corresponding to the bits with non-zero
  985. * coefficient. Array is terminated with -1. Up to max elements of the array
  986. * will be filled. Return value is total number of array elements that would
  987. * be filled if array was large enough.
  988. */
  989. int BN_GF2m_poly2arr(const BIGNUM *a, int p[], int max)
  990. {
  991. int i, j, k = 0;
  992. BN_ULONG mask;
  993. if (BN_is_zero(a))
  994. return 0;
  995. for (i = a->top - 1; i >= 0; i--) {
  996. if (!a->d[i])
  997. /* skip word if a->d[i] == 0 */
  998. continue;
  999. mask = BN_TBIT;
  1000. for (j = BN_BITS2 - 1; j >= 0; j--) {
  1001. if (a->d[i] & mask) {
  1002. if (k < max)
  1003. p[k] = BN_BITS2 * i + j;
  1004. k++;
  1005. }
  1006. mask >>= 1;
  1007. }
  1008. }
  1009. if (k < max) {
  1010. p[k] = -1;
  1011. k++;
  1012. }
  1013. return k;
  1014. }
  1015. /*
  1016. * Convert the coefficient array representation of a polynomial to a
  1017. * bit-string. The array must be terminated by -1.
  1018. */
  1019. int BN_GF2m_arr2poly(const int p[], BIGNUM *a)
  1020. {
  1021. int i;
  1022. bn_check_top(a);
  1023. BN_zero(a);
  1024. for (i = 0; p[i] != -1; i++) {
  1025. if (BN_set_bit(a, p[i]) == 0)
  1026. return 0;
  1027. }
  1028. bn_check_top(a);
  1029. return 1;
  1030. }
  1031. #endif