bn_gf2m.c 29 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181
  1. /*
  2. * Copyright 2002-2021 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. ERR_raise(ERR_LIB_BN, 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;
  418. bn_check_top(a);
  419. bn_check_top(b);
  420. bn_check_top(p);
  421. arr = OPENSSL_malloc(sizeof(*arr) * max);
  422. if (arr == NULL)
  423. return 0;
  424. ret = BN_GF2m_poly2arr(p, arr, max);
  425. if (!ret || ret > max) {
  426. ERR_raise(ERR_LIB_BN, BN_R_INVALID_LENGTH);
  427. goto err;
  428. }
  429. ret = BN_GF2m_mod_mul_arr(r, a, b, arr, ctx);
  430. bn_check_top(r);
  431. err:
  432. OPENSSL_free(arr);
  433. return ret;
  434. }
  435. /* Square a, reduce the result mod p, and store it in a. r could be a. */
  436. int BN_GF2m_mod_sqr_arr(BIGNUM *r, const BIGNUM *a, const int p[],
  437. BN_CTX *ctx)
  438. {
  439. int i, ret = 0;
  440. BIGNUM *s;
  441. bn_check_top(a);
  442. BN_CTX_start(ctx);
  443. if ((s = BN_CTX_get(ctx)) == NULL)
  444. goto err;
  445. if (!bn_wexpand(s, 2 * a->top))
  446. goto err;
  447. for (i = a->top - 1; i >= 0; i--) {
  448. s->d[2 * i + 1] = SQR1(a->d[i]);
  449. s->d[2 * i] = SQR0(a->d[i]);
  450. }
  451. s->top = 2 * a->top;
  452. bn_correct_top(s);
  453. if (!BN_GF2m_mod_arr(r, s, p))
  454. goto err;
  455. bn_check_top(r);
  456. ret = 1;
  457. err:
  458. BN_CTX_end(ctx);
  459. return ret;
  460. }
  461. /*
  462. * Square a, reduce the result mod p, and store it in a. r could be a. This
  463. * function calls down to the BN_GF2m_mod_sqr_arr implementation; this
  464. * wrapper function is only provided for convenience; for best performance,
  465. * use the BN_GF2m_mod_sqr_arr function.
  466. */
  467. int BN_GF2m_mod_sqr(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  468. {
  469. int ret = 0;
  470. const int max = BN_num_bits(p) + 1;
  471. int *arr;
  472. bn_check_top(a);
  473. bn_check_top(p);
  474. arr = OPENSSL_malloc(sizeof(*arr) * max);
  475. if (arr == NULL)
  476. return 0;
  477. ret = BN_GF2m_poly2arr(p, arr, max);
  478. if (!ret || ret > max) {
  479. ERR_raise(ERR_LIB_BN, BN_R_INVALID_LENGTH);
  480. goto err;
  481. }
  482. ret = BN_GF2m_mod_sqr_arr(r, a, arr, ctx);
  483. bn_check_top(r);
  484. err:
  485. OPENSSL_free(arr);
  486. return ret;
  487. }
  488. /*
  489. * Invert a, reduce modulo p, and store the result in r. r could be a. Uses
  490. * Modified Almost Inverse Algorithm (Algorithm 10) from Hankerson, D.,
  491. * Hernandez, J.L., and Menezes, A. "Software Implementation of Elliptic
  492. * Curve Cryptography Over Binary Fields".
  493. */
  494. static int BN_GF2m_mod_inv_vartime(BIGNUM *r, const BIGNUM *a,
  495. const BIGNUM *p, BN_CTX *ctx)
  496. {
  497. BIGNUM *b, *c = NULL, *u = NULL, *v = NULL, *tmp;
  498. int ret = 0;
  499. bn_check_top(a);
  500. bn_check_top(p);
  501. BN_CTX_start(ctx);
  502. b = BN_CTX_get(ctx);
  503. c = BN_CTX_get(ctx);
  504. u = BN_CTX_get(ctx);
  505. v = BN_CTX_get(ctx);
  506. if (v == NULL)
  507. goto err;
  508. if (!BN_GF2m_mod(u, a, p))
  509. goto err;
  510. if (BN_is_zero(u))
  511. goto err;
  512. if (!BN_copy(v, p))
  513. goto err;
  514. # if 0
  515. if (!BN_one(b))
  516. goto err;
  517. while (1) {
  518. while (!BN_is_odd(u)) {
  519. if (BN_is_zero(u))
  520. goto err;
  521. if (!BN_rshift1(u, u))
  522. goto err;
  523. if (BN_is_odd(b)) {
  524. if (!BN_GF2m_add(b, b, p))
  525. goto err;
  526. }
  527. if (!BN_rshift1(b, b))
  528. goto err;
  529. }
  530. if (BN_abs_is_word(u, 1))
  531. break;
  532. if (BN_num_bits(u) < BN_num_bits(v)) {
  533. tmp = u;
  534. u = v;
  535. v = tmp;
  536. tmp = b;
  537. b = c;
  538. c = tmp;
  539. }
  540. if (!BN_GF2m_add(u, u, v))
  541. goto err;
  542. if (!BN_GF2m_add(b, b, c))
  543. goto err;
  544. }
  545. # else
  546. {
  547. int i;
  548. int ubits = BN_num_bits(u);
  549. int vbits = BN_num_bits(v); /* v is copy of p */
  550. int top = p->top;
  551. BN_ULONG *udp, *bdp, *vdp, *cdp;
  552. if (!bn_wexpand(u, top))
  553. goto err;
  554. udp = u->d;
  555. for (i = u->top; i < top; i++)
  556. udp[i] = 0;
  557. u->top = top;
  558. if (!bn_wexpand(b, top))
  559. goto err;
  560. bdp = b->d;
  561. bdp[0] = 1;
  562. for (i = 1; i < top; i++)
  563. bdp[i] = 0;
  564. b->top = top;
  565. if (!bn_wexpand(c, top))
  566. goto err;
  567. cdp = c->d;
  568. for (i = 0; i < top; i++)
  569. cdp[i] = 0;
  570. c->top = top;
  571. vdp = v->d; /* It pays off to "cache" *->d pointers,
  572. * because it allows optimizer to be more
  573. * aggressive. But we don't have to "cache"
  574. * p->d, because *p is declared 'const'... */
  575. while (1) {
  576. while (ubits && !(udp[0] & 1)) {
  577. BN_ULONG u0, u1, b0, b1, mask;
  578. u0 = udp[0];
  579. b0 = bdp[0];
  580. mask = (BN_ULONG)0 - (b0 & 1);
  581. b0 ^= p->d[0] & mask;
  582. for (i = 0; i < top - 1; i++) {
  583. u1 = udp[i + 1];
  584. udp[i] = ((u0 >> 1) | (u1 << (BN_BITS2 - 1))) & BN_MASK2;
  585. u0 = u1;
  586. b1 = bdp[i + 1] ^ (p->d[i + 1] & mask);
  587. bdp[i] = ((b0 >> 1) | (b1 << (BN_BITS2 - 1))) & BN_MASK2;
  588. b0 = b1;
  589. }
  590. udp[i] = u0 >> 1;
  591. bdp[i] = b0 >> 1;
  592. ubits--;
  593. }
  594. if (ubits <= BN_BITS2) {
  595. if (udp[0] == 0) /* poly was reducible */
  596. goto err;
  597. if (udp[0] == 1)
  598. break;
  599. }
  600. if (ubits < vbits) {
  601. i = ubits;
  602. ubits = vbits;
  603. vbits = i;
  604. tmp = u;
  605. u = v;
  606. v = tmp;
  607. tmp = b;
  608. b = c;
  609. c = tmp;
  610. udp = vdp;
  611. vdp = v->d;
  612. bdp = cdp;
  613. cdp = c->d;
  614. }
  615. for (i = 0; i < top; i++) {
  616. udp[i] ^= vdp[i];
  617. bdp[i] ^= cdp[i];
  618. }
  619. if (ubits == vbits) {
  620. BN_ULONG ul;
  621. int utop = (ubits - 1) / BN_BITS2;
  622. while ((ul = udp[utop]) == 0 && utop)
  623. utop--;
  624. ubits = utop * BN_BITS2 + BN_num_bits_word(ul);
  625. }
  626. }
  627. bn_correct_top(b);
  628. }
  629. # endif
  630. if (!BN_copy(r, b))
  631. goto err;
  632. bn_check_top(r);
  633. ret = 1;
  634. err:
  635. # ifdef BN_DEBUG
  636. /* BN_CTX_end would complain about the expanded form */
  637. bn_correct_top(c);
  638. bn_correct_top(u);
  639. bn_correct_top(v);
  640. # endif
  641. BN_CTX_end(ctx);
  642. return ret;
  643. }
  644. /*-
  645. * Wrapper for BN_GF2m_mod_inv_vartime that blinds the input before calling.
  646. * This is not constant time.
  647. * But it does eliminate first order deduction on the input.
  648. */
  649. int BN_GF2m_mod_inv(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  650. {
  651. BIGNUM *b = NULL;
  652. int ret = 0;
  653. BN_CTX_start(ctx);
  654. if ((b = BN_CTX_get(ctx)) == NULL)
  655. goto err;
  656. /* generate blinding value */
  657. do {
  658. if (!BN_priv_rand_ex(b, BN_num_bits(p) - 1,
  659. BN_RAND_TOP_ANY, BN_RAND_BOTTOM_ANY, 0, ctx))
  660. goto err;
  661. } while (BN_is_zero(b));
  662. /* r := a * b */
  663. if (!BN_GF2m_mod_mul(r, a, b, p, ctx))
  664. goto err;
  665. /* r := 1/(a * b) */
  666. if (!BN_GF2m_mod_inv_vartime(r, r, p, ctx))
  667. goto err;
  668. /* r := b/(a * b) = 1/a */
  669. if (!BN_GF2m_mod_mul(r, r, b, p, ctx))
  670. goto err;
  671. ret = 1;
  672. err:
  673. BN_CTX_end(ctx);
  674. return ret;
  675. }
  676. /*
  677. * Invert xx, reduce modulo p, and store the result in r. r could be xx.
  678. * This function calls down to the BN_GF2m_mod_inv implementation; this
  679. * wrapper function is only provided for convenience; for best performance,
  680. * use the BN_GF2m_mod_inv function.
  681. */
  682. int BN_GF2m_mod_inv_arr(BIGNUM *r, const BIGNUM *xx, const int p[],
  683. BN_CTX *ctx)
  684. {
  685. BIGNUM *field;
  686. int ret = 0;
  687. bn_check_top(xx);
  688. BN_CTX_start(ctx);
  689. if ((field = BN_CTX_get(ctx)) == NULL)
  690. goto err;
  691. if (!BN_GF2m_arr2poly(p, field))
  692. goto err;
  693. ret = BN_GF2m_mod_inv(r, xx, field, ctx);
  694. bn_check_top(r);
  695. err:
  696. BN_CTX_end(ctx);
  697. return ret;
  698. }
  699. /*
  700. * Divide y by x, reduce modulo p, and store the result in r. r could be x
  701. * or y, x could equal y.
  702. */
  703. int BN_GF2m_mod_div(BIGNUM *r, const BIGNUM *y, const BIGNUM *x,
  704. const BIGNUM *p, BN_CTX *ctx)
  705. {
  706. BIGNUM *xinv = NULL;
  707. int ret = 0;
  708. bn_check_top(y);
  709. bn_check_top(x);
  710. bn_check_top(p);
  711. BN_CTX_start(ctx);
  712. xinv = BN_CTX_get(ctx);
  713. if (xinv == NULL)
  714. goto err;
  715. if (!BN_GF2m_mod_inv(xinv, x, p, ctx))
  716. goto err;
  717. if (!BN_GF2m_mod_mul(r, y, xinv, p, ctx))
  718. goto err;
  719. bn_check_top(r);
  720. ret = 1;
  721. err:
  722. BN_CTX_end(ctx);
  723. return ret;
  724. }
  725. /*
  726. * Divide yy by xx, reduce modulo p, and store the result in r. r could be xx
  727. * * or yy, xx could equal yy. This function calls down to the
  728. * BN_GF2m_mod_div implementation; this wrapper function is only provided for
  729. * convenience; for best performance, use the BN_GF2m_mod_div function.
  730. */
  731. int BN_GF2m_mod_div_arr(BIGNUM *r, const BIGNUM *yy, const BIGNUM *xx,
  732. const int p[], BN_CTX *ctx)
  733. {
  734. BIGNUM *field;
  735. int ret = 0;
  736. bn_check_top(yy);
  737. bn_check_top(xx);
  738. BN_CTX_start(ctx);
  739. if ((field = BN_CTX_get(ctx)) == NULL)
  740. goto err;
  741. if (!BN_GF2m_arr2poly(p, field))
  742. goto err;
  743. ret = BN_GF2m_mod_div(r, yy, xx, field, ctx);
  744. bn_check_top(r);
  745. err:
  746. BN_CTX_end(ctx);
  747. return ret;
  748. }
  749. /*
  750. * Compute the bth power of a, reduce modulo p, and store the result in r. r
  751. * could be a. Uses simple square-and-multiply algorithm A.5.1 from IEEE
  752. * P1363.
  753. */
  754. int BN_GF2m_mod_exp_arr(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  755. const int p[], BN_CTX *ctx)
  756. {
  757. int ret = 0, i, n;
  758. BIGNUM *u;
  759. bn_check_top(a);
  760. bn_check_top(b);
  761. if (BN_is_zero(b))
  762. return BN_one(r);
  763. if (BN_abs_is_word(b, 1))
  764. return (BN_copy(r, a) != NULL);
  765. BN_CTX_start(ctx);
  766. if ((u = BN_CTX_get(ctx)) == NULL)
  767. goto err;
  768. if (!BN_GF2m_mod_arr(u, a, p))
  769. goto err;
  770. n = BN_num_bits(b) - 1;
  771. for (i = n - 1; i >= 0; i--) {
  772. if (!BN_GF2m_mod_sqr_arr(u, u, p, ctx))
  773. goto err;
  774. if (BN_is_bit_set(b, i)) {
  775. if (!BN_GF2m_mod_mul_arr(u, u, a, p, ctx))
  776. goto err;
  777. }
  778. }
  779. if (!BN_copy(r, u))
  780. goto err;
  781. bn_check_top(r);
  782. ret = 1;
  783. err:
  784. BN_CTX_end(ctx);
  785. return ret;
  786. }
  787. /*
  788. * Compute the bth power of a, reduce modulo p, and store the result in r. r
  789. * could be a. This function calls down to the BN_GF2m_mod_exp_arr
  790. * implementation; this wrapper function is only provided for convenience;
  791. * for best performance, use the BN_GF2m_mod_exp_arr function.
  792. */
  793. int BN_GF2m_mod_exp(BIGNUM *r, const BIGNUM *a, const BIGNUM *b,
  794. const BIGNUM *p, BN_CTX *ctx)
  795. {
  796. int ret = 0;
  797. const int max = BN_num_bits(p) + 1;
  798. int *arr;
  799. bn_check_top(a);
  800. bn_check_top(b);
  801. bn_check_top(p);
  802. arr = OPENSSL_malloc(sizeof(*arr) * max);
  803. if (arr == NULL)
  804. return 0;
  805. ret = BN_GF2m_poly2arr(p, arr, max);
  806. if (!ret || ret > max) {
  807. ERR_raise(ERR_LIB_BN, BN_R_INVALID_LENGTH);
  808. goto err;
  809. }
  810. ret = BN_GF2m_mod_exp_arr(r, a, b, arr, ctx);
  811. bn_check_top(r);
  812. err:
  813. OPENSSL_free(arr);
  814. return ret;
  815. }
  816. /*
  817. * Compute the square root of a, reduce modulo p, and store the result in r.
  818. * r could be a. Uses exponentiation as in algorithm A.4.1 from IEEE P1363.
  819. */
  820. int BN_GF2m_mod_sqrt_arr(BIGNUM *r, const BIGNUM *a, const int p[],
  821. BN_CTX *ctx)
  822. {
  823. int ret = 0;
  824. BIGNUM *u;
  825. bn_check_top(a);
  826. if (p[0] == 0) {
  827. /* reduction mod 1 => return 0 */
  828. BN_zero(r);
  829. return 1;
  830. }
  831. BN_CTX_start(ctx);
  832. if ((u = BN_CTX_get(ctx)) == NULL)
  833. goto err;
  834. if (!BN_set_bit(u, p[0] - 1))
  835. goto err;
  836. ret = BN_GF2m_mod_exp_arr(r, a, u, p, ctx);
  837. bn_check_top(r);
  838. err:
  839. BN_CTX_end(ctx);
  840. return ret;
  841. }
  842. /*
  843. * Compute the square root of a, reduce modulo p, and store the result in r.
  844. * r could be a. This function calls down to the BN_GF2m_mod_sqrt_arr
  845. * implementation; this wrapper function is only provided for convenience;
  846. * for best performance, use the BN_GF2m_mod_sqrt_arr function.
  847. */
  848. int BN_GF2m_mod_sqrt(BIGNUM *r, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx)
  849. {
  850. int ret = 0;
  851. const int max = BN_num_bits(p) + 1;
  852. int *arr;
  853. bn_check_top(a);
  854. bn_check_top(p);
  855. arr = OPENSSL_malloc(sizeof(*arr) * max);
  856. if (arr == NULL)
  857. return 0;
  858. ret = BN_GF2m_poly2arr(p, arr, max);
  859. if (!ret || ret > max) {
  860. ERR_raise(ERR_LIB_BN, BN_R_INVALID_LENGTH);
  861. goto err;
  862. }
  863. ret = BN_GF2m_mod_sqrt_arr(r, a, arr, ctx);
  864. bn_check_top(r);
  865. err:
  866. OPENSSL_free(arr);
  867. return ret;
  868. }
  869. /*
  870. * Find r such that r^2 + r = a mod p. r could be a. If no r exists returns
  871. * 0. Uses algorithms A.4.7 and A.4.6 from IEEE P1363.
  872. */
  873. int BN_GF2m_mod_solve_quad_arr(BIGNUM *r, const BIGNUM *a_, const int p[],
  874. BN_CTX *ctx)
  875. {
  876. int ret = 0, count = 0, j;
  877. BIGNUM *a, *z, *rho, *w, *w2, *tmp;
  878. bn_check_top(a_);
  879. if (p[0] == 0) {
  880. /* reduction mod 1 => return 0 */
  881. BN_zero(r);
  882. return 1;
  883. }
  884. BN_CTX_start(ctx);
  885. a = BN_CTX_get(ctx);
  886. z = BN_CTX_get(ctx);
  887. w = BN_CTX_get(ctx);
  888. if (w == NULL)
  889. goto err;
  890. if (!BN_GF2m_mod_arr(a, a_, p))
  891. goto err;
  892. if (BN_is_zero(a)) {
  893. BN_zero(r);
  894. ret = 1;
  895. goto err;
  896. }
  897. if (p[0] & 0x1) { /* m is odd */
  898. /* compute half-trace of a */
  899. if (!BN_copy(z, a))
  900. goto err;
  901. for (j = 1; j <= (p[0] - 1) / 2; j++) {
  902. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  903. goto err;
  904. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  905. goto err;
  906. if (!BN_GF2m_add(z, z, a))
  907. goto err;
  908. }
  909. } else { /* m is even */
  910. rho = BN_CTX_get(ctx);
  911. w2 = BN_CTX_get(ctx);
  912. tmp = BN_CTX_get(ctx);
  913. if (tmp == NULL)
  914. goto err;
  915. do {
  916. if (!BN_priv_rand_ex(rho, p[0], BN_RAND_TOP_ONE, BN_RAND_BOTTOM_ANY,
  917. 0, ctx))
  918. goto err;
  919. if (!BN_GF2m_mod_arr(rho, rho, p))
  920. goto err;
  921. BN_zero(z);
  922. if (!BN_copy(w, rho))
  923. goto err;
  924. for (j = 1; j <= p[0] - 1; j++) {
  925. if (!BN_GF2m_mod_sqr_arr(z, z, p, ctx))
  926. goto err;
  927. if (!BN_GF2m_mod_sqr_arr(w2, w, p, ctx))
  928. goto err;
  929. if (!BN_GF2m_mod_mul_arr(tmp, w2, a, p, ctx))
  930. goto err;
  931. if (!BN_GF2m_add(z, z, tmp))
  932. goto err;
  933. if (!BN_GF2m_add(w, w2, rho))
  934. goto err;
  935. }
  936. count++;
  937. } while (BN_is_zero(w) && (count < MAX_ITERATIONS));
  938. if (BN_is_zero(w)) {
  939. ERR_raise(ERR_LIB_BN, BN_R_TOO_MANY_ITERATIONS);
  940. goto err;
  941. }
  942. }
  943. if (!BN_GF2m_mod_sqr_arr(w, z, p, ctx))
  944. goto err;
  945. if (!BN_GF2m_add(w, z, w))
  946. goto err;
  947. if (BN_GF2m_cmp(w, a)) {
  948. ERR_raise(ERR_LIB_BN, BN_R_NO_SOLUTION);
  949. goto err;
  950. }
  951. if (!BN_copy(r, z))
  952. goto err;
  953. bn_check_top(r);
  954. ret = 1;
  955. err:
  956. BN_CTX_end(ctx);
  957. return ret;
  958. }
  959. /*
  960. * Find r such that r^2 + r = a mod p. r could be a. If no r exists returns
  961. * 0. This function calls down to the BN_GF2m_mod_solve_quad_arr
  962. * implementation; this wrapper function is only provided for convenience;
  963. * for best performance, use the BN_GF2m_mod_solve_quad_arr function.
  964. */
  965. int BN_GF2m_mod_solve_quad(BIGNUM *r, const BIGNUM *a, const BIGNUM *p,
  966. BN_CTX *ctx)
  967. {
  968. int ret = 0;
  969. const int max = BN_num_bits(p) + 1;
  970. int *arr;
  971. bn_check_top(a);
  972. bn_check_top(p);
  973. arr = OPENSSL_malloc(sizeof(*arr) * max);
  974. if (arr == NULL)
  975. goto err;
  976. ret = BN_GF2m_poly2arr(p, arr, max);
  977. if (!ret || ret > max) {
  978. ERR_raise(ERR_LIB_BN, BN_R_INVALID_LENGTH);
  979. goto err;
  980. }
  981. ret = BN_GF2m_mod_solve_quad_arr(r, a, arr, ctx);
  982. bn_check_top(r);
  983. err:
  984. OPENSSL_free(arr);
  985. return ret;
  986. }
  987. /*
  988. * Convert the bit-string representation of a polynomial ( \sum_{i=0}^n a_i *
  989. * x^i) into an array of integers corresponding to the bits with non-zero
  990. * coefficient. Array is terminated with -1. Up to max elements of the array
  991. * will be filled. Return value is total number of array elements that would
  992. * be filled if array was large enough.
  993. */
  994. int BN_GF2m_poly2arr(const BIGNUM *a, int p[], int max)
  995. {
  996. int i, j, k = 0;
  997. BN_ULONG mask;
  998. if (BN_is_zero(a))
  999. return 0;
  1000. for (i = a->top - 1; i >= 0; i--) {
  1001. if (!a->d[i])
  1002. /* skip word if a->d[i] == 0 */
  1003. continue;
  1004. mask = BN_TBIT;
  1005. for (j = BN_BITS2 - 1; j >= 0; j--) {
  1006. if (a->d[i] & mask) {
  1007. if (k < max)
  1008. p[k] = BN_BITS2 * i + j;
  1009. k++;
  1010. }
  1011. mask >>= 1;
  1012. }
  1013. }
  1014. if (k < max) {
  1015. p[k] = -1;
  1016. k++;
  1017. }
  1018. return k;
  1019. }
  1020. /*
  1021. * Convert the coefficient array representation of a polynomial to a
  1022. * bit-string. The array must be terminated by -1.
  1023. */
  1024. int BN_GF2m_arr2poly(const int p[], BIGNUM *a)
  1025. {
  1026. int i;
  1027. bn_check_top(a);
  1028. BN_zero(a);
  1029. for (i = 0; p[i] != -1; i++) {
  1030. if (BN_set_bit(a, p[i]) == 0)
  1031. return 0;
  1032. }
  1033. bn_check_top(a);
  1034. return 1;
  1035. }
  1036. #endif