3 This file is part of the ARM-Crypto-Lib.
4 Copyright (C) 2008 Daniel Otte (daniel.otte@rub.de)
6 This program is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
16 You should have received a copy of the GNU General Public License
17 along with this program. If not, see <http://www.gnu.org/licenses/>.
24 * \license GPLv3 or later
30 #define STRING(x) STRING2(x)
31 #define STR_LINE STRING(__LINE__)
37 #define PREFERE_HEAP_SPACE 1
39 #if PREFERE_HEAP_SPACE
42 #define ALLOC_BIGINT_WORDS(var,words) \
43 bigint_word_t *(var) = malloc((words) * sizeof(bigint_word_t)); \
45 puts_P(PSTR("\n\nDBG: OOM ERROR (in arithmeics)!\n")); \
51 #define FREE(x) free(x)
55 #define ALLOC_BIGINT_WORDS(var,words) bigint_word_t var[words]
65 #include "bigint_io.h"
69 #define MAX(a,b) (((a)>(b))?(a):(b))
73 #define MIN(a,b) (((a) < (b)) ? (a) : (b))
76 #define SET_FBS(a, v) do {(a)->info &= ~BIGINT_FBS_MASK; (a)->info |= (v);} while(0)
77 #define GET_FBS(a) ((a)->info & BIGINT_FBS_MASK)
78 #define SET_NEG(a) (a)->info |= BIGINT_NEG_MASK
79 #define SET_POS(a) (a)->info &= ~BIGINT_NEG_MASK
80 #define XCHG(a,b) do{(a) ^= (b); (b) ^= (a); (a) ^= (b);} while(0)
81 #define XCHG_PTR(a,b) do{ a = (void*)(((intptr_t)(a)) ^ ((intptr_t)(b))); \
82 b = (void*)(((intptr_t)(a)) ^ ((intptr_t)(b))); \
83 a = (void*)(((intptr_t)(a)) ^ ((intptr_t)(b)));} while(0)
85 #define GET_SIGN(a) ((a)->info & BIGINT_NEG_MASK)
87 /******************************************************************************/
88 void bigint_adjust(bigint_t *a){
89 while (a->length_W != 0 && a->wordv[a->length_W - 1] == 0) {
92 if (a->length_W == 0) {
97 uint8_t i = BIGINT_WORD_SIZE - 1;
98 t = a->wordv[a->length_W - 1];
99 while ((t & (((bigint_word_t)1) << (BIGINT_WORD_SIZE - 1))) == 0 && i) {
106 /******************************************************************************/
108 bigint_length_t bigint_length_b(const bigint_t *a){
109 if(!a->length_W || a->length_W==0){
112 return (a->length_W-1) * BIGINT_WORD_SIZE + GET_FBS(a);
115 /******************************************************************************/
117 bigint_length_t bigint_length_B(const bigint_t *a){
118 return a->length_W * sizeof(bigint_word_t);
121 /******************************************************************************/
123 int32_t bigint_get_first_set_bit(const bigint_t *a){
124 if(a->length_W == 0) {
127 return (a->length_W-1) * sizeof(bigint_word_t) * CHAR_BIT + GET_FBS(a);
131 /******************************************************************************/
133 int32_t bigint_get_last_set_bit(const bigint_t *a){
137 if (a->length_W == 0) {
140 while (a->wordv[r] == 0 && r < a->length_W) {
143 if (a->wordv[r] == 0) {
144 return (uint32_t)(-1);
146 while ((x&a->wordv[r])==0) {
150 return r * BIGINT_WORD_SIZE + b;
153 /******************************************************************************/
155 void bigint_copy(bigint_t *dest, const bigint_t *src){
156 if(dest->wordv != src->wordv){
157 memcpy(dest->wordv, src->wordv, src->length_W * sizeof(bigint_word_t));
159 dest->length_W = src->length_W;
160 dest->info = src->info;
163 /******************************************************************************/
165 /* this should be implemented in assembly */
166 void bigint_add_u(bigint_t *dest, const bigint_t *a, const bigint_t *b){
168 bigint_wordplus_t t = 0LL;
169 if (a->length_W < b->length_W) {
172 for(i = 0; i < b->length_W; ++i) {
175 dest->wordv[i] = (bigint_word_t)t;
176 t >>= BIGINT_WORD_SIZE;
178 for(; i < a->length_W; ++i){
180 dest->wordv[i] = (bigint_word_t)t;
181 t >>= BIGINT_WORD_SIZE;
184 dest->wordv[i++] = (bigint_word_t)t;
191 /******************************************************************************/
193 /* this should be implemented in assembly */
194 void bigint_add_scale_u(bigint_t *dest, const bigint_t *a, bigint_length_t scale){
195 if(a->length_W == 0){
199 bigint_add_u(dest, dest, a);
203 #if BIGINT_WORD_SIZE == 8
204 memset(dest->wordv + dest->length_W, 0, MAX(dest->length_W, a->length_W + scale) - dest->length_W);
205 x.wordv = dest->wordv + scale;
206 x.length_W = dest->length_W - scale;
207 if((int16_t)x.length_W < 0){
213 bigint_add_u(&x, &x, a);
214 dest->length_W = x.length_W + scale;
219 bigint_length_t word_shift = scale / sizeof(bigint_word_t), byte_shift = scale % sizeof(bigint_word_t);
220 bigint_word_t bv[a->length_W + 1];
222 bv[0] = bv[a->length_W] = 0;
223 memcpy((uint8_t*)bv + byte_shift, a->wordv, a->length_W * sizeof(bigint_word_t));
224 s.length_W = a->length_W + 1;
226 memset(dest->wordv + dest->length_W, 0, (MAX(dest->length_W, s.length_W + word_shift) - dest->length_W) * sizeof(bigint_word_t));
227 x.wordv = dest->wordv + word_shift;
228 x.length_W = dest->length_W - word_shift;
229 if((int16_t)x.length_W < 0){
235 bigint_add_u(&x, &x, &s);
236 dest->length_W = x.length_W + word_shift;
242 /******************************************************************************/
244 /* this should be implemented in assembly */
245 void bigint_sub_u(bigint_t *dest, const bigint_t *a, const bigint_t *b){
248 bigint_wordplus_signed_t t = 0;
250 if(b->length_W == 0){
251 bigint_copy(dest, a);
255 if(a->length_W == 0){
256 bigint_copy(dest, b);
260 r = bigint_cmp_u(a,b);
262 bigint_set_zero(dest);
266 bigint_sub_u(dest, b, a);
270 for(i = 0; i < a->length_W; ++i){
276 dest->wordv[i] = (bigint_word_t)t;
277 borrow = t < 0 ? 1 : 0;
284 /******************************************************************************/
286 int8_t bigint_cmp_u(const bigint_t *a, const bigint_t *b){
287 if(a->length_W > b->length_W){
290 if(a->length_W < b->length_W){
293 if(a->length_W == 0){
299 if(a->wordv[i] != b->wordv[i]){
300 if(a->wordv[i] > b->wordv[i]){
310 /******************************************************************************/
312 void bigint_add_s(bigint_t *dest, const bigint_t *a, const bigint_t *b){
315 s |= GET_SIGN(b)?1:0;
317 case 0: /* both positive */
318 bigint_add_u(dest, a,b);
321 case 1: /* a positive, b negative */
322 bigint_sub_u(dest, a, b);
324 case 2: /* a negative, b positive */
325 bigint_sub_u(dest, b, a);
327 case 3: /* both negative */
328 bigint_add_u(dest, a, b);
331 default: /* how can this happen?*/
336 /******************************************************************************/
338 void bigint_sub_s(bigint_t *dest, const bigint_t *a, const bigint_t *b){
341 s |= GET_SIGN(b)?1:0;
343 case 0: /* both positive */
344 bigint_sub_u(dest, a,b);
346 case 1: /* a positive, b negative */
347 bigint_add_u(dest, a, b);
350 case 2: /* a negative, b positive */
351 bigint_add_u(dest, a, b);
354 case 3: /* both negative */
355 bigint_sub_u(dest, b, a);
357 default: /* how can this happen?*/
363 /******************************************************************************/
365 int8_t bigint_cmp_s(const bigint_t *a, const bigint_t *b){
367 if(a->length_W==0 && b->length_W==0){
371 s |= GET_SIGN(b)?1:0;
373 case 0: /* both positive */
374 return bigint_cmp_u(a, b);
376 case 1: /* a positive, b negative */
379 case 2: /* a negative, b positive */
382 case 3: /* both negative */
383 return bigint_cmp_u(b, a);
385 default: /* how can this happen?*/
388 return 0; /* just to satisfy the compiler */
391 /******************************************************************************/
393 void bigint_shiftleft_bits(bigint_t *a, uint8_t shift) {
395 bigint_wordplus_t t = 0;
397 for (i = 0; i < a->length_W; ++i) {
398 t |= ((bigint_wordplus_t)a->wordv[i]) << shift;
399 a->wordv[i] = (bigint_word_t)t;
400 t >>= BIGINT_WORD_SIZE;
403 a->wordv[a->length_W++] = (bigint_word_t)t;
408 /******************************************************************************/
410 void bigint_shiftleft(bigint_t *a, bigint_length_t shift){
411 bigint_length_t byteshift;
414 if (a->length_W == 0 || shift == 0) {
417 byteshift = shift / 8;
418 bitshift = shift & 7;
419 if (byteshift % sizeof(bigint_word_t)) {
420 a->wordv[a->length_W + byteshift / sizeof(bigint_t)] = 0;
423 memmove(((uint8_t*)a->wordv) + byteshift, a->wordv, a->length_W * sizeof(bigint_word_t));
424 memset(a->wordv, 0, byteshift);
425 a->length_W += (byteshift + sizeof(bigint_word_t) - 1) / sizeof(bigint_word_t);
431 bigint_shiftleft_bits(a, bitshift);
435 /******************************************************************************/
437 void bigint_shiftright_bits(bigint_t *a, uint8_t shift){
439 bigint_wordplus_t t = 0;
443 t |= ((bigint_wordplus_t)(a->wordv[i])) << (BIGINT_WORD_SIZE - shift);
444 a->wordv[i] = (bigint_word_t)(t >> BIGINT_WORD_SIZE);
445 t <<= BIGINT_WORD_SIZE;
451 /******************************************************************************/
453 void bigint_shiftright_1bit(bigint_t *a){
455 bigint_word_t t1 = 0, t2;
459 t2 = a->wordv[i] & 1;
462 t1 = t2 << (BIGINT_WORD_SIZE - 1);
467 /******************************************************************************/
469 void bigint_shiftright_1word(bigint_t *a){
470 if (a->length_W == 0) {
473 memmove(a->wordv, &a->wordv[1], (a->length_W - 1) * sizeof(bigint_word_t));
477 /******************************************************************************/
479 void bigint_shiftright(bigint_t *a, bigint_length_t shift){
480 bigint_length_t byteshift = shift / 8;
481 uint8_t bitshift = shift & 7;
483 if (a->length_W == 0) {
487 if(bigint_get_first_set_bit(a) < shift){ /* we would shift out more than we have */
493 memmove(a->wordv, (uint8_t*)a->wordv + byteshift, a->length_W * sizeof(bigint_word_t) - byteshift);
494 memset((uint8_t*)&a->wordv[a->length_W] - byteshift, 0, byteshift);
495 a->length_W -= byteshift / sizeof(bigint_word_t);
500 bigint_shiftright_bits(a, bitshift);
504 /******************************************************************************/
506 void bigint_xor(bigint_t *dest, const bigint_t *a){
508 for(i=0; i<a->length_W; ++i){
509 dest->wordv[i] ^= a->wordv[i];
514 /******************************************************************************/
516 void bigint_set_zero(bigint_t *a){
520 /******************************************************************************/
522 /* using the Karatsuba-Algorithm */
523 /* x*y = (xh*yh)*b**2n + ((xh+xl)*(yh+yl) - xh*yh - xl*yl)*b**n + yh*yl */
524 void bigint_mul_u(bigint_t *dest, const bigint_t *a, const bigint_t *b){
525 if (a->length_W == 0 || b->length_W == 0) {
526 bigint_set_zero(dest);
529 if (dest == a || dest == b) {
531 bigint_word_t d_b[a->length_W + b->length_W];
533 bigint_mul_u(&d, a, b);
534 bigint_copy(dest, &d);
537 if (a->length_W == 1 || b->length_W == 1) {
538 if (a->length_W != 1) {
541 bigint_wordplus_t t = 0;
543 bigint_word_t x = a->wordv[0];
544 for (i = 0; i < b->length_W; ++i) {
545 t += ((bigint_wordplus_t)b->wordv[i]) * ((bigint_wordplus_t)x);
546 dest->wordv[i] = (bigint_word_t)t;
547 t >>= BIGINT_WORD_SIZE;
551 dest->wordv[i] = (bigint_word_t)t;
558 if (a->length_W * sizeof(bigint_word_t) <= 4 && b->length_W * sizeof(bigint_word_t) <= 4) {
559 uint32_t p = 0, q = 0;
561 memcpy(&p, a->wordv, a->length_W * sizeof(bigint_word_t));
562 memcpy(&q, b->wordv, b->length_W * sizeof(bigint_word_t));
563 r = (uint64_t)p * (uint64_t)q;
564 memcpy(dest->wordv, &r, (dest->length_W = a->length_W + b->length_W) * sizeof(bigint_word_t));
568 /* split a in xh & xl; split b in yh & yl */
569 const bigint_length_t n = (MAX(a->length_W, b->length_W) + 1) / 2;
570 bigint_t xl, xh, yl, yh;
573 if (a->length_W <= n) {
574 bigint_set_zero(&xh);
575 xl.length_W = a->length_W;
581 xh.wordv = &(a->wordv[n]);
582 xh.length_W = a->length_W-n;
585 if (b->length_W <= n) {
586 bigint_set_zero(&yh);
587 yl.length_W = b->length_W;
593 yh.wordv = &(b->wordv[n]);
594 yh.length_W = b->length_W-n;
597 /* now we have split up a and b */
598 /* remember we want to do:
599 * x*y = (xh * b**n + xl) * (yh * b**n + yl)
600 * = (xh * yh) * b**2n + xh * b**n * yl + yh * b**n * xl + xl * yl
601 * = (xh * yh) * b**2n + (xh * yl + yh * xl) * b**n + xl *yl
602 * // xh * yl + yh * xl = (xh + yh) * (xl + yl) - xh * yh - xl * yl
603 * x*y = (xh * yh) * b**2n + ((xh+xl)*(yh+yl) - xh*yh - xl*yl)*b**n + xl*yl
604 * 5 9 2 4 3 7 5 6 1 8 1
606 ALLOC_BIGINT_WORDS(tmp_w, 2 * n + 2);
607 ALLOC_BIGINT_WORDS(m_w, 2 * n + 2);
608 bigint_t tmp, tmp2, m;
610 tmp2.wordv = &(tmp_w[n + 1]);
613 bigint_mul_u(dest, &xl, &yl); /* 1: dest <= xl*yl */
614 bigint_add_u(&tmp2, &xh, &xl); /* 2: tmp2 <= xh+xl */
615 bigint_add_u(&tmp, &yh, &yl); /* 3: tmp <= yh+yl */
616 bigint_mul_u(&m, &tmp2, &tmp); /* 4: m <= tmp2*tmp */
617 bigint_mul_u(&tmp, &xh, &yh); /* 5: tmp <= xh*yh */
618 bigint_sub_u(&m, &m, dest); /* 6: m <= m-dest */
619 bigint_sub_u(&m, &m, &tmp); /* 7: m <= m-tmp */
620 bigint_add_scale_u(dest, &m, n * sizeof(bigint_word_t)); /* 8: dest <= dest+m**n*/
621 bigint_add_scale_u(dest, &tmp, 2 * n * sizeof(bigint_word_t)); /* 9: dest <= dest+tmp**(2*n) */
626 /******************************************************************************/
628 void bigint_mul_s(bigint_t *dest, const bigint_t *a, const bigint_t *b){
631 s |= GET_SIGN(b)?1:0;
633 case 0: /* both positive */
634 bigint_mul_u(dest, a,b);
637 case 1: /* a positive, b negative */
638 bigint_mul_u(dest, a,b);
641 case 2: /* a negative, b positive */
642 bigint_mul_u(dest, a,b);
645 case 3: /* both negative */
646 bigint_mul_u(dest, a,b);
649 default: /* how can this happen?*/
654 /******************************************************************************/
656 void bigint_square(bigint_t *dest, const bigint_t *a) {
657 union __attribute__((packed)) {
659 bigint_wordplus_t uv;
661 bigint_word_t q, c1, c2;
662 bigint_length_t i, j;
664 if (a->length_W * sizeof(bigint_word_t) <= 4) {
666 memcpy(&r, a->wordv, a->length_W * sizeof(bigint_word_t));
668 memcpy(dest->wordv, &r, 2 * a->length_W * sizeof(bigint_word_t));
670 dest->length_W = 2 * a->length_W;
675 if (dest->wordv == a->wordv) {
677 ALLOC_BIGINT_WORDS(d_w, a->length_W * 2);
679 bigint_square(&d, a);
680 bigint_copy(dest, &d);
685 memset(dest->wordv, 0, a->length_W * 2 * sizeof(bigint_word_t));
687 for (i = 0; i < a->length_W; ++i) {
688 acc.uv = (bigint_wordplus_t)a->wordv[i] * (bigint_wordplus_t)a->wordv[i] + (bigint_wordplus_t)dest->wordv[2 * i];
689 dest->wordv[2 * i] = acc.u[0];
692 for (j = i + 1; j < a->length_W; ++j) {
693 acc.uv = (bigint_wordplus_t)a->wordv[i] * (bigint_wordplus_t)a->wordv[j];
694 q = acc.u[1] >> (BIGINT_WORD_SIZE - 1);
696 acc.uv += dest->wordv[i + j];
697 q += (acc.uv < dest->wordv[i + j]);
700 dest->wordv[i + j] = acc.u[0];
701 c1 = (bigint_wordplus_t)acc.u[1] + c2;
704 dest->wordv[i + a->length_W] += c1;
705 if (i < a->length_W - 1) {
706 dest->wordv[i + a->length_W + 1] = c2 + (dest->wordv[i + a->length_W] < c1);
710 dest->length_W = 2 * a->length_W;
714 /******************************************************************************/
716 void bigint_sub_u_bitscale(bigint_t *a, const bigint_t *b, bigint_length_t bitscale){
718 const bigint_length_t word_shift = bitscale / BIGINT_WORD_SIZE;
720 if (a->length_W < b->length_W + word_shift) {
722 cli_putstr("\r\nDBG: *bang*\r\n");
727 ALLOC_BIGINT_WORDS(tmp_w, b->length_W + 1);
729 bigint_copy(&tmp, b);
730 bigint_shiftleft_bits(&tmp, bitscale % BIGINT_WORD_SIZE);
733 x.wordv = &(a->wordv[word_shift]);
734 x.length_W = a->length_W - word_shift;
736 bigint_sub_u(&x, &x, &tmp);
742 /******************************************************************************/
744 void bigint_reduce(bigint_t *a, const bigint_t *r){
745 uint8_t rfbs = GET_FBS(r);
746 if (r->length_W == 0 || a->length_W == 0) {
750 if (bigint_length_b(a) + 3 > bigint_length_b(r)) {
751 if ((r->length_W * sizeof(bigint_word_t) <= 4) && (a->length_W * sizeof(bigint_word_t) <= 4)) {
752 uint32_t p = 0, q = 0;
753 memcpy(&p, a->wordv, a->length_W * sizeof(bigint_word_t));
754 memcpy(&q, r->wordv, r->length_W * sizeof(bigint_word_t));
756 memcpy(a->wordv, &p, a->length_W * sizeof(bigint_word_t));
757 a->length_W = r->length_W;
762 while (a->length_W > r->length_W) {
763 shift = (a->length_W - r->length_W) * CHAR_BIT * sizeof(bigint_word_t) + GET_FBS(a) - rfbs - 1;
764 bigint_sub_u_bitscale(a, r, shift);
766 while ((GET_FBS(a) > rfbs) && (a->length_W == r->length_W)) {
767 shift = GET_FBS(a) - rfbs - 1;
768 bigint_sub_u_bitscale(a, r, shift);
771 while (bigint_cmp_u(a, r) >= 0) {
772 bigint_sub_u(a, a, r);
777 /******************************************************************************/
779 /* calculate dest = a**exp % r */
780 /* using square&multiply */
781 void bigint_expmod_u_sam(bigint_t *dest, const bigint_t *a, const bigint_t *exp, const bigint_t *r){
782 if (a->length_W == 0) {
783 bigint_set_zero(dest);
787 if(exp->length_W == 0){
799 ALLOC_BIGINT_WORDS(base_w, MAX(a->length_W, r->length_W));
800 ALLOC_BIGINT_WORDS(res_w, r->length_W * 2);
804 bigint_copy(&base, a);
805 bigint_reduce(&base, r);
810 bigint_copy(&res, &base);
812 t = exp->wordv[exp->length_W - 1];
813 for (i = exp->length_W; i > 0; --i) {
814 t = exp->wordv[i - 1];
815 for (j = BIGINT_WORD_SIZE; j > 0; --j) {
817 if (t & (((bigint_word_t)1) << (BIGINT_WORD_SIZE - 1))) {
821 bigint_square(&res, &res);
822 bigint_reduce(&res, r);
823 if(t & (((bigint_word_t)1) << (BIGINT_WORD_SIZE - 1))){
824 bigint_mul_u(&res, &res, &base);
825 bigint_reduce(&res, r);
833 bigint_copy(dest, &res);
838 /******************************************************************************/
839 /* gcd <-- gcd(x,y) a*x+b*y=gcd */
840 void bigint_gcdext(bigint_t *gcd, bigint_t *a, bigint_t *b, const bigint_t *x, const bigint_t *y){
841 bigint_length_t i = 0;
842 if(x->length_W == 0 || y->length_W == 0){
844 bigint_set_zero(gcd);
854 if(x->length_W == 1 && x->wordv[0] == 1){
871 if(y->length_W == 1 && y->wordv[0] == 1){
889 while(x->wordv[i] == 0 && y->wordv[i] == 0){
893 ALLOC_BIGINT_WORDS(g_w, i + 2);
894 ALLOC_BIGINT_WORDS(x_w, x->length_W - i);
895 ALLOC_BIGINT_WORDS(y_w, y->length_W - i);
896 ALLOC_BIGINT_WORDS(u_w, x->length_W - i);
897 ALLOC_BIGINT_WORDS(v_w, y->length_W - i);
898 ALLOC_BIGINT_WORDS(a_w, y->length_W + 2);
899 ALLOC_BIGINT_WORDS(c_w, y->length_W + 2);
900 ALLOC_BIGINT_WORDS(b_w, x->length_W + 2);
901 ALLOC_BIGINT_WORDS(d_w, x->length_W + 2);
903 bigint_t g, x_, y_, u, v, a_, b_, c_, d_;
908 memset(g_w, 0, i * sizeof(bigint_word_t));
912 x_.info = y_.info = 0;
913 x_.length_W = x->length_W - i;
914 y_.length_W = y->length_W - i;
915 memcpy(x_.wordv, x->wordv + i, x_.length_W * sizeof(bigint_word_t));
916 memcpy(y_.wordv, y->wordv + i, y_.length_W * sizeof(bigint_word_t));
918 for(i = 0; (x_.wordv[0] & ((bigint_word_t)1 << i)) == 0 && (y_.wordv[0] & ((bigint_word_t)1 << i)) == 0; ++i)
925 bigint_shiftleft_bits(&g, i);
926 bigint_shiftright_bits(&x_, i);
927 bigint_shiftright_bits(&y_, i);
937 bigint_copy(&u, &x_);
938 bigint_copy(&v, &y_);
945 bigint_set_zero(&b_);
946 bigint_set_zero(&c_);
948 while ((u.wordv[0] & 1) == 0) {
949 bigint_shiftright_1bit(&u);
950 if((a_.wordv[0] & 1) || (b_.wordv[0] & 1)){
951 bigint_add_s(&a_, &a_, &y_);
952 bigint_sub_s(&b_, &b_, &x_);
954 bigint_shiftright_1bit(&a_);
955 bigint_shiftright_1bit(&b_);
957 while ((v.wordv[0] & 1) == 0) {
958 bigint_shiftright_1bit(&v);
959 if((c_.wordv[0] & 1) || (d_.wordv[0] & 1)){
960 bigint_add_s(&c_, &c_, &y_);
961 bigint_sub_s(&d_, &d_, &x_);
963 bigint_shiftright_1bit(&c_);
964 bigint_shiftright_1bit(&d_);
966 if(bigint_cmp_u(&u, &v) >= 0){
967 bigint_sub_u(&u, &u, &v);
968 bigint_sub_s(&a_, &a_, &c_);
969 bigint_sub_s(&b_, &b_, &d_);
971 bigint_sub_u(&v, &v, &u);
972 bigint_sub_s(&c_, &c_, &a_);
973 bigint_sub_s(&d_, &d_, &b_);
977 bigint_mul_s(gcd, &v, &g);
997 /******************************************************************************/
999 void bigint_inverse(bigint_t *dest, const bigint_t *a, const bigint_t *m){
1000 bigint_gcdext(NULL, dest, NULL, a, m);
1001 while(dest->info&BIGINT_NEG_MASK){
1002 bigint_add_s(dest, dest, m);
1006 /******************************************************************************/
1008 void bigint_changeendianess(bigint_t *a){
1010 p = (uint8_t*)a->wordv;
1011 q = p + a->length_W * sizeof(bigint_word_t) - 1;
1020 /******************************************************************************/
1022 void bigint_mul_word_u(bigint_t *a, bigint_word_t b){
1023 bigint_wordplus_t c0 = 0, c1 = 0;
1031 for(i = 0; i < a->length_W; ++i){
1032 c1 = ((bigint_wordplus_t)(a->wordv[i])) * ((bigint_wordplus_t)b);
1034 a->wordv[i] = (bigint_word_t)c1;
1035 c0 = c1 >> BIGINT_WORD_SIZE;
1038 a->wordv[a->length_W] = (bigint_word_t)c0;
1044 /******************************************************************************/
1046 void bigint_clip(bigint_t *dest, bigint_length_t length_W){
1047 if(dest->length_W > length_W){
1048 dest->length_W = length_W;
1050 bigint_adjust(dest);
1053 /******************************************************************************/
1056 * dest = (a * b) % m (?)
1059 void bigint_mont_mul(bigint_t *dest, const bigint_t *a, const bigint_t *b, const bigint_t *m, const bigint_t *m_){
1060 const bigint_length_t s = MAX(MAX(a->length_W, b->length_W), m->length_W);
1065 if (a->length_W == 0 || b->length_W == 0) {
1066 bigint_set_zero(dest);
1069 ALLOC_BIGINT_WORDS(u_w, s + 2);
1070 ALLOC_BIGINT_WORDS(t_w, s + 2);
1075 for (i = 0; i < a->length_W; ++i) {
1077 bigint_mul_word_u(&t, a->wordv[i]);
1078 bigint_add_u(&u, &u, &t);
1079 bigint_copy(&t, m_);
1080 if (u.length_W != 0) {
1081 bigint_mul_word_u(&t, u.wordv[0]);
1082 bigint_add_u(&u, &u, &t);
1084 bigint_shiftright_1word(&u);
1086 for (; i < s; ++i) {
1087 bigint_copy(&t, m_);
1088 if (u.length_W != 0) {
1089 bigint_mul_word_u(&t, u.wordv[0]);
1090 bigint_add_u(&u, &u, &t);
1092 bigint_shiftright_1word(&u);
1094 bigint_reduce(&u, m);
1095 bigint_copy(dest, &u);
1100 /******************************************************************************/
1102 void bigint_mont_red(bigint_t *dest, const bigint_t *a, const bigint_t *m, const bigint_t *m_){
1104 bigint_length_t i, s = MAX(a->length_W, MAX(m->length_W, m_->length_W));
1106 if (a->length_W == 0) {
1107 bigint_set_zero(dest);
1111 ALLOC_BIGINT_WORDS(u_w, s + 2);
1112 ALLOC_BIGINT_WORDS(t_w, s + 2);
1116 for (i = 0; i < m->length_W; ++i) {
1117 bigint_copy(&t, m_);
1118 if (u.length_W != 0) {
1119 bigint_mul_word_u(&t, u.wordv[0]);
1120 bigint_add_u(&u, &u, &t);
1122 bigint_shiftright_1word(&u);
1124 bigint_reduce(&u, m);
1125 bigint_copy(dest, &u);
1130 /******************************************************************************/
1132 * m_ = m * (- m0^-1 (mod 2^W))
1134 void bigint_mont_gen_m_(bigint_t* dest, const bigint_t* m){
1135 bigint_word_t x_w[2], m_w_0[1];
1137 if (m->length_W == 0) {
1138 bigint_set_zero(dest);
1141 if ((m->wordv[0] & 1) == 0) {
1143 printf_P(PSTR("ERROR: m must not be even, m = "));
1144 bigint_print_hex(m);
1158 m_0.wordv[0] = m->wordv[0];
1160 bigint_adjust(&m_0);
1161 bigint_inverse(dest, &m_0, &x);
1162 bigint_sub_s(&x, &x, dest);
1163 bigint_copy(dest, m);
1164 bigint_mul_word_u(dest, x.wordv[0]);
1168 /******************************************************************************/
1171 * dest = a * R mod m
1173 void bigint_mont_trans(bigint_t *dest, const bigint_t *a, const bigint_t *m){
1176 ALLOC_BIGINT_WORDS(t_w, a->length_W + m->length_W);
1178 memset(t_w, 0, m->length_W * sizeof(bigint_word_t));
1179 memcpy(&t_w[m->length_W], a->wordv, a->length_W * sizeof(bigint_word_t));
1181 t.length_W = a->length_W + m->length_W;
1182 bigint_reduce(&t, m);
1183 bigint_copy(dest, &t);
1187 /******************************************************************************/
1189 /* calculate dest = a**exp % r */
1190 /* using square&multiply */
1191 void bigint_expmod_u_mont_accel(bigint_t *dest, const bigint_t *a, const bigint_t *exp, const bigint_t *r, const bigint_t *m_){
1192 if(r->length_W == 0) {
1201 if (exp->length_W == 0) {
1208 ALLOC_BIGINT_WORDS(res_w, r->length_W * 2);
1209 ALLOC_BIGINT_WORDS(ax_w, MAX(r->length_W, a->length_W));
1218 bigint_copy(&ax, a);
1219 bigint_reduce(&ax, r);
1221 bigint_mont_trans(&ax, &ax, r);
1222 bigint_mont_trans(&res, &res, r);
1225 t = exp->wordv[exp->length_W - 1];
1226 for (i = exp->length_W; i > 0; --i) {
1227 t = exp->wordv[i - 1];
1228 for(j = BIGINT_WORD_SIZE; j > 0; --j){
1230 if(t & (((bigint_word_t)1) << (BIGINT_WORD_SIZE - 1))){
1235 bigint_square(&res, &res);
1236 bigint_mont_red(&res, &res, r, m_);
1237 if (t & (((bigint_word_t)1) << (BIGINT_WORD_SIZE - 1))) {
1238 bigint_mont_mul(&res, &res, &ax, r, m_);
1245 bigint_mont_red(dest, &res, r, m_);
1250 /******************************************************************************/
1252 void bigint_expmod_u_mont_sam(bigint_t *dest, const bigint_t *a, const bigint_t *exp, const bigint_t *r){
1253 if(r->length_W == 0) {
1256 if(a->length_W == 0) {
1257 bigint_set_zero(dest);
1261 bigint_word_t m_w_[r->length_W + 1];
1263 bigint_mont_gen_m_(&m_, r);
1264 bigint_expmod_u_mont_accel(dest, a, exp, r,&m_);
1267 /******************************************************************************/
1269 void bigint_expmod_u(bigint_t *dest, const bigint_t *a, const bigint_t *exp, const bigint_t *r){
1271 printf("\nDBG: expmod_u (a ** e %% m) <%s %s %d>\n\ta: ", __FILE__, __func__, __LINE__);
1272 bigint_print_hex(a);
1274 bigint_print_hex(exp);
1276 bigint_print_hex(r);
1278 if (r->wordv[0] & 1) {
1279 bigint_expmod_u_mont_sam(dest, a, exp, r);
1281 bigint_expmod_u_sam(dest, a, exp, r);