PolarSSL v1.2.5
bignum.c
Go to the documentation of this file.
1 /*
2  * Multi-precision integer library
3  *
4  * Copyright (C) 2006-2010, Brainspark B.V.
5  *
6  * This file is part of PolarSSL (http://www.polarssl.org)
7  * Lead Maintainer: Paul Bakker <polarssl_maintainer at polarssl.org>
8  *
9  * All rights reserved.
10  *
11  * This program is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License as published by
13  * the Free Software Foundation; either version 2 of the License, or
14  * (at your option) any later version.
15  *
16  * This program is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19  * GNU General Public License for more details.
20  *
21  * You should have received a copy of the GNU General Public License along
22  * with this program; if not, write to the Free Software Foundation, Inc.,
23  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24  */
25 /*
26  * This MPI implementation is based on:
27  *
28  * http://www.cacr.math.uwaterloo.ca/hac/about/chap14.pdf
29  * http://www.stillhq.com/extracted/gnupg-api/mpi/
30  * http://math.libtomcrypt.com/files/tommath.pdf
31  */
32 
33 #include "polarssl/config.h"
34 
35 #if defined(POLARSSL_BIGNUM_C)
36 
37 #include "polarssl/bignum.h"
38 #include "polarssl/bn_mul.h"
39 
40 #include <stdlib.h>
41 
42 #define ciL (sizeof(t_uint)) /* chars in limb */
43 #define biL (ciL << 3) /* bits in limb */
44 #define biH (ciL << 2) /* half limb size */
45 
46 /*
47  * Convert between bits/chars and number of limbs
48  */
49 #define BITS_TO_LIMBS(i) (((i) + biL - 1) / biL)
50 #define CHARS_TO_LIMBS(i) (((i) + ciL - 1) / ciL)
51 
52 /*
53  * Initialize one MPI
54  */
55 void mpi_init( mpi *X )
56 {
57  if( X == NULL )
58  return;
59 
60  X->s = 1;
61  X->n = 0;
62  X->p = NULL;
63 }
64 
65 /*
66  * Unallocate one MPI
67  */
68 void mpi_free( mpi *X )
69 {
70  if( X == NULL )
71  return;
72 
73  if( X->p != NULL )
74  {
75  memset( X->p, 0, X->n * ciL );
76  free( X->p );
77  }
78 
79  X->s = 1;
80  X->n = 0;
81  X->p = NULL;
82 }
83 
84 /*
85  * Enlarge to the specified number of limbs
86  */
87 int mpi_grow( mpi *X, size_t nblimbs )
88 {
89  t_uint *p;
90 
91  if( nblimbs > POLARSSL_MPI_MAX_LIMBS )
93 
94  if( X->n < nblimbs )
95  {
96  if( ( p = (t_uint *) malloc( nblimbs * ciL ) ) == NULL )
98 
99  memset( p, 0, nblimbs * ciL );
100 
101  if( X->p != NULL )
102  {
103  memcpy( p, X->p, X->n * ciL );
104  memset( X->p, 0, X->n * ciL );
105  free( X->p );
106  }
107 
108  X->n = nblimbs;
109  X->p = p;
110  }
111 
112  return( 0 );
113 }
114 
115 /*
116  * Copy the contents of Y into X
117  */
118 int mpi_copy( mpi *X, const mpi *Y )
119 {
120  int ret;
121  size_t i;
122 
123  if( X == Y )
124  return( 0 );
125 
126  for( i = Y->n - 1; i > 0; i-- )
127  if( Y->p[i] != 0 )
128  break;
129  i++;
130 
131  X->s = Y->s;
132 
133  MPI_CHK( mpi_grow( X, i ) );
134 
135  memset( X->p, 0, X->n * ciL );
136  memcpy( X->p, Y->p, i * ciL );
137 
138 cleanup:
139 
140  return( ret );
141 }
142 
143 /*
144  * Swap the contents of X and Y
145  */
146 void mpi_swap( mpi *X, mpi *Y )
147 {
148  mpi T;
149 
150  memcpy( &T, X, sizeof( mpi ) );
151  memcpy( X, Y, sizeof( mpi ) );
152  memcpy( Y, &T, sizeof( mpi ) );
153 }
154 
155 /*
156  * Set value from integer
157  */
158 int mpi_lset( mpi *X, t_sint z )
159 {
160  int ret;
161 
162  MPI_CHK( mpi_grow( X, 1 ) );
163  memset( X->p, 0, X->n * ciL );
164 
165  X->p[0] = ( z < 0 ) ? -z : z;
166  X->s = ( z < 0 ) ? -1 : 1;
167 
168 cleanup:
169 
170  return( ret );
171 }
172 
173 /*
174  * Get a specific bit
175  */
176 int mpi_get_bit( const mpi *X, size_t pos )
177 {
178  if( X->n * biL <= pos )
179  return( 0 );
180 
181  return ( X->p[pos / biL] >> ( pos % biL ) ) & 0x01;
182 }
183 
184 /*
185  * Set a bit to a specific value of 0 or 1
186  */
187 int mpi_set_bit( mpi *X, size_t pos, unsigned char val )
188 {
189  int ret = 0;
190  size_t off = pos / biL;
191  size_t idx = pos % biL;
192 
193  if( val != 0 && val != 1 )
195 
196  if( X->n * biL <= pos )
197  {
198  if( val == 0 )
199  return ( 0 );
200 
201  MPI_CHK( mpi_grow( X, off + 1 ) );
202  }
203 
204  X->p[off] = ( X->p[off] & ~( 0x01 << idx ) ) | ( val << idx );
205 
206 cleanup:
207 
208  return( ret );
209 }
210 
211 /*
212  * Return the number of least significant bits
213  */
214 size_t mpi_lsb( const mpi *X )
215 {
216  size_t i, j, count = 0;
217 
218  for( i = 0; i < X->n; i++ )
219  for( j = 0; j < biL; j++, count++ )
220  if( ( ( X->p[i] >> j ) & 1 ) != 0 )
221  return( count );
222 
223  return( 0 );
224 }
225 
226 /*
227  * Return the number of most significant bits
228  */
229 size_t mpi_msb( const mpi *X )
230 {
231  size_t i, j;
232 
233  for( i = X->n - 1; i > 0; i-- )
234  if( X->p[i] != 0 )
235  break;
236 
237  for( j = biL; j > 0; j-- )
238  if( ( ( X->p[i] >> ( j - 1 ) ) & 1 ) != 0 )
239  break;
240 
241  return( ( i * biL ) + j );
242 }
243 
244 /*
245  * Return the total size in bytes
246  */
247 size_t mpi_size( const mpi *X )
248 {
249  return( ( mpi_msb( X ) + 7 ) >> 3 );
250 }
251 
252 /*
253  * Convert an ASCII character to digit value
254  */
255 static int mpi_get_digit( t_uint *d, int radix, char c )
256 {
257  *d = 255;
258 
259  if( c >= 0x30 && c <= 0x39 ) *d = c - 0x30;
260  if( c >= 0x41 && c <= 0x46 ) *d = c - 0x37;
261  if( c >= 0x61 && c <= 0x66 ) *d = c - 0x57;
262 
263  if( *d >= (t_uint) radix )
265 
266  return( 0 );
267 }
268 
269 /*
270  * Import from an ASCII string
271  */
272 int mpi_read_string( mpi *X, int radix, const char *s )
273 {
274  int ret;
275  size_t i, j, slen, n;
276  t_uint d;
277  mpi T;
278 
279  if( radix < 2 || radix > 16 )
281 
282  mpi_init( &T );
283 
284  slen = strlen( s );
285 
286  if( radix == 16 )
287  {
288  n = BITS_TO_LIMBS( slen << 2 );
289 
290  MPI_CHK( mpi_grow( X, n ) );
291  MPI_CHK( mpi_lset( X, 0 ) );
292 
293  for( i = slen, j = 0; i > 0; i--, j++ )
294  {
295  if( i == 1 && s[i - 1] == '-' )
296  {
297  X->s = -1;
298  break;
299  }
300 
301  MPI_CHK( mpi_get_digit( &d, radix, s[i - 1] ) );
302  X->p[j / (2 * ciL)] |= d << ( (j % (2 * ciL)) << 2 );
303  }
304  }
305  else
306  {
307  MPI_CHK( mpi_lset( X, 0 ) );
308 
309  for( i = 0; i < slen; i++ )
310  {
311  if( i == 0 && s[i] == '-' )
312  {
313  X->s = -1;
314  continue;
315  }
316 
317  MPI_CHK( mpi_get_digit( &d, radix, s[i] ) );
318  MPI_CHK( mpi_mul_int( &T, X, radix ) );
319 
320  if( X->s == 1 )
321  {
322  MPI_CHK( mpi_add_int( X, &T, d ) );
323  }
324  else
325  {
326  MPI_CHK( mpi_sub_int( X, &T, d ) );
327  }
328  }
329  }
330 
331 cleanup:
332 
333  mpi_free( &T );
334 
335  return( ret );
336 }
337 
338 /*
339  * Helper to write the digits high-order first
340  */
341 static int mpi_write_hlp( mpi *X, int radix, char **p )
342 {
343  int ret;
344  t_uint r;
345 
346  if( radix < 2 || radix > 16 )
348 
349  MPI_CHK( mpi_mod_int( &r, X, radix ) );
350  MPI_CHK( mpi_div_int( X, NULL, X, radix ) );
351 
352  if( mpi_cmp_int( X, 0 ) != 0 )
353  MPI_CHK( mpi_write_hlp( X, radix, p ) );
354 
355  if( r < 10 )
356  *(*p)++ = (char)( r + 0x30 );
357  else
358  *(*p)++ = (char)( r + 0x37 );
359 
360 cleanup:
361 
362  return( ret );
363 }
364 
365 /*
366  * Export into an ASCII string
367  */
368 int mpi_write_string( const mpi *X, int radix, char *s, size_t *slen )
369 {
370  int ret = 0;
371  size_t n;
372  char *p;
373  mpi T;
374 
375  if( radix < 2 || radix > 16 )
377 
378  n = mpi_msb( X );
379  if( radix >= 4 ) n >>= 1;
380  if( radix >= 16 ) n >>= 1;
381  n += 3;
382 
383  if( *slen < n )
384  {
385  *slen = n;
387  }
388 
389  p = s;
390  mpi_init( &T );
391 
392  if( X->s == -1 )
393  *p++ = '-';
394 
395  if( radix == 16 )
396  {
397  int c;
398  size_t i, j, k;
399 
400  for( i = X->n, k = 0; i > 0; i-- )
401  {
402  for( j = ciL; j > 0; j-- )
403  {
404  c = ( X->p[i - 1] >> ( ( j - 1 ) << 3) ) & 0xFF;
405 
406  if( c == 0 && k == 0 && ( i + j + 3 ) != 0 )
407  continue;
408 
409  *(p++) = "0123456789ABCDEF" [c / 16];
410  *(p++) = "0123456789ABCDEF" [c % 16];
411  k = 1;
412  }
413  }
414  }
415  else
416  {
417  MPI_CHK( mpi_copy( &T, X ) );
418 
419  if( T.s == -1 )
420  T.s = 1;
421 
422  MPI_CHK( mpi_write_hlp( &T, radix, &p ) );
423  }
424 
425  *p++ = '\0';
426  *slen = p - s;
427 
428 cleanup:
429 
430  mpi_free( &T );
431 
432  return( ret );
433 }
434 
435 #if defined(POLARSSL_FS_IO)
436 /*
437  * Read X from an opened file
438  */
439 int mpi_read_file( mpi *X, int radix, FILE *fin )
440 {
441  t_uint d;
442  size_t slen;
443  char *p;
444  /*
445  * Buffer should have space for (short) label and decimal formatted MPI,
446  * newline characters and '\0'
447  */
448  char s[ POLARSSL_MPI_RW_BUFFER_SIZE ];
449 
450  memset( s, 0, sizeof( s ) );
451  if( fgets( s, sizeof( s ) - 1, fin ) == NULL )
453 
454  slen = strlen( s );
455  if( slen == sizeof( s ) - 2 )
457 
458  if( s[slen - 1] == '\n' ) { slen--; s[slen] = '\0'; }
459  if( s[slen - 1] == '\r' ) { slen--; s[slen] = '\0'; }
460 
461  p = s + slen;
462  while( --p >= s )
463  if( mpi_get_digit( &d, radix, *p ) != 0 )
464  break;
465 
466  return( mpi_read_string( X, radix, p + 1 ) );
467 }
468 
469 /*
470  * Write X into an opened file (or stdout if fout == NULL)
471  */
472 int mpi_write_file( const char *p, const mpi *X, int radix, FILE *fout )
473 {
474  int ret;
475  size_t n, slen, plen;
476  /*
477  * Buffer should have space for (short) label and decimal formatted MPI,
478  * newline characters and '\0'
479  */
480  char s[ POLARSSL_MPI_RW_BUFFER_SIZE ];
481 
482  n = sizeof( s );
483  memset( s, 0, n );
484  n -= 2;
485 
486  MPI_CHK( mpi_write_string( X, radix, s, (size_t *) &n ) );
487 
488  if( p == NULL ) p = "";
489 
490  plen = strlen( p );
491  slen = strlen( s );
492  s[slen++] = '\r';
493  s[slen++] = '\n';
494 
495  if( fout != NULL )
496  {
497  if( fwrite( p, 1, plen, fout ) != plen ||
498  fwrite( s, 1, slen, fout ) != slen )
500  }
501  else
502  printf( "%s%s", p, s );
503 
504 cleanup:
505 
506  return( ret );
507 }
508 #endif /* POLARSSL_FS_IO */
509 
510 /*
511  * Import X from unsigned binary data, big endian
512  */
513 int mpi_read_binary( mpi *X, const unsigned char *buf, size_t buflen )
514 {
515  int ret;
516  size_t i, j, n;
517 
518  for( n = 0; n < buflen; n++ )
519  if( buf[n] != 0 )
520  break;
521 
522  MPI_CHK( mpi_grow( X, CHARS_TO_LIMBS( buflen - n ) ) );
523  MPI_CHK( mpi_lset( X, 0 ) );
524 
525  for( i = buflen, j = 0; i > n; i--, j++ )
526  X->p[j / ciL] |= ((t_uint) buf[i - 1]) << ((j % ciL) << 3);
527 
528 cleanup:
529 
530  return( ret );
531 }
532 
533 /*
534  * Export X into unsigned binary data, big endian
535  */
536 int mpi_write_binary( const mpi *X, unsigned char *buf, size_t buflen )
537 {
538  size_t i, j, n;
539 
540  n = mpi_size( X );
541 
542  if( buflen < n )
544 
545  memset( buf, 0, buflen );
546 
547  for( i = buflen - 1, j = 0; n > 0; i--, j++, n-- )
548  buf[i] = (unsigned char)( X->p[j / ciL] >> ((j % ciL) << 3) );
549 
550  return( 0 );
551 }
552 
553 /*
554  * Left-shift: X <<= count
555  */
556 int mpi_shift_l( mpi *X, size_t count )
557 {
558  int ret;
559  size_t i, v0, t1;
560  t_uint r0 = 0, r1;
561 
562  v0 = count / (biL );
563  t1 = count & (biL - 1);
564 
565  i = mpi_msb( X ) + count;
566 
567  if( X->n * biL < i )
568  MPI_CHK( mpi_grow( X, BITS_TO_LIMBS( i ) ) );
569 
570  ret = 0;
571 
572  /*
573  * shift by count / limb_size
574  */
575  if( v0 > 0 )
576  {
577  for( i = X->n; i > v0; i-- )
578  X->p[i - 1] = X->p[i - v0 - 1];
579 
580  for( ; i > 0; i-- )
581  X->p[i - 1] = 0;
582  }
583 
584  /*
585  * shift by count % limb_size
586  */
587  if( t1 > 0 )
588  {
589  for( i = v0; i < X->n; i++ )
590  {
591  r1 = X->p[i] >> (biL - t1);
592  X->p[i] <<= t1;
593  X->p[i] |= r0;
594  r0 = r1;
595  }
596  }
597 
598 cleanup:
599 
600  return( ret );
601 }
602 
603 /*
604  * Right-shift: X >>= count
605  */
606 int mpi_shift_r( mpi *X, size_t count )
607 {
608  size_t i, v0, v1;
609  t_uint r0 = 0, r1;
610 
611  v0 = count / biL;
612  v1 = count & (biL - 1);
613 
614  if( v0 > X->n || ( v0 == X->n && v1 > 0 ) )
615  return mpi_lset( X, 0 );
616 
617  /*
618  * shift by count / limb_size
619  */
620  if( v0 > 0 )
621  {
622  for( i = 0; i < X->n - v0; i++ )
623  X->p[i] = X->p[i + v0];
624 
625  for( ; i < X->n; i++ )
626  X->p[i] = 0;
627  }
628 
629  /*
630  * shift by count % limb_size
631  */
632  if( v1 > 0 )
633  {
634  for( i = X->n; i > 0; i-- )
635  {
636  r1 = X->p[i - 1] << (biL - v1);
637  X->p[i - 1] >>= v1;
638  X->p[i - 1] |= r0;
639  r0 = r1;
640  }
641  }
642 
643  return( 0 );
644 }
645 
646 /*
647  * Compare unsigned values
648  */
649 int mpi_cmp_abs( const mpi *X, const mpi *Y )
650 {
651  size_t i, j;
652 
653  for( i = X->n; i > 0; i-- )
654  if( X->p[i - 1] != 0 )
655  break;
656 
657  for( j = Y->n; j > 0; j-- )
658  if( Y->p[j - 1] != 0 )
659  break;
660 
661  if( i == 0 && j == 0 )
662  return( 0 );
663 
664  if( i > j ) return( 1 );
665  if( j > i ) return( -1 );
666 
667  for( ; i > 0; i-- )
668  {
669  if( X->p[i - 1] > Y->p[i - 1] ) return( 1 );
670  if( X->p[i - 1] < Y->p[i - 1] ) return( -1 );
671  }
672 
673  return( 0 );
674 }
675 
676 /*
677  * Compare signed values
678  */
679 int mpi_cmp_mpi( const mpi *X, const mpi *Y )
680 {
681  size_t i, j;
682 
683  for( i = X->n; i > 0; i-- )
684  if( X->p[i - 1] != 0 )
685  break;
686 
687  for( j = Y->n; j > 0; j-- )
688  if( Y->p[j - 1] != 0 )
689  break;
690 
691  if( i == 0 && j == 0 )
692  return( 0 );
693 
694  if( i > j ) return( X->s );
695  if( j > i ) return( -Y->s );
696 
697  if( X->s > 0 && Y->s < 0 ) return( 1 );
698  if( Y->s > 0 && X->s < 0 ) return( -1 );
699 
700  for( ; i > 0; i-- )
701  {
702  if( X->p[i - 1] > Y->p[i - 1] ) return( X->s );
703  if( X->p[i - 1] < Y->p[i - 1] ) return( -X->s );
704  }
705 
706  return( 0 );
707 }
708 
709 /*
710  * Compare signed values
711  */
712 int mpi_cmp_int( const mpi *X, t_sint z )
713 {
714  mpi Y;
715  t_uint p[1];
716 
717  *p = ( z < 0 ) ? -z : z;
718  Y.s = ( z < 0 ) ? -1 : 1;
719  Y.n = 1;
720  Y.p = p;
721 
722  return( mpi_cmp_mpi( X, &Y ) );
723 }
724 
725 /*
726  * Unsigned addition: X = |A| + |B| (HAC 14.7)
727  */
728 int mpi_add_abs( mpi *X, const mpi *A, const mpi *B )
729 {
730  int ret;
731  size_t i, j;
732  t_uint *o, *p, c;
733 
734  if( X == B )
735  {
736  const mpi *T = A; A = X; B = T;
737  }
738 
739  if( X != A )
740  MPI_CHK( mpi_copy( X, A ) );
741 
742  /*
743  * X should always be positive as a result of unsigned additions.
744  */
745  X->s = 1;
746 
747  for( j = B->n; j > 0; j-- )
748  if( B->p[j - 1] != 0 )
749  break;
750 
751  MPI_CHK( mpi_grow( X, j ) );
752 
753  o = B->p; p = X->p; c = 0;
754 
755  for( i = 0; i < j; i++, o++, p++ )
756  {
757  *p += c; c = ( *p < c );
758  *p += *o; c += ( *p < *o );
759  }
760 
761  while( c != 0 )
762  {
763  if( i >= X->n )
764  {
765  MPI_CHK( mpi_grow( X, i + 1 ) );
766  p = X->p + i;
767  }
768 
769  *p += c; c = ( *p < c ); i++; p++;
770  }
771 
772 cleanup:
773 
774  return( ret );
775 }
776 
777 /*
778  * Helper for mpi substraction
779  */
780 static void mpi_sub_hlp( size_t n, t_uint *s, t_uint *d )
781 {
782  size_t i;
783  t_uint c, z;
784 
785  for( i = c = 0; i < n; i++, s++, d++ )
786  {
787  z = ( *d < c ); *d -= c;
788  c = ( *d < *s ) + z; *d -= *s;
789  }
790 
791  while( c != 0 )
792  {
793  z = ( *d < c ); *d -= c;
794  c = z; i++; d++;
795  }
796 }
797 
798 /*
799  * Unsigned substraction: X = |A| - |B| (HAC 14.9)
800  */
801 int mpi_sub_abs( mpi *X, const mpi *A, const mpi *B )
802 {
803  mpi TB;
804  int ret;
805  size_t n;
806 
807  if( mpi_cmp_abs( A, B ) < 0 )
809 
810  mpi_init( &TB );
811 
812  if( X == B )
813  {
814  MPI_CHK( mpi_copy( &TB, B ) );
815  B = &TB;
816  }
817 
818  if( X != A )
819  MPI_CHK( mpi_copy( X, A ) );
820 
821  /*
822  * X should always be positive as a result of unsigned substractions.
823  */
824  X->s = 1;
825 
826  ret = 0;
827 
828  for( n = B->n; n > 0; n-- )
829  if( B->p[n - 1] != 0 )
830  break;
831 
832  mpi_sub_hlp( n, B->p, X->p );
833 
834 cleanup:
835 
836  mpi_free( &TB );
837 
838  return( ret );
839 }
840 
841 /*
842  * Signed addition: X = A + B
843  */
844 int mpi_add_mpi( mpi *X, const mpi *A, const mpi *B )
845 {
846  int ret, s = A->s;
847 
848  if( A->s * B->s < 0 )
849  {
850  if( mpi_cmp_abs( A, B ) >= 0 )
851  {
852  MPI_CHK( mpi_sub_abs( X, A, B ) );
853  X->s = s;
854  }
855  else
856  {
857  MPI_CHK( mpi_sub_abs( X, B, A ) );
858  X->s = -s;
859  }
860  }
861  else
862  {
863  MPI_CHK( mpi_add_abs( X, A, B ) );
864  X->s = s;
865  }
866 
867 cleanup:
868 
869  return( ret );
870 }
871 
872 /*
873  * Signed substraction: X = A - B
874  */
875 int mpi_sub_mpi( mpi *X, const mpi *A, const mpi *B )
876 {
877  int ret, s = A->s;
878 
879  if( A->s * B->s > 0 )
880  {
881  if( mpi_cmp_abs( A, B ) >= 0 )
882  {
883  MPI_CHK( mpi_sub_abs( X, A, B ) );
884  X->s = s;
885  }
886  else
887  {
888  MPI_CHK( mpi_sub_abs( X, B, A ) );
889  X->s = -s;
890  }
891  }
892  else
893  {
894  MPI_CHK( mpi_add_abs( X, A, B ) );
895  X->s = s;
896  }
897 
898 cleanup:
899 
900  return( ret );
901 }
902 
903 /*
904  * Signed addition: X = A + b
905  */
906 int mpi_add_int( mpi *X, const mpi *A, t_sint b )
907 {
908  mpi _B;
909  t_uint p[1];
910 
911  p[0] = ( b < 0 ) ? -b : b;
912  _B.s = ( b < 0 ) ? -1 : 1;
913  _B.n = 1;
914  _B.p = p;
915 
916  return( mpi_add_mpi( X, A, &_B ) );
917 }
918 
919 /*
920  * Signed substraction: X = A - b
921  */
922 int mpi_sub_int( mpi *X, const mpi *A, t_sint b )
923 {
924  mpi _B;
925  t_uint p[1];
926 
927  p[0] = ( b < 0 ) ? -b : b;
928  _B.s = ( b < 0 ) ? -1 : 1;
929  _B.n = 1;
930  _B.p = p;
931 
932  return( mpi_sub_mpi( X, A, &_B ) );
933 }
934 
935 /*
936  * Helper for mpi multiplication
937  */
938 static void mpi_mul_hlp( size_t i, t_uint *s, t_uint *d, t_uint b )
939 {
940  t_uint c = 0, t = 0;
941 
942 #if defined(MULADDC_HUIT)
943  for( ; i >= 8; i -= 8 )
944  {
946  MULADDC_HUIT
948  }
949 
950  for( ; i > 0; i-- )
951  {
955  }
956 #else
957  for( ; i >= 16; i -= 16 )
958  {
964 
970  }
971 
972  for( ; i >= 8; i -= 8 )
973  {
977 
981  }
982 
983  for( ; i > 0; i-- )
984  {
988  }
989 #endif
990 
991  t++;
992 
993  do {
994  *d += c; c = ( *d < c ); d++;
995  }
996  while( c != 0 );
997 }
998 
999 /*
1000  * Baseline multiplication: X = A * B (HAC 14.12)
1001  */
1002 int mpi_mul_mpi( mpi *X, const mpi *A, const mpi *B )
1003 {
1004  int ret;
1005  size_t i, j;
1006  mpi TA, TB;
1007 
1008  mpi_init( &TA ); mpi_init( &TB );
1009 
1010  if( X == A ) { MPI_CHK( mpi_copy( &TA, A ) ); A = &TA; }
1011  if( X == B ) { MPI_CHK( mpi_copy( &TB, B ) ); B = &TB; }
1012 
1013  for( i = A->n; i > 0; i-- )
1014  if( A->p[i - 1] != 0 )
1015  break;
1016 
1017  for( j = B->n; j > 0; j-- )
1018  if( B->p[j - 1] != 0 )
1019  break;
1020 
1021  MPI_CHK( mpi_grow( X, i + j ) );
1022  MPI_CHK( mpi_lset( X, 0 ) );
1023 
1024  for( i++; j > 0; j-- )
1025  mpi_mul_hlp( i - 1, A->p, X->p + j - 1, B->p[j - 1] );
1026 
1027  X->s = A->s * B->s;
1028 
1029 cleanup:
1030 
1031  mpi_free( &TB ); mpi_free( &TA );
1032 
1033  return( ret );
1034 }
1035 
1036 /*
1037  * Baseline multiplication: X = A * b
1038  */
1039 int mpi_mul_int( mpi *X, const mpi *A, t_sint b )
1040 {
1041  mpi _B;
1042  t_uint p[1];
1043 
1044  _B.s = 1;
1045  _B.n = 1;
1046  _B.p = p;
1047  p[0] = b;
1048 
1049  return( mpi_mul_mpi( X, A, &_B ) );
1050 }
1051 
1052 /*
1053  * Division by mpi: A = Q * B + R (HAC 14.20)
1054  */
1055 int mpi_div_mpi( mpi *Q, mpi *R, const mpi *A, const mpi *B )
1056 {
1057  int ret;
1058  size_t i, n, t, k;
1059  mpi X, Y, Z, T1, T2;
1060 
1061  if( mpi_cmp_int( B, 0 ) == 0 )
1063 
1064  mpi_init( &X ); mpi_init( &Y ); mpi_init( &Z );
1065  mpi_init( &T1 ); mpi_init( &T2 );
1066 
1067  if( mpi_cmp_abs( A, B ) < 0 )
1068  {
1069  if( Q != NULL ) MPI_CHK( mpi_lset( Q, 0 ) );
1070  if( R != NULL ) MPI_CHK( mpi_copy( R, A ) );
1071  return( 0 );
1072  }
1073 
1074  MPI_CHK( mpi_copy( &X, A ) );
1075  MPI_CHK( mpi_copy( &Y, B ) );
1076  X.s = Y.s = 1;
1077 
1078  MPI_CHK( mpi_grow( &Z, A->n + 2 ) );
1079  MPI_CHK( mpi_lset( &Z, 0 ) );
1080  MPI_CHK( mpi_grow( &T1, 2 ) );
1081  MPI_CHK( mpi_grow( &T2, 3 ) );
1082 
1083  k = mpi_msb( &Y ) % biL;
1084  if( k < biL - 1 )
1085  {
1086  k = biL - 1 - k;
1087  MPI_CHK( mpi_shift_l( &X, k ) );
1088  MPI_CHK( mpi_shift_l( &Y, k ) );
1089  }
1090  else k = 0;
1091 
1092  n = X.n - 1;
1093  t = Y.n - 1;
1094  MPI_CHK( mpi_shift_l( &Y, biL * (n - t) ) );
1095 
1096  while( mpi_cmp_mpi( &X, &Y ) >= 0 )
1097  {
1098  Z.p[n - t]++;
1099  mpi_sub_mpi( &X, &X, &Y );
1100  }
1101  mpi_shift_r( &Y, biL * (n - t) );
1102 
1103  for( i = n; i > t ; i-- )
1104  {
1105  if( X.p[i] >= Y.p[t] )
1106  Z.p[i - t - 1] = ~0;
1107  else
1108  {
1109 #if defined(POLARSSL_HAVE_UDBL)
1110  t_udbl r;
1111 
1112  r = (t_udbl) X.p[i] << biL;
1113  r |= (t_udbl) X.p[i - 1];
1114  r /= Y.p[t];
1115  if( r > ((t_udbl) 1 << biL) - 1)
1116  r = ((t_udbl) 1 << biL) - 1;
1117 
1118  Z.p[i - t - 1] = (t_uint) r;
1119 #else
1120  /*
1121  * __udiv_qrnnd_c, from gmp/longlong.h
1122  */
1123  t_uint q0, q1, r0, r1;
1124  t_uint d0, d1, d, m;
1125 
1126  d = Y.p[t];
1127  d0 = ( d << biH ) >> biH;
1128  d1 = ( d >> biH );
1129 
1130  q1 = X.p[i] / d1;
1131  r1 = X.p[i] - d1 * q1;
1132  r1 <<= biH;
1133  r1 |= ( X.p[i - 1] >> biH );
1134 
1135  m = q1 * d0;
1136  if( r1 < m )
1137  {
1138  q1--, r1 += d;
1139  while( r1 >= d && r1 < m )
1140  q1--, r1 += d;
1141  }
1142  r1 -= m;
1143 
1144  q0 = r1 / d1;
1145  r0 = r1 - d1 * q0;
1146  r0 <<= biH;
1147  r0 |= ( X.p[i - 1] << biH ) >> biH;
1148 
1149  m = q0 * d0;
1150  if( r0 < m )
1151  {
1152  q0--, r0 += d;
1153  while( r0 >= d && r0 < m )
1154  q0--, r0 += d;
1155  }
1156  r0 -= m;
1157 
1158  Z.p[i - t - 1] = ( q1 << biH ) | q0;
1159 #endif
1160  }
1161 
1162  Z.p[i - t - 1]++;
1163  do
1164  {
1165  Z.p[i - t - 1]--;
1166 
1167  MPI_CHK( mpi_lset( &T1, 0 ) );
1168  T1.p[0] = (t < 1) ? 0 : Y.p[t - 1];
1169  T1.p[1] = Y.p[t];
1170  MPI_CHK( mpi_mul_int( &T1, &T1, Z.p[i - t - 1] ) );
1171 
1172  MPI_CHK( mpi_lset( &T2, 0 ) );
1173  T2.p[0] = (i < 2) ? 0 : X.p[i - 2];
1174  T2.p[1] = (i < 1) ? 0 : X.p[i - 1];
1175  T2.p[2] = X.p[i];
1176  }
1177  while( mpi_cmp_mpi( &T1, &T2 ) > 0 );
1178 
1179  MPI_CHK( mpi_mul_int( &T1, &Y, Z.p[i - t - 1] ) );
1180  MPI_CHK( mpi_shift_l( &T1, biL * (i - t - 1) ) );
1181  MPI_CHK( mpi_sub_mpi( &X, &X, &T1 ) );
1182 
1183  if( mpi_cmp_int( &X, 0 ) < 0 )
1184  {
1185  MPI_CHK( mpi_copy( &T1, &Y ) );
1186  MPI_CHK( mpi_shift_l( &T1, biL * (i - t - 1) ) );
1187  MPI_CHK( mpi_add_mpi( &X, &X, &T1 ) );
1188  Z.p[i - t - 1]--;
1189  }
1190  }
1191 
1192  if( Q != NULL )
1193  {
1194  mpi_copy( Q, &Z );
1195  Q->s = A->s * B->s;
1196  }
1197 
1198  if( R != NULL )
1199  {
1200  mpi_shift_r( &X, k );
1201  X.s = A->s;
1202  mpi_copy( R, &X );
1203 
1204  if( mpi_cmp_int( R, 0 ) == 0 )
1205  R->s = 1;
1206  }
1207 
1208 cleanup:
1209 
1210  mpi_free( &X ); mpi_free( &Y ); mpi_free( &Z );
1211  mpi_free( &T1 ); mpi_free( &T2 );
1212 
1213  return( ret );
1214 }
1215 
1216 /*
1217  * Division by int: A = Q * b + R
1218  */
1219 int mpi_div_int( mpi *Q, mpi *R, const mpi *A, t_sint b )
1220 {
1221  mpi _B;
1222  t_uint p[1];
1223 
1224  p[0] = ( b < 0 ) ? -b : b;
1225  _B.s = ( b < 0 ) ? -1 : 1;
1226  _B.n = 1;
1227  _B.p = p;
1228 
1229  return( mpi_div_mpi( Q, R, A, &_B ) );
1230 }
1231 
1232 /*
1233  * Modulo: R = A mod B
1234  */
1235 int mpi_mod_mpi( mpi *R, const mpi *A, const mpi *B )
1236 {
1237  int ret;
1238 
1239  if( mpi_cmp_int( B, 0 ) < 0 )
1241 
1242  MPI_CHK( mpi_div_mpi( NULL, R, A, B ) );
1243 
1244  while( mpi_cmp_int( R, 0 ) < 0 )
1245  MPI_CHK( mpi_add_mpi( R, R, B ) );
1246 
1247  while( mpi_cmp_mpi( R, B ) >= 0 )
1248  MPI_CHK( mpi_sub_mpi( R, R, B ) );
1249 
1250 cleanup:
1251 
1252  return( ret );
1253 }
1254 
1255 /*
1256  * Modulo: r = A mod b
1257  */
1258 int mpi_mod_int( t_uint *r, const mpi *A, t_sint b )
1259 {
1260  size_t i;
1261  t_uint x, y, z;
1262 
1263  if( b == 0 )
1265 
1266  if( b < 0 )
1268 
1269  /*
1270  * handle trivial cases
1271  */
1272  if( b == 1 )
1273  {
1274  *r = 0;
1275  return( 0 );
1276  }
1277 
1278  if( b == 2 )
1279  {
1280  *r = A->p[0] & 1;
1281  return( 0 );
1282  }
1283 
1284  /*
1285  * general case
1286  */
1287  for( i = A->n, y = 0; i > 0; i-- )
1288  {
1289  x = A->p[i - 1];
1290  y = ( y << biH ) | ( x >> biH );
1291  z = y / b;
1292  y -= z * b;
1293 
1294  x <<= biH;
1295  y = ( y << biH ) | ( x >> biH );
1296  z = y / b;
1297  y -= z * b;
1298  }
1299 
1300  /*
1301  * If A is negative, then the current y represents a negative value.
1302  * Flipping it to the positive side.
1303  */
1304  if( A->s < 0 && y != 0 )
1305  y = b - y;
1306 
1307  *r = y;
1308 
1309  return( 0 );
1310 }
1311 
1312 /*
1313  * Fast Montgomery initialization (thanks to Tom St Denis)
1314  */
1315 static void mpi_montg_init( t_uint *mm, const mpi *N )
1316 {
1317  t_uint x, m0 = N->p[0];
1318 
1319  x = m0;
1320  x += ( ( m0 + 2 ) & 4 ) << 1;
1321  x *= ( 2 - ( m0 * x ) );
1322 
1323  if( biL >= 16 ) x *= ( 2 - ( m0 * x ) );
1324  if( biL >= 32 ) x *= ( 2 - ( m0 * x ) );
1325  if( biL >= 64 ) x *= ( 2 - ( m0 * x ) );
1326 
1327  *mm = ~x + 1;
1328 }
1329 
1330 /*
1331  * Montgomery multiplication: A = A * B * R^-1 mod N (HAC 14.36)
1332  */
1333 static void mpi_montmul( mpi *A, const mpi *B, const mpi *N, t_uint mm, const mpi *T )
1334 {
1335  size_t i, n, m;
1336  t_uint u0, u1, *d;
1337 
1338  memset( T->p, 0, T->n * ciL );
1339 
1340  d = T->p;
1341  n = N->n;
1342  m = ( B->n < n ) ? B->n : n;
1343 
1344  for( i = 0; i < n; i++ )
1345  {
1346  /*
1347  * T = (T + u0*B + u1*N) / 2^biL
1348  */
1349  u0 = A->p[i];
1350  u1 = ( d[0] + u0 * B->p[0] ) * mm;
1351 
1352  mpi_mul_hlp( m, B->p, d, u0 );
1353  mpi_mul_hlp( n, N->p, d, u1 );
1354 
1355  *d++ = u0; d[n + 1] = 0;
1356  }
1357 
1358  memcpy( A->p, d, (n + 1) * ciL );
1359 
1360  if( mpi_cmp_abs( A, N ) >= 0 )
1361  mpi_sub_hlp( n, N->p, A->p );
1362  else
1363  /* prevent timing attacks */
1364  mpi_sub_hlp( n, A->p, T->p );
1365 }
1366 
1367 /*
1368  * Montgomery reduction: A = A * R^-1 mod N
1369  */
1370 static void mpi_montred( mpi *A, const mpi *N, t_uint mm, const mpi *T )
1371 {
1372  t_uint z = 1;
1373  mpi U;
1374 
1375  U.n = U.s = z;
1376  U.p = &z;
1377 
1378  mpi_montmul( A, &U, N, mm, T );
1379 }
1380 
1381 /*
1382  * Sliding-window exponentiation: X = A^E mod N (HAC 14.85)
1383  */
1384 int mpi_exp_mod( mpi *X, const mpi *A, const mpi *E, const mpi *N, mpi *_RR )
1385 {
1386  int ret;
1387  size_t wbits, wsize, one = 1;
1388  size_t i, j, nblimbs;
1389  size_t bufsize, nbits;
1390  t_uint ei, mm, state;
1391  mpi RR, T, W[ 2 << POLARSSL_MPI_WINDOW_SIZE ], Apos;
1392  int neg;
1393 
1394  if( mpi_cmp_int( N, 0 ) < 0 || ( N->p[0] & 1 ) == 0 )
1396 
1397  if( mpi_cmp_int( E, 0 ) < 0 )
1399 
1400  /*
1401  * Init temps and window size
1402  */
1403  mpi_montg_init( &mm, N );
1404  mpi_init( &RR ); mpi_init( &T );
1405  memset( W, 0, sizeof( W ) );
1406 
1407  i = mpi_msb( E );
1408 
1409  wsize = ( i > 671 ) ? 6 : ( i > 239 ) ? 5 :
1410  ( i > 79 ) ? 4 : ( i > 23 ) ? 3 : 1;
1411 
1412  if( wsize > POLARSSL_MPI_WINDOW_SIZE )
1413  wsize = POLARSSL_MPI_WINDOW_SIZE;
1414 
1415  j = N->n + 1;
1416  MPI_CHK( mpi_grow( X, j ) );
1417  MPI_CHK( mpi_grow( &W[1], j ) );
1418  MPI_CHK( mpi_grow( &T, j * 2 ) );
1419 
1420  /*
1421  * Compensate for negative A (and correct at the end)
1422  */
1423  neg = ( A->s == -1 );
1424 
1425  mpi_init( &Apos );
1426  if( neg )
1427  {
1428  MPI_CHK( mpi_copy( &Apos, A ) );
1429  Apos.s = 1;
1430  A = &Apos;
1431  }
1432 
1433  /*
1434  * If 1st call, pre-compute R^2 mod N
1435  */
1436  if( _RR == NULL || _RR->p == NULL )
1437  {
1438  MPI_CHK( mpi_lset( &RR, 1 ) );
1439  MPI_CHK( mpi_shift_l( &RR, N->n * 2 * biL ) );
1440  MPI_CHK( mpi_mod_mpi( &RR, &RR, N ) );
1441 
1442  if( _RR != NULL )
1443  memcpy( _RR, &RR, sizeof( mpi ) );
1444  }
1445  else
1446  memcpy( &RR, _RR, sizeof( mpi ) );
1447 
1448  /*
1449  * W[1] = A * R^2 * R^-1 mod N = A * R mod N
1450  */
1451  if( mpi_cmp_mpi( A, N ) >= 0 )
1452  mpi_mod_mpi( &W[1], A, N );
1453  else mpi_copy( &W[1], A );
1454 
1455  mpi_montmul( &W[1], &RR, N, mm, &T );
1456 
1457  /*
1458  * X = R^2 * R^-1 mod N = R mod N
1459  */
1460  MPI_CHK( mpi_copy( X, &RR ) );
1461  mpi_montred( X, N, mm, &T );
1462 
1463  if( wsize > 1 )
1464  {
1465  /*
1466  * W[1 << (wsize - 1)] = W[1] ^ (wsize - 1)
1467  */
1468  j = one << (wsize - 1);
1469 
1470  MPI_CHK( mpi_grow( &W[j], N->n + 1 ) );
1471  MPI_CHK( mpi_copy( &W[j], &W[1] ) );
1472 
1473  for( i = 0; i < wsize - 1; i++ )
1474  mpi_montmul( &W[j], &W[j], N, mm, &T );
1475 
1476  /*
1477  * W[i] = W[i - 1] * W[1]
1478  */
1479  for( i = j + 1; i < (one << wsize); i++ )
1480  {
1481  MPI_CHK( mpi_grow( &W[i], N->n + 1 ) );
1482  MPI_CHK( mpi_copy( &W[i], &W[i - 1] ) );
1483 
1484  mpi_montmul( &W[i], &W[1], N, mm, &T );
1485  }
1486  }
1487 
1488  nblimbs = E->n;
1489  bufsize = 0;
1490  nbits = 0;
1491  wbits = 0;
1492  state = 0;
1493 
1494  while( 1 )
1495  {
1496  if( bufsize == 0 )
1497  {
1498  if( nblimbs-- == 0 )
1499  break;
1500 
1501  bufsize = sizeof( t_uint ) << 3;
1502  }
1503 
1504  bufsize--;
1505 
1506  ei = (E->p[nblimbs] >> bufsize) & 1;
1507 
1508  /*
1509  * skip leading 0s
1510  */
1511  if( ei == 0 && state == 0 )
1512  continue;
1513 
1514  if( ei == 0 && state == 1 )
1515  {
1516  /*
1517  * out of window, square X
1518  */
1519  mpi_montmul( X, X, N, mm, &T );
1520  continue;
1521  }
1522 
1523  /*
1524  * add ei to current window
1525  */
1526  state = 2;
1527 
1528  nbits++;
1529  wbits |= (ei << (wsize - nbits));
1530 
1531  if( nbits == wsize )
1532  {
1533  /*
1534  * X = X^wsize R^-1 mod N
1535  */
1536  for( i = 0; i < wsize; i++ )
1537  mpi_montmul( X, X, N, mm, &T );
1538 
1539  /*
1540  * X = X * W[wbits] R^-1 mod N
1541  */
1542  mpi_montmul( X, &W[wbits], N, mm, &T );
1543 
1544  state--;
1545  nbits = 0;
1546  wbits = 0;
1547  }
1548  }
1549 
1550  /*
1551  * process the remaining bits
1552  */
1553  for( i = 0; i < nbits; i++ )
1554  {
1555  mpi_montmul( X, X, N, mm, &T );
1556 
1557  wbits <<= 1;
1558 
1559  if( (wbits & (one << wsize)) != 0 )
1560  mpi_montmul( X, &W[1], N, mm, &T );
1561  }
1562 
1563  /*
1564  * X = A^E * R * R^-1 mod N = A^E mod N
1565  */
1566  mpi_montred( X, N, mm, &T );
1567 
1568  if( neg )
1569  {
1570  X->s = -1;
1571  mpi_add_mpi( X, N, X );
1572  }
1573 
1574 cleanup:
1575 
1576  for( i = (one << (wsize - 1)); i < (one << wsize); i++ )
1577  mpi_free( &W[i] );
1578 
1579  mpi_free( &W[1] ); mpi_free( &T ); mpi_free( &Apos );
1580 
1581  if( _RR == NULL )
1582  mpi_free( &RR );
1583 
1584  return( ret );
1585 }
1586 
1587 /*
1588  * Greatest common divisor: G = gcd(A, B) (HAC 14.54)
1589  */
1590 int mpi_gcd( mpi *G, const mpi *A, const mpi *B )
1591 {
1592  int ret;
1593  size_t lz, lzt;
1594  mpi TG, TA, TB;
1595 
1596  mpi_init( &TG ); mpi_init( &TA ); mpi_init( &TB );
1597 
1598  MPI_CHK( mpi_copy( &TA, A ) );
1599  MPI_CHK( mpi_copy( &TB, B ) );
1600 
1601  lz = mpi_lsb( &TA );
1602  lzt = mpi_lsb( &TB );
1603 
1604  if ( lzt < lz )
1605  lz = lzt;
1606 
1607  MPI_CHK( mpi_shift_r( &TA, lz ) );
1608  MPI_CHK( mpi_shift_r( &TB, lz ) );
1609 
1610  TA.s = TB.s = 1;
1611 
1612  while( mpi_cmp_int( &TA, 0 ) != 0 )
1613  {
1614  MPI_CHK( mpi_shift_r( &TA, mpi_lsb( &TA ) ) );
1615  MPI_CHK( mpi_shift_r( &TB, mpi_lsb( &TB ) ) );
1616 
1617  if( mpi_cmp_mpi( &TA, &TB ) >= 0 )
1618  {
1619  MPI_CHK( mpi_sub_abs( &TA, &TA, &TB ) );
1620  MPI_CHK( mpi_shift_r( &TA, 1 ) );
1621  }
1622  else
1623  {
1624  MPI_CHK( mpi_sub_abs( &TB, &TB, &TA ) );
1625  MPI_CHK( mpi_shift_r( &TB, 1 ) );
1626  }
1627  }
1628 
1629  MPI_CHK( mpi_shift_l( &TB, lz ) );
1630  MPI_CHK( mpi_copy( G, &TB ) );
1631 
1632 cleanup:
1633 
1634  mpi_free( &TG ); mpi_free( &TA ); mpi_free( &TB );
1635 
1636  return( ret );
1637 }
1638 
1639 int mpi_fill_random( mpi *X, size_t size,
1640  int (*f_rng)(void *, unsigned char *, size_t),
1641  void *p_rng )
1642 {
1643  int ret;
1644 
1645  MPI_CHK( mpi_grow( X, CHARS_TO_LIMBS( size ) ) );
1646  MPI_CHK( mpi_lset( X, 0 ) );
1647 
1648  MPI_CHK( f_rng( p_rng, (unsigned char *) X->p, size ) );
1649 
1650 cleanup:
1651  return( ret );
1652 }
1653 
1654 /*
1655  * Modular inverse: X = A^-1 mod N (HAC 14.61 / 14.64)
1656  */
1657 int mpi_inv_mod( mpi *X, const mpi *A, const mpi *N )
1658 {
1659  int ret;
1660  mpi G, TA, TU, U1, U2, TB, TV, V1, V2;
1661 
1662  if( mpi_cmp_int( N, 0 ) <= 0 )
1664 
1665  mpi_init( &TA ); mpi_init( &TU ); mpi_init( &U1 ); mpi_init( &U2 );
1666  mpi_init( &G ); mpi_init( &TB ); mpi_init( &TV );
1667  mpi_init( &V1 ); mpi_init( &V2 );
1668 
1669  MPI_CHK( mpi_gcd( &G, A, N ) );
1670 
1671  if( mpi_cmp_int( &G, 1 ) != 0 )
1672  {
1674  goto cleanup;
1675  }
1676 
1677  MPI_CHK( mpi_mod_mpi( &TA, A, N ) );
1678  MPI_CHK( mpi_copy( &TU, &TA ) );
1679  MPI_CHK( mpi_copy( &TB, N ) );
1680  MPI_CHK( mpi_copy( &TV, N ) );
1681 
1682  MPI_CHK( mpi_lset( &U1, 1 ) );
1683  MPI_CHK( mpi_lset( &U2, 0 ) );
1684  MPI_CHK( mpi_lset( &V1, 0 ) );
1685  MPI_CHK( mpi_lset( &V2, 1 ) );
1686 
1687  do
1688  {
1689  while( ( TU.p[0] & 1 ) == 0 )
1690  {
1691  MPI_CHK( mpi_shift_r( &TU, 1 ) );
1692 
1693  if( ( U1.p[0] & 1 ) != 0 || ( U2.p[0] & 1 ) != 0 )
1694  {
1695  MPI_CHK( mpi_add_mpi( &U1, &U1, &TB ) );
1696  MPI_CHK( mpi_sub_mpi( &U2, &U2, &TA ) );
1697  }
1698 
1699  MPI_CHK( mpi_shift_r( &U1, 1 ) );
1700  MPI_CHK( mpi_shift_r( &U2, 1 ) );
1701  }
1702 
1703  while( ( TV.p[0] & 1 ) == 0 )
1704  {
1705  MPI_CHK( mpi_shift_r( &TV, 1 ) );
1706 
1707  if( ( V1.p[0] & 1 ) != 0 || ( V2.p[0] & 1 ) != 0 )
1708  {
1709  MPI_CHK( mpi_add_mpi( &V1, &V1, &TB ) );
1710  MPI_CHK( mpi_sub_mpi( &V2, &V2, &TA ) );
1711  }
1712 
1713  MPI_CHK( mpi_shift_r( &V1, 1 ) );
1714  MPI_CHK( mpi_shift_r( &V2, 1 ) );
1715  }
1716 
1717  if( mpi_cmp_mpi( &TU, &TV ) >= 0 )
1718  {
1719  MPI_CHK( mpi_sub_mpi( &TU, &TU, &TV ) );
1720  MPI_CHK( mpi_sub_mpi( &U1, &U1, &V1 ) );
1721  MPI_CHK( mpi_sub_mpi( &U2, &U2, &V2 ) );
1722  }
1723  else
1724  {
1725  MPI_CHK( mpi_sub_mpi( &TV, &TV, &TU ) );
1726  MPI_CHK( mpi_sub_mpi( &V1, &V1, &U1 ) );
1727  MPI_CHK( mpi_sub_mpi( &V2, &V2, &U2 ) );
1728  }
1729  }
1730  while( mpi_cmp_int( &TU, 0 ) != 0 );
1731 
1732  while( mpi_cmp_int( &V1, 0 ) < 0 )
1733  MPI_CHK( mpi_add_mpi( &V1, &V1, N ) );
1734 
1735  while( mpi_cmp_mpi( &V1, N ) >= 0 )
1736  MPI_CHK( mpi_sub_mpi( &V1, &V1, N ) );
1737 
1738  MPI_CHK( mpi_copy( X, &V1 ) );
1739 
1740 cleanup:
1741 
1742  mpi_free( &TA ); mpi_free( &TU ); mpi_free( &U1 ); mpi_free( &U2 );
1743  mpi_free( &G ); mpi_free( &TB ); mpi_free( &TV );
1744  mpi_free( &V1 ); mpi_free( &V2 );
1745 
1746  return( ret );
1747 }
1748 
1749 #if defined(POLARSSL_GENPRIME)
1750 
1751 static const int small_prime[] =
1752 {
1753  3, 5, 7, 11, 13, 17, 19, 23,
1754  29, 31, 37, 41, 43, 47, 53, 59,
1755  61, 67, 71, 73, 79, 83, 89, 97,
1756  101, 103, 107, 109, 113, 127, 131, 137,
1757  139, 149, 151, 157, 163, 167, 173, 179,
1758  181, 191, 193, 197, 199, 211, 223, 227,
1759  229, 233, 239, 241, 251, 257, 263, 269,
1760  271, 277, 281, 283, 293, 307, 311, 313,
1761  317, 331, 337, 347, 349, 353, 359, 367,
1762  373, 379, 383, 389, 397, 401, 409, 419,
1763  421, 431, 433, 439, 443, 449, 457, 461,
1764  463, 467, 479, 487, 491, 499, 503, 509,
1765  521, 523, 541, 547, 557, 563, 569, 571,
1766  577, 587, 593, 599, 601, 607, 613, 617,
1767  619, 631, 641, 643, 647, 653, 659, 661,
1768  673, 677, 683, 691, 701, 709, 719, 727,
1769  733, 739, 743, 751, 757, 761, 769, 773,
1770  787, 797, 809, 811, 821, 823, 827, 829,
1771  839, 853, 857, 859, 863, 877, 881, 883,
1772  887, 907, 911, 919, 929, 937, 941, 947,
1773  953, 967, 971, 977, 983, 991, 997, -103
1774 };
1775 
1776 /*
1777  * Miller-Rabin primality test (HAC 4.24)
1778  */
1779 int mpi_is_prime( mpi *X,
1780  int (*f_rng)(void *, unsigned char *, size_t),
1781  void *p_rng )
1782 {
1783  int ret, xs;
1784  size_t i, j, n, s;
1785  mpi W, R, T, A, RR;
1786 
1787  if( mpi_cmp_int( X, 0 ) == 0 ||
1788  mpi_cmp_int( X, 1 ) == 0 )
1790 
1791  if( mpi_cmp_int( X, 2 ) == 0 )
1792  return( 0 );
1793 
1794  mpi_init( &W ); mpi_init( &R ); mpi_init( &T ); mpi_init( &A );
1795  mpi_init( &RR );
1796 
1797  xs = X->s; X->s = 1;
1798 
1799  /*
1800  * test trivial factors first
1801  */
1802  if( ( X->p[0] & 1 ) == 0 )
1804 
1805  for( i = 0; small_prime[i] > 0; i++ )
1806  {
1807  t_uint r;
1808 
1809  if( mpi_cmp_int( X, small_prime[i] ) <= 0 )
1810  return( 0 );
1811 
1812  MPI_CHK( mpi_mod_int( &r, X, small_prime[i] ) );
1813 
1814  if( r == 0 )
1816  }
1817 
1818  /*
1819  * W = |X| - 1
1820  * R = W >> lsb( W )
1821  */
1822  MPI_CHK( mpi_sub_int( &W, X, 1 ) );
1823  s = mpi_lsb( &W );
1824  MPI_CHK( mpi_copy( &R, &W ) );
1825  MPI_CHK( mpi_shift_r( &R, s ) );
1826 
1827  i = mpi_msb( X );
1828  /*
1829  * HAC, table 4.4
1830  */
1831  n = ( ( i >= 1300 ) ? 2 : ( i >= 850 ) ? 3 :
1832  ( i >= 650 ) ? 4 : ( i >= 350 ) ? 8 :
1833  ( i >= 250 ) ? 12 : ( i >= 150 ) ? 18 : 27 );
1834 
1835  for( i = 0; i < n; i++ )
1836  {
1837  /*
1838  * pick a random A, 1 < A < |X| - 1
1839  */
1840  MPI_CHK( mpi_fill_random( &A, X->n * ciL, f_rng, p_rng ) );
1841 
1842  if( mpi_cmp_mpi( &A, &W ) >= 0 )
1843  {
1844  j = mpi_msb( &A ) - mpi_msb( &W );
1845  MPI_CHK( mpi_shift_r( &A, j + 1 ) );
1846  }
1847  A.p[0] |= 3;
1848 
1849  /*
1850  * A = A^R mod |X|
1851  */
1852  MPI_CHK( mpi_exp_mod( &A, &A, &R, X, &RR ) );
1853 
1854  if( mpi_cmp_mpi( &A, &W ) == 0 ||
1855  mpi_cmp_int( &A, 1 ) == 0 )
1856  continue;
1857 
1858  j = 1;
1859  while( j < s && mpi_cmp_mpi( &A, &W ) != 0 )
1860  {
1861  /*
1862  * A = A * A mod |X|
1863  */
1864  MPI_CHK( mpi_mul_mpi( &T, &A, &A ) );
1865  MPI_CHK( mpi_mod_mpi( &A, &T, X ) );
1866 
1867  if( mpi_cmp_int( &A, 1 ) == 0 )
1868  break;
1869 
1870  j++;
1871  }
1872 
1873  /*
1874  * not prime if A != |X| - 1 or A == 1
1875  */
1876  if( mpi_cmp_mpi( &A, &W ) != 0 ||
1877  mpi_cmp_int( &A, 1 ) == 0 )
1878  {
1880  break;
1881  }
1882  }
1883 
1884 cleanup:
1885 
1886  X->s = xs;
1887 
1888  mpi_free( &W ); mpi_free( &R ); mpi_free( &T ); mpi_free( &A );
1889  mpi_free( &RR );
1890 
1891  return( ret );
1892 }
1893 
1894 /*
1895  * Prime number generation
1896  */
1897 int mpi_gen_prime( mpi *X, size_t nbits, int dh_flag,
1898  int (*f_rng)(void *, unsigned char *, size_t),
1899  void *p_rng )
1900 {
1901  int ret;
1902  size_t k, n;
1903  mpi Y;
1904 
1905  if( nbits < 3 || nbits > POLARSSL_MPI_MAX_BITS )
1907 
1908  mpi_init( &Y );
1909 
1910  n = BITS_TO_LIMBS( nbits );
1911 
1912  MPI_CHK( mpi_fill_random( X, n * ciL, f_rng, p_rng ) );
1913 
1914  k = mpi_msb( X );
1915  if( k < nbits ) MPI_CHK( mpi_shift_l( X, nbits - k ) );
1916  if( k > nbits ) MPI_CHK( mpi_shift_r( X, k - nbits ) );
1917 
1918  X->p[0] |= 3;
1919 
1920  if( dh_flag == 0 )
1921  {
1922  while( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) != 0 )
1923  {
1924  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1925  goto cleanup;
1926 
1927  MPI_CHK( mpi_add_int( X, X, 2 ) );
1928  }
1929  }
1930  else
1931  {
1932  MPI_CHK( mpi_sub_int( &Y, X, 1 ) );
1933  MPI_CHK( mpi_shift_r( &Y, 1 ) );
1934 
1935  while( 1 )
1936  {
1937  if( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) == 0 )
1938  {
1939  if( ( ret = mpi_is_prime( &Y, f_rng, p_rng ) ) == 0 )
1940  break;
1941 
1942  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1943  goto cleanup;
1944  }
1945 
1946  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1947  goto cleanup;
1948 
1949  MPI_CHK( mpi_add_int( &Y, X, 1 ) );
1950  MPI_CHK( mpi_add_int( X, X, 2 ) );
1951  MPI_CHK( mpi_shift_r( &Y, 1 ) );
1952  }
1953  }
1954 
1955 cleanup:
1956 
1957  mpi_free( &Y );
1958 
1959  return( ret );
1960 }
1961 
1962 #endif
1963 
1964 #if defined(POLARSSL_SELF_TEST)
1965 
1966 #define GCD_PAIR_COUNT 3
1967 
1968 static const int gcd_pairs[GCD_PAIR_COUNT][3] =
1969 {
1970  { 693, 609, 21 },
1971  { 1764, 868, 28 },
1972  { 768454923, 542167814, 1 }
1973 };
1974 
1975 /*
1976  * Checkup routine
1977  */
1978 int mpi_self_test( int verbose )
1979 {
1980  int ret, i;
1981  mpi A, E, N, X, Y, U, V;
1982 
1983  mpi_init( &A ); mpi_init( &E ); mpi_init( &N ); mpi_init( &X );
1984  mpi_init( &Y ); mpi_init( &U ); mpi_init( &V );
1985 
1986  MPI_CHK( mpi_read_string( &A, 16,
1987  "EFE021C2645FD1DC586E69184AF4A31E" \
1988  "D5F53E93B5F123FA41680867BA110131" \
1989  "944FE7952E2517337780CB0DB80E61AA" \
1990  "E7C8DDC6C5C6AADEB34EB38A2F40D5E6" ) );
1991 
1992  MPI_CHK( mpi_read_string( &E, 16,
1993  "B2E7EFD37075B9F03FF989C7C5051C20" \
1994  "34D2A323810251127E7BF8625A4F49A5" \
1995  "F3E27F4DA8BD59C47D6DAABA4C8127BD" \
1996  "5B5C25763222FEFCCFC38B832366C29E" ) );
1997 
1998  MPI_CHK( mpi_read_string( &N, 16,
1999  "0066A198186C18C10B2F5ED9B522752A" \
2000  "9830B69916E535C8F047518A889A43A5" \
2001  "94B6BED27A168D31D4A52F88925AA8F5" ) );
2002 
2003  MPI_CHK( mpi_mul_mpi( &X, &A, &N ) );
2004 
2005  MPI_CHK( mpi_read_string( &U, 16,
2006  "602AB7ECA597A3D6B56FF9829A5E8B85" \
2007  "9E857EA95A03512E2BAE7391688D264A" \
2008  "A5663B0341DB9CCFD2C4C5F421FEC814" \
2009  "8001B72E848A38CAE1C65F78E56ABDEF" \
2010  "E12D3C039B8A02D6BE593F0BBBDA56F1" \
2011  "ECF677152EF804370C1A305CAF3B5BF1" \
2012  "30879B56C61DE584A0F53A2447A51E" ) );
2013 
2014  if( verbose != 0 )
2015  printf( " MPI test #1 (mul_mpi): " );
2016 
2017  if( mpi_cmp_mpi( &X, &U ) != 0 )
2018  {
2019  if( verbose != 0 )
2020  printf( "failed\n" );
2021 
2022  return( 1 );
2023  }
2024 
2025  if( verbose != 0 )
2026  printf( "passed\n" );
2027 
2028  MPI_CHK( mpi_div_mpi( &X, &Y, &A, &N ) );
2029 
2030  MPI_CHK( mpi_read_string( &U, 16,
2031  "256567336059E52CAE22925474705F39A94" ) );
2032 
2033  MPI_CHK( mpi_read_string( &V, 16,
2034  "6613F26162223DF488E9CD48CC132C7A" \
2035  "0AC93C701B001B092E4E5B9F73BCD27B" \
2036  "9EE50D0657C77F374E903CDFA4C642" ) );
2037 
2038  if( verbose != 0 )
2039  printf( " MPI test #2 (div_mpi): " );
2040 
2041  if( mpi_cmp_mpi( &X, &U ) != 0 ||
2042  mpi_cmp_mpi( &Y, &V ) != 0 )
2043  {
2044  if( verbose != 0 )
2045  printf( "failed\n" );
2046 
2047  return( 1 );
2048  }
2049 
2050  if( verbose != 0 )
2051  printf( "passed\n" );
2052 
2053  MPI_CHK( mpi_exp_mod( &X, &A, &E, &N, NULL ) );
2054 
2055  MPI_CHK( mpi_read_string( &U, 16,
2056  "36E139AEA55215609D2816998ED020BB" \
2057  "BD96C37890F65171D948E9BC7CBAA4D9" \
2058  "325D24D6A3C12710F10A09FA08AB87" ) );
2059 
2060  if( verbose != 0 )
2061  printf( " MPI test #3 (exp_mod): " );
2062 
2063  if( mpi_cmp_mpi( &X, &U ) != 0 )
2064  {
2065  if( verbose != 0 )
2066  printf( "failed\n" );
2067 
2068  return( 1 );
2069  }
2070 
2071  if( verbose != 0 )
2072  printf( "passed\n" );
2073 
2074 #if defined(POLARSSL_GENPRIME)
2075  MPI_CHK( mpi_inv_mod( &X, &A, &N ) );
2076 
2077  MPI_CHK( mpi_read_string( &U, 16,
2078  "003A0AAEDD7E784FC07D8F9EC6E3BFD5" \
2079  "C3DBA76456363A10869622EAC2DD84EC" \
2080  "C5B8A74DAC4D09E03B5E0BE779F2DF61" ) );
2081 
2082  if( verbose != 0 )
2083  printf( " MPI test #4 (inv_mod): " );
2084 
2085  if( mpi_cmp_mpi( &X, &U ) != 0 )
2086  {
2087  if( verbose != 0 )
2088  printf( "failed\n" );
2089 
2090  return( 1 );
2091  }
2092 
2093  if( verbose != 0 )
2094  printf( "passed\n" );
2095 #endif
2096 
2097  if( verbose != 0 )
2098  printf( " MPI test #5 (simple gcd): " );
2099 
2100  for ( i = 0; i < GCD_PAIR_COUNT; i++)
2101  {
2102  MPI_CHK( mpi_lset( &X, gcd_pairs[i][0] ) );
2103  MPI_CHK( mpi_lset( &Y, gcd_pairs[i][1] ) );
2104 
2105  MPI_CHK( mpi_gcd( &A, &X, &Y ) );
2106 
2107  if( mpi_cmp_int( &A, gcd_pairs[i][2] ) != 0 )
2108  {
2109  if( verbose != 0 )
2110  printf( "failed at %d\n", i );
2111 
2112  return( 1 );
2113  }
2114  }
2115 
2116  if( verbose != 0 )
2117  printf( "passed\n" );
2118 
2119 cleanup:
2120 
2121  if( ret != 0 && verbose != 0 )
2122  printf( "Unexpected error, return code = %08X\n", ret );
2123 
2124  mpi_free( &A ); mpi_free( &E ); mpi_free( &N ); mpi_free( &X );
2125  mpi_free( &Y ); mpi_free( &U ); mpi_free( &V );
2126 
2127  if( verbose != 0 )
2128  printf( "\n" );
2129 
2130  return( ret );
2131 }
2132 
2133 #endif
2134 
2135 #endif