2
0

bn_gcd.c 17 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616
  1. /*
  2. * Copyright 1995-2017 The OpenSSL Project Authors. All Rights Reserved.
  3. *
  4. * Licensed under the OpenSSL license (the "License"). You may not use
  5. * this file except in compliance with the License. You can obtain a copy
  6. * in the file LICENSE in the source distribution or at
  7. * https://www.openssl.org/source/license.html
  8. */
  9. #include "internal/cryptlib.h"
  10. #include "bn_lcl.h"
  11. static BIGNUM *euclid(BIGNUM *a, BIGNUM *b);
  12. int BN_gcd(BIGNUM *r, const BIGNUM *in_a, const BIGNUM *in_b, BN_CTX *ctx)
  13. {
  14. BIGNUM *a, *b, *t;
  15. int ret = 0;
  16. bn_check_top(in_a);
  17. bn_check_top(in_b);
  18. BN_CTX_start(ctx);
  19. a = BN_CTX_get(ctx);
  20. b = BN_CTX_get(ctx);
  21. if (b == NULL)
  22. goto err;
  23. if (BN_copy(a, in_a) == NULL)
  24. goto err;
  25. if (BN_copy(b, in_b) == NULL)
  26. goto err;
  27. a->neg = 0;
  28. b->neg = 0;
  29. if (BN_cmp(a, b) < 0) {
  30. t = a;
  31. a = b;
  32. b = t;
  33. }
  34. t = euclid(a, b);
  35. if (t == NULL)
  36. goto err;
  37. if (BN_copy(r, t) == NULL)
  38. goto err;
  39. ret = 1;
  40. err:
  41. BN_CTX_end(ctx);
  42. bn_check_top(r);
  43. return ret;
  44. }
  45. static BIGNUM *euclid(BIGNUM *a, BIGNUM *b)
  46. {
  47. BIGNUM *t;
  48. int shifts = 0;
  49. bn_check_top(a);
  50. bn_check_top(b);
  51. /* 0 <= b <= a */
  52. while (!BN_is_zero(b)) {
  53. /* 0 < b <= a */
  54. if (BN_is_odd(a)) {
  55. if (BN_is_odd(b)) {
  56. if (!BN_sub(a, a, b))
  57. goto err;
  58. if (!BN_rshift1(a, a))
  59. goto err;
  60. if (BN_cmp(a, b) < 0) {
  61. t = a;
  62. a = b;
  63. b = t;
  64. }
  65. } else { /* a odd - b even */
  66. if (!BN_rshift1(b, b))
  67. goto err;
  68. if (BN_cmp(a, b) < 0) {
  69. t = a;
  70. a = b;
  71. b = t;
  72. }
  73. }
  74. } else { /* a is even */
  75. if (BN_is_odd(b)) {
  76. if (!BN_rshift1(a, a))
  77. goto err;
  78. if (BN_cmp(a, b) < 0) {
  79. t = a;
  80. a = b;
  81. b = t;
  82. }
  83. } else { /* a even - b even */
  84. if (!BN_rshift1(a, a))
  85. goto err;
  86. if (!BN_rshift1(b, b))
  87. goto err;
  88. shifts++;
  89. }
  90. }
  91. /* 0 <= b <= a */
  92. }
  93. if (shifts) {
  94. if (!BN_lshift(a, a, shifts))
  95. goto err;
  96. }
  97. bn_check_top(a);
  98. return a;
  99. err:
  100. return NULL;
  101. }
  102. /* solves ax == 1 (mod n) */
  103. static BIGNUM *BN_mod_inverse_no_branch(BIGNUM *in,
  104. const BIGNUM *a, const BIGNUM *n,
  105. BN_CTX *ctx);
  106. BIGNUM *BN_mod_inverse(BIGNUM *in,
  107. const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx)
  108. {
  109. BIGNUM *rv;
  110. int noinv;
  111. rv = int_bn_mod_inverse(in, a, n, ctx, &noinv);
  112. if (noinv)
  113. BNerr(BN_F_BN_MOD_INVERSE, BN_R_NO_INVERSE);
  114. return rv;
  115. }
  116. BIGNUM *int_bn_mod_inverse(BIGNUM *in,
  117. const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx,
  118. int *pnoinv)
  119. {
  120. BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
  121. BIGNUM *ret = NULL;
  122. int sign;
  123. if (pnoinv)
  124. *pnoinv = 0;
  125. if ((BN_get_flags(a, BN_FLG_CONSTTIME) != 0)
  126. || (BN_get_flags(n, BN_FLG_CONSTTIME) != 0)) {
  127. return BN_mod_inverse_no_branch(in, a, n, ctx);
  128. }
  129. bn_check_top(a);
  130. bn_check_top(n);
  131. BN_CTX_start(ctx);
  132. A = BN_CTX_get(ctx);
  133. B = BN_CTX_get(ctx);
  134. X = BN_CTX_get(ctx);
  135. D = BN_CTX_get(ctx);
  136. M = BN_CTX_get(ctx);
  137. Y = BN_CTX_get(ctx);
  138. T = BN_CTX_get(ctx);
  139. if (T == NULL)
  140. goto err;
  141. if (in == NULL)
  142. R = BN_new();
  143. else
  144. R = in;
  145. if (R == NULL)
  146. goto err;
  147. BN_one(X);
  148. BN_zero(Y);
  149. if (BN_copy(B, a) == NULL)
  150. goto err;
  151. if (BN_copy(A, n) == NULL)
  152. goto err;
  153. A->neg = 0;
  154. if (B->neg || (BN_ucmp(B, A) >= 0)) {
  155. if (!BN_nnmod(B, B, A, ctx))
  156. goto err;
  157. }
  158. sign = -1;
  159. /*-
  160. * From B = a mod |n|, A = |n| it follows that
  161. *
  162. * 0 <= B < A,
  163. * -sign*X*a == B (mod |n|),
  164. * sign*Y*a == A (mod |n|).
  165. */
  166. if (BN_is_odd(n) && (BN_num_bits(n) <= 2048)) {
  167. /*
  168. * Binary inversion algorithm; requires odd modulus. This is faster
  169. * than the general algorithm if the modulus is sufficiently small
  170. * (about 400 .. 500 bits on 32-bit systems, but much more on 64-bit
  171. * systems)
  172. */
  173. int shift;
  174. while (!BN_is_zero(B)) {
  175. /*-
  176. * 0 < B < |n|,
  177. * 0 < A <= |n|,
  178. * (1) -sign*X*a == B (mod |n|),
  179. * (2) sign*Y*a == A (mod |n|)
  180. */
  181. /*
  182. * Now divide B by the maximum possible power of two in the
  183. * integers, and divide X by the same value mod |n|. When we're
  184. * done, (1) still holds.
  185. */
  186. shift = 0;
  187. while (!BN_is_bit_set(B, shift)) { /* note that 0 < B */
  188. shift++;
  189. if (BN_is_odd(X)) {
  190. if (!BN_uadd(X, X, n))
  191. goto err;
  192. }
  193. /*
  194. * now X is even, so we can easily divide it by two
  195. */
  196. if (!BN_rshift1(X, X))
  197. goto err;
  198. }
  199. if (shift > 0) {
  200. if (!BN_rshift(B, B, shift))
  201. goto err;
  202. }
  203. /*
  204. * Same for A and Y. Afterwards, (2) still holds.
  205. */
  206. shift = 0;
  207. while (!BN_is_bit_set(A, shift)) { /* note that 0 < A */
  208. shift++;
  209. if (BN_is_odd(Y)) {
  210. if (!BN_uadd(Y, Y, n))
  211. goto err;
  212. }
  213. /* now Y is even */
  214. if (!BN_rshift1(Y, Y))
  215. goto err;
  216. }
  217. if (shift > 0) {
  218. if (!BN_rshift(A, A, shift))
  219. goto err;
  220. }
  221. /*-
  222. * We still have (1) and (2).
  223. * Both A and B are odd.
  224. * The following computations ensure that
  225. *
  226. * 0 <= B < |n|,
  227. * 0 < A < |n|,
  228. * (1) -sign*X*a == B (mod |n|),
  229. * (2) sign*Y*a == A (mod |n|),
  230. *
  231. * and that either A or B is even in the next iteration.
  232. */
  233. if (BN_ucmp(B, A) >= 0) {
  234. /* -sign*(X + Y)*a == B - A (mod |n|) */
  235. if (!BN_uadd(X, X, Y))
  236. goto err;
  237. /*
  238. * NB: we could use BN_mod_add_quick(X, X, Y, n), but that
  239. * actually makes the algorithm slower
  240. */
  241. if (!BN_usub(B, B, A))
  242. goto err;
  243. } else {
  244. /* sign*(X + Y)*a == A - B (mod |n|) */
  245. if (!BN_uadd(Y, Y, X))
  246. goto err;
  247. /*
  248. * as above, BN_mod_add_quick(Y, Y, X, n) would slow things down
  249. */
  250. if (!BN_usub(A, A, B))
  251. goto err;
  252. }
  253. }
  254. } else {
  255. /* general inversion algorithm */
  256. while (!BN_is_zero(B)) {
  257. BIGNUM *tmp;
  258. /*-
  259. * 0 < B < A,
  260. * (*) -sign*X*a == B (mod |n|),
  261. * sign*Y*a == A (mod |n|)
  262. */
  263. /* (D, M) := (A/B, A%B) ... */
  264. if (BN_num_bits(A) == BN_num_bits(B)) {
  265. if (!BN_one(D))
  266. goto err;
  267. if (!BN_sub(M, A, B))
  268. goto err;
  269. } else if (BN_num_bits(A) == BN_num_bits(B) + 1) {
  270. /* A/B is 1, 2, or 3 */
  271. if (!BN_lshift1(T, B))
  272. goto err;
  273. if (BN_ucmp(A, T) < 0) {
  274. /* A < 2*B, so D=1 */
  275. if (!BN_one(D))
  276. goto err;
  277. if (!BN_sub(M, A, B))
  278. goto err;
  279. } else {
  280. /* A >= 2*B, so D=2 or D=3 */
  281. if (!BN_sub(M, A, T))
  282. goto err;
  283. if (!BN_add(D, T, B))
  284. goto err; /* use D (:= 3*B) as temp */
  285. if (BN_ucmp(A, D) < 0) {
  286. /* A < 3*B, so D=2 */
  287. if (!BN_set_word(D, 2))
  288. goto err;
  289. /*
  290. * M (= A - 2*B) already has the correct value
  291. */
  292. } else {
  293. /* only D=3 remains */
  294. if (!BN_set_word(D, 3))
  295. goto err;
  296. /*
  297. * currently M = A - 2*B, but we need M = A - 3*B
  298. */
  299. if (!BN_sub(M, M, B))
  300. goto err;
  301. }
  302. }
  303. } else {
  304. if (!BN_div(D, M, A, B, ctx))
  305. goto err;
  306. }
  307. /*-
  308. * Now
  309. * A = D*B + M;
  310. * thus we have
  311. * (**) sign*Y*a == D*B + M (mod |n|).
  312. */
  313. tmp = A; /* keep the BIGNUM object, the value does not matter */
  314. /* (A, B) := (B, A mod B) ... */
  315. A = B;
  316. B = M;
  317. /* ... so we have 0 <= B < A again */
  318. /*-
  319. * Since the former M is now B and the former B is now A,
  320. * (**) translates into
  321. * sign*Y*a == D*A + B (mod |n|),
  322. * i.e.
  323. * sign*Y*a - D*A == B (mod |n|).
  324. * Similarly, (*) translates into
  325. * -sign*X*a == A (mod |n|).
  326. *
  327. * Thus,
  328. * sign*Y*a + D*sign*X*a == B (mod |n|),
  329. * i.e.
  330. * sign*(Y + D*X)*a == B (mod |n|).
  331. *
  332. * So if we set (X, Y, sign) := (Y + D*X, X, -sign), we arrive back at
  333. * -sign*X*a == B (mod |n|),
  334. * sign*Y*a == A (mod |n|).
  335. * Note that X and Y stay non-negative all the time.
  336. */
  337. /*
  338. * most of the time D is very small, so we can optimize tmp := D*X+Y
  339. */
  340. if (BN_is_one(D)) {
  341. if (!BN_add(tmp, X, Y))
  342. goto err;
  343. } else {
  344. if (BN_is_word(D, 2)) {
  345. if (!BN_lshift1(tmp, X))
  346. goto err;
  347. } else if (BN_is_word(D, 4)) {
  348. if (!BN_lshift(tmp, X, 2))
  349. goto err;
  350. } else if (D->top == 1) {
  351. if (!BN_copy(tmp, X))
  352. goto err;
  353. if (!BN_mul_word(tmp, D->d[0]))
  354. goto err;
  355. } else {
  356. if (!BN_mul(tmp, D, X, ctx))
  357. goto err;
  358. }
  359. if (!BN_add(tmp, tmp, Y))
  360. goto err;
  361. }
  362. M = Y; /* keep the BIGNUM object, the value does not matter */
  363. Y = X;
  364. X = tmp;
  365. sign = -sign;
  366. }
  367. }
  368. /*-
  369. * The while loop (Euclid's algorithm) ends when
  370. * A == gcd(a,n);
  371. * we have
  372. * sign*Y*a == A (mod |n|),
  373. * where Y is non-negative.
  374. */
  375. if (sign < 0) {
  376. if (!BN_sub(Y, n, Y))
  377. goto err;
  378. }
  379. /* Now Y*a == A (mod |n|). */
  380. if (BN_is_one(A)) {
  381. /* Y*a == 1 (mod |n|) */
  382. if (!Y->neg && BN_ucmp(Y, n) < 0) {
  383. if (!BN_copy(R, Y))
  384. goto err;
  385. } else {
  386. if (!BN_nnmod(R, Y, n, ctx))
  387. goto err;
  388. }
  389. } else {
  390. if (pnoinv)
  391. *pnoinv = 1;
  392. goto err;
  393. }
  394. ret = R;
  395. err:
  396. if ((ret == NULL) && (in == NULL))
  397. BN_free(R);
  398. BN_CTX_end(ctx);
  399. bn_check_top(ret);
  400. return ret;
  401. }
  402. /*
  403. * BN_mod_inverse_no_branch is a special version of BN_mod_inverse. It does
  404. * not contain branches that may leak sensitive information.
  405. */
  406. static BIGNUM *BN_mod_inverse_no_branch(BIGNUM *in,
  407. const BIGNUM *a, const BIGNUM *n,
  408. BN_CTX *ctx)
  409. {
  410. BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
  411. BIGNUM *ret = NULL;
  412. int sign;
  413. bn_check_top(a);
  414. bn_check_top(n);
  415. BN_CTX_start(ctx);
  416. A = BN_CTX_get(ctx);
  417. B = BN_CTX_get(ctx);
  418. X = BN_CTX_get(ctx);
  419. D = BN_CTX_get(ctx);
  420. M = BN_CTX_get(ctx);
  421. Y = BN_CTX_get(ctx);
  422. T = BN_CTX_get(ctx);
  423. if (T == NULL)
  424. goto err;
  425. if (in == NULL)
  426. R = BN_new();
  427. else
  428. R = in;
  429. if (R == NULL)
  430. goto err;
  431. BN_one(X);
  432. BN_zero(Y);
  433. if (BN_copy(B, a) == NULL)
  434. goto err;
  435. if (BN_copy(A, n) == NULL)
  436. goto err;
  437. A->neg = 0;
  438. if (B->neg || (BN_ucmp(B, A) >= 0)) {
  439. /*
  440. * Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
  441. * BN_div_no_branch will be called eventually.
  442. */
  443. {
  444. BIGNUM local_B;
  445. bn_init(&local_B);
  446. BN_with_flags(&local_B, B, BN_FLG_CONSTTIME);
  447. if (!BN_nnmod(B, &local_B, A, ctx))
  448. goto err;
  449. /* Ensure local_B goes out of scope before any further use of B */
  450. }
  451. }
  452. sign = -1;
  453. /*-
  454. * From B = a mod |n|, A = |n| it follows that
  455. *
  456. * 0 <= B < A,
  457. * -sign*X*a == B (mod |n|),
  458. * sign*Y*a == A (mod |n|).
  459. */
  460. while (!BN_is_zero(B)) {
  461. BIGNUM *tmp;
  462. /*-
  463. * 0 < B < A,
  464. * (*) -sign*X*a == B (mod |n|),
  465. * sign*Y*a == A (mod |n|)
  466. */
  467. /*
  468. * Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
  469. * BN_div_no_branch will be called eventually.
  470. */
  471. {
  472. BIGNUM local_A;
  473. bn_init(&local_A);
  474. BN_with_flags(&local_A, A, BN_FLG_CONSTTIME);
  475. /* (D, M) := (A/B, A%B) ... */
  476. if (!BN_div(D, M, &local_A, B, ctx))
  477. goto err;
  478. /* Ensure local_A goes out of scope before any further use of A */
  479. }
  480. /*-
  481. * Now
  482. * A = D*B + M;
  483. * thus we have
  484. * (**) sign*Y*a == D*B + M (mod |n|).
  485. */
  486. tmp = A; /* keep the BIGNUM object, the value does not
  487. * matter */
  488. /* (A, B) := (B, A mod B) ... */
  489. A = B;
  490. B = M;
  491. /* ... so we have 0 <= B < A again */
  492. /*-
  493. * Since the former M is now B and the former B is now A,
  494. * (**) translates into
  495. * sign*Y*a == D*A + B (mod |n|),
  496. * i.e.
  497. * sign*Y*a - D*A == B (mod |n|).
  498. * Similarly, (*) translates into
  499. * -sign*X*a == A (mod |n|).
  500. *
  501. * Thus,
  502. * sign*Y*a + D*sign*X*a == B (mod |n|),
  503. * i.e.
  504. * sign*(Y + D*X)*a == B (mod |n|).
  505. *
  506. * So if we set (X, Y, sign) := (Y + D*X, X, -sign), we arrive back at
  507. * -sign*X*a == B (mod |n|),
  508. * sign*Y*a == A (mod |n|).
  509. * Note that X and Y stay non-negative all the time.
  510. */
  511. if (!BN_mul(tmp, D, X, ctx))
  512. goto err;
  513. if (!BN_add(tmp, tmp, Y))
  514. goto err;
  515. M = Y; /* keep the BIGNUM object, the value does not
  516. * matter */
  517. Y = X;
  518. X = tmp;
  519. sign = -sign;
  520. }
  521. /*-
  522. * The while loop (Euclid's algorithm) ends when
  523. * A == gcd(a,n);
  524. * we have
  525. * sign*Y*a == A (mod |n|),
  526. * where Y is non-negative.
  527. */
  528. if (sign < 0) {
  529. if (!BN_sub(Y, n, Y))
  530. goto err;
  531. }
  532. /* Now Y*a == A (mod |n|). */
  533. if (BN_is_one(A)) {
  534. /* Y*a == 1 (mod |n|) */
  535. if (!Y->neg && BN_ucmp(Y, n) < 0) {
  536. if (!BN_copy(R, Y))
  537. goto err;
  538. } else {
  539. if (!BN_nnmod(R, Y, n, ctx))
  540. goto err;
  541. }
  542. } else {
  543. BNerr(BN_F_BN_MOD_INVERSE_NO_BRANCH, BN_R_NO_INVERSE);
  544. goto err;
  545. }
  546. ret = R;
  547. err:
  548. if ((ret == NULL) && (in == NULL))
  549. BN_free(R);
  550. BN_CTX_end(ctx);
  551. bn_check_top(ret);
  552. return ret;
  553. }