d7ee059473407006498537c407a0579c62312c73
[arm-crypto-lib.git] / bigint / bigint.c
1 /* bigint.c */
2 /*
3     This file is part of the AVR-Crypto-Lib.
4     Copyright (C) 2008  Daniel Otte (daniel.otte@rub.de)
5
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.
10
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.
15
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/>.
18 */
19 /**
20  * \file                bigint.c
21  * \author              Daniel Otte
22  * \date                2010-02-22
23  * 
24  * \license         GPLv3 or later
25  * 
26  */
27  
28
29 #define STRING2(x) #x
30 #define STRING(x) STRING2(x)
31 #define STR_LINE STRING(__LINE__)
32
33 #include "bigint.h"
34 #include <string.h>
35
36 #define DEBUG 1
37
38 #if DEBUG
39 #include "cli.h"
40 #include "uart_lowlevel.h"
41 #include "bigint_io.h"
42 #endif
43
44 #ifndef MAX
45  #define MAX(a,b) (((a)>(b))?(a):(b))
46 #endif
47
48 #ifndef MIN
49  #define MIN(a,b) (((a)<(b))?(a):(b))
50 #endif
51
52 #define SET_FBS(a, v) do{(a)->info &=~BIGINT_FBS_MASK; (a)->info |= (v);}while(0)
53 #define GET_FBS(a)   ((a)->info&BIGINT_FBS_MASK)
54 #define SET_NEG(a)   (a)->info |= BIGINT_NEG_MASK
55 #define SET_POS(a)   (a)->info &= ~BIGINT_NEG_MASK
56 #define XCHG(a,b)    do{(a)^=(b); (b)^=(a); (a)^=(b);}while(0)
57 #define XCHG_PTR(a,b)    do{ a = (void*)(((uint32_t)(a)) ^ ((uint32_t)(b))); \
58                                  b = (void*)(((uint32_t)(a)) ^ ((uint32_t)(b))); \
59                                  a = (void*)(((uint32_t)(a)) ^ ((uint32_t)(b)));}while(0)
60
61 #define GET_SIGN(a) ((a)->info&BIGINT_NEG_MASK)
62
63 /******************************************************************************/
64 void bigint_adjust(bigint_t* a){
65         while(a->length_B!=0 && a->wordv[a->length_B-1]==0){
66                 a->length_B--;
67         }
68         if(a->length_B==0){
69                 a->info=0;
70                 return;
71         }
72         bigint_word_t t;
73         uint8_t i = BIGINT_WORD_SIZE-1;
74         t = a->wordv[a->length_B-1];
75         while((t&(1<<(BIGINT_WORD_SIZE-1)))==0 && i){
76                 t<<=1;
77                 i--;
78         }
79         SET_FBS(a, i);
80 }
81
82 /******************************************************************************/
83
84 void bigint_copy(bigint_t* dest, const bigint_t* src){
85         memcpy(dest->wordv, src->wordv, src->length_B*sizeof(bigint_word_t));
86         dest->length_B = src->length_B;
87         dest->info = src->info;
88 }
89
90 /******************************************************************************/
91
92 /* this should be implemented in assembly */
93 void bigint_add_u(bigint_t* dest, const bigint_t* a, const bigint_t* b){
94         uint16_t i;
95         bigint_wordplus_t t=0LL;
96         if(a->length_B < b->length_B){
97                 XCHG_PTR(a,b);
98         }
99         for(i=0; i<b->length_B; ++i){
100 //              t = (bigint_wordplus_t)(a->wordv[i]) + (bigint_wordplus_t)(b->wordv[i]) + t;
101                 t += a->wordv[i];
102                 t += b->wordv[i];
103                 dest->wordv[i] = (bigint_word_t)t;
104                 t>>=BIGINT_WORD_SIZE;
105         }
106         for(; i<a->length_B; ++i){
107                 t += a->wordv[i];
108                 dest->wordv[i] = (bigint_word_t)t;
109                 t>>=BIGINT_WORD_SIZE;
110         }
111         dest->wordv[i++] = (bigint_word_t)t;
112         dest->length_B = i;
113         bigint_adjust(dest);
114 }
115
116 /******************************************************************************/
117
118 /* this should be implemented in assembly */
119 void bigint_add_scale_u(bigint_t* dest, const bigint_t* a, uint16_t scale){
120         uint16_t i,j=0;
121         uint16_t scale_w;
122         uint32_t *dst;
123         bigint_wordplus_t t=0;
124         scale_w = (scale+sizeof(bigint_word_t)-1)/sizeof(bigint_word_t);
125         if(scale>dest->length_B*sizeof(bigint_word_t)){
126                 memset(((uint8_t*)dest->wordv)+dest->length_B*sizeof(bigint_word_t), 0, scale-dest->length_B*sizeof(bigint_word_t));
127         }
128         // a->wordv = (const uint32_t*)(((uint8_t*)a->wordv)+(scale&3));
129         dst  = dest->wordv + (scale&(sizeof(bigint_word_t)-1));
130         for(i=scale/sizeof(bigint_word_t); i<a->length_B+scale_w; ++i,++j){
131                 t += a->wordv[j];
132                 if(dest->length_B>i){
133                         t += dst[i];
134                 }
135                 dst[i] = (bigint_word_t)t;
136                 t>>=BIGINT_WORD_SIZE;
137         }
138         while(t){
139                 if(dest->length_B>i){
140                         t += dst[i];
141                 }
142                 dst[i] = (bigint_word_t)t;
143                 t>>=BIGINT_WORD_SIZE;
144                 ++i;
145         }
146         if(dest->length_B < i){
147                 dest->length_B = i;
148         }
149         bigint_adjust(dest);
150 }
151
152 /******************************************************************************/
153
154 /* this should be implemented in assembly */
155 void bigint_sub_u(bigint_t* dest, const bigint_t* a, const bigint_t* b){
156         int8_t borrow=0;
157         int8_t  r;
158         bigint_wordplus_signed_t t=0LL;
159         uint16_t i, min, max;
160         min = MIN(a->length_B, b->length_B);
161         max = MAX(a->length_B, b->length_B);
162         r = bigint_cmp_u(a,b);
163         if(r==0){
164                 bigint_set_zero(dest);
165                 return;
166         }
167         if(b->length_B==0){
168                 bigint_copy(dest, a);
169                 SET_POS(dest);
170                 return;
171         }
172         if(a->length_B==0){
173                 bigint_copy(dest, b);
174                 SET_NEG(dest);
175                 return;
176         }
177         if(r<0){
178                 bigint_sub_u(dest, b, a);
179                 SET_NEG(dest);
180         }else{
181                 for(i=0; i<min; ++i){
182                         t = a->wordv[i];
183                         t -= b->wordv[i];
184                         t -= borrow;
185                         if(t<0){
186                                 borrow = 1;
187                                 dest->wordv[i]=(bigint_word_t)t;
188                         }else{
189                                 borrow = 0;
190                                 dest->wordv[i]=(bigint_word_t)t;
191                         }
192                 }
193                 for(;i<max; ++i){
194                         t = a->wordv[i] - borrow;
195                         if(t<0){
196                                 borrow = 1;
197                                 dest->wordv[i]=(bigint_word_t)t;
198                         }else{
199                                 borrow = 0;
200                                 dest->wordv[i]=(bigint_word_t)t;
201                         }
202
203                 }
204                 SET_POS(dest);
205                 dest->length_B = i;
206                 bigint_adjust(dest);
207         }
208 }
209
210 /******************************************************************************/
211
212 int8_t bigint_cmp_u(const bigint_t* a, const bigint_t* b){
213         if(a->length_B > b->length_B){
214                 return 1;
215         }
216         if(a->length_B < b->length_B){
217                 return -1;
218         }
219         if(a->length_B==0){
220                 return 0;
221         }
222         uint16_t i;
223         i = a->length_B-1;
224         do{
225                 if(a->wordv[i]!=b->wordv[i]){
226                         if(a->wordv[i]>b->wordv[i]){
227                                 return 1;
228                         }else{
229                                 return -1;
230                         }
231                 }
232         }while(i--);
233         return 0;
234 }
235
236 /******************************************************************************/
237
238 void bigint_add_s(bigint_t* dest, const bigint_t* a, const bigint_t* b){
239         uint8_t s;
240         s  = GET_SIGN(a)?2:0;
241         s |= GET_SIGN(b)?1:0;
242         switch(s){
243                 case 0: /* both positive */
244                         bigint_add_u(dest, a,b);
245                         SET_POS(dest);
246                         break;
247                 case 1: /* a positive, b negative */
248                         bigint_sub_u(dest, a, b);
249                         break;
250                 case 2: /* a negative, b positive */
251                         bigint_sub_u(dest, b, a);
252                         break;
253                 case 3: /* both negative */
254                         bigint_add_u(dest, a, b);
255                         SET_NEG(dest);
256                         break;
257                 default: /* how can this happen?*/
258                         break;
259         }
260 }
261
262 /******************************************************************************/
263
264 void bigint_sub_s(bigint_t* dest, const bigint_t* a, const bigint_t* b){
265         uint8_t s;
266         s  = GET_SIGN(a)?2:0;
267         s |= GET_SIGN(b)?1:0;
268         switch(s){
269                 case 0: /* both positive */
270                         bigint_sub_u(dest, a,b);
271                         break;
272                 case 1: /* a positive, b negative */
273                         bigint_add_u(dest, a, b);
274                         SET_POS(dest);
275                         break;
276                 case 2: /* a negative, b positive */
277                         bigint_add_u(dest, a, b);
278                         SET_NEG(dest);
279                         break;
280                 case 3: /* both negative */
281                         bigint_sub_u(dest, b, a);
282                         break;
283                 default: /* how can this happen?*/
284                                         break;
285         }
286
287 }
288
289 /******************************************************************************/
290
291 int8_t bigint_cmp_s(const bigint_t* a, const bigint_t* b){
292         uint8_t s;
293         if(a->length_B==0 && b->length_B==0){
294                 return 0;
295         }
296         s  = GET_SIGN(a)?2:0;
297         s |= GET_SIGN(b)?1:0;
298         switch(s){
299                 case 0: /* both positive */
300                         return bigint_cmp_u(a, b);
301                         break;
302                 case 1: /* a positive, b negative */
303                         return 1;
304                         break;
305                 case 2: /* a negative, b positive */
306                         return -1;
307                         break;
308                 case 3: /* both negative */
309                         return bigint_cmp_u(b, a);
310                         break;
311                 default: /* how can this happen?*/
312                                         break;
313         }
314         return 0; /* just to satisfy the compiler */
315 }
316
317 /******************************************************************************/
318
319 void bigint_shiftleft(bigint_t* a, uint16_t shift){
320         uint16_t byteshift, word_alloc;
321         int16_t i;
322         uint8_t bitshift;
323         bigint_word_t *p;
324         bigint_wordplus_t t=0;
325         if(shift==0){
326                 return;
327         }
328         byteshift = shift/8;
329         bitshift = shift&7;
330         for(i=0;i<=byteshift/sizeof(bigint_word_t); ++i){
331                 a->wordv[a->length_B+i] = 0;
332         }
333         if(byteshift){
334                 memmove(((uint8_t*)a->wordv)+byteshift, a->wordv, a->length_B*sizeof(bigint_word_t));
335                 memset(a->wordv, 0, byteshift);
336         }
337         p = (bigint_word_t*)(((uint8_t*)a->wordv)+byteshift);
338         word_alloc = a->length_B+(byteshift+sizeof(bigint_word_t)-1)/sizeof(bigint_word_t)+1;
339         a->wordv[word_alloc-1]=0;
340         if(bitshift!=0){
341                 for(i=0; i<a->length_B; ++i){
342                         t |= ((bigint_wordplus_t)p[i])<<bitshift;
343                         p[i] = (bigint_word_t)t;
344                         t >>= BIGINT_WORD_SIZE;
345                 }
346                 p[i] = (bigint_word_t)t;
347         }
348         a->length_B = word_alloc;
349         bigint_adjust(a);
350 }
351
352 /******************************************************************************/
353
354 void bigint_shiftright(bigint_t* a, uint16_t shift){
355         uint16_t byteshift;
356         uint16_t i;
357         uint8_t bitshift;
358         bigint_wordplus_t t=0;
359         byteshift = shift/8;
360         bitshift = shift&7;
361         if(byteshift >= a->length_B*sizeof(bigint_word_t)){ /* we would shift out more than we have */
362                 bigint_set_zero(a);
363                 return;
364         }
365         if(byteshift == a->length_B*sizeof(bigint_word_t)-1 && bitshift>GET_FBS(a)){
366                 bigint_set_zero(a);
367                 return;
368         }
369         if(byteshift){
370                 memmove(a->wordv, (uint8_t*)a->wordv+byteshift, a->length_B-byteshift);
371                 memset((uint8_t*)a->wordv+a->length_B-byteshift, 0,  byteshift);
372         }
373         byteshift /= sizeof(bigint_word_t);
374         if(bitshift!=0){
375          /* shift to the right */
376                 for(i=a->length_B-byteshift-1; i>0; --i){
377                         t |= ((bigint_wordplus_t)(a->wordv[i]))<<(BIGINT_WORD_SIZE-bitshift);
378                         a->wordv[i] = (bigint_word_t)(t>>BIGINT_WORD_SIZE);
379                         t <<= BIGINT_WORD_SIZE;
380                 }
381                 t |= ((bigint_wordplus_t)(a->wordv[0]))<<(BIGINT_WORD_SIZE-bitshift);
382                 a->wordv[0] = (bigint_word_t)(t>>BIGINT_WORD_SIZE);
383         }
384         a->length_B -= ((shift/8)+sizeof(bigint_word_t)-1)/sizeof(bigint_word_t);
385         bigint_adjust(a);
386 }
387
388 /******************************************************************************/
389
390 void bigint_xor(bigint_t* dest, const bigint_t* a){
391         uint16_t i;
392         for(i=0; i<a->length_B; ++i){
393                 dest->wordv[i] ^= a->wordv[i];
394         }
395         bigint_adjust(dest);
396 }
397
398 /******************************************************************************/
399
400 void bigint_set_zero(bigint_t* a){
401         a->length_B=0;
402 }
403
404 /******************************************************************************/
405
406 /* using the Karatsuba-Algorithm */
407 /* x*y = (xh*yh)*b**2n + ((xh+xl)*(yh+yl) - xh*yh - xl*yl)*b**n + yh*yl */
408 void bigint_mul_u(bigint_t* dest, const bigint_t* a, const bigint_t* b){
409         if(a->length_B==0 || b->length_B==0){
410                 bigint_set_zero(dest);
411                 return;
412         }
413         if(dest==a || dest==b){
414                 bigint_t d;
415                 bigint_word_t d_b[a->length_B+b->length_B];
416                 d.wordv = d_b;
417                 bigint_mul_u(&d, a, b);
418                 bigint_copy(dest, &d);
419                 return;
420         }
421         if(a->length_B==1 || b->length_B==1){
422                 if(a->length_B!=1){
423                         XCHG_PTR(a,b);
424                 }
425                 bigint_wordplus_t i, t=0;
426                 bigint_word_t x = a->wordv[0];
427                 for(i=0; i<b->length_B; ++i){
428                         t += ((bigint_wordplus_t)b->wordv[i])*((bigint_wordplus_t)x);
429                         dest->wordv[i] = (bigint_word_t)t;
430                         t>>=BIGINT_WORD_SIZE;
431                 }
432                 dest->wordv[i] = (bigint_word_t)t;
433                 dest->length_B=i+1;
434                 bigint_adjust(dest);
435                 return;
436         }
437         if(a->length_B<=4/sizeof(bigint_word_t) && b->length_B<=4/sizeof(bigint_word_t)){
438                 uint32_t p=0, q=0;
439                 uint64_t r;
440                 memcpy(&p, a->wordv, a->length_B*sizeof(bigint_word_t));
441                 memcpy(&q, b->wordv, b->length_B*sizeof(bigint_word_t));
442                 r = (uint64_t)p*(uint64_t)q;
443                 memcpy(dest->wordv, &r, (a->length_B+b->length_B)*sizeof(bigint_word_t));
444                 dest->length_B =  a->length_B+b->length_B;
445                 bigint_adjust(dest);
446                 return;
447         }
448         bigint_set_zero(dest);
449         /* split a in xh & xl; split b in yh & yl */
450         uint16_t n;
451         n=(MAX(a->length_B, b->length_B)+1)/2;
452         bigint_t xl, xh, yl, yh;
453         xl.wordv = a->wordv;
454         yl.wordv = b->wordv;
455         if(a->length_B<=n){
456                 xh.info=0;
457                 xh.length_B = 0;
458                 xl.length_B = a->length_B;
459                 xl.info = 0;
460         }else{
461                 xl.length_B=n;
462                 xl.info = 0;
463                 bigint_adjust(&xl);
464                 xh.wordv = a->wordv+n;
465                 xh.length_B = a->length_B-n;
466                 xh.info = 0;
467         }
468         if(b->length_B<=n){
469                 yh.info=0;
470                 yh.length_B = 0;
471                 yl.length_B = b->length_B;
472                 yl.info = b->info;
473         }else{
474                 yl.length_B=n;
475                 yl.info = 0;
476                 bigint_adjust(&yl);
477                 yh.wordv = b->wordv+n;
478                 yh.length_B = b->length_B-n;
479                 yh.info = 0;
480         }
481         /* now we have split up a and b */
482         bigint_word_t  tmp_b[2*n+2], m_b[2*(n+1)];
483         bigint_t tmp, tmp2, m;
484         tmp.wordv = tmp_b;
485         tmp2.wordv = tmp_b+n+1;
486         m.wordv = m_b;
487
488         bigint_mul_u(dest, &xl, &yl);  /* dest <= xl*yl     */
489         bigint_add_u(&tmp2, &xh, &xl); /* tmp2 <= xh+xl     */
490         bigint_add_u(&tmp, &yh, &yl);  /* tmp  <= yh+yl     */
491         bigint_mul_u(&m, &tmp2, &tmp); /* m    <= tmp2*tmp  */
492         bigint_mul_u(&tmp, &xh, &yh);  /* h    <= xh*yh     */
493         bigint_sub_u(&m, &m, dest);    /* m    <= m-dest    */
494     bigint_sub_u(&m, &m, &tmp);    /* m    <= m-h       */
495         bigint_add_scale_u(dest, &m, n*sizeof(bigint_word_t));
496         bigint_add_scale_u(dest, &tmp, 2*n*sizeof(bigint_word_t));
497 }
498
499 /******************************************************************************/
500
501 void bigint_mul_s(bigint_t* dest, const bigint_t* a, const bigint_t* b){
502         uint8_t s;
503         s  = GET_SIGN(a)?2:0;
504         s |= GET_SIGN(b)?1:0;
505         switch(s){
506                 case 0: /* both positive */
507                         bigint_mul_u(dest, a,b);
508                         SET_POS(dest);
509                         break;
510                 case 1: /* a positive, b negative */
511                         bigint_mul_u(dest, a,b);
512                         SET_NEG(dest);
513                         break;
514                 case 2: /* a negative, b positive */
515                         bigint_mul_u(dest, a,b);
516                         SET_NEG(dest);
517                         break;
518                 case 3: /* both negative */
519                         bigint_mul_u(dest, a,b);
520                         SET_POS(dest);
521                         break;
522                 default: /* how can this happen?*/
523                         break;
524         }
525 }
526
527 /******************************************************************************/
528
529 /* square */
530 /* (xh*b^n+xl)^2 = xh^2*b^2n + 2*xh*xl*b^n + xl^2 */
531 void bigint_square(bigint_t* dest, const bigint_t* a){
532         if(a->length_B*sizeof(bigint_word_t)<=4){
533                 uint64_t r=0;
534                 memcpy(&r, a->wordv, a->length_B*sizeof(bigint_word_t));
535                 r = r*r;
536                 memcpy(dest->wordv, &r, 2*a->length_B*sizeof(bigint_word_t));
537                 SET_POS(dest);
538                 dest->length_B=2*a->length_B;
539                 bigint_adjust(dest);
540                 return;
541         }
542         if(dest==a){
543                 bigint_t d;
544                 bigint_word_t d_b[a->length_B*2];
545                 d.wordv = d_b;
546                 bigint_square(&d, a);
547                 bigint_copy(dest, &d);
548                 return;
549         }
550         uint16_t n;
551         n=(a->length_B+1)/2;
552         bigint_t xh, xl, tmp; /* x-high, x-low, temp */
553         bigint_word_t buffer[2*n+1];
554         xl.wordv = a->wordv;
555         xl.length_B = n;
556         xh.wordv = &(a->wordv[n]);
557         xh.length_B = a->length_B-n;
558         tmp.wordv = buffer;
559 //      cli_putstr("\r\nDBG (a): xl: "); bigint_print_hex(&xl);
560 //      cli_putstr("\r\nDBG (b): xh: "); bigint_print_hex(&xh);
561         bigint_square(dest, &xl);
562 //      cli_putstr("\r\nDBG (1): xl**2: "); bigint_print_hex(dest);
563         bigint_square(&tmp, &xh);
564 //      cli_putstr("\r\nDBG (2): xh**2: "); bigint_print_hex(&tmp);
565         bigint_add_scale_u(dest, &tmp, 2*n*sizeof(bigint_word_t));
566 //      cli_putstr("\r\nDBG (3): xl**2 + xh**2*n**2: "); bigint_print_hex(dest);
567         bigint_mul_u(&tmp, &xl, &xh);
568 //      cli_putstr("\r\nDBG (4): xl*xh: "); bigint_print_hex(&tmp);
569         bigint_shiftleft(&tmp, 1);
570 //      cli_putstr("\r\nDBG (5): xl*xh*2: "); bigint_print_hex(&tmp);
571         bigint_add_scale_u(dest, &tmp, n*sizeof(bigint_word_t));
572 //      cli_putstr("\r\nDBG (6): x**2: "); bigint_print_hex(dest);
573 //      cli_putstr("\r\n");
574 }
575
576 /******************************************************************************/
577
578 #define cli_putstr(a)
579 #define bigint_print_hex(a)
580 #define cli_hexdump_rev(a,b)
581 #define uart_flush(a)
582
583 void bigint_sub_u_bitscale(bigint_t* a, const bigint_t* b, uint16_t bitscale){
584         bigint_t tmp;
585         bigint_word_t tmp_b[b->length_B+4];
586         uint16_t i,j,word_shift=bitscale/(8*sizeof(bigint_word_t));
587         uint8_t borrow=0;
588         bigint_wordplus_signed_t t;
589
590         if(a->length_B < b->length_B+word_shift){
591                 cli_putstr("\r\nDBG: *bang*\r\n");
592                 bigint_set_zero(a);
593                 return;
594         }
595         tmp.wordv = tmp_b;
596         bigint_copy(&tmp, b);
597         bigint_shiftleft(&tmp, bitscale&(BIGINT_WORD_SIZE-1));
598         cli_putstr("\r\nDBG(sub_ub.0) tmp_shift    = "); bigint_print_hex(&tmp);
599         for(j=0,i=word_shift; i<tmp.length_B+word_shift; ++i, ++j){
600                 t = a->wordv[i];
601                 t -= tmp.wordv[j];
602                 t -= borrow;
603                 a->wordv[i] = (bigint_word_t)t;
604                 if(t<0){
605                         borrow = 1;
606                 }else{
607                         borrow = 0;
608                 }
609         }
610         while(borrow){
611                 if(i+1 > a->length_B){
612                         cli_putstr("\r\nDBG: *boom*\r\n");
613                         bigint_set_zero(a);
614                         return;
615                 }
616                 a->wordv[i] -= borrow;
617                 if(a->wordv[i]!=0xff){
618                         borrow=0;
619                 }
620                 ++i;
621         }
622         bigint_adjust(a);
623 }
624
625 /******************************************************************************/
626
627 void bigint_reduce(bigint_t* a, const bigint_t* r){
628 //      bigint_adjust(r);
629         uint8_t rfbs = GET_FBS(r);
630
631         cli_putstr("\r\nDBG: (a) = "); bigint_print_hex(a);
632         if(r->length_B==0 || a->length_B==0){
633                 return;
634         }
635         if((r->length_B*sizeof(bigint_word_t)<=4) && (a->length_B*sizeof(bigint_word_t)<=4)){
636                 uint32_t p=0, q=0;
637                 memcpy(&p, a->wordv, a->length_B*sizeof(bigint_word_t));
638                 memcpy(&q, r->wordv, r->length_B*sizeof(bigint_word_t));
639                 p %= q;
640                 memcpy(a->wordv, &p, a->length_B*sizeof(bigint_word_t));
641                 bigint_adjust(a);
642                 cli_putstr("\r\nDBG: (0) = "); bigint_print_hex(a);
643                 return;
644         }
645         uint16_t shift;
646         while(a->length_B > r->length_B){
647                 shift = (a->length_B-r->length_B)*8*sizeof(bigint_word_t)+GET_FBS(a)-rfbs-1;
648                 cli_putstr("\r\nDBG: (p) shift = "); cli_hexdump_rev(&shift, 2);
649                 uart_flush(0);
650                 bigint_sub_u_bitscale(a, r, shift);
651                 cli_putstr("\r\nDBG: (1) = "); bigint_print_hex(a);
652         }
653         while((GET_FBS(a) > rfbs+1) && (a->length_B == r->length_B)){
654                 shift = GET_FBS(a)-rfbs-1;
655                 cli_putstr("\r\nDBG: (q) shift = "); cli_hexdump_rev(&shift, 2);
656                 bigint_sub_u_bitscale(a, r, GET_FBS(a)-rfbs-1);
657                 cli_putstr("\r\nDBG: (2) = "); bigint_print_hex(a);
658         }
659         while(bigint_cmp_u(a,r)>=0){
660                 bigint_sub_u(a,a,r);
661                 cli_putstr("\r\nDBG: (3) = "); bigint_print_hex(a);
662         }
663         bigint_adjust(a);
664         cli_putstr("\r\nDBG: (a) = "); bigint_print_hex(a);
665         cli_putstr("\r\n");
666 }
667
668 /******************************************************************************/
669
670 /* calculate dest = a**exp % r */
671 /* using square&multiply */
672 void bigint_expmod_u(bigint_t* dest, const bigint_t* a, const bigint_t* exp, const bigint_t* r){
673         if(a->length_B==0 || r->length_B==0){
674                 return;
675         }
676
677         bigint_t res, base;
678         bigint_word_t base_b[MAX(a->length_B,r->length_B*2)], res_b[r->length_B*2];
679         uint16_t i;
680         uint8_t j, t;
681         res.wordv = res_b;
682         base.wordv = base_b;
683         bigint_copy(&base, a);
684         bigint_reduce(&base, r);
685         res.wordv[0]=1;
686         res.length_B=1;
687         res.info = 0;
688         bigint_adjust(&res);
689         for(i=0; i+1<exp->length_B; ++i){
690                 t=exp->wordv[i];
691                 for(j=0; j<8; ++j){
692                         if(t&1){
693                                 bigint_mul_u(&res, &res, &base);
694                                 bigint_reduce(&res, r);
695                         }
696                         bigint_square(&base, &base);
697                         bigint_reduce(&base, r);
698                         t>>=1;
699                 }
700         }
701         t=exp->wordv[i];
702         while(t){
703                 if(t&1){
704                         bigint_mul_u(&res, &res, &base);
705                         bigint_reduce(&res, r);
706                 }
707                 bigint_square(&base, &base);
708                 bigint_reduce(&base, r);
709                 t>>=1;
710         }
711         SET_POS(&res);
712         bigint_copy(dest, &res);
713 }
714
715 /******************************************************************************/
716 /* gcd <-- gcd(x,y) a*x+b*y=gcd */
717 void bigint_gcdext(bigint_t* gcd, bigint_t* a, bigint_t* b, const bigint_t* x, const bigint_t* y){
718          bigint_t g, x_, y_, u, v, a_, b_, c_, d_;
719          volatile uint16_t i=0;
720          if(x->length_B==0 || y->length_B==0){
721                  return;
722          }
723          while(x->wordv[i]==0 && y->wordv[i]==0){
724                  ++i;
725          }
726          bigint_word_t g_b[i+2], x_b[x->length_B-i], y_b[y->length_B-i];
727          bigint_word_t u_b[x->length_B-i], v_b[y->length_B-i];
728          bigint_word_t a_b[y->length_B+2], c_b[y->length_B+2];
729          bigint_word_t b_b[x->length_B+2], d_b[x->length_B+2];
730
731          g.wordv = g_b;
732          x_.wordv = x_b;
733          y_.wordv = y_b;
734          memset(g_b, 0, i);
735          g_b[i]=1;
736          g.length_B = i+1;
737          g.info=0;
738          x_.info = y_.info = 0;
739          x_.length_B = x->length_B-i;
740          y_.length_B = y->length_B-i;
741          memcpy(x_.wordv, x->wordv+i, x_.length_B);
742          memcpy(y_.wordv, y->wordv+i, y_.length_B);
743          for(i=0; (x_.wordv[0]&(1<<i))==0 && (y_.wordv[0]&(1<<i))==0; ++i){
744          }
745
746          bigint_adjust(&x_);
747          bigint_adjust(&y_);
748
749          if(i){
750                  bigint_shiftleft(&g, i);
751                  bigint_shiftright(&x_, i);
752                  bigint_shiftright(&y_, i);
753          }
754          u.wordv = u_b;
755          v.wordv = v_b;
756          a_.wordv = a_b;
757          b_.wordv = b_b;
758          c_.wordv = c_b;
759          d_.wordv = d_b;
760
761          bigint_copy(&u, &x_);
762          bigint_copy(&v, &y_);
763          a_.wordv[0] = 1;
764          a_.length_B = 1;
765          a_.info = 0;
766          d_.wordv[0] = 1;
767          d_.length_B = 1;
768          d_.info = 0;
769          bigint_set_zero(&b_);
770          bigint_set_zero(&c_);
771          do{
772                  while((u.wordv[0]&1)==0){
773                          bigint_shiftright(&u, 1);
774                          if((a_.wordv[0]&1) || (b_.wordv[0]&1)){
775                                  bigint_add_s(&a_, &a_, &y_);
776                                  bigint_sub_s(&b_, &b_, &x_);
777                          }
778                          bigint_shiftright(&a_, 1);
779                          bigint_shiftright(&b_, 1);
780                  }
781                  while((v.wordv[0]&1)==0){
782                          bigint_shiftright(&v, 1);
783                          if((c_.wordv[0]&1) || (d_.wordv[0]&1)){
784                                  bigint_add_s(&c_, &c_, &y_);
785                                  bigint_sub_s(&d_, &d_, &x_);
786                          }
787                          bigint_shiftright(&c_, 1);
788                          bigint_shiftright(&d_, 1);
789
790                  }
791                  if(bigint_cmp_u(&u, &v)>=0){
792                         bigint_sub_u(&u, &u, &v);
793                         bigint_sub_s(&a_, &a_, &c_);
794                         bigint_sub_s(&b_, &b_, &d_);
795                  }else{
796                         bigint_sub_u(&v, &v, &u);
797                         bigint_sub_s(&c_, &c_, &a_);
798                         bigint_sub_s(&d_, &d_, &b_);
799                  }
800          }while(u.length_B);
801          if(gcd){
802                  bigint_mul_s(gcd, &v, &g);
803          }
804          if(a){
805                 bigint_copy(a, &c_);
806          }
807          if(b){
808                  bigint_copy(b, &d_);
809          }
810 }
811
812 /******************************************************************************/
813
814 void bigint_inverse(bigint_t* dest, const bigint_t* a, const bigint_t* m){
815         bigint_gcdext(NULL, dest, NULL, a, m);
816         while(dest->info&BIGINT_NEG_MASK){
817                 bigint_add_s(dest, dest, m);
818         }
819 }
820
821 /******************************************************************************/
822
823 void bigint_changeendianess(bigint_t* a){
824         uint8_t t, *p, *q;
825         p = (uint8_t*)(a->wordv);
826         q = ((uint8_t*)p)+a->length_B*sizeof(bigint_word_t)-1;
827         while(p<q){
828                 t = *p;
829                 *p = *q;
830                 *q = t;
831                 ++p; --q;
832         }
833 }
834
835 /******************************************************************************/
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858