Mercurial > dropbear
comparison pre_gen/mpi.c @ 190:d8254fc979e9 libtommath-orig LTM_0.35
Initial import of libtommath 0.35
author | Matt Johnston <matt@ucc.asn.au> |
---|---|
date | Fri, 06 May 2005 08:59:30 +0000 |
parents | d29b64170cf0 |
children |
comparison
equal
deleted
inserted
replaced
142:d29b64170cf0 | 190:d8254fc979e9 |
---|---|
67 * that is c = 1/a mod b | 67 * that is c = 1/a mod b |
68 * | 68 * |
69 * Based on slow invmod except this is optimized for the case where b is | 69 * Based on slow invmod except this is optimized for the case where b is |
70 * odd as per HAC Note 14.64 on pp. 610 | 70 * odd as per HAC Note 14.64 on pp. 610 |
71 */ | 71 */ |
72 int | 72 int fast_mp_invmod (mp_int * a, mp_int * b, mp_int * c) |
73 fast_mp_invmod (mp_int * a, mp_int * b, mp_int * c) | |
74 { | 73 { |
75 mp_int x, y, u, v, B, D; | 74 mp_int x, y, u, v, B, D; |
76 int res, neg; | 75 int res, neg; |
77 | 76 |
78 /* 2. [modified] b must be odd */ | 77 /* 2. [modified] b must be odd */ |
85 return res; | 84 return res; |
86 } | 85 } |
87 | 86 |
88 /* x == modulus, y == value to invert */ | 87 /* x == modulus, y == value to invert */ |
89 if ((res = mp_copy (b, &x)) != MP_OKAY) { | 88 if ((res = mp_copy (b, &x)) != MP_OKAY) { |
90 goto __ERR; | 89 goto LBL_ERR; |
91 } | 90 } |
92 | 91 |
93 /* we need y = |a| */ | 92 /* we need y = |a| */ |
94 if ((res = mp_abs (a, &y)) != MP_OKAY) { | 93 if ((res = mp_mod (a, b, &y)) != MP_OKAY) { |
95 goto __ERR; | 94 goto LBL_ERR; |
96 } | 95 } |
97 | 96 |
98 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */ | 97 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */ |
99 if ((res = mp_copy (&x, &u)) != MP_OKAY) { | 98 if ((res = mp_copy (&x, &u)) != MP_OKAY) { |
100 goto __ERR; | 99 goto LBL_ERR; |
101 } | 100 } |
102 if ((res = mp_copy (&y, &v)) != MP_OKAY) { | 101 if ((res = mp_copy (&y, &v)) != MP_OKAY) { |
103 goto __ERR; | 102 goto LBL_ERR; |
104 } | 103 } |
105 mp_set (&D, 1); | 104 mp_set (&D, 1); |
106 | 105 |
107 top: | 106 top: |
108 /* 4. while u is even do */ | 107 /* 4. while u is even do */ |
109 while (mp_iseven (&u) == 1) { | 108 while (mp_iseven (&u) == 1) { |
110 /* 4.1 u = u/2 */ | 109 /* 4.1 u = u/2 */ |
111 if ((res = mp_div_2 (&u, &u)) != MP_OKAY) { | 110 if ((res = mp_div_2 (&u, &u)) != MP_OKAY) { |
112 goto __ERR; | 111 goto LBL_ERR; |
113 } | 112 } |
114 /* 4.2 if B is odd then */ | 113 /* 4.2 if B is odd then */ |
115 if (mp_isodd (&B) == 1) { | 114 if (mp_isodd (&B) == 1) { |
116 if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) { | 115 if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) { |
117 goto __ERR; | 116 goto LBL_ERR; |
118 } | 117 } |
119 } | 118 } |
120 /* B = B/2 */ | 119 /* B = B/2 */ |
121 if ((res = mp_div_2 (&B, &B)) != MP_OKAY) { | 120 if ((res = mp_div_2 (&B, &B)) != MP_OKAY) { |
122 goto __ERR; | 121 goto LBL_ERR; |
123 } | 122 } |
124 } | 123 } |
125 | 124 |
126 /* 5. while v is even do */ | 125 /* 5. while v is even do */ |
127 while (mp_iseven (&v) == 1) { | 126 while (mp_iseven (&v) == 1) { |
128 /* 5.1 v = v/2 */ | 127 /* 5.1 v = v/2 */ |
129 if ((res = mp_div_2 (&v, &v)) != MP_OKAY) { | 128 if ((res = mp_div_2 (&v, &v)) != MP_OKAY) { |
130 goto __ERR; | 129 goto LBL_ERR; |
131 } | 130 } |
132 /* 5.2 if D is odd then */ | 131 /* 5.2 if D is odd then */ |
133 if (mp_isodd (&D) == 1) { | 132 if (mp_isodd (&D) == 1) { |
134 /* D = (D-x)/2 */ | 133 /* D = (D-x)/2 */ |
135 if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) { | 134 if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) { |
136 goto __ERR; | 135 goto LBL_ERR; |
137 } | 136 } |
138 } | 137 } |
139 /* D = D/2 */ | 138 /* D = D/2 */ |
140 if ((res = mp_div_2 (&D, &D)) != MP_OKAY) { | 139 if ((res = mp_div_2 (&D, &D)) != MP_OKAY) { |
141 goto __ERR; | 140 goto LBL_ERR; |
142 } | 141 } |
143 } | 142 } |
144 | 143 |
145 /* 6. if u >= v then */ | 144 /* 6. if u >= v then */ |
146 if (mp_cmp (&u, &v) != MP_LT) { | 145 if (mp_cmp (&u, &v) != MP_LT) { |
147 /* u = u - v, B = B - D */ | 146 /* u = u - v, B = B - D */ |
148 if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) { | 147 if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) { |
149 goto __ERR; | 148 goto LBL_ERR; |
150 } | 149 } |
151 | 150 |
152 if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) { | 151 if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) { |
153 goto __ERR; | 152 goto LBL_ERR; |
154 } | 153 } |
155 } else { | 154 } else { |
156 /* v - v - u, D = D - B */ | 155 /* v - v - u, D = D - B */ |
157 if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) { | 156 if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) { |
158 goto __ERR; | 157 goto LBL_ERR; |
159 } | 158 } |
160 | 159 |
161 if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) { | 160 if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) { |
162 goto __ERR; | 161 goto LBL_ERR; |
163 } | 162 } |
164 } | 163 } |
165 | 164 |
166 /* if not zero goto step 4 */ | 165 /* if not zero goto step 4 */ |
167 if (mp_iszero (&u) == 0) { | 166 if (mp_iszero (&u) == 0) { |
171 /* now a = C, b = D, gcd == g*v */ | 170 /* now a = C, b = D, gcd == g*v */ |
172 | 171 |
173 /* if v != 1 then there is no inverse */ | 172 /* if v != 1 then there is no inverse */ |
174 if (mp_cmp_d (&v, 1) != MP_EQ) { | 173 if (mp_cmp_d (&v, 1) != MP_EQ) { |
175 res = MP_VAL; | 174 res = MP_VAL; |
176 goto __ERR; | 175 goto LBL_ERR; |
177 } | 176 } |
178 | 177 |
179 /* b is now the inverse */ | 178 /* b is now the inverse */ |
180 neg = a->sign; | 179 neg = a->sign; |
181 while (D.sign == MP_NEG) { | 180 while (D.sign == MP_NEG) { |
182 if ((res = mp_add (&D, b, &D)) != MP_OKAY) { | 181 if ((res = mp_add (&D, b, &D)) != MP_OKAY) { |
183 goto __ERR; | 182 goto LBL_ERR; |
184 } | 183 } |
185 } | 184 } |
186 mp_exch (&D, c); | 185 mp_exch (&D, c); |
187 c->sign = neg; | 186 c->sign = neg; |
188 res = MP_OKAY; | 187 res = MP_OKAY; |
189 | 188 |
190 __ERR:mp_clear_multi (&x, &y, &u, &v, &B, &D, NULL); | 189 LBL_ERR:mp_clear_multi (&x, &y, &u, &v, &B, &D, NULL); |
191 return res; | 190 return res; |
192 } | 191 } |
193 #endif | 192 #endif |
194 | 193 |
195 /* End: bn_fast_mp_invmod.c */ | 194 /* End: bn_fast_mp_invmod.c */ |
218 * which uses the comba method to quickly calculate the columns of the | 217 * which uses the comba method to quickly calculate the columns of the |
219 * reduction. | 218 * reduction. |
220 * | 219 * |
221 * Based on Algorithm 14.32 on pp.601 of HAC. | 220 * Based on Algorithm 14.32 on pp.601 of HAC. |
222 */ | 221 */ |
223 int | 222 int fast_mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho) |
224 fast_mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho) | |
225 { | 223 { |
226 int ix, res, olduse; | 224 int ix, res, olduse; |
227 mp_word W[MP_WARRAY]; | 225 mp_word W[MP_WARRAY]; |
228 | 226 |
229 /* get old used count */ | 227 /* get old used count */ |
399 * required for fast Barrett reduction). | 397 * required for fast Barrett reduction). |
400 * | 398 * |
401 * Based on Algorithm 14.12 on pp.595 of HAC. | 399 * Based on Algorithm 14.12 on pp.595 of HAC. |
402 * | 400 * |
403 */ | 401 */ |
404 int | 402 int fast_s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs) |
405 fast_s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs) | |
406 { | 403 { |
407 int olduse, res, pa, ix, iz; | 404 int olduse, res, pa, ix, iz; |
408 mp_digit W[MP_WARRAY]; | 405 mp_digit W[MP_WARRAY]; |
409 register mp_word _W; | 406 register mp_word _W; |
410 | 407 |
418 /* number of output digits to produce */ | 415 /* number of output digits to produce */ |
419 pa = MIN(digs, a->used + b->used); | 416 pa = MIN(digs, a->used + b->used); |
420 | 417 |
421 /* clear the carry */ | 418 /* clear the carry */ |
422 _W = 0; | 419 _W = 0; |
423 for (ix = 0; ix <= pa; ix++) { | 420 for (ix = 0; ix < pa; ix++) { |
424 int tx, ty; | 421 int tx, ty; |
425 int iy; | 422 int iy; |
426 mp_digit *tmpx, *tmpy; | 423 mp_digit *tmpx, *tmpy; |
427 | 424 |
428 /* get offsets into the two bignums */ | 425 /* get offsets into the two bignums */ |
431 | 428 |
432 /* setup temp aliases */ | 429 /* setup temp aliases */ |
433 tmpx = a->dp + tx; | 430 tmpx = a->dp + tx; |
434 tmpy = b->dp + ty; | 431 tmpy = b->dp + ty; |
435 | 432 |
436 /* this is the number of times the loop will iterrate, essentially its | 433 /* this is the number of times the loop will iterrate, essentially |
437 while (tx++ < a->used && ty-- >= 0) { ... } | 434 while (tx++ < a->used && ty-- >= 0) { ... } |
438 */ | 435 */ |
439 iy = MIN(a->used-tx, ty+1); | 436 iy = MIN(a->used-tx, ty+1); |
440 | 437 |
441 /* execute loop */ | 438 /* execute loop */ |
448 | 445 |
449 /* make next carry */ | 446 /* make next carry */ |
450 _W = _W >> ((mp_word)DIGIT_BIT); | 447 _W = _W >> ((mp_word)DIGIT_BIT); |
451 } | 448 } |
452 | 449 |
450 /* store final carry */ | |
451 W[ix] = (mp_digit)(_W & MP_MASK); | |
452 | |
453 /* setup dest */ | 453 /* setup dest */ |
454 olduse = c->used; | 454 olduse = c->used; |
455 c->used = digs; | 455 c->used = pa; |
456 | 456 |
457 { | 457 { |
458 register mp_digit *tmpc; | 458 register mp_digit *tmpc; |
459 tmpc = c->dp; | 459 tmpc = c->dp; |
460 for (ix = 0; ix < digs; ix++) { | 460 for (ix = 0; ix < pa+1; ix++) { |
461 /* now extract the previous digit [below the carry] */ | 461 /* now extract the previous digit [below the carry] */ |
462 *tmpc++ = W[ix]; | 462 *tmpc++ = W[ix]; |
463 } | 463 } |
464 | 464 |
465 /* clear unused digits [that existed in the old copy of c] */ | 465 /* clear unused digits [that existed in the old copy of c] */ |
499 * This is used in the Barrett reduction since for one of the multiplications | 499 * This is used in the Barrett reduction since for one of the multiplications |
500 * only the higher digits were needed. This essentially halves the work. | 500 * only the higher digits were needed. This essentially halves the work. |
501 * | 501 * |
502 * Based on Algorithm 14.12 on pp.595 of HAC. | 502 * Based on Algorithm 14.12 on pp.595 of HAC. |
503 */ | 503 */ |
504 int | 504 int fast_s_mp_mul_high_digs (mp_int * a, mp_int * b, mp_int * c, int digs) |
505 fast_s_mp_mul_high_digs (mp_int * a, mp_int * b, mp_int * c, int digs) | |
506 { | 505 { |
507 int olduse, res, pa, ix, iz; | 506 int olduse, res, pa, ix, iz; |
508 mp_digit W[MP_WARRAY]; | 507 mp_digit W[MP_WARRAY]; |
509 mp_word _W; | 508 mp_word _W; |
510 | 509 |
517 } | 516 } |
518 | 517 |
519 /* number of output digits to produce */ | 518 /* number of output digits to produce */ |
520 pa = a->used + b->used; | 519 pa = a->used + b->used; |
521 _W = 0; | 520 _W = 0; |
522 for (ix = digs; ix <= pa; ix++) { | 521 for (ix = digs; ix < pa; ix++) { |
523 int tx, ty, iy; | 522 int tx, ty, iy; |
524 mp_digit *tmpx, *tmpy; | 523 mp_digit *tmpx, *tmpy; |
525 | 524 |
526 /* get offsets into the two bignums */ | 525 /* get offsets into the two bignums */ |
527 ty = MIN(b->used-1, ix); | 526 ty = MIN(b->used-1, ix); |
545 W[ix] = ((mp_digit)_W) & MP_MASK; | 544 W[ix] = ((mp_digit)_W) & MP_MASK; |
546 | 545 |
547 /* make next carry */ | 546 /* make next carry */ |
548 _W = _W >> ((mp_word)DIGIT_BIT); | 547 _W = _W >> ((mp_word)DIGIT_BIT); |
549 } | 548 } |
549 | |
550 /* store final carry */ | |
551 W[ix] = (mp_digit)(_W & MP_MASK); | |
550 | 552 |
551 /* setup dest */ | 553 /* setup dest */ |
552 olduse = c->used; | 554 olduse = c->used; |
553 c->used = pa; | 555 c->used = pa; |
554 | 556 |
589 * guarantee it works. | 591 * guarantee it works. |
590 * | 592 * |
591 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 593 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
592 */ | 594 */ |
593 | 595 |
594 /* fast squaring | |
595 * | |
596 * This is the comba method where the columns of the product | |
597 * are computed first then the carries are computed. This | |
598 * has the effect of making a very simple inner loop that | |
599 * is executed the most | |
600 * | |
601 * W2 represents the outer products and W the inner. | |
602 * | |
603 * A further optimizations is made because the inner | |
604 * products are of the form "A * B * 2". The *2 part does | |
605 * not need to be computed until the end which is good | |
606 * because 64-bit shifts are slow! | |
607 * | |
608 * Based on Algorithm 14.16 on pp.597 of HAC. | |
609 * | |
610 */ | |
611 /* the jist of squaring... | 596 /* the jist of squaring... |
612 | 597 * you do like mult except the offset of the tmpx [one that |
613 you do like mult except the offset of the tmpx [one that starts closer to zero] | 598 * starts closer to zero] can't equal the offset of tmpy. |
614 can't equal the offset of tmpy. So basically you set up iy like before then you min it with | 599 * So basically you set up iy like before then you min it with |
615 (ty-tx) so that it never happens. You double all those you add in the inner loop | 600 * (ty-tx) so that it never happens. You double all those |
601 * you add in the inner loop | |
616 | 602 |
617 After that loop you do the squares and add them in. | 603 After that loop you do the squares and add them in. |
618 | |
619 Remove W2 and don't memset W | |
620 | |
621 */ | 604 */ |
622 | 605 |
623 int fast_s_mp_sqr (mp_int * a, mp_int * b) | 606 int fast_s_mp_sqr (mp_int * a, mp_int * b) |
624 { | 607 { |
625 int olduse, res, pa, ix, iz; | 608 int olduse, res, pa, ix, iz; |
634 } | 617 } |
635 } | 618 } |
636 | 619 |
637 /* number of output digits to produce */ | 620 /* number of output digits to produce */ |
638 W1 = 0; | 621 W1 = 0; |
639 for (ix = 0; ix <= pa; ix++) { | 622 for (ix = 0; ix < pa; ix++) { |
640 int tx, ty, iy; | 623 int tx, ty, iy; |
641 mp_word _W; | 624 mp_word _W; |
642 mp_digit *tmpy; | 625 mp_digit *tmpy; |
643 | 626 |
644 /* clear counter */ | 627 /* clear counter */ |
650 | 633 |
651 /* setup temp aliases */ | 634 /* setup temp aliases */ |
652 tmpx = a->dp + tx; | 635 tmpx = a->dp + tx; |
653 tmpy = a->dp + ty; | 636 tmpy = a->dp + ty; |
654 | 637 |
655 /* this is the number of times the loop will iterrate, essentially its | 638 /* this is the number of times the loop will iterrate, essentially |
656 while (tx++ < a->used && ty-- >= 0) { ... } | 639 while (tx++ < a->used && ty-- >= 0) { ... } |
657 */ | 640 */ |
658 iy = MIN(a->used-tx, ty+1); | 641 iy = MIN(a->used-tx, ty+1); |
659 | 642 |
660 /* now for squaring tx can never equal ty | 643 /* now for squaring tx can never equal ty |
675 if ((ix&1) == 0) { | 658 if ((ix&1) == 0) { |
676 _W += ((mp_word)a->dp[ix>>1])*((mp_word)a->dp[ix>>1]); | 659 _W += ((mp_word)a->dp[ix>>1])*((mp_word)a->dp[ix>>1]); |
677 } | 660 } |
678 | 661 |
679 /* store it */ | 662 /* store it */ |
680 W[ix] = _W; | 663 W[ix] = (mp_digit)(_W & MP_MASK); |
681 | 664 |
682 /* make next carry */ | 665 /* make next carry */ |
683 W1 = _W >> ((mp_word)DIGIT_BIT); | 666 W1 = _W >> ((mp_word)DIGIT_BIT); |
684 } | 667 } |
685 | 668 |
1537 } | 1520 } |
1538 | 1521 |
1539 | 1522 |
1540 mp_set(&tq, 1); | 1523 mp_set(&tq, 1); |
1541 n = mp_count_bits(a) - mp_count_bits(b); | 1524 n = mp_count_bits(a) - mp_count_bits(b); |
1542 if (((res = mp_copy(a, &ta)) != MP_OKAY) || | 1525 if (((res = mp_abs(a, &ta)) != MP_OKAY) || |
1543 ((res = mp_copy(b, &tb)) != MP_OKAY) || | 1526 ((res = mp_abs(b, &tb)) != MP_OKAY) || |
1544 ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) || | 1527 ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) || |
1545 ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) { | 1528 ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) { |
1546 goto __ERR; | 1529 goto LBL_ERR; |
1547 } | 1530 } |
1548 | 1531 |
1549 while (n-- >= 0) { | 1532 while (n-- >= 0) { |
1550 if (mp_cmp(&tb, &ta) != MP_GT) { | 1533 if (mp_cmp(&tb, &ta) != MP_GT) { |
1551 if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) || | 1534 if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) || |
1552 ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) { | 1535 ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) { |
1553 goto __ERR; | 1536 goto LBL_ERR; |
1554 } | 1537 } |
1555 } | 1538 } |
1556 if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) || | 1539 if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) || |
1557 ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) { | 1540 ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) { |
1558 goto __ERR; | 1541 goto LBL_ERR; |
1559 } | 1542 } |
1560 } | 1543 } |
1561 | 1544 |
1562 /* now q == quotient and ta == remainder */ | 1545 /* now q == quotient and ta == remainder */ |
1563 n = a->sign; | 1546 n = a->sign; |
1564 n2 = (a->sign == b->sign ? MP_ZPOS : MP_NEG); | 1547 n2 = (a->sign == b->sign ? MP_ZPOS : MP_NEG); |
1565 if (c != NULL) { | 1548 if (c != NULL) { |
1566 mp_exch(c, &q); | 1549 mp_exch(c, &q); |
1567 c->sign = n2; | 1550 c->sign = (mp_iszero(c) == MP_YES) ? MP_ZPOS : n2; |
1568 } | 1551 } |
1569 if (d != NULL) { | 1552 if (d != NULL) { |
1570 mp_exch(d, &ta); | 1553 mp_exch(d, &ta); |
1571 d->sign = n; | 1554 d->sign = (mp_iszero(d) == MP_YES) ? MP_ZPOS : n; |
1572 } | 1555 } |
1573 __ERR: | 1556 LBL_ERR: |
1574 mp_clear_multi(&ta, &tb, &tq, &q, NULL); | 1557 mp_clear_multi(&ta, &tb, &tq, &q, NULL); |
1575 return res; | 1558 return res; |
1576 } | 1559 } |
1577 | 1560 |
1578 #else | 1561 #else |
1617 return res; | 1600 return res; |
1618 } | 1601 } |
1619 q.used = a->used + 2; | 1602 q.used = a->used + 2; |
1620 | 1603 |
1621 if ((res = mp_init (&t1)) != MP_OKAY) { | 1604 if ((res = mp_init (&t1)) != MP_OKAY) { |
1622 goto __Q; | 1605 goto LBL_Q; |
1623 } | 1606 } |
1624 | 1607 |
1625 if ((res = mp_init (&t2)) != MP_OKAY) { | 1608 if ((res = mp_init (&t2)) != MP_OKAY) { |
1626 goto __T1; | 1609 goto LBL_T1; |
1627 } | 1610 } |
1628 | 1611 |
1629 if ((res = mp_init_copy (&x, a)) != MP_OKAY) { | 1612 if ((res = mp_init_copy (&x, a)) != MP_OKAY) { |
1630 goto __T2; | 1613 goto LBL_T2; |
1631 } | 1614 } |
1632 | 1615 |
1633 if ((res = mp_init_copy (&y, b)) != MP_OKAY) { | 1616 if ((res = mp_init_copy (&y, b)) != MP_OKAY) { |
1634 goto __X; | 1617 goto LBL_X; |
1635 } | 1618 } |
1636 | 1619 |
1637 /* fix the sign */ | 1620 /* fix the sign */ |
1638 neg = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; | 1621 neg = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; |
1639 x.sign = y.sign = MP_ZPOS; | 1622 x.sign = y.sign = MP_ZPOS; |
1641 /* normalize both x and y, ensure that y >= b/2, [b == 2**DIGIT_BIT] */ | 1624 /* normalize both x and y, ensure that y >= b/2, [b == 2**DIGIT_BIT] */ |
1642 norm = mp_count_bits(&y) % DIGIT_BIT; | 1625 norm = mp_count_bits(&y) % DIGIT_BIT; |
1643 if (norm < (int)(DIGIT_BIT-1)) { | 1626 if (norm < (int)(DIGIT_BIT-1)) { |
1644 norm = (DIGIT_BIT-1) - norm; | 1627 norm = (DIGIT_BIT-1) - norm; |
1645 if ((res = mp_mul_2d (&x, norm, &x)) != MP_OKAY) { | 1628 if ((res = mp_mul_2d (&x, norm, &x)) != MP_OKAY) { |
1646 goto __Y; | 1629 goto LBL_Y; |
1647 } | 1630 } |
1648 if ((res = mp_mul_2d (&y, norm, &y)) != MP_OKAY) { | 1631 if ((res = mp_mul_2d (&y, norm, &y)) != MP_OKAY) { |
1649 goto __Y; | 1632 goto LBL_Y; |
1650 } | 1633 } |
1651 } else { | 1634 } else { |
1652 norm = 0; | 1635 norm = 0; |
1653 } | 1636 } |
1654 | 1637 |
1656 n = x.used - 1; | 1639 n = x.used - 1; |
1657 t = y.used - 1; | 1640 t = y.used - 1; |
1658 | 1641 |
1659 /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */ | 1642 /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */ |
1660 if ((res = mp_lshd (&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */ | 1643 if ((res = mp_lshd (&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */ |
1661 goto __Y; | 1644 goto LBL_Y; |
1662 } | 1645 } |
1663 | 1646 |
1664 while (mp_cmp (&x, &y) != MP_LT) { | 1647 while (mp_cmp (&x, &y) != MP_LT) { |
1665 ++(q.dp[n - t]); | 1648 ++(q.dp[n - t]); |
1666 if ((res = mp_sub (&x, &y, &x)) != MP_OKAY) { | 1649 if ((res = mp_sub (&x, &y, &x)) != MP_OKAY) { |
1667 goto __Y; | 1650 goto LBL_Y; |
1668 } | 1651 } |
1669 } | 1652 } |
1670 | 1653 |
1671 /* reset y by shifting it back down */ | 1654 /* reset y by shifting it back down */ |
1672 mp_rshd (&y, n - t); | 1655 mp_rshd (&y, n - t); |
1704 mp_zero (&t1); | 1687 mp_zero (&t1); |
1705 t1.dp[0] = (t - 1 < 0) ? 0 : y.dp[t - 1]; | 1688 t1.dp[0] = (t - 1 < 0) ? 0 : y.dp[t - 1]; |
1706 t1.dp[1] = y.dp[t]; | 1689 t1.dp[1] = y.dp[t]; |
1707 t1.used = 2; | 1690 t1.used = 2; |
1708 if ((res = mp_mul_d (&t1, q.dp[i - t - 1], &t1)) != MP_OKAY) { | 1691 if ((res = mp_mul_d (&t1, q.dp[i - t - 1], &t1)) != MP_OKAY) { |
1709 goto __Y; | 1692 goto LBL_Y; |
1710 } | 1693 } |
1711 | 1694 |
1712 /* find right hand */ | 1695 /* find right hand */ |
1713 t2.dp[0] = (i - 2 < 0) ? 0 : x.dp[i - 2]; | 1696 t2.dp[0] = (i - 2 < 0) ? 0 : x.dp[i - 2]; |
1714 t2.dp[1] = (i - 1 < 0) ? 0 : x.dp[i - 1]; | 1697 t2.dp[1] = (i - 1 < 0) ? 0 : x.dp[i - 1]; |
1716 t2.used = 3; | 1699 t2.used = 3; |
1717 } while (mp_cmp_mag(&t1, &t2) == MP_GT); | 1700 } while (mp_cmp_mag(&t1, &t2) == MP_GT); |
1718 | 1701 |
1719 /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */ | 1702 /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */ |
1720 if ((res = mp_mul_d (&y, q.dp[i - t - 1], &t1)) != MP_OKAY) { | 1703 if ((res = mp_mul_d (&y, q.dp[i - t - 1], &t1)) != MP_OKAY) { |
1721 goto __Y; | 1704 goto LBL_Y; |
1722 } | 1705 } |
1723 | 1706 |
1724 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) { | 1707 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) { |
1725 goto __Y; | 1708 goto LBL_Y; |
1726 } | 1709 } |
1727 | 1710 |
1728 if ((res = mp_sub (&x, &t1, &x)) != MP_OKAY) { | 1711 if ((res = mp_sub (&x, &t1, &x)) != MP_OKAY) { |
1729 goto __Y; | 1712 goto LBL_Y; |
1730 } | 1713 } |
1731 | 1714 |
1732 /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */ | 1715 /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */ |
1733 if (x.sign == MP_NEG) { | 1716 if (x.sign == MP_NEG) { |
1734 if ((res = mp_copy (&y, &t1)) != MP_OKAY) { | 1717 if ((res = mp_copy (&y, &t1)) != MP_OKAY) { |
1735 goto __Y; | 1718 goto LBL_Y; |
1736 } | 1719 } |
1737 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) { | 1720 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) { |
1738 goto __Y; | 1721 goto LBL_Y; |
1739 } | 1722 } |
1740 if ((res = mp_add (&x, &t1, &x)) != MP_OKAY) { | 1723 if ((res = mp_add (&x, &t1, &x)) != MP_OKAY) { |
1741 goto __Y; | 1724 goto LBL_Y; |
1742 } | 1725 } |
1743 | 1726 |
1744 q.dp[i - t - 1] = (q.dp[i - t - 1] - 1UL) & MP_MASK; | 1727 q.dp[i - t - 1] = (q.dp[i - t - 1] - 1UL) & MP_MASK; |
1745 } | 1728 } |
1746 } | 1729 } |
1763 mp_exch (&x, d); | 1746 mp_exch (&x, d); |
1764 } | 1747 } |
1765 | 1748 |
1766 res = MP_OKAY; | 1749 res = MP_OKAY; |
1767 | 1750 |
1768 __Y:mp_clear (&y); | 1751 LBL_Y:mp_clear (&y); |
1769 __X:mp_clear (&x); | 1752 LBL_X:mp_clear (&x); |
1770 __T2:mp_clear (&t2); | 1753 LBL_T2:mp_clear (&t2); |
1771 __T1:mp_clear (&t1); | 1754 LBL_T1:mp_clear (&t1); |
1772 __Q:mp_clear (&q); | 1755 LBL_Q:mp_clear (&q); |
1773 return res; | 1756 return res; |
1774 } | 1757 } |
1775 | 1758 |
1776 #endif | 1759 #endif |
1777 | 1760 |
2197 /* reduce "x" in place modulo "n" using the Diminished Radix algorithm. | 2180 /* reduce "x" in place modulo "n" using the Diminished Radix algorithm. |
2198 * | 2181 * |
2199 * Based on algorithm from the paper | 2182 * Based on algorithm from the paper |
2200 * | 2183 * |
2201 * "Generating Efficient Primes for Discrete Log Cryptosystems" | 2184 * "Generating Efficient Primes for Discrete Log Cryptosystems" |
2202 * Chae Hoon Lim, Pil Loong Lee, | 2185 * Chae Hoon Lim, Pil Joong Lee, |
2203 * POSTECH Information Research Laboratories | 2186 * POSTECH Information Research Laboratories |
2204 * | 2187 * |
2205 * The modulus must be of a special format [see manual] | 2188 * The modulus must be of a special format [see manual] |
2206 * | 2189 * |
2207 * Has been modified to use algorithm 7.10 from the LTM book instead | 2190 * Has been modified to use algorithm 7.10 from the LTM book instead |
2455 err = mp_exptmod(&tmpG, &tmpX, P, Y); | 2438 err = mp_exptmod(&tmpG, &tmpX, P, Y); |
2456 mp_clear_multi(&tmpG, &tmpX, NULL); | 2439 mp_clear_multi(&tmpG, &tmpX, NULL); |
2457 return err; | 2440 return err; |
2458 #else | 2441 #else |
2459 /* no invmod */ | 2442 /* no invmod */ |
2460 return MP_VAL | 2443 return MP_VAL; |
2461 #endif | 2444 #endif |
2462 } | 2445 } |
2446 | |
2447 /* modified diminished radix reduction */ | |
2448 #if defined(BN_MP_REDUCE_IS_2K_L_C) && defined(BN_MP_REDUCE_2K_L_C) | |
2449 if (mp_reduce_is_2k_l(P) == MP_YES) { | |
2450 return s_mp_exptmod(G, X, P, Y, 1); | |
2451 } | |
2452 #endif | |
2463 | 2453 |
2464 #ifdef BN_MP_DR_IS_MODULUS_C | 2454 #ifdef BN_MP_DR_IS_MODULUS_C |
2465 /* is it a DR modulus? */ | 2455 /* is it a DR modulus? */ |
2466 dr = mp_dr_is_modulus(P); | 2456 dr = mp_dr_is_modulus(P); |
2467 #else | 2457 #else |
2458 /* default to no */ | |
2468 dr = 0; | 2459 dr = 0; |
2469 #endif | 2460 #endif |
2470 | 2461 |
2471 #ifdef BN_MP_REDUCE_IS_2K_C | 2462 #ifdef BN_MP_REDUCE_IS_2K_C |
2472 /* if not, is it a uDR modulus? */ | 2463 /* if not, is it a unrestricted DR modulus? */ |
2473 if (dr == 0) { | 2464 if (dr == 0) { |
2474 dr = mp_reduce_is_2k(P) << 1; | 2465 dr = mp_reduce_is_2k(P) << 1; |
2475 } | 2466 } |
2476 #endif | 2467 #endif |
2477 | 2468 |
2478 /* if the modulus is odd or dr != 0 use the fast method */ | 2469 /* if the modulus is odd or dr != 0 use the montgomery method */ |
2479 #ifdef BN_MP_EXPTMOD_FAST_C | 2470 #ifdef BN_MP_EXPTMOD_FAST_C |
2480 if (mp_isodd (P) == 1 || dr != 0) { | 2471 if (mp_isodd (P) == 1 || dr != 0) { |
2481 return mp_exptmod_fast (G, X, P, Y, dr); | 2472 return mp_exptmod_fast (G, X, P, Y, dr); |
2482 } else { | 2473 } else { |
2483 #endif | 2474 #endif |
2484 #ifdef BN_S_MP_EXPTMOD_C | 2475 #ifdef BN_S_MP_EXPTMOD_C |
2485 /* otherwise use the generic Barrett reduction technique */ | 2476 /* otherwise use the generic Barrett reduction technique */ |
2486 return s_mp_exptmod (G, X, P, Y); | 2477 return s_mp_exptmod (G, X, P, Y, 0); |
2487 #else | 2478 #else |
2488 /* no exptmod for evens */ | 2479 /* no exptmod for evens */ |
2489 return MP_VAL; | 2480 return MP_VAL; |
2490 #endif | 2481 #endif |
2491 #ifdef BN_MP_EXPTMOD_FAST_C | 2482 #ifdef BN_MP_EXPTMOD_FAST_C |
2527 #define TAB_SIZE 32 | 2518 #define TAB_SIZE 32 |
2528 #else | 2519 #else |
2529 #define TAB_SIZE 256 | 2520 #define TAB_SIZE 256 |
2530 #endif | 2521 #endif |
2531 | 2522 |
2532 int | 2523 int mp_exptmod_fast (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode) |
2533 mp_exptmod_fast (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode) | |
2534 { | 2524 { |
2535 mp_int M[TAB_SIZE], res; | 2525 mp_int M[TAB_SIZE], res; |
2536 mp_digit buf, mp; | 2526 mp_digit buf, mp; |
2537 int err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize; | 2527 int err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize; |
2538 | 2528 |
2586 /* determine and setup reduction code */ | 2576 /* determine and setup reduction code */ |
2587 if (redmode == 0) { | 2577 if (redmode == 0) { |
2588 #ifdef BN_MP_MONTGOMERY_SETUP_C | 2578 #ifdef BN_MP_MONTGOMERY_SETUP_C |
2589 /* now setup montgomery */ | 2579 /* now setup montgomery */ |
2590 if ((err = mp_montgomery_setup (P, &mp)) != MP_OKAY) { | 2580 if ((err = mp_montgomery_setup (P, &mp)) != MP_OKAY) { |
2591 goto __M; | 2581 goto LBL_M; |
2592 } | 2582 } |
2593 #else | 2583 #else |
2594 err = MP_VAL; | 2584 err = MP_VAL; |
2595 goto __M; | 2585 goto LBL_M; |
2596 #endif | 2586 #endif |
2597 | 2587 |
2598 /* automatically pick the comba one if available (saves quite a few calls/ifs) */ | 2588 /* automatically pick the comba one if available (saves quite a few calls/ifs) */ |
2599 #ifdef BN_FAST_MP_MONTGOMERY_REDUCE_C | 2589 #ifdef BN_FAST_MP_MONTGOMERY_REDUCE_C |
2600 if (((P->used * 2 + 1) < MP_WARRAY) && | 2590 if (((P->used * 2 + 1) < MP_WARRAY) && |
2606 #ifdef BN_MP_MONTGOMERY_REDUCE_C | 2596 #ifdef BN_MP_MONTGOMERY_REDUCE_C |
2607 /* use slower baseline Montgomery method */ | 2597 /* use slower baseline Montgomery method */ |
2608 redux = mp_montgomery_reduce; | 2598 redux = mp_montgomery_reduce; |
2609 #else | 2599 #else |
2610 err = MP_VAL; | 2600 err = MP_VAL; |
2611 goto __M; | 2601 goto LBL_M; |
2612 #endif | 2602 #endif |
2613 } | 2603 } |
2614 } else if (redmode == 1) { | 2604 } else if (redmode == 1) { |
2615 #if defined(BN_MP_DR_SETUP_C) && defined(BN_MP_DR_REDUCE_C) | 2605 #if defined(BN_MP_DR_SETUP_C) && defined(BN_MP_DR_REDUCE_C) |
2616 /* setup DR reduction for moduli of the form B**k - b */ | 2606 /* setup DR reduction for moduli of the form B**k - b */ |
2617 mp_dr_setup(P, &mp); | 2607 mp_dr_setup(P, &mp); |
2618 redux = mp_dr_reduce; | 2608 redux = mp_dr_reduce; |
2619 #else | 2609 #else |
2620 err = MP_VAL; | 2610 err = MP_VAL; |
2621 goto __M; | 2611 goto LBL_M; |
2622 #endif | 2612 #endif |
2623 } else { | 2613 } else { |
2624 #if defined(BN_MP_REDUCE_2K_SETUP_C) && defined(BN_MP_REDUCE_2K_C) | 2614 #if defined(BN_MP_REDUCE_2K_SETUP_C) && defined(BN_MP_REDUCE_2K_C) |
2625 /* setup DR reduction for moduli of the form 2**k - b */ | 2615 /* setup DR reduction for moduli of the form 2**k - b */ |
2626 if ((err = mp_reduce_2k_setup(P, &mp)) != MP_OKAY) { | 2616 if ((err = mp_reduce_2k_setup(P, &mp)) != MP_OKAY) { |
2627 goto __M; | 2617 goto LBL_M; |
2628 } | 2618 } |
2629 redux = mp_reduce_2k; | 2619 redux = mp_reduce_2k; |
2630 #else | 2620 #else |
2631 err = MP_VAL; | 2621 err = MP_VAL; |
2632 goto __M; | 2622 goto LBL_M; |
2633 #endif | 2623 #endif |
2634 } | 2624 } |
2635 | 2625 |
2636 /* setup result */ | 2626 /* setup result */ |
2637 if ((err = mp_init (&res)) != MP_OKAY) { | 2627 if ((err = mp_init (&res)) != MP_OKAY) { |
2638 goto __M; | 2628 goto LBL_M; |
2639 } | 2629 } |
2640 | 2630 |
2641 /* create M table | 2631 /* create M table |
2642 * | 2632 * |
2643 | 2633 |
2647 | 2637 |
2648 if (redmode == 0) { | 2638 if (redmode == 0) { |
2649 #ifdef BN_MP_MONTGOMERY_CALC_NORMALIZATION_C | 2639 #ifdef BN_MP_MONTGOMERY_CALC_NORMALIZATION_C |
2650 /* now we need R mod m */ | 2640 /* now we need R mod m */ |
2651 if ((err = mp_montgomery_calc_normalization (&res, P)) != MP_OKAY) { | 2641 if ((err = mp_montgomery_calc_normalization (&res, P)) != MP_OKAY) { |
2652 goto __RES; | 2642 goto LBL_RES; |
2653 } | 2643 } |
2654 #else | 2644 #else |
2655 err = MP_VAL; | 2645 err = MP_VAL; |
2656 goto __RES; | 2646 goto LBL_RES; |
2657 #endif | 2647 #endif |
2658 | 2648 |
2659 /* now set M[1] to G * R mod m */ | 2649 /* now set M[1] to G * R mod m */ |
2660 if ((err = mp_mulmod (G, &res, P, &M[1])) != MP_OKAY) { | 2650 if ((err = mp_mulmod (G, &res, P, &M[1])) != MP_OKAY) { |
2661 goto __RES; | 2651 goto LBL_RES; |
2662 } | 2652 } |
2663 } else { | 2653 } else { |
2664 mp_set(&res, 1); | 2654 mp_set(&res, 1); |
2665 if ((err = mp_mod(G, P, &M[1])) != MP_OKAY) { | 2655 if ((err = mp_mod(G, P, &M[1])) != MP_OKAY) { |
2666 goto __RES; | 2656 goto LBL_RES; |
2667 } | 2657 } |
2668 } | 2658 } |
2669 | 2659 |
2670 /* compute the value at M[1<<(winsize-1)] by squaring M[1] (winsize-1) times */ | 2660 /* compute the value at M[1<<(winsize-1)] by squaring M[1] (winsize-1) times */ |
2671 if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) { | 2661 if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) { |
2672 goto __RES; | 2662 goto LBL_RES; |
2673 } | 2663 } |
2674 | 2664 |
2675 for (x = 0; x < (winsize - 1); x++) { | 2665 for (x = 0; x < (winsize - 1); x++) { |
2676 if ((err = mp_sqr (&M[1 << (winsize - 1)], &M[1 << (winsize - 1)])) != MP_OKAY) { | 2666 if ((err = mp_sqr (&M[1 << (winsize - 1)], &M[1 << (winsize - 1)])) != MP_OKAY) { |
2677 goto __RES; | 2667 goto LBL_RES; |
2678 } | 2668 } |
2679 if ((err = redux (&M[1 << (winsize - 1)], P, mp)) != MP_OKAY) { | 2669 if ((err = redux (&M[1 << (winsize - 1)], P, mp)) != MP_OKAY) { |
2680 goto __RES; | 2670 goto LBL_RES; |
2681 } | 2671 } |
2682 } | 2672 } |
2683 | 2673 |
2684 /* create upper table */ | 2674 /* create upper table */ |
2685 for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) { | 2675 for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) { |
2686 if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) { | 2676 if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) { |
2687 goto __RES; | 2677 goto LBL_RES; |
2688 } | 2678 } |
2689 if ((err = redux (&M[x], P, mp)) != MP_OKAY) { | 2679 if ((err = redux (&M[x], P, mp)) != MP_OKAY) { |
2690 goto __RES; | 2680 goto LBL_RES; |
2691 } | 2681 } |
2692 } | 2682 } |
2693 | 2683 |
2694 /* set initial mode and bit cnt */ | 2684 /* set initial mode and bit cnt */ |
2695 mode = 0; | 2685 mode = 0; |
2725 } | 2715 } |
2726 | 2716 |
2727 /* if the bit is zero and mode == 1 then we square */ | 2717 /* if the bit is zero and mode == 1 then we square */ |
2728 if (mode == 1 && y == 0) { | 2718 if (mode == 1 && y == 0) { |
2729 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 2719 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
2730 goto __RES; | 2720 goto LBL_RES; |
2731 } | 2721 } |
2732 if ((err = redux (&res, P, mp)) != MP_OKAY) { | 2722 if ((err = redux (&res, P, mp)) != MP_OKAY) { |
2733 goto __RES; | 2723 goto LBL_RES; |
2734 } | 2724 } |
2735 continue; | 2725 continue; |
2736 } | 2726 } |
2737 | 2727 |
2738 /* else we add it to the window */ | 2728 /* else we add it to the window */ |
2742 if (bitcpy == winsize) { | 2732 if (bitcpy == winsize) { |
2743 /* ok window is filled so square as required and multiply */ | 2733 /* ok window is filled so square as required and multiply */ |
2744 /* square first */ | 2734 /* square first */ |
2745 for (x = 0; x < winsize; x++) { | 2735 for (x = 0; x < winsize; x++) { |
2746 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 2736 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
2747 goto __RES; | 2737 goto LBL_RES; |
2748 } | 2738 } |
2749 if ((err = redux (&res, P, mp)) != MP_OKAY) { | 2739 if ((err = redux (&res, P, mp)) != MP_OKAY) { |
2750 goto __RES; | 2740 goto LBL_RES; |
2751 } | 2741 } |
2752 } | 2742 } |
2753 | 2743 |
2754 /* then multiply */ | 2744 /* then multiply */ |
2755 if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) { | 2745 if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) { |
2756 goto __RES; | 2746 goto LBL_RES; |
2757 } | 2747 } |
2758 if ((err = redux (&res, P, mp)) != MP_OKAY) { | 2748 if ((err = redux (&res, P, mp)) != MP_OKAY) { |
2759 goto __RES; | 2749 goto LBL_RES; |
2760 } | 2750 } |
2761 | 2751 |
2762 /* empty window and reset */ | 2752 /* empty window and reset */ |
2763 bitcpy = 0; | 2753 bitcpy = 0; |
2764 bitbuf = 0; | 2754 bitbuf = 0; |
2769 /* if bits remain then square/multiply */ | 2759 /* if bits remain then square/multiply */ |
2770 if (mode == 2 && bitcpy > 0) { | 2760 if (mode == 2 && bitcpy > 0) { |
2771 /* square then multiply if the bit is set */ | 2761 /* square then multiply if the bit is set */ |
2772 for (x = 0; x < bitcpy; x++) { | 2762 for (x = 0; x < bitcpy; x++) { |
2773 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 2763 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
2774 goto __RES; | 2764 goto LBL_RES; |
2775 } | 2765 } |
2776 if ((err = redux (&res, P, mp)) != MP_OKAY) { | 2766 if ((err = redux (&res, P, mp)) != MP_OKAY) { |
2777 goto __RES; | 2767 goto LBL_RES; |
2778 } | 2768 } |
2779 | 2769 |
2780 /* get next bit of the window */ | 2770 /* get next bit of the window */ |
2781 bitbuf <<= 1; | 2771 bitbuf <<= 1; |
2782 if ((bitbuf & (1 << winsize)) != 0) { | 2772 if ((bitbuf & (1 << winsize)) != 0) { |
2783 /* then multiply */ | 2773 /* then multiply */ |
2784 if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) { | 2774 if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) { |
2785 goto __RES; | 2775 goto LBL_RES; |
2786 } | 2776 } |
2787 if ((err = redux (&res, P, mp)) != MP_OKAY) { | 2777 if ((err = redux (&res, P, mp)) != MP_OKAY) { |
2788 goto __RES; | 2778 goto LBL_RES; |
2789 } | 2779 } |
2790 } | 2780 } |
2791 } | 2781 } |
2792 } | 2782 } |
2793 | 2783 |
2797 * actually multiplied by R mod n. So we have | 2787 * actually multiplied by R mod n. So we have |
2798 * to reduce one more time to cancel out the factor | 2788 * to reduce one more time to cancel out the factor |
2799 * of R. | 2789 * of R. |
2800 */ | 2790 */ |
2801 if ((err = redux(&res, P, mp)) != MP_OKAY) { | 2791 if ((err = redux(&res, P, mp)) != MP_OKAY) { |
2802 goto __RES; | 2792 goto LBL_RES; |
2803 } | 2793 } |
2804 } | 2794 } |
2805 | 2795 |
2806 /* swap res with Y */ | 2796 /* swap res with Y */ |
2807 mp_exch (&res, Y); | 2797 mp_exch (&res, Y); |
2808 err = MP_OKAY; | 2798 err = MP_OKAY; |
2809 __RES:mp_clear (&res); | 2799 LBL_RES:mp_clear (&res); |
2810 __M: | 2800 LBL_M: |
2811 mp_clear(&M[1]); | 2801 mp_clear(&M[1]); |
2812 for (x = 1<<(winsize-1); x < (1 << winsize); x++) { | 2802 for (x = 1<<(winsize-1); x < (1 << winsize); x++) { |
2813 mp_clear (&M[x]); | 2803 mp_clear (&M[x]); |
2814 } | 2804 } |
2815 return err; | 2805 return err; |
2877 | 2867 |
2878 /* (v1,v2,v3) = (t1,t2,t3) */ | 2868 /* (v1,v2,v3) = (t1,t2,t3) */ |
2879 if ((err = mp_copy(&t1, &v1)) != MP_OKAY) { goto _ERR; } | 2869 if ((err = mp_copy(&t1, &v1)) != MP_OKAY) { goto _ERR; } |
2880 if ((err = mp_copy(&t2, &v2)) != MP_OKAY) { goto _ERR; } | 2870 if ((err = mp_copy(&t2, &v2)) != MP_OKAY) { goto _ERR; } |
2881 if ((err = mp_copy(&t3, &v3)) != MP_OKAY) { goto _ERR; } | 2871 if ((err = mp_copy(&t3, &v3)) != MP_OKAY) { goto _ERR; } |
2872 } | |
2873 | |
2874 /* make sure U3 >= 0 */ | |
2875 if (u3.sign == MP_NEG) { | |
2876 mp_neg(&u1, &u1); | |
2877 mp_neg(&u2, &u2); | |
2878 mp_neg(&u3, &u3); | |
2882 } | 2879 } |
2883 | 2880 |
2884 /* copy result out */ | 2881 /* copy result out */ |
2885 if (U1 != NULL) { mp_exch(U1, &u1); } | 2882 if (U1 != NULL) { mp_exch(U1, &u1); } |
2886 if (U2 != NULL) { mp_exch(U2, &u2); } | 2883 if (U2 != NULL) { mp_exch(U2, &u2); } |
3057 if ((res = mp_init_copy (&u, a)) != MP_OKAY) { | 3054 if ((res = mp_init_copy (&u, a)) != MP_OKAY) { |
3058 return res; | 3055 return res; |
3059 } | 3056 } |
3060 | 3057 |
3061 if ((res = mp_init_copy (&v, b)) != MP_OKAY) { | 3058 if ((res = mp_init_copy (&v, b)) != MP_OKAY) { |
3062 goto __U; | 3059 goto LBL_U; |
3063 } | 3060 } |
3064 | 3061 |
3065 /* must be positive for the remainder of the algorithm */ | 3062 /* must be positive for the remainder of the algorithm */ |
3066 u.sign = v.sign = MP_ZPOS; | 3063 u.sign = v.sign = MP_ZPOS; |
3067 | 3064 |
3071 k = MIN(u_lsb, v_lsb); | 3068 k = MIN(u_lsb, v_lsb); |
3072 | 3069 |
3073 if (k > 0) { | 3070 if (k > 0) { |
3074 /* divide the power of two out */ | 3071 /* divide the power of two out */ |
3075 if ((res = mp_div_2d(&u, k, &u, NULL)) != MP_OKAY) { | 3072 if ((res = mp_div_2d(&u, k, &u, NULL)) != MP_OKAY) { |
3076 goto __V; | 3073 goto LBL_V; |
3077 } | 3074 } |
3078 | 3075 |
3079 if ((res = mp_div_2d(&v, k, &v, NULL)) != MP_OKAY) { | 3076 if ((res = mp_div_2d(&v, k, &v, NULL)) != MP_OKAY) { |
3080 goto __V; | 3077 goto LBL_V; |
3081 } | 3078 } |
3082 } | 3079 } |
3083 | 3080 |
3084 /* divide any remaining factors of two out */ | 3081 /* divide any remaining factors of two out */ |
3085 if (u_lsb != k) { | 3082 if (u_lsb != k) { |
3086 if ((res = mp_div_2d(&u, u_lsb - k, &u, NULL)) != MP_OKAY) { | 3083 if ((res = mp_div_2d(&u, u_lsb - k, &u, NULL)) != MP_OKAY) { |
3087 goto __V; | 3084 goto LBL_V; |
3088 } | 3085 } |
3089 } | 3086 } |
3090 | 3087 |
3091 if (v_lsb != k) { | 3088 if (v_lsb != k) { |
3092 if ((res = mp_div_2d(&v, v_lsb - k, &v, NULL)) != MP_OKAY) { | 3089 if ((res = mp_div_2d(&v, v_lsb - k, &v, NULL)) != MP_OKAY) { |
3093 goto __V; | 3090 goto LBL_V; |
3094 } | 3091 } |
3095 } | 3092 } |
3096 | 3093 |
3097 while (mp_iszero(&v) == 0) { | 3094 while (mp_iszero(&v) == 0) { |
3098 /* make sure v is the largest */ | 3095 /* make sure v is the largest */ |
3101 mp_exch(&u, &v); | 3098 mp_exch(&u, &v); |
3102 } | 3099 } |
3103 | 3100 |
3104 /* subtract smallest from largest */ | 3101 /* subtract smallest from largest */ |
3105 if ((res = s_mp_sub(&v, &u, &v)) != MP_OKAY) { | 3102 if ((res = s_mp_sub(&v, &u, &v)) != MP_OKAY) { |
3106 goto __V; | 3103 goto LBL_V; |
3107 } | 3104 } |
3108 | 3105 |
3109 /* Divide out all factors of two */ | 3106 /* Divide out all factors of two */ |
3110 if ((res = mp_div_2d(&v, mp_cnt_lsb(&v), &v, NULL)) != MP_OKAY) { | 3107 if ((res = mp_div_2d(&v, mp_cnt_lsb(&v), &v, NULL)) != MP_OKAY) { |
3111 goto __V; | 3108 goto LBL_V; |
3112 } | 3109 } |
3113 } | 3110 } |
3114 | 3111 |
3115 /* multiply by 2**k which we divided out at the beginning */ | 3112 /* multiply by 2**k which we divided out at the beginning */ |
3116 if ((res = mp_mul_2d (&u, k, c)) != MP_OKAY) { | 3113 if ((res = mp_mul_2d (&u, k, c)) != MP_OKAY) { |
3117 goto __V; | 3114 goto LBL_V; |
3118 } | 3115 } |
3119 c->sign = MP_ZPOS; | 3116 c->sign = MP_ZPOS; |
3120 res = MP_OKAY; | 3117 res = MP_OKAY; |
3121 __V:mp_clear (&u); | 3118 LBL_V:mp_clear (&u); |
3122 __U:mp_clear (&v); | 3119 LBL_U:mp_clear (&v); |
3123 return res; | 3120 return res; |
3124 } | 3121 } |
3125 #endif | 3122 #endif |
3126 | 3123 |
3127 /* End: bn_mp_gcd.c */ | 3124 /* End: bn_mp_gcd.c */ |
3553 &A, &B, &C, &D, NULL)) != MP_OKAY) { | 3550 &A, &B, &C, &D, NULL)) != MP_OKAY) { |
3554 return res; | 3551 return res; |
3555 } | 3552 } |
3556 | 3553 |
3557 /* x = a, y = b */ | 3554 /* x = a, y = b */ |
3558 if ((res = mp_copy (a, &x)) != MP_OKAY) { | 3555 if ((res = mp_mod(a, b, &x)) != MP_OKAY) { |
3559 goto __ERR; | 3556 goto LBL_ERR; |
3560 } | 3557 } |
3561 if ((res = mp_copy (b, &y)) != MP_OKAY) { | 3558 if ((res = mp_copy (b, &y)) != MP_OKAY) { |
3562 goto __ERR; | 3559 goto LBL_ERR; |
3563 } | 3560 } |
3564 | 3561 |
3565 /* 2. [modified] if x,y are both even then return an error! */ | 3562 /* 2. [modified] if x,y are both even then return an error! */ |
3566 if (mp_iseven (&x) == 1 && mp_iseven (&y) == 1) { | 3563 if (mp_iseven (&x) == 1 && mp_iseven (&y) == 1) { |
3567 res = MP_VAL; | 3564 res = MP_VAL; |
3568 goto __ERR; | 3565 goto LBL_ERR; |
3569 } | 3566 } |
3570 | 3567 |
3571 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */ | 3568 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */ |
3572 if ((res = mp_copy (&x, &u)) != MP_OKAY) { | 3569 if ((res = mp_copy (&x, &u)) != MP_OKAY) { |
3573 goto __ERR; | 3570 goto LBL_ERR; |
3574 } | 3571 } |
3575 if ((res = mp_copy (&y, &v)) != MP_OKAY) { | 3572 if ((res = mp_copy (&y, &v)) != MP_OKAY) { |
3576 goto __ERR; | 3573 goto LBL_ERR; |
3577 } | 3574 } |
3578 mp_set (&A, 1); | 3575 mp_set (&A, 1); |
3579 mp_set (&D, 1); | 3576 mp_set (&D, 1); |
3580 | 3577 |
3581 top: | 3578 top: |
3582 /* 4. while u is even do */ | 3579 /* 4. while u is even do */ |
3583 while (mp_iseven (&u) == 1) { | 3580 while (mp_iseven (&u) == 1) { |
3584 /* 4.1 u = u/2 */ | 3581 /* 4.1 u = u/2 */ |
3585 if ((res = mp_div_2 (&u, &u)) != MP_OKAY) { | 3582 if ((res = mp_div_2 (&u, &u)) != MP_OKAY) { |
3586 goto __ERR; | 3583 goto LBL_ERR; |
3587 } | 3584 } |
3588 /* 4.2 if A or B is odd then */ | 3585 /* 4.2 if A or B is odd then */ |
3589 if (mp_isodd (&A) == 1 || mp_isodd (&B) == 1) { | 3586 if (mp_isodd (&A) == 1 || mp_isodd (&B) == 1) { |
3590 /* A = (A+y)/2, B = (B-x)/2 */ | 3587 /* A = (A+y)/2, B = (B-x)/2 */ |
3591 if ((res = mp_add (&A, &y, &A)) != MP_OKAY) { | 3588 if ((res = mp_add (&A, &y, &A)) != MP_OKAY) { |
3592 goto __ERR; | 3589 goto LBL_ERR; |
3593 } | 3590 } |
3594 if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) { | 3591 if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) { |
3595 goto __ERR; | 3592 goto LBL_ERR; |
3596 } | 3593 } |
3597 } | 3594 } |
3598 /* A = A/2, B = B/2 */ | 3595 /* A = A/2, B = B/2 */ |
3599 if ((res = mp_div_2 (&A, &A)) != MP_OKAY) { | 3596 if ((res = mp_div_2 (&A, &A)) != MP_OKAY) { |
3600 goto __ERR; | 3597 goto LBL_ERR; |
3601 } | 3598 } |
3602 if ((res = mp_div_2 (&B, &B)) != MP_OKAY) { | 3599 if ((res = mp_div_2 (&B, &B)) != MP_OKAY) { |
3603 goto __ERR; | 3600 goto LBL_ERR; |
3604 } | 3601 } |
3605 } | 3602 } |
3606 | 3603 |
3607 /* 5. while v is even do */ | 3604 /* 5. while v is even do */ |
3608 while (mp_iseven (&v) == 1) { | 3605 while (mp_iseven (&v) == 1) { |
3609 /* 5.1 v = v/2 */ | 3606 /* 5.1 v = v/2 */ |
3610 if ((res = mp_div_2 (&v, &v)) != MP_OKAY) { | 3607 if ((res = mp_div_2 (&v, &v)) != MP_OKAY) { |
3611 goto __ERR; | 3608 goto LBL_ERR; |
3612 } | 3609 } |
3613 /* 5.2 if C or D is odd then */ | 3610 /* 5.2 if C or D is odd then */ |
3614 if (mp_isodd (&C) == 1 || mp_isodd (&D) == 1) { | 3611 if (mp_isodd (&C) == 1 || mp_isodd (&D) == 1) { |
3615 /* C = (C+y)/2, D = (D-x)/2 */ | 3612 /* C = (C+y)/2, D = (D-x)/2 */ |
3616 if ((res = mp_add (&C, &y, &C)) != MP_OKAY) { | 3613 if ((res = mp_add (&C, &y, &C)) != MP_OKAY) { |
3617 goto __ERR; | 3614 goto LBL_ERR; |
3618 } | 3615 } |
3619 if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) { | 3616 if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) { |
3620 goto __ERR; | 3617 goto LBL_ERR; |
3621 } | 3618 } |
3622 } | 3619 } |
3623 /* C = C/2, D = D/2 */ | 3620 /* C = C/2, D = D/2 */ |
3624 if ((res = mp_div_2 (&C, &C)) != MP_OKAY) { | 3621 if ((res = mp_div_2 (&C, &C)) != MP_OKAY) { |
3625 goto __ERR; | 3622 goto LBL_ERR; |
3626 } | 3623 } |
3627 if ((res = mp_div_2 (&D, &D)) != MP_OKAY) { | 3624 if ((res = mp_div_2 (&D, &D)) != MP_OKAY) { |
3628 goto __ERR; | 3625 goto LBL_ERR; |
3629 } | 3626 } |
3630 } | 3627 } |
3631 | 3628 |
3632 /* 6. if u >= v then */ | 3629 /* 6. if u >= v then */ |
3633 if (mp_cmp (&u, &v) != MP_LT) { | 3630 if (mp_cmp (&u, &v) != MP_LT) { |
3634 /* u = u - v, A = A - C, B = B - D */ | 3631 /* u = u - v, A = A - C, B = B - D */ |
3635 if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) { | 3632 if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) { |
3636 goto __ERR; | 3633 goto LBL_ERR; |
3637 } | 3634 } |
3638 | 3635 |
3639 if ((res = mp_sub (&A, &C, &A)) != MP_OKAY) { | 3636 if ((res = mp_sub (&A, &C, &A)) != MP_OKAY) { |
3640 goto __ERR; | 3637 goto LBL_ERR; |
3641 } | 3638 } |
3642 | 3639 |
3643 if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) { | 3640 if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) { |
3644 goto __ERR; | 3641 goto LBL_ERR; |
3645 } | 3642 } |
3646 } else { | 3643 } else { |
3647 /* v - v - u, C = C - A, D = D - B */ | 3644 /* v - v - u, C = C - A, D = D - B */ |
3648 if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) { | 3645 if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) { |
3649 goto __ERR; | 3646 goto LBL_ERR; |
3650 } | 3647 } |
3651 | 3648 |
3652 if ((res = mp_sub (&C, &A, &C)) != MP_OKAY) { | 3649 if ((res = mp_sub (&C, &A, &C)) != MP_OKAY) { |
3653 goto __ERR; | 3650 goto LBL_ERR; |
3654 } | 3651 } |
3655 | 3652 |
3656 if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) { | 3653 if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) { |
3657 goto __ERR; | 3654 goto LBL_ERR; |
3658 } | 3655 } |
3659 } | 3656 } |
3660 | 3657 |
3661 /* if not zero goto step 4 */ | 3658 /* if not zero goto step 4 */ |
3662 if (mp_iszero (&u) == 0) | 3659 if (mp_iszero (&u) == 0) |
3665 /* now a = C, b = D, gcd == g*v */ | 3662 /* now a = C, b = D, gcd == g*v */ |
3666 | 3663 |
3667 /* if v != 1 then there is no inverse */ | 3664 /* if v != 1 then there is no inverse */ |
3668 if (mp_cmp_d (&v, 1) != MP_EQ) { | 3665 if (mp_cmp_d (&v, 1) != MP_EQ) { |
3669 res = MP_VAL; | 3666 res = MP_VAL; |
3670 goto __ERR; | 3667 goto LBL_ERR; |
3671 } | 3668 } |
3672 | 3669 |
3673 /* if its too low */ | 3670 /* if its too low */ |
3674 while (mp_cmp_d(&C, 0) == MP_LT) { | 3671 while (mp_cmp_d(&C, 0) == MP_LT) { |
3675 if ((res = mp_add(&C, b, &C)) != MP_OKAY) { | 3672 if ((res = mp_add(&C, b, &C)) != MP_OKAY) { |
3676 goto __ERR; | 3673 goto LBL_ERR; |
3677 } | 3674 } |
3678 } | 3675 } |
3679 | 3676 |
3680 /* too big */ | 3677 /* too big */ |
3681 while (mp_cmp_mag(&C, b) != MP_LT) { | 3678 while (mp_cmp_mag(&C, b) != MP_LT) { |
3682 if ((res = mp_sub(&C, b, &C)) != MP_OKAY) { | 3679 if ((res = mp_sub(&C, b, &C)) != MP_OKAY) { |
3683 goto __ERR; | 3680 goto LBL_ERR; |
3684 } | 3681 } |
3685 } | 3682 } |
3686 | 3683 |
3687 /* C is now the inverse */ | 3684 /* C is now the inverse */ |
3688 mp_exch (&C, c); | 3685 mp_exch (&C, c); |
3689 res = MP_OKAY; | 3686 res = MP_OKAY; |
3690 __ERR:mp_clear_multi (&x, &y, &u, &v, &A, &B, &C, &D, NULL); | 3687 LBL_ERR:mp_clear_multi (&x, &y, &u, &v, &A, &B, &C, &D, NULL); |
3691 return res; | 3688 return res; |
3692 } | 3689 } |
3693 #endif | 3690 #endif |
3694 | 3691 |
3695 /* End: bn_mp_invmod_slow.c */ | 3692 /* End: bn_mp_invmod_slow.c */ |
3854 if ((res = mp_init_copy (&a1, a)) != MP_OKAY) { | 3851 if ((res = mp_init_copy (&a1, a)) != MP_OKAY) { |
3855 return res; | 3852 return res; |
3856 } | 3853 } |
3857 | 3854 |
3858 if ((res = mp_init (&p1)) != MP_OKAY) { | 3855 if ((res = mp_init (&p1)) != MP_OKAY) { |
3859 goto __A1; | 3856 goto LBL_A1; |
3860 } | 3857 } |
3861 | 3858 |
3862 /* divide out larger power of two */ | 3859 /* divide out larger power of two */ |
3863 k = mp_cnt_lsb(&a1); | 3860 k = mp_cnt_lsb(&a1); |
3864 if ((res = mp_div_2d(&a1, k, &a1, NULL)) != MP_OKAY) { | 3861 if ((res = mp_div_2d(&a1, k, &a1, NULL)) != MP_OKAY) { |
3865 goto __P1; | 3862 goto LBL_P1; |
3866 } | 3863 } |
3867 | 3864 |
3868 /* step 4. if e is even set s=1 */ | 3865 /* step 4. if e is even set s=1 */ |
3869 if ((k & 1) == 0) { | 3866 if ((k & 1) == 0) { |
3870 s = 1; | 3867 s = 1; |
3888 if (mp_cmp_d (&a1, 1) == MP_EQ) { | 3885 if (mp_cmp_d (&a1, 1) == MP_EQ) { |
3889 *c = s; | 3886 *c = s; |
3890 } else { | 3887 } else { |
3891 /* n1 = n mod a1 */ | 3888 /* n1 = n mod a1 */ |
3892 if ((res = mp_mod (p, &a1, &p1)) != MP_OKAY) { | 3889 if ((res = mp_mod (p, &a1, &p1)) != MP_OKAY) { |
3893 goto __P1; | 3890 goto LBL_P1; |
3894 } | 3891 } |
3895 if ((res = mp_jacobi (&p1, &a1, &r)) != MP_OKAY) { | 3892 if ((res = mp_jacobi (&p1, &a1, &r)) != MP_OKAY) { |
3896 goto __P1; | 3893 goto LBL_P1; |
3897 } | 3894 } |
3898 *c = s * r; | 3895 *c = s * r; |
3899 } | 3896 } |
3900 | 3897 |
3901 /* done */ | 3898 /* done */ |
3902 res = MP_OKAY; | 3899 res = MP_OKAY; |
3903 __P1:mp_clear (&p1); | 3900 LBL_P1:mp_clear (&p1); |
3904 __A1:mp_clear (&a1); | 3901 LBL_A1:mp_clear (&a1); |
3905 return res; | 3902 return res; |
3906 } | 3903 } |
3907 #endif | 3904 #endif |
3908 | 3905 |
3909 /* End: bn_mp_jacobi.c */ | 3906 /* End: bn_mp_jacobi.c */ |
4225 return res; | 4222 return res; |
4226 } | 4223 } |
4227 | 4224 |
4228 /* t1 = get the GCD of the two inputs */ | 4225 /* t1 = get the GCD of the two inputs */ |
4229 if ((res = mp_gcd (a, b, &t1)) != MP_OKAY) { | 4226 if ((res = mp_gcd (a, b, &t1)) != MP_OKAY) { |
4230 goto __T; | 4227 goto LBL_T; |
4231 } | 4228 } |
4232 | 4229 |
4233 /* divide the smallest by the GCD */ | 4230 /* divide the smallest by the GCD */ |
4234 if (mp_cmp_mag(a, b) == MP_LT) { | 4231 if (mp_cmp_mag(a, b) == MP_LT) { |
4235 /* store quotient in t2 such that t2 * b is the LCM */ | 4232 /* store quotient in t2 such that t2 * b is the LCM */ |
4236 if ((res = mp_div(a, &t1, &t2, NULL)) != MP_OKAY) { | 4233 if ((res = mp_div(a, &t1, &t2, NULL)) != MP_OKAY) { |
4237 goto __T; | 4234 goto LBL_T; |
4238 } | 4235 } |
4239 res = mp_mul(b, &t2, c); | 4236 res = mp_mul(b, &t2, c); |
4240 } else { | 4237 } else { |
4241 /* store quotient in t2 such that t2 * a is the LCM */ | 4238 /* store quotient in t2 such that t2 * a is the LCM */ |
4242 if ((res = mp_div(b, &t1, &t2, NULL)) != MP_OKAY) { | 4239 if ((res = mp_div(b, &t1, &t2, NULL)) != MP_OKAY) { |
4243 goto __T; | 4240 goto LBL_T; |
4244 } | 4241 } |
4245 res = mp_mul(a, &t2, c); | 4242 res = mp_mul(a, &t2, c); |
4246 } | 4243 } |
4247 | 4244 |
4248 /* fix the sign to positive */ | 4245 /* fix the sign to positive */ |
4249 c->sign = MP_ZPOS; | 4246 c->sign = MP_ZPOS; |
4250 | 4247 |
4251 __T: | 4248 LBL_T: |
4252 mp_clear_multi (&t1, &t2, NULL); | 4249 mp_clear_multi (&t1, &t2, NULL); |
4253 return res; | 4250 return res; |
4254 } | 4251 } |
4255 #endif | 4252 #endif |
4256 | 4253 |
4400 mp_zero (c); | 4397 mp_zero (c); |
4401 return MP_OKAY; | 4398 return MP_OKAY; |
4402 } | 4399 } |
4403 | 4400 |
4404 /* if the modulus is larger than the value than return */ | 4401 /* if the modulus is larger than the value than return */ |
4405 if (b > (int) (a->used * DIGIT_BIT)) { | 4402 if (b >= (int) (a->used * DIGIT_BIT)) { |
4406 res = mp_copy (a, c); | 4403 res = mp_copy (a, c); |
4407 return res; | 4404 return res; |
4408 } | 4405 } |
4409 | 4406 |
4410 /* copy */ | 4407 /* copy */ |
4481 { | 4478 { |
4482 int x, bits, res; | 4479 int x, bits, res; |
4483 | 4480 |
4484 /* how many bits of last digit does b use */ | 4481 /* how many bits of last digit does b use */ |
4485 bits = mp_count_bits (b) % DIGIT_BIT; | 4482 bits = mp_count_bits (b) % DIGIT_BIT; |
4486 | |
4487 | 4483 |
4488 if (b->used > 1) { | 4484 if (b->used > 1) { |
4489 if ((res = mp_2expt (a, (b->used - 1) * DIGIT_BIT + bits - 1)) != MP_OKAY) { | 4485 if ((res = mp_2expt (a, (b->used - 1) * DIGIT_BIT + bits - 1)) != MP_OKAY) { |
4490 return res; | 4486 return res; |
4491 } | 4487 } |
4981 | 4977 |
4982 /* send carry into next iteration */ | 4978 /* send carry into next iteration */ |
4983 u = (mp_digit) (r >> ((mp_word) DIGIT_BIT)); | 4979 u = (mp_digit) (r >> ((mp_word) DIGIT_BIT)); |
4984 } | 4980 } |
4985 | 4981 |
4986 /* store final carry [if any] */ | 4982 /* store final carry [if any] and increment ix offset */ |
4987 *tmpc++ = u; | 4983 *tmpc++ = u; |
4984 ++ix; | |
4988 | 4985 |
4989 /* now zero digits above the top */ | 4986 /* now zero digits above the top */ |
4990 while (ix++ < olduse) { | 4987 while (ix++ < olduse) { |
4991 *tmpc++ = 0; | 4988 *tmpc++ = 0; |
4992 } | 4989 } |
5083 if ((res = mp_init (&t1)) != MP_OKAY) { | 5080 if ((res = mp_init (&t1)) != MP_OKAY) { |
5084 return res; | 5081 return res; |
5085 } | 5082 } |
5086 | 5083 |
5087 if ((res = mp_init (&t2)) != MP_OKAY) { | 5084 if ((res = mp_init (&t2)) != MP_OKAY) { |
5088 goto __T1; | 5085 goto LBL_T1; |
5089 } | 5086 } |
5090 | 5087 |
5091 if ((res = mp_init (&t3)) != MP_OKAY) { | 5088 if ((res = mp_init (&t3)) != MP_OKAY) { |
5092 goto __T2; | 5089 goto LBL_T2; |
5093 } | 5090 } |
5094 | 5091 |
5095 /* if a is negative fudge the sign but keep track */ | 5092 /* if a is negative fudge the sign but keep track */ |
5096 neg = a->sign; | 5093 neg = a->sign; |
5097 a->sign = MP_ZPOS; | 5094 a->sign = MP_ZPOS; |
5100 mp_set (&t2, 2); | 5097 mp_set (&t2, 2); |
5101 | 5098 |
5102 do { | 5099 do { |
5103 /* t1 = t2 */ | 5100 /* t1 = t2 */ |
5104 if ((res = mp_copy (&t2, &t1)) != MP_OKAY) { | 5101 if ((res = mp_copy (&t2, &t1)) != MP_OKAY) { |
5105 goto __T3; | 5102 goto LBL_T3; |
5106 } | 5103 } |
5107 | 5104 |
5108 /* t2 = t1 - ((t1**b - a) / (b * t1**(b-1))) */ | 5105 /* t2 = t1 - ((t1**b - a) / (b * t1**(b-1))) */ |
5109 | 5106 |
5110 /* t3 = t1**(b-1) */ | 5107 /* t3 = t1**(b-1) */ |
5111 if ((res = mp_expt_d (&t1, b - 1, &t3)) != MP_OKAY) { | 5108 if ((res = mp_expt_d (&t1, b - 1, &t3)) != MP_OKAY) { |
5112 goto __T3; | 5109 goto LBL_T3; |
5113 } | 5110 } |
5114 | 5111 |
5115 /* numerator */ | 5112 /* numerator */ |
5116 /* t2 = t1**b */ | 5113 /* t2 = t1**b */ |
5117 if ((res = mp_mul (&t3, &t1, &t2)) != MP_OKAY) { | 5114 if ((res = mp_mul (&t3, &t1, &t2)) != MP_OKAY) { |
5118 goto __T3; | 5115 goto LBL_T3; |
5119 } | 5116 } |
5120 | 5117 |
5121 /* t2 = t1**b - a */ | 5118 /* t2 = t1**b - a */ |
5122 if ((res = mp_sub (&t2, a, &t2)) != MP_OKAY) { | 5119 if ((res = mp_sub (&t2, a, &t2)) != MP_OKAY) { |
5123 goto __T3; | 5120 goto LBL_T3; |
5124 } | 5121 } |
5125 | 5122 |
5126 /* denominator */ | 5123 /* denominator */ |
5127 /* t3 = t1**(b-1) * b */ | 5124 /* t3 = t1**(b-1) * b */ |
5128 if ((res = mp_mul_d (&t3, b, &t3)) != MP_OKAY) { | 5125 if ((res = mp_mul_d (&t3, b, &t3)) != MP_OKAY) { |
5129 goto __T3; | 5126 goto LBL_T3; |
5130 } | 5127 } |
5131 | 5128 |
5132 /* t3 = (t1**b - a)/(b * t1**(b-1)) */ | 5129 /* t3 = (t1**b - a)/(b * t1**(b-1)) */ |
5133 if ((res = mp_div (&t2, &t3, &t3, NULL)) != MP_OKAY) { | 5130 if ((res = mp_div (&t2, &t3, &t3, NULL)) != MP_OKAY) { |
5134 goto __T3; | 5131 goto LBL_T3; |
5135 } | 5132 } |
5136 | 5133 |
5137 if ((res = mp_sub (&t1, &t3, &t2)) != MP_OKAY) { | 5134 if ((res = mp_sub (&t1, &t3, &t2)) != MP_OKAY) { |
5138 goto __T3; | 5135 goto LBL_T3; |
5139 } | 5136 } |
5140 } while (mp_cmp (&t1, &t2) != MP_EQ); | 5137 } while (mp_cmp (&t1, &t2) != MP_EQ); |
5141 | 5138 |
5142 /* result can be off by a few so check */ | 5139 /* result can be off by a few so check */ |
5143 for (;;) { | 5140 for (;;) { |
5144 if ((res = mp_expt_d (&t1, b, &t2)) != MP_OKAY) { | 5141 if ((res = mp_expt_d (&t1, b, &t2)) != MP_OKAY) { |
5145 goto __T3; | 5142 goto LBL_T3; |
5146 } | 5143 } |
5147 | 5144 |
5148 if (mp_cmp (&t2, a) == MP_GT) { | 5145 if (mp_cmp (&t2, a) == MP_GT) { |
5149 if ((res = mp_sub_d (&t1, 1, &t1)) != MP_OKAY) { | 5146 if ((res = mp_sub_d (&t1, 1, &t1)) != MP_OKAY) { |
5150 goto __T3; | 5147 goto LBL_T3; |
5151 } | 5148 } |
5152 } else { | 5149 } else { |
5153 break; | 5150 break; |
5154 } | 5151 } |
5155 } | 5152 } |
5163 /* set the sign of the result */ | 5160 /* set the sign of the result */ |
5164 c->sign = neg; | 5161 c->sign = neg; |
5165 | 5162 |
5166 res = MP_OKAY; | 5163 res = MP_OKAY; |
5167 | 5164 |
5168 __T3:mp_clear (&t3); | 5165 LBL_T3:mp_clear (&t3); |
5169 __T2:mp_clear (&t2); | 5166 LBL_T2:mp_clear (&t2); |
5170 __T1:mp_clear (&t1); | 5167 LBL_T1:mp_clear (&t1); |
5171 return res; | 5168 return res; |
5172 } | 5169 } |
5173 #endif | 5170 #endif |
5174 | 5171 |
5175 /* End: bn_mp_n_root.c */ | 5172 /* End: bn_mp_n_root.c */ |
5194 | 5191 |
5195 /* b = -a */ | 5192 /* b = -a */ |
5196 int mp_neg (mp_int * a, mp_int * b) | 5193 int mp_neg (mp_int * a, mp_int * b) |
5197 { | 5194 { |
5198 int res; | 5195 int res; |
5199 if ((res = mp_copy (a, b)) != MP_OKAY) { | 5196 if (a != b) { |
5200 return res; | 5197 if ((res = mp_copy (a, b)) != MP_OKAY) { |
5201 } | 5198 return res; |
5199 } | |
5200 } | |
5201 | |
5202 if (mp_iszero(b) != MP_YES) { | 5202 if (mp_iszero(b) != MP_YES) { |
5203 b->sign = (a->sign == MP_ZPOS) ? MP_NEG : MP_ZPOS; | 5203 b->sign = (a->sign == MP_ZPOS) ? MP_NEG : MP_ZPOS; |
5204 } | 5204 } else { |
5205 b->sign = MP_ZPOS; | |
5206 } | |
5207 | |
5205 return MP_OKAY; | 5208 return MP_OKAY; |
5206 } | 5209 } |
5207 #endif | 5210 #endif |
5208 | 5211 |
5209 /* End: bn_mp_neg.c */ | 5212 /* End: bn_mp_neg.c */ |
5302 return err; | 5305 return err; |
5303 } | 5306 } |
5304 | 5307 |
5305 /* compute t = b**a mod a */ | 5308 /* compute t = b**a mod a */ |
5306 if ((err = mp_exptmod (b, a, a, &t)) != MP_OKAY) { | 5309 if ((err = mp_exptmod (b, a, a, &t)) != MP_OKAY) { |
5307 goto __T; | 5310 goto LBL_T; |
5308 } | 5311 } |
5309 | 5312 |
5310 /* is it equal to b? */ | 5313 /* is it equal to b? */ |
5311 if (mp_cmp (&t, b) == MP_EQ) { | 5314 if (mp_cmp (&t, b) == MP_EQ) { |
5312 *result = MP_YES; | 5315 *result = MP_YES; |
5313 } | 5316 } |
5314 | 5317 |
5315 err = MP_OKAY; | 5318 err = MP_OKAY; |
5316 __T:mp_clear (&t); | 5319 LBL_T:mp_clear (&t); |
5317 return err; | 5320 return err; |
5318 } | 5321 } |
5319 #endif | 5322 #endif |
5320 | 5323 |
5321 /* End: bn_mp_prime_fermat.c */ | 5324 /* End: bn_mp_prime_fermat.c */ |
5350 | 5353 |
5351 /* default to not */ | 5354 /* default to not */ |
5352 *result = MP_NO; | 5355 *result = MP_NO; |
5353 | 5356 |
5354 for (ix = 0; ix < PRIME_SIZE; ix++) { | 5357 for (ix = 0; ix < PRIME_SIZE; ix++) { |
5355 /* what is a mod __prime_tab[ix] */ | 5358 /* what is a mod LBL_prime_tab[ix] */ |
5356 if ((err = mp_mod_d (a, __prime_tab[ix], &res)) != MP_OKAY) { | 5359 if ((err = mp_mod_d (a, ltm_prime_tab[ix], &res)) != MP_OKAY) { |
5357 return err; | 5360 return err; |
5358 } | 5361 } |
5359 | 5362 |
5360 /* is the residue zero? */ | 5363 /* is the residue zero? */ |
5361 if (res == 0) { | 5364 if (res == 0) { |
5408 return MP_VAL; | 5411 return MP_VAL; |
5409 } | 5412 } |
5410 | 5413 |
5411 /* is the input equal to one of the primes in the table? */ | 5414 /* is the input equal to one of the primes in the table? */ |
5412 for (ix = 0; ix < PRIME_SIZE; ix++) { | 5415 for (ix = 0; ix < PRIME_SIZE; ix++) { |
5413 if (mp_cmp_d(a, __prime_tab[ix]) == MP_EQ) { | 5416 if (mp_cmp_d(a, ltm_prime_tab[ix]) == MP_EQ) { |
5414 *result = 1; | 5417 *result = 1; |
5415 return MP_OKAY; | 5418 return MP_OKAY; |
5416 } | 5419 } |
5417 } | 5420 } |
5418 | 5421 |
5431 return err; | 5434 return err; |
5432 } | 5435 } |
5433 | 5436 |
5434 for (ix = 0; ix < t; ix++) { | 5437 for (ix = 0; ix < t; ix++) { |
5435 /* set the prime */ | 5438 /* set the prime */ |
5436 mp_set (&b, __prime_tab[ix]); | 5439 mp_set (&b, ltm_prime_tab[ix]); |
5437 | 5440 |
5438 if ((err = mp_prime_miller_rabin (a, &b, &res)) != MP_OKAY) { | 5441 if ((err = mp_prime_miller_rabin (a, &b, &res)) != MP_OKAY) { |
5439 goto __B; | 5442 goto LBL_B; |
5440 } | 5443 } |
5441 | 5444 |
5442 if (res == MP_NO) { | 5445 if (res == MP_NO) { |
5443 goto __B; | 5446 goto LBL_B; |
5444 } | 5447 } |
5445 } | 5448 } |
5446 | 5449 |
5447 /* passed the test */ | 5450 /* passed the test */ |
5448 *result = MP_YES; | 5451 *result = MP_YES; |
5449 __B:mp_clear (&b); | 5452 LBL_B:mp_clear (&b); |
5450 return err; | 5453 return err; |
5451 } | 5454 } |
5452 #endif | 5455 #endif |
5453 | 5456 |
5454 /* End: bn_mp_prime_is_prime.c */ | 5457 /* End: bn_mp_prime_is_prime.c */ |
5494 /* get n1 = a - 1 */ | 5497 /* get n1 = a - 1 */ |
5495 if ((err = mp_init_copy (&n1, a)) != MP_OKAY) { | 5498 if ((err = mp_init_copy (&n1, a)) != MP_OKAY) { |
5496 return err; | 5499 return err; |
5497 } | 5500 } |
5498 if ((err = mp_sub_d (&n1, 1, &n1)) != MP_OKAY) { | 5501 if ((err = mp_sub_d (&n1, 1, &n1)) != MP_OKAY) { |
5499 goto __N1; | 5502 goto LBL_N1; |
5500 } | 5503 } |
5501 | 5504 |
5502 /* set 2**s * r = n1 */ | 5505 /* set 2**s * r = n1 */ |
5503 if ((err = mp_init_copy (&r, &n1)) != MP_OKAY) { | 5506 if ((err = mp_init_copy (&r, &n1)) != MP_OKAY) { |
5504 goto __N1; | 5507 goto LBL_N1; |
5505 } | 5508 } |
5506 | 5509 |
5507 /* count the number of least significant bits | 5510 /* count the number of least significant bits |
5508 * which are zero | 5511 * which are zero |
5509 */ | 5512 */ |
5510 s = mp_cnt_lsb(&r); | 5513 s = mp_cnt_lsb(&r); |
5511 | 5514 |
5512 /* now divide n - 1 by 2**s */ | 5515 /* now divide n - 1 by 2**s */ |
5513 if ((err = mp_div_2d (&r, s, &r, NULL)) != MP_OKAY) { | 5516 if ((err = mp_div_2d (&r, s, &r, NULL)) != MP_OKAY) { |
5514 goto __R; | 5517 goto LBL_R; |
5515 } | 5518 } |
5516 | 5519 |
5517 /* compute y = b**r mod a */ | 5520 /* compute y = b**r mod a */ |
5518 if ((err = mp_init (&y)) != MP_OKAY) { | 5521 if ((err = mp_init (&y)) != MP_OKAY) { |
5519 goto __R; | 5522 goto LBL_R; |
5520 } | 5523 } |
5521 if ((err = mp_exptmod (b, &r, a, &y)) != MP_OKAY) { | 5524 if ((err = mp_exptmod (b, &r, a, &y)) != MP_OKAY) { |
5522 goto __Y; | 5525 goto LBL_Y; |
5523 } | 5526 } |
5524 | 5527 |
5525 /* if y != 1 and y != n1 do */ | 5528 /* if y != 1 and y != n1 do */ |
5526 if (mp_cmp_d (&y, 1) != MP_EQ && mp_cmp (&y, &n1) != MP_EQ) { | 5529 if (mp_cmp_d (&y, 1) != MP_EQ && mp_cmp (&y, &n1) != MP_EQ) { |
5527 j = 1; | 5530 j = 1; |
5528 /* while j <= s-1 and y != n1 */ | 5531 /* while j <= s-1 and y != n1 */ |
5529 while ((j <= (s - 1)) && mp_cmp (&y, &n1) != MP_EQ) { | 5532 while ((j <= (s - 1)) && mp_cmp (&y, &n1) != MP_EQ) { |
5530 if ((err = mp_sqrmod (&y, a, &y)) != MP_OKAY) { | 5533 if ((err = mp_sqrmod (&y, a, &y)) != MP_OKAY) { |
5531 goto __Y; | 5534 goto LBL_Y; |
5532 } | 5535 } |
5533 | 5536 |
5534 /* if y == 1 then composite */ | 5537 /* if y == 1 then composite */ |
5535 if (mp_cmp_d (&y, 1) == MP_EQ) { | 5538 if (mp_cmp_d (&y, 1) == MP_EQ) { |
5536 goto __Y; | 5539 goto LBL_Y; |
5537 } | 5540 } |
5538 | 5541 |
5539 ++j; | 5542 ++j; |
5540 } | 5543 } |
5541 | 5544 |
5542 /* if y != n1 then composite */ | 5545 /* if y != n1 then composite */ |
5543 if (mp_cmp (&y, &n1) != MP_EQ) { | 5546 if (mp_cmp (&y, &n1) != MP_EQ) { |
5544 goto __Y; | 5547 goto LBL_Y; |
5545 } | 5548 } |
5546 } | 5549 } |
5547 | 5550 |
5548 /* probably prime now */ | 5551 /* probably prime now */ |
5549 *result = MP_YES; | 5552 *result = MP_YES; |
5550 __Y:mp_clear (&y); | 5553 LBL_Y:mp_clear (&y); |
5551 __R:mp_clear (&r); | 5554 LBL_R:mp_clear (&r); |
5552 __N1:mp_clear (&n1); | 5555 LBL_N1:mp_clear (&n1); |
5553 return err; | 5556 return err; |
5554 } | 5557 } |
5555 #endif | 5558 #endif |
5556 | 5559 |
5557 /* End: bn_mp_prime_miller_rabin.c */ | 5560 /* End: bn_mp_prime_miller_rabin.c */ |
5592 | 5595 |
5593 /* force positive */ | 5596 /* force positive */ |
5594 a->sign = MP_ZPOS; | 5597 a->sign = MP_ZPOS; |
5595 | 5598 |
5596 /* simple algo if a is less than the largest prime in the table */ | 5599 /* simple algo if a is less than the largest prime in the table */ |
5597 if (mp_cmp_d(a, __prime_tab[PRIME_SIZE-1]) == MP_LT) { | 5600 if (mp_cmp_d(a, ltm_prime_tab[PRIME_SIZE-1]) == MP_LT) { |
5598 /* find which prime it is bigger than */ | 5601 /* find which prime it is bigger than */ |
5599 for (x = PRIME_SIZE - 2; x >= 0; x--) { | 5602 for (x = PRIME_SIZE - 2; x >= 0; x--) { |
5600 if (mp_cmp_d(a, __prime_tab[x]) != MP_LT) { | 5603 if (mp_cmp_d(a, ltm_prime_tab[x]) != MP_LT) { |
5601 if (bbs_style == 1) { | 5604 if (bbs_style == 1) { |
5602 /* ok we found a prime smaller or | 5605 /* ok we found a prime smaller or |
5603 * equal [so the next is larger] | 5606 * equal [so the next is larger] |
5604 * | 5607 * |
5605 * however, the prime must be | 5608 * however, the prime must be |
5606 * congruent to 3 mod 4 | 5609 * congruent to 3 mod 4 |
5607 */ | 5610 */ |
5608 if ((__prime_tab[x + 1] & 3) != 3) { | 5611 if ((ltm_prime_tab[x + 1] & 3) != 3) { |
5609 /* scan upwards for a prime congruent to 3 mod 4 */ | 5612 /* scan upwards for a prime congruent to 3 mod 4 */ |
5610 for (y = x + 1; y < PRIME_SIZE; y++) { | 5613 for (y = x + 1; y < PRIME_SIZE; y++) { |
5611 if ((__prime_tab[y] & 3) == 3) { | 5614 if ((ltm_prime_tab[y] & 3) == 3) { |
5612 mp_set(a, __prime_tab[y]); | 5615 mp_set(a, ltm_prime_tab[y]); |
5613 return MP_OKAY; | 5616 return MP_OKAY; |
5614 } | 5617 } |
5615 } | 5618 } |
5616 } | 5619 } |
5617 } else { | 5620 } else { |
5618 mp_set(a, __prime_tab[x + 1]); | 5621 mp_set(a, ltm_prime_tab[x + 1]); |
5619 return MP_OKAY; | 5622 return MP_OKAY; |
5620 } | 5623 } |
5621 } | 5624 } |
5622 } | 5625 } |
5623 /* at this point a maybe 1 */ | 5626 /* at this point a maybe 1 */ |
5651 } | 5654 } |
5652 } | 5655 } |
5653 | 5656 |
5654 /* generate the restable */ | 5657 /* generate the restable */ |
5655 for (x = 1; x < PRIME_SIZE; x++) { | 5658 for (x = 1; x < PRIME_SIZE; x++) { |
5656 if ((err = mp_mod_d(a, __prime_tab[x], res_tab + x)) != MP_OKAY) { | 5659 if ((err = mp_mod_d(a, ltm_prime_tab[x], res_tab + x)) != MP_OKAY) { |
5657 return err; | 5660 return err; |
5658 } | 5661 } |
5659 } | 5662 } |
5660 | 5663 |
5661 /* init temp used for Miller-Rabin Testing */ | 5664 /* init temp used for Miller-Rabin Testing */ |
5677 for (x = 1; x < PRIME_SIZE; x++) { | 5680 for (x = 1; x < PRIME_SIZE; x++) { |
5678 /* add the step to each residue */ | 5681 /* add the step to each residue */ |
5679 res_tab[x] += kstep; | 5682 res_tab[x] += kstep; |
5680 | 5683 |
5681 /* subtract the modulus [instead of using division] */ | 5684 /* subtract the modulus [instead of using division] */ |
5682 if (res_tab[x] >= __prime_tab[x]) { | 5685 if (res_tab[x] >= ltm_prime_tab[x]) { |
5683 res_tab[x] -= __prime_tab[x]; | 5686 res_tab[x] -= ltm_prime_tab[x]; |
5684 } | 5687 } |
5685 | 5688 |
5686 /* set flag if zero */ | 5689 /* set flag if zero */ |
5687 if (res_tab[x] == 0) { | 5690 if (res_tab[x] == 0) { |
5688 y = 1; | 5691 y = 1; |
5690 } | 5693 } |
5691 } while (y == 1 && step < ((((mp_digit)1)<<DIGIT_BIT) - kstep)); | 5694 } while (y == 1 && step < ((((mp_digit)1)<<DIGIT_BIT) - kstep)); |
5692 | 5695 |
5693 /* add the step */ | 5696 /* add the step */ |
5694 if ((err = mp_add_d(a, step, a)) != MP_OKAY) { | 5697 if ((err = mp_add_d(a, step, a)) != MP_OKAY) { |
5695 goto __ERR; | 5698 goto LBL_ERR; |
5696 } | 5699 } |
5697 | 5700 |
5698 /* if didn't pass sieve and step == MAX then skip test */ | 5701 /* if didn't pass sieve and step == MAX then skip test */ |
5699 if (y == 1 && step >= ((((mp_digit)1)<<DIGIT_BIT) - kstep)) { | 5702 if (y == 1 && step >= ((((mp_digit)1)<<DIGIT_BIT) - kstep)) { |
5700 continue; | 5703 continue; |
5701 } | 5704 } |
5702 | 5705 |
5703 /* is this prime? */ | 5706 /* is this prime? */ |
5704 for (x = 0; x < t; x++) { | 5707 for (x = 0; x < t; x++) { |
5705 mp_set(&b, __prime_tab[t]); | 5708 mp_set(&b, ltm_prime_tab[t]); |
5706 if ((err = mp_prime_miller_rabin(a, &b, &res)) != MP_OKAY) { | 5709 if ((err = mp_prime_miller_rabin(a, &b, &res)) != MP_OKAY) { |
5707 goto __ERR; | 5710 goto LBL_ERR; |
5708 } | 5711 } |
5709 if (res == MP_NO) { | 5712 if (res == MP_NO) { |
5710 break; | 5713 break; |
5711 } | 5714 } |
5712 } | 5715 } |
5715 break; | 5718 break; |
5716 } | 5719 } |
5717 } | 5720 } |
5718 | 5721 |
5719 err = MP_OKAY; | 5722 err = MP_OKAY; |
5720 __ERR: | 5723 LBL_ERR: |
5721 mp_clear(&b); | 5724 mp_clear(&b); |
5722 return err; | 5725 return err; |
5723 } | 5726 } |
5724 | 5727 |
5725 #endif | 5728 #endif |
5826 if (flags & LTM_PRIME_SAFE) { | 5829 if (flags & LTM_PRIME_SAFE) { |
5827 flags |= LTM_PRIME_BBS; | 5830 flags |= LTM_PRIME_BBS; |
5828 } | 5831 } |
5829 | 5832 |
5830 /* calc the byte size */ | 5833 /* calc the byte size */ |
5831 bsize = (size>>3)+(size&7?1:0); | 5834 bsize = (size>>3) + ((size&7)?1:0); |
5832 | 5835 |
5833 /* we need a buffer of bsize bytes */ | 5836 /* we need a buffer of bsize bytes */ |
5834 tmp = OPT_CAST(unsigned char) XMALLOC(bsize); | 5837 tmp = OPT_CAST(unsigned char) XMALLOC(bsize); |
5835 if (tmp == NULL) { | 5838 if (tmp == NULL) { |
5836 return MP_MEM; | 5839 return MP_MEM; |
5837 } | 5840 } |
5838 | 5841 |
5839 /* calc the maskAND value for the MSbyte*/ | 5842 /* calc the maskAND value for the MSbyte*/ |
5840 maskAND = 0xFF >> (8 - (size & 7)); | 5843 maskAND = ((size&7) == 0) ? 0xFF : (0xFF >> (8 - (size & 7))); |
5841 | 5844 |
5842 /* calc the maskOR_msb */ | 5845 /* calc the maskOR_msb */ |
5843 maskOR_msb = 0; | 5846 maskOR_msb = 0; |
5844 maskOR_msb_offset = (size - 2) >> 3; | 5847 maskOR_msb_offset = ((size & 7) == 1) ? 1 : 0; |
5845 if (flags & LTM_PRIME_2MSB_ON) { | 5848 if (flags & LTM_PRIME_2MSB_ON) { |
5846 maskOR_msb |= 1 << ((size - 2) & 7); | 5849 maskOR_msb |= 1 << ((size - 2) & 7); |
5847 } else if (flags & LTM_PRIME_2MSB_OFF) { | 5850 } else if (flags & LTM_PRIME_2MSB_OFF) { |
5848 maskAND &= ~(1 << ((size - 2) & 7)); | 5851 maskAND &= ~(1 << ((size - 2) & 7)); |
5849 } | 5852 } |
5850 | 5853 |
5851 /* get the maskOR_lsb */ | 5854 /* get the maskOR_lsb */ |
5852 maskOR_lsb = 0; | 5855 maskOR_lsb = 1; |
5853 if (flags & LTM_PRIME_BBS) { | 5856 if (flags & LTM_PRIME_BBS) { |
5854 maskOR_lsb |= 3; | 5857 maskOR_lsb |= 3; |
5855 } | 5858 } |
5856 | 5859 |
5857 do { | 5860 do { |
5941 /* make sure the radix is in range */ | 5944 /* make sure the radix is in range */ |
5942 if (radix < 2 || radix > 64) { | 5945 if (radix < 2 || radix > 64) { |
5943 return MP_VAL; | 5946 return MP_VAL; |
5944 } | 5947 } |
5945 | 5948 |
5949 if (mp_iszero(a) == MP_YES) { | |
5950 *size = 2; | |
5951 return MP_OKAY; | |
5952 } | |
5953 | |
5954 /* digs is the digit count */ | |
5955 digs = 0; | |
5956 | |
5957 /* if it's negative add one for the sign */ | |
5958 if (a->sign == MP_NEG) { | |
5959 ++digs; | |
5960 } | |
5961 | |
5946 /* init a copy of the input */ | 5962 /* init a copy of the input */ |
5947 if ((res = mp_init_copy (&t, a)) != MP_OKAY) { | 5963 if ((res = mp_init_copy (&t, a)) != MP_OKAY) { |
5948 return res; | 5964 return res; |
5949 } | 5965 } |
5950 | 5966 |
5951 /* digs is the digit count */ | 5967 /* force temp to positive */ |
5952 digs = 0; | 5968 t.sign = MP_ZPOS; |
5953 | |
5954 /* if it's negative add one for the sign */ | |
5955 if (t.sign == MP_NEG) { | |
5956 ++digs; | |
5957 t.sign = MP_ZPOS; | |
5958 } | |
5959 | 5969 |
5960 /* fetch out all of the digits */ | 5970 /* fetch out all of the digits */ |
5961 while (mp_iszero (&t) == 0) { | 5971 while (mp_iszero (&t) == MP_NO) { |
5962 if ((res = mp_div_d (&t, (mp_digit) radix, &t, &d)) != MP_OKAY) { | 5972 if ((res = mp_div_d (&t, (mp_digit) radix, &t, &d)) != MP_OKAY) { |
5963 mp_clear (&t); | 5973 mp_clear (&t); |
5964 return res; | 5974 return res; |
5965 } | 5975 } |
5966 ++digs; | 5976 ++digs; |
6030 return MP_OKAY; | 6040 return MP_OKAY; |
6031 } | 6041 } |
6032 | 6042 |
6033 /* first place a random non-zero digit */ | 6043 /* first place a random non-zero digit */ |
6034 do { | 6044 do { |
6035 d = ((mp_digit) abs (rand ())); | 6045 d = ((mp_digit) abs (rand ())) & MP_MASK; |
6036 } while (d == 0); | 6046 } while (d == 0); |
6037 | 6047 |
6038 if ((res = mp_add_d (a, d, a)) != MP_OKAY) { | 6048 if ((res = mp_add_d (a, d, a)) != MP_OKAY) { |
6039 return res; | 6049 return res; |
6040 } | 6050 } |
6041 | 6051 |
6042 while (digits-- > 0) { | 6052 while (--digits > 0) { |
6043 if ((res = mp_lshd (a, 1)) != MP_OKAY) { | 6053 if ((res = mp_lshd (a, 1)) != MP_OKAY) { |
6044 return res; | 6054 return res; |
6045 } | 6055 } |
6046 | 6056 |
6047 if ((res = mp_add_d (a, ((mp_digit) abs (rand ())), a)) != MP_OKAY) { | 6057 if ((res = mp_add_d (a, ((mp_digit) abs (rand ())), a)) != MP_OKAY) { |
6072 * | 6082 * |
6073 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 6083 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
6074 */ | 6084 */ |
6075 | 6085 |
6076 /* read a string [ASCII] in a given radix */ | 6086 /* read a string [ASCII] in a given radix */ |
6077 int mp_read_radix (mp_int * a, char *str, int radix) | 6087 int mp_read_radix (mp_int * a, const char *str, int radix) |
6078 { | 6088 { |
6079 int y, res, neg; | 6089 int y, res, neg; |
6080 char ch; | 6090 char ch; |
6081 | 6091 |
6082 /* make sure the radix is ok */ | 6092 /* make sure the radix is ok */ |
6255 | 6265 |
6256 /* reduces x mod m, assumes 0 < x < m**2, mu is | 6266 /* reduces x mod m, assumes 0 < x < m**2, mu is |
6257 * precomputed via mp_reduce_setup. | 6267 * precomputed via mp_reduce_setup. |
6258 * From HAC pp.604 Algorithm 14.42 | 6268 * From HAC pp.604 Algorithm 14.42 |
6259 */ | 6269 */ |
6260 int | 6270 int mp_reduce (mp_int * x, mp_int * m, mp_int * mu) |
6261 mp_reduce (mp_int * x, mp_int * m, mp_int * mu) | |
6262 { | 6271 { |
6263 mp_int q; | 6272 mp_int q; |
6264 int res, um = m->used; | 6273 int res, um = m->used; |
6265 | 6274 |
6266 /* q = x */ | 6275 /* q = x */ |
6276 if ((res = mp_mul (&q, mu, &q)) != MP_OKAY) { | 6285 if ((res = mp_mul (&q, mu, &q)) != MP_OKAY) { |
6277 goto CLEANUP; | 6286 goto CLEANUP; |
6278 } | 6287 } |
6279 } else { | 6288 } else { |
6280 #ifdef BN_S_MP_MUL_HIGH_DIGS_C | 6289 #ifdef BN_S_MP_MUL_HIGH_DIGS_C |
6281 if ((res = s_mp_mul_high_digs (&q, mu, &q, um - 1)) != MP_OKAY) { | 6290 if ((res = s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) { |
6282 goto CLEANUP; | 6291 goto CLEANUP; |
6283 } | 6292 } |
6284 #elif defined(BN_FAST_S_MP_MUL_HIGH_DIGS_C) | 6293 #elif defined(BN_FAST_S_MP_MUL_HIGH_DIGS_C) |
6285 if ((res = fast_s_mp_mul_high_digs (&q, mu, &q, um - 1)) != MP_OKAY) { | 6294 if ((res = fast_s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) { |
6286 goto CLEANUP; | 6295 goto CLEANUP; |
6287 } | 6296 } |
6288 #else | 6297 #else |
6289 { | 6298 { |
6290 res = MP_VAL; | 6299 res = MP_VAL; |
6353 * | 6362 * |
6354 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 6363 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
6355 */ | 6364 */ |
6356 | 6365 |
6357 /* reduces a modulo n where n is of the form 2**p - d */ | 6366 /* reduces a modulo n where n is of the form 2**p - d */ |
6358 int | 6367 int mp_reduce_2k(mp_int *a, mp_int *n, mp_digit d) |
6359 mp_reduce_2k(mp_int *a, mp_int *n, mp_digit d) | |
6360 { | 6368 { |
6361 mp_int q; | 6369 mp_int q; |
6362 int p, res; | 6370 int p, res; |
6363 | 6371 |
6364 if ((res = mp_init(&q)) != MP_OKAY) { | 6372 if ((res = mp_init(&q)) != MP_OKAY) { |
6396 | 6404 |
6397 #endif | 6405 #endif |
6398 | 6406 |
6399 /* End: bn_mp_reduce_2k.c */ | 6407 /* End: bn_mp_reduce_2k.c */ |
6400 | 6408 |
6409 /* Start: bn_mp_reduce_2k_l.c */ | |
6410 #include <tommath.h> | |
6411 #ifdef BN_MP_REDUCE_2K_L_C | |
6412 /* LibTomMath, multiple-precision integer library -- Tom St Denis | |
6413 * | |
6414 * LibTomMath is a library that provides multiple-precision | |
6415 * integer arithmetic as well as number theoretic functionality. | |
6416 * | |
6417 * The library was designed directly after the MPI library by | |
6418 * Michael Fromberger but has been written from scratch with | |
6419 * additional optimizations in place. | |
6420 * | |
6421 * The library is free for all purposes without any express | |
6422 * guarantee it works. | |
6423 * | |
6424 * Tom St Denis, [email protected], http://math.libtomcrypt.org | |
6425 */ | |
6426 | |
6427 /* reduces a modulo n where n is of the form 2**p - d | |
6428 This differs from reduce_2k since "d" can be larger | |
6429 than a single digit. | |
6430 */ | |
6431 int mp_reduce_2k_l(mp_int *a, mp_int *n, mp_int *d) | |
6432 { | |
6433 mp_int q; | |
6434 int p, res; | |
6435 | |
6436 if ((res = mp_init(&q)) != MP_OKAY) { | |
6437 return res; | |
6438 } | |
6439 | |
6440 p = mp_count_bits(n); | |
6441 top: | |
6442 /* q = a/2**p, a = a mod 2**p */ | |
6443 if ((res = mp_div_2d(a, p, &q, a)) != MP_OKAY) { | |
6444 goto ERR; | |
6445 } | |
6446 | |
6447 /* q = q * d */ | |
6448 if ((res = mp_mul(&q, d, &q)) != MP_OKAY) { | |
6449 goto ERR; | |
6450 } | |
6451 | |
6452 /* a = a + q */ | |
6453 if ((res = s_mp_add(a, &q, a)) != MP_OKAY) { | |
6454 goto ERR; | |
6455 } | |
6456 | |
6457 if (mp_cmp_mag(a, n) != MP_LT) { | |
6458 s_mp_sub(a, n, a); | |
6459 goto top; | |
6460 } | |
6461 | |
6462 ERR: | |
6463 mp_clear(&q); | |
6464 return res; | |
6465 } | |
6466 | |
6467 #endif | |
6468 | |
6469 /* End: bn_mp_reduce_2k_l.c */ | |
6470 | |
6401 /* Start: bn_mp_reduce_2k_setup.c */ | 6471 /* Start: bn_mp_reduce_2k_setup.c */ |
6402 #include <tommath.h> | 6472 #include <tommath.h> |
6403 #ifdef BN_MP_REDUCE_2K_SETUP_C | 6473 #ifdef BN_MP_REDUCE_2K_SETUP_C |
6404 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 6474 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
6405 * | 6475 * |
6415 * | 6485 * |
6416 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 6486 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
6417 */ | 6487 */ |
6418 | 6488 |
6419 /* determines the setup value */ | 6489 /* determines the setup value */ |
6420 int | 6490 int mp_reduce_2k_setup(mp_int *a, mp_digit *d) |
6421 mp_reduce_2k_setup(mp_int *a, mp_digit *d) | |
6422 { | 6491 { |
6423 int res, p; | 6492 int res, p; |
6424 mp_int tmp; | 6493 mp_int tmp; |
6425 | 6494 |
6426 if ((res = mp_init(&tmp)) != MP_OKAY) { | 6495 if ((res = mp_init(&tmp)) != MP_OKAY) { |
6444 } | 6513 } |
6445 #endif | 6514 #endif |
6446 | 6515 |
6447 /* End: bn_mp_reduce_2k_setup.c */ | 6516 /* End: bn_mp_reduce_2k_setup.c */ |
6448 | 6517 |
6518 /* Start: bn_mp_reduce_2k_setup_l.c */ | |
6519 #include <tommath.h> | |
6520 #ifdef BN_MP_REDUCE_2K_SETUP_L_C | |
6521 /* LibTomMath, multiple-precision integer library -- Tom St Denis | |
6522 * | |
6523 * LibTomMath is a library that provides multiple-precision | |
6524 * integer arithmetic as well as number theoretic functionality. | |
6525 * | |
6526 * The library was designed directly after the MPI library by | |
6527 * Michael Fromberger but has been written from scratch with | |
6528 * additional optimizations in place. | |
6529 * | |
6530 * The library is free for all purposes without any express | |
6531 * guarantee it works. | |
6532 * | |
6533 * Tom St Denis, [email protected], http://math.libtomcrypt.org | |
6534 */ | |
6535 | |
6536 /* determines the setup value */ | |
6537 int mp_reduce_2k_setup_l(mp_int *a, mp_int *d) | |
6538 { | |
6539 int res; | |
6540 mp_int tmp; | |
6541 | |
6542 if ((res = mp_init(&tmp)) != MP_OKAY) { | |
6543 return res; | |
6544 } | |
6545 | |
6546 if ((res = mp_2expt(&tmp, mp_count_bits(a))) != MP_OKAY) { | |
6547 goto ERR; | |
6548 } | |
6549 | |
6550 if ((res = s_mp_sub(&tmp, a, d)) != MP_OKAY) { | |
6551 goto ERR; | |
6552 } | |
6553 | |
6554 ERR: | |
6555 mp_clear(&tmp); | |
6556 return res; | |
6557 } | |
6558 #endif | |
6559 | |
6560 /* End: bn_mp_reduce_2k_setup_l.c */ | |
6561 | |
6449 /* Start: bn_mp_reduce_is_2k.c */ | 6562 /* Start: bn_mp_reduce_is_2k.c */ |
6450 #include <tommath.h> | 6563 #include <tommath.h> |
6451 #ifdef BN_MP_REDUCE_IS_2K_C | 6564 #ifdef BN_MP_REDUCE_IS_2K_C |
6452 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 6565 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
6453 * | 6566 * |
6469 { | 6582 { |
6470 int ix, iy, iw; | 6583 int ix, iy, iw; |
6471 mp_digit iz; | 6584 mp_digit iz; |
6472 | 6585 |
6473 if (a->used == 0) { | 6586 if (a->used == 0) { |
6474 return 0; | 6587 return MP_NO; |
6475 } else if (a->used == 1) { | 6588 } else if (a->used == 1) { |
6476 return 1; | 6589 return MP_YES; |
6477 } else if (a->used > 1) { | 6590 } else if (a->used > 1) { |
6478 iy = mp_count_bits(a); | 6591 iy = mp_count_bits(a); |
6479 iz = 1; | 6592 iz = 1; |
6480 iw = 1; | 6593 iw = 1; |
6481 | 6594 |
6482 /* Test every bit from the second digit up, must be 1 */ | 6595 /* Test every bit from the second digit up, must be 1 */ |
6483 for (ix = DIGIT_BIT; ix < iy; ix++) { | 6596 for (ix = DIGIT_BIT; ix < iy; ix++) { |
6484 if ((a->dp[iw] & iz) == 0) { | 6597 if ((a->dp[iw] & iz) == 0) { |
6485 return 0; | 6598 return MP_NO; |
6486 } | 6599 } |
6487 iz <<= 1; | 6600 iz <<= 1; |
6488 if (iz > (mp_digit)MP_MASK) { | 6601 if (iz > (mp_digit)MP_MASK) { |
6489 ++iw; | 6602 ++iw; |
6490 iz = 1; | 6603 iz = 1; |
6491 } | 6604 } |
6492 } | 6605 } |
6493 } | 6606 } |
6494 return 1; | 6607 return MP_YES; |
6495 } | 6608 } |
6496 | 6609 |
6497 #endif | 6610 #endif |
6498 | 6611 |
6499 /* End: bn_mp_reduce_is_2k.c */ | 6612 /* End: bn_mp_reduce_is_2k.c */ |
6613 | |
6614 /* Start: bn_mp_reduce_is_2k_l.c */ | |
6615 #include <tommath.h> | |
6616 #ifdef BN_MP_REDUCE_IS_2K_L_C | |
6617 /* LibTomMath, multiple-precision integer library -- Tom St Denis | |
6618 * | |
6619 * LibTomMath is a library that provides multiple-precision | |
6620 * integer arithmetic as well as number theoretic functionality. | |
6621 * | |
6622 * The library was designed directly after the MPI library by | |
6623 * Michael Fromberger but has been written from scratch with | |
6624 * additional optimizations in place. | |
6625 * | |
6626 * The library is free for all purposes without any express | |
6627 * guarantee it works. | |
6628 * | |
6629 * Tom St Denis, [email protected], http://math.libtomcrypt.org | |
6630 */ | |
6631 | |
6632 /* determines if reduce_2k_l can be used */ | |
6633 int mp_reduce_is_2k_l(mp_int *a) | |
6634 { | |
6635 int ix, iy; | |
6636 | |
6637 if (a->used == 0) { | |
6638 return MP_NO; | |
6639 } else if (a->used == 1) { | |
6640 return MP_YES; | |
6641 } else if (a->used > 1) { | |
6642 /* if more than half of the digits are -1 we're sold */ | |
6643 for (iy = ix = 0; ix < a->used; ix++) { | |
6644 if (a->dp[ix] == MP_MASK) { | |
6645 ++iy; | |
6646 } | |
6647 } | |
6648 return (iy >= (a->used/2)) ? MP_YES : MP_NO; | |
6649 | |
6650 } | |
6651 return MP_NO; | |
6652 } | |
6653 | |
6654 #endif | |
6655 | |
6656 /* End: bn_mp_reduce_is_2k_l.c */ | |
6500 | 6657 |
6501 /* Start: bn_mp_reduce_setup.c */ | 6658 /* Start: bn_mp_reduce_setup.c */ |
6502 #include <tommath.h> | 6659 #include <tommath.h> |
6503 #ifdef BN_MP_REDUCE_SETUP_C | 6660 #ifdef BN_MP_REDUCE_SETUP_C |
6504 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 6661 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
7130 * | 7287 * |
7131 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 7288 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7132 */ | 7289 */ |
7133 | 7290 |
7134 /* store in signed [big endian] format */ | 7291 /* store in signed [big endian] format */ |
7135 int | 7292 int mp_to_signed_bin (mp_int * a, unsigned char *b) |
7136 mp_to_signed_bin (mp_int * a, unsigned char *b) | |
7137 { | 7293 { |
7138 int res; | 7294 int res; |
7139 | 7295 |
7140 if ((res = mp_to_unsigned_bin (a, b + 1)) != MP_OKAY) { | 7296 if ((res = mp_to_unsigned_bin (a, b + 1)) != MP_OKAY) { |
7141 return res; | 7297 return res; |
7145 } | 7301 } |
7146 #endif | 7302 #endif |
7147 | 7303 |
7148 /* End: bn_mp_to_signed_bin.c */ | 7304 /* End: bn_mp_to_signed_bin.c */ |
7149 | 7305 |
7306 /* Start: bn_mp_to_signed_bin_n.c */ | |
7307 #include <tommath.h> | |
7308 #ifdef BN_MP_TO_SIGNED_BIN_N_C | |
7309 /* LibTomMath, multiple-precision integer library -- Tom St Denis | |
7310 * | |
7311 * LibTomMath is a library that provides multiple-precision | |
7312 * integer arithmetic as well as number theoretic functionality. | |
7313 * | |
7314 * The library was designed directly after the MPI library by | |
7315 * Michael Fromberger but has been written from scratch with | |
7316 * additional optimizations in place. | |
7317 * | |
7318 * The library is free for all purposes without any express | |
7319 * guarantee it works. | |
7320 * | |
7321 * Tom St Denis, [email protected], http://math.libtomcrypt.org | |
7322 */ | |
7323 | |
7324 /* store in signed [big endian] format */ | |
7325 int mp_to_signed_bin_n (mp_int * a, unsigned char *b, unsigned long *outlen) | |
7326 { | |
7327 if (*outlen < (unsigned long)mp_signed_bin_size(a)) { | |
7328 return MP_VAL; | |
7329 } | |
7330 *outlen = mp_signed_bin_size(a); | |
7331 return mp_to_signed_bin(a, b); | |
7332 } | |
7333 #endif | |
7334 | |
7335 /* End: bn_mp_to_signed_bin_n.c */ | |
7336 | |
7150 /* Start: bn_mp_to_unsigned_bin.c */ | 7337 /* Start: bn_mp_to_unsigned_bin.c */ |
7151 #include <tommath.h> | 7338 #include <tommath.h> |
7152 #ifdef BN_MP_TO_UNSIGNED_BIN_C | 7339 #ifdef BN_MP_TO_UNSIGNED_BIN_C |
7153 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 7340 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
7154 * | 7341 * |
7164 * | 7351 * |
7165 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 7352 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7166 */ | 7353 */ |
7167 | 7354 |
7168 /* store in unsigned [big endian] format */ | 7355 /* store in unsigned [big endian] format */ |
7169 int | 7356 int mp_to_unsigned_bin (mp_int * a, unsigned char *b) |
7170 mp_to_unsigned_bin (mp_int * a, unsigned char *b) | |
7171 { | 7357 { |
7172 int x, res; | 7358 int x, res; |
7173 mp_int t; | 7359 mp_int t; |
7174 | 7360 |
7175 if ((res = mp_init_copy (&t, a)) != MP_OKAY) { | 7361 if ((res = mp_init_copy (&t, a)) != MP_OKAY) { |
7194 } | 7380 } |
7195 #endif | 7381 #endif |
7196 | 7382 |
7197 /* End: bn_mp_to_unsigned_bin.c */ | 7383 /* End: bn_mp_to_unsigned_bin.c */ |
7198 | 7384 |
7385 /* Start: bn_mp_to_unsigned_bin_n.c */ | |
7386 #include <tommath.h> | |
7387 #ifdef BN_MP_TO_UNSIGNED_BIN_N_C | |
7388 /* LibTomMath, multiple-precision integer library -- Tom St Denis | |
7389 * | |
7390 * LibTomMath is a library that provides multiple-precision | |
7391 * integer arithmetic as well as number theoretic functionality. | |
7392 * | |
7393 * The library was designed directly after the MPI library by | |
7394 * Michael Fromberger but has been written from scratch with | |
7395 * additional optimizations in place. | |
7396 * | |
7397 * The library is free for all purposes without any express | |
7398 * guarantee it works. | |
7399 * | |
7400 * Tom St Denis, [email protected], http://math.libtomcrypt.org | |
7401 */ | |
7402 | |
7403 /* store in unsigned [big endian] format */ | |
7404 int mp_to_unsigned_bin_n (mp_int * a, unsigned char *b, unsigned long *outlen) | |
7405 { | |
7406 if (*outlen < (unsigned long)mp_unsigned_bin_size(a)) { | |
7407 return MP_VAL; | |
7408 } | |
7409 *outlen = mp_unsigned_bin_size(a); | |
7410 return mp_to_unsigned_bin(a, b); | |
7411 } | |
7412 #endif | |
7413 | |
7414 /* End: bn_mp_to_unsigned_bin_n.c */ | |
7415 | |
7199 /* Start: bn_mp_toom_mul.c */ | 7416 /* Start: bn_mp_toom_mul.c */ |
7200 #include <tommath.h> | 7417 #include <tommath.h> |
7201 #ifdef BN_MP_TOOM_MUL_C | 7418 #ifdef BN_MP_TOOM_MUL_C |
7202 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 7419 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
7203 * | 7420 * |
7214 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 7431 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7215 */ | 7432 */ |
7216 | 7433 |
7217 /* multiplication using the Toom-Cook 3-way algorithm | 7434 /* multiplication using the Toom-Cook 3-way algorithm |
7218 * | 7435 * |
7219 * Much more complicated than Karatsuba but has a lower asymptotic running time of | 7436 * Much more complicated than Karatsuba but has a lower |
7220 * O(N**1.464). This algorithm is only particularly useful on VERY large | 7437 * asymptotic running time of O(N**1.464). This algorithm is |
7221 * inputs (we're talking 1000s of digits here...). | 7438 * only particularly useful on VERY large inputs |
7439 * (we're talking 1000s of digits here...). | |
7222 */ | 7440 */ |
7223 int mp_toom_mul(mp_int *a, mp_int *b, mp_int *c) | 7441 int mp_toom_mul(mp_int *a, mp_int *b, mp_int *c) |
7224 { | 7442 { |
7225 mp_int w0, w1, w2, w3, w4, tmp1, tmp2, a0, a1, a2, b0, b1, b2; | 7443 mp_int w0, w1, w2, w3, w4, tmp1, tmp2, a0, a1, a2, b0, b1, b2; |
7226 int res, B; | 7444 int res, B; |
7886 * | 8104 * |
7887 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 8105 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7888 */ | 8106 */ |
7889 | 8107 |
7890 /* get the size for an unsigned equivalent */ | 8108 /* get the size for an unsigned equivalent */ |
7891 int | 8109 int mp_unsigned_bin_size (mp_int * a) |
7892 mp_unsigned_bin_size (mp_int * a) | |
7893 { | 8110 { |
7894 int size = mp_count_bits (a); | 8111 int size = mp_count_bits (a); |
7895 return (size / 8 + ((size & 7) != 0 ? 1 : 0)); | 8112 return (size / 8 + ((size & 7) != 0 ? 1 : 0)); |
7896 } | 8113 } |
7897 #endif | 8114 #endif |
7936 px = a->used; | 8153 px = a->used; |
7937 x = a; | 8154 x = a; |
7938 } | 8155 } |
7939 | 8156 |
7940 for (ix = 0; ix < px; ix++) { | 8157 for (ix = 0; ix < px; ix++) { |
7941 | 8158 t.dp[ix] ^= x->dp[ix]; |
7942 } | 8159 } |
7943 mp_clamp (&t); | 8160 mp_clamp (&t); |
7944 mp_exch (c, &t); | 8161 mp_exch (c, &t); |
7945 mp_clear (&t); | 8162 mp_clear (&t); |
7946 return MP_OKAY; | 8163 return MP_OKAY; |
7966 * | 8183 * |
7967 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 8184 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7968 */ | 8185 */ |
7969 | 8186 |
7970 /* set to zero */ | 8187 /* set to zero */ |
7971 void | 8188 void mp_zero (mp_int * a) |
7972 mp_zero (mp_int * a) | 8189 { |
7973 { | 8190 int n; |
8191 mp_digit *tmp; | |
8192 | |
7974 a->sign = MP_ZPOS; | 8193 a->sign = MP_ZPOS; |
7975 a->used = 0; | 8194 a->used = 0; |
7976 memset (a->dp, 0, sizeof (mp_digit) * a->alloc); | 8195 |
8196 tmp = a->dp; | |
8197 for (n = 0; n < a->alloc; n++) { | |
8198 *tmp++ = 0; | |
8199 } | |
7977 } | 8200 } |
7978 #endif | 8201 #endif |
7979 | 8202 |
7980 /* End: bn_mp_zero.c */ | 8203 /* End: bn_mp_zero.c */ |
7981 | 8204 |
7994 * The library is free for all purposes without any express | 8217 * The library is free for all purposes without any express |
7995 * guarantee it works. | 8218 * guarantee it works. |
7996 * | 8219 * |
7997 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 8220 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
7998 */ | 8221 */ |
7999 const mp_digit __prime_tab[] = { | 8222 const mp_digit ltm_prime_tab[] = { |
8000 0x0002, 0x0003, 0x0005, 0x0007, 0x000B, 0x000D, 0x0011, 0x0013, | 8223 0x0002, 0x0003, 0x0005, 0x0007, 0x000B, 0x000D, 0x0011, 0x0013, |
8001 0x0017, 0x001D, 0x001F, 0x0025, 0x0029, 0x002B, 0x002F, 0x0035, | 8224 0x0017, 0x001D, 0x001F, 0x0025, 0x0029, 0x002B, 0x002F, 0x0035, |
8002 0x003B, 0x003D, 0x0043, 0x0047, 0x0049, 0x004F, 0x0053, 0x0059, | 8225 0x003B, 0x003D, 0x0043, 0x0047, 0x0049, 0x004F, 0x0053, 0x0059, |
8003 0x0061, 0x0065, 0x0067, 0x006B, 0x006D, 0x0071, 0x007F, | 8226 0x0061, 0x0065, 0x0067, 0x006B, 0x006D, 0x0071, 0x007F, |
8004 #ifndef MP_8BIT | 8227 #ifndef MP_8BIT |
8210 #define TAB_SIZE 32 | 8433 #define TAB_SIZE 32 |
8211 #else | 8434 #else |
8212 #define TAB_SIZE 256 | 8435 #define TAB_SIZE 256 |
8213 #endif | 8436 #endif |
8214 | 8437 |
8215 int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y) | 8438 int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode) |
8216 { | 8439 { |
8217 mp_int M[TAB_SIZE], res, mu; | 8440 mp_int M[TAB_SIZE], res, mu; |
8218 mp_digit buf; | 8441 mp_digit buf; |
8219 int err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize; | 8442 int err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize; |
8443 int (*redux)(mp_int*,mp_int*,mp_int*); | |
8220 | 8444 |
8221 /* find window size */ | 8445 /* find window size */ |
8222 x = mp_count_bits (X); | 8446 x = mp_count_bits (X); |
8223 if (x <= 7) { | 8447 if (x <= 7) { |
8224 winsize = 2; | 8448 winsize = 2; |
8259 } | 8483 } |
8260 } | 8484 } |
8261 | 8485 |
8262 /* create mu, used for Barrett reduction */ | 8486 /* create mu, used for Barrett reduction */ |
8263 if ((err = mp_init (&mu)) != MP_OKAY) { | 8487 if ((err = mp_init (&mu)) != MP_OKAY) { |
8264 goto __M; | 8488 goto LBL_M; |
8265 } | 8489 } |
8266 if ((err = mp_reduce_setup (&mu, P)) != MP_OKAY) { | 8490 |
8267 goto __MU; | 8491 if (redmode == 0) { |
8268 } | 8492 if ((err = mp_reduce_setup (&mu, P)) != MP_OKAY) { |
8493 goto LBL_MU; | |
8494 } | |
8495 redux = mp_reduce; | |
8496 } else { | |
8497 if ((err = mp_reduce_2k_setup_l (P, &mu)) != MP_OKAY) { | |
8498 goto LBL_MU; | |
8499 } | |
8500 redux = mp_reduce_2k_l; | |
8501 } | |
8269 | 8502 |
8270 /* create M table | 8503 /* create M table |
8271 * | 8504 * |
8272 * The M table contains powers of the base, | 8505 * The M table contains powers of the base, |
8273 * e.g. M[x] = G**x mod P | 8506 * e.g. M[x] = G**x mod P |
8274 * | 8507 * |
8275 * The first half of the table is not | 8508 * The first half of the table is not |
8276 * computed though accept for M[0] and M[1] | 8509 * computed though accept for M[0] and M[1] |
8277 */ | 8510 */ |
8278 if ((err = mp_mod (G, P, &M[1])) != MP_OKAY) { | 8511 if ((err = mp_mod (G, P, &M[1])) != MP_OKAY) { |
8279 goto __MU; | 8512 goto LBL_MU; |
8280 } | 8513 } |
8281 | 8514 |
8282 /* compute the value at M[1<<(winsize-1)] by squaring | 8515 /* compute the value at M[1<<(winsize-1)] by squaring |
8283 * M[1] (winsize-1) times | 8516 * M[1] (winsize-1) times |
8284 */ | 8517 */ |
8285 if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) { | 8518 if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) { |
8286 goto __MU; | 8519 goto LBL_MU; |
8287 } | 8520 } |
8288 | 8521 |
8289 for (x = 0; x < (winsize - 1); x++) { | 8522 for (x = 0; x < (winsize - 1); x++) { |
8523 /* square it */ | |
8290 if ((err = mp_sqr (&M[1 << (winsize - 1)], | 8524 if ((err = mp_sqr (&M[1 << (winsize - 1)], |
8291 &M[1 << (winsize - 1)])) != MP_OKAY) { | 8525 &M[1 << (winsize - 1)])) != MP_OKAY) { |
8292 goto __MU; | 8526 goto LBL_MU; |
8293 } | 8527 } |
8294 if ((err = mp_reduce (&M[1 << (winsize - 1)], P, &mu)) != MP_OKAY) { | 8528 |
8295 goto __MU; | 8529 /* reduce modulo P */ |
8530 if ((err = redux (&M[1 << (winsize - 1)], P, &mu)) != MP_OKAY) { | |
8531 goto LBL_MU; | |
8296 } | 8532 } |
8297 } | 8533 } |
8298 | 8534 |
8299 /* create upper table, that is M[x] = M[x-1] * M[1] (mod P) | 8535 /* create upper table, that is M[x] = M[x-1] * M[1] (mod P) |
8300 * for x = (2**(winsize - 1) + 1) to (2**winsize - 1) | 8536 * for x = (2**(winsize - 1) + 1) to (2**winsize - 1) |
8301 */ | 8537 */ |
8302 for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) { | 8538 for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) { |
8303 if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) { | 8539 if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) { |
8304 goto __MU; | 8540 goto LBL_MU; |
8305 } | 8541 } |
8306 if ((err = mp_reduce (&M[x], P, &mu)) != MP_OKAY) { | 8542 if ((err = redux (&M[x], P, &mu)) != MP_OKAY) { |
8307 goto __MU; | 8543 goto LBL_MU; |
8308 } | 8544 } |
8309 } | 8545 } |
8310 | 8546 |
8311 /* setup result */ | 8547 /* setup result */ |
8312 if ((err = mp_init (&res)) != MP_OKAY) { | 8548 if ((err = mp_init (&res)) != MP_OKAY) { |
8313 goto __MU; | 8549 goto LBL_MU; |
8314 } | 8550 } |
8315 mp_set (&res, 1); | 8551 mp_set (&res, 1); |
8316 | 8552 |
8317 /* set initial mode and bit cnt */ | 8553 /* set initial mode and bit cnt */ |
8318 mode = 0; | 8554 mode = 0; |
8348 } | 8584 } |
8349 | 8585 |
8350 /* if the bit is zero and mode == 1 then we square */ | 8586 /* if the bit is zero and mode == 1 then we square */ |
8351 if (mode == 1 && y == 0) { | 8587 if (mode == 1 && y == 0) { |
8352 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 8588 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
8353 goto __RES; | 8589 goto LBL_RES; |
8354 } | 8590 } |
8355 if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) { | 8591 if ((err = redux (&res, P, &mu)) != MP_OKAY) { |
8356 goto __RES; | 8592 goto LBL_RES; |
8357 } | 8593 } |
8358 continue; | 8594 continue; |
8359 } | 8595 } |
8360 | 8596 |
8361 /* else we add it to the window */ | 8597 /* else we add it to the window */ |
8365 if (bitcpy == winsize) { | 8601 if (bitcpy == winsize) { |
8366 /* ok window is filled so square as required and multiply */ | 8602 /* ok window is filled so square as required and multiply */ |
8367 /* square first */ | 8603 /* square first */ |
8368 for (x = 0; x < winsize; x++) { | 8604 for (x = 0; x < winsize; x++) { |
8369 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 8605 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
8370 goto __RES; | 8606 goto LBL_RES; |
8371 } | 8607 } |
8372 if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) { | 8608 if ((err = redux (&res, P, &mu)) != MP_OKAY) { |
8373 goto __RES; | 8609 goto LBL_RES; |
8374 } | 8610 } |
8375 } | 8611 } |
8376 | 8612 |
8377 /* then multiply */ | 8613 /* then multiply */ |
8378 if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) { | 8614 if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) { |
8379 goto __RES; | 8615 goto LBL_RES; |
8380 } | 8616 } |
8381 if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) { | 8617 if ((err = redux (&res, P, &mu)) != MP_OKAY) { |
8382 goto __RES; | 8618 goto LBL_RES; |
8383 } | 8619 } |
8384 | 8620 |
8385 /* empty window and reset */ | 8621 /* empty window and reset */ |
8386 bitcpy = 0; | 8622 bitcpy = 0; |
8387 bitbuf = 0; | 8623 bitbuf = 0; |
8392 /* if bits remain then square/multiply */ | 8628 /* if bits remain then square/multiply */ |
8393 if (mode == 2 && bitcpy > 0) { | 8629 if (mode == 2 && bitcpy > 0) { |
8394 /* square then multiply if the bit is set */ | 8630 /* square then multiply if the bit is set */ |
8395 for (x = 0; x < bitcpy; x++) { | 8631 for (x = 0; x < bitcpy; x++) { |
8396 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { | 8632 if ((err = mp_sqr (&res, &res)) != MP_OKAY) { |
8397 goto __RES; | 8633 goto LBL_RES; |
8398 } | 8634 } |
8399 if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) { | 8635 if ((err = redux (&res, P, &mu)) != MP_OKAY) { |
8400 goto __RES; | 8636 goto LBL_RES; |
8401 } | 8637 } |
8402 | 8638 |
8403 bitbuf <<= 1; | 8639 bitbuf <<= 1; |
8404 if ((bitbuf & (1 << winsize)) != 0) { | 8640 if ((bitbuf & (1 << winsize)) != 0) { |
8405 /* then multiply */ | 8641 /* then multiply */ |
8406 if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) { | 8642 if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) { |
8407 goto __RES; | 8643 goto LBL_RES; |
8408 } | 8644 } |
8409 if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) { | 8645 if ((err = redux (&res, P, &mu)) != MP_OKAY) { |
8410 goto __RES; | 8646 goto LBL_RES; |
8411 } | 8647 } |
8412 } | 8648 } |
8413 } | 8649 } |
8414 } | 8650 } |
8415 | 8651 |
8416 mp_exch (&res, Y); | 8652 mp_exch (&res, Y); |
8417 err = MP_OKAY; | 8653 err = MP_OKAY; |
8418 __RES:mp_clear (&res); | 8654 LBL_RES:mp_clear (&res); |
8419 __MU:mp_clear (&mu); | 8655 LBL_MU:mp_clear (&mu); |
8420 __M: | 8656 LBL_M: |
8421 mp_clear(&M[1]); | 8657 mp_clear(&M[1]); |
8422 for (x = 1<<(winsize-1); x < (1 << winsize); x++) { | 8658 for (x = 1<<(winsize-1); x < (1 << winsize); x++) { |
8423 mp_clear (&M[x]); | 8659 mp_clear (&M[x]); |
8424 } | 8660 } |
8425 return err; | 8661 return err; |
8448 | 8684 |
8449 /* multiplies |a| * |b| and only computes upto digs digits of result | 8685 /* multiplies |a| * |b| and only computes upto digs digits of result |
8450 * HAC pp. 595, Algorithm 14.12 Modified so you can control how | 8686 * HAC pp. 595, Algorithm 14.12 Modified so you can control how |
8451 * many digits of output are created. | 8687 * many digits of output are created. |
8452 */ | 8688 */ |
8453 int | 8689 int s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs) |
8454 s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs) | |
8455 { | 8690 { |
8456 mp_int t; | 8691 mp_int t; |
8457 int res, pa, pb, ix, iy; | 8692 int res, pa, pb, ix, iy; |
8458 mp_digit u; | 8693 mp_digit u; |
8459 mp_word r; | 8694 mp_word r; |
8617 * | 8852 * |
8618 * Tom St Denis, [email protected], http://math.libtomcrypt.org | 8853 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
8619 */ | 8854 */ |
8620 | 8855 |
8621 /* low level squaring, b = a*a, HAC pp.596-597, Algorithm 14.16 */ | 8856 /* low level squaring, b = a*a, HAC pp.596-597, Algorithm 14.16 */ |
8622 int | 8857 int s_mp_sqr (mp_int * a, mp_int * b) |
8623 s_mp_sqr (mp_int * a, mp_int * b) | |
8624 { | 8858 { |
8625 mp_int t; | 8859 mp_int t; |
8626 int res, ix, iy, pa; | 8860 int res, ix, iy, pa; |
8627 mp_word r; | 8861 mp_word r; |
8628 mp_digit u, tmpx, *tmpt; | 8862 mp_digit u, tmpx, *tmpt; |
8795 /* Known optimal configurations | 9029 /* Known optimal configurations |
8796 | 9030 |
8797 CPU /Compiler /MUL CUTOFF/SQR CUTOFF | 9031 CPU /Compiler /MUL CUTOFF/SQR CUTOFF |
8798 ------------------------------------------------------------- | 9032 ------------------------------------------------------------- |
8799 Intel P4 Northwood /GCC v3.4.1 / 88/ 128/LTM 0.32 ;-) | 9033 Intel P4 Northwood /GCC v3.4.1 / 88/ 128/LTM 0.32 ;-) |
9034 AMD Athlon64 /GCC v3.4.4 / 74/ 124/LTM 0.34 | |
8800 | 9035 |
8801 */ | 9036 */ |
8802 | 9037 |
8803 int KARATSUBA_MUL_CUTOFF = 88, /* Min. number of digits before Karatsuba multiplication is used. */ | 9038 int KARATSUBA_MUL_CUTOFF = 74, /* Min. number of digits before Karatsuba multiplication is used. */ |
8804 KARATSUBA_SQR_CUTOFF = 128, /* Min. number of digits before Karatsuba squaring is used. */ | 9039 KARATSUBA_SQR_CUTOFF = 124, /* Min. number of digits before Karatsuba squaring is used. */ |
8805 | 9040 |
8806 TOOM_MUL_CUTOFF = 350, /* no optimal values of these are known yet so set em high */ | 9041 TOOM_MUL_CUTOFF = 350, /* no optimal values of these are known yet so set em high */ |
8807 TOOM_SQR_CUTOFF = 400; | 9042 TOOM_SQR_CUTOFF = 400; |
8808 #endif | 9043 #endif |
8809 | 9044 |