PolarSSL v1.1.4
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( 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 += sprintf( p, "%02X", c );
410  k = 1;
411  }
412  }
413  }
414  else
415  {
416  MPI_CHK( mpi_copy( &T, X ) );
417 
418  if( T.s == -1 )
419  T.s = 1;
420 
421  MPI_CHK( mpi_write_hlp( &T, radix, &p ) );
422  }
423 
424  *p++ = '\0';
425  *slen = p - s;
426 
427 cleanup:
428 
429  mpi_free( &T );
430 
431  return( ret );
432 }
433 
434 #if defined(POLARSSL_FS_IO)
435 /*
436  * Read X from an opened file
437  */
438 int mpi_read_file( mpi *X, int radix, FILE *fin )
439 {
440  t_uint d;
441  size_t slen;
442  char *p;
443  /*
444  * Buffer should have space for (short) label and decimal formatted MPI,
445  * newline characters and '\0'
446  */
448 
449  memset( s, 0, sizeof( s ) );
450  if( fgets( s, sizeof( s ) - 1, fin ) == NULL )
452 
453  slen = strlen( s );
454  if( slen == sizeof( s ) - 2 )
456 
457  if( s[slen - 1] == '\n' ) { slen--; s[slen] = '\0'; }
458  if( s[slen - 1] == '\r' ) { slen--; s[slen] = '\0'; }
459 
460  p = s + slen;
461  while( --p >= s )
462  if( mpi_get_digit( &d, radix, *p ) != 0 )
463  break;
464 
465  return( mpi_read_string( X, radix, p + 1 ) );
466 }
467 
468 /*
469  * Write X into an opened file (or stdout if fout == NULL)
470  */
471 int mpi_write_file( const char *p, const mpi *X, int radix, FILE *fout )
472 {
473  int ret;
474  size_t n, slen, plen;
475  /*
476  * Buffer should have space for minus sign, hexified MPI and '\0'
477  */
478  char s[ 2 * POLARSSL_MPI_MAX_SIZE + 2 ];
479 
480  n = sizeof( s );
481  memset( s, 0, n );
482  n -= 2;
483 
484  MPI_CHK( mpi_write_string( X, radix, s, (size_t *) &n ) );
485 
486  if( p == NULL ) p = "";
487 
488  plen = strlen( p );
489  slen = strlen( s );
490  s[slen++] = '\r';
491  s[slen++] = '\n';
492 
493  if( fout != NULL )
494  {
495  if( fwrite( p, 1, plen, fout ) != plen ||
496  fwrite( s, 1, slen, fout ) != slen )
498  }
499  else
500  printf( "%s%s", p, s );
501 
502 cleanup:
503 
504  return( ret );
505 }
506 #endif /* POLARSSL_FS_IO */
507 
508 /*
509  * Import X from unsigned binary data, big endian
510  */
511 int mpi_read_binary( mpi *X, const unsigned char *buf, size_t buflen )
512 {
513  int ret;
514  size_t i, j, n;
515 
516  for( n = 0; n < buflen; n++ )
517  if( buf[n] != 0 )
518  break;
519 
520  MPI_CHK( mpi_grow( X, CHARS_TO_LIMBS( buflen - n ) ) );
521  MPI_CHK( mpi_lset( X, 0 ) );
522 
523  for( i = buflen, j = 0; i > n; i--, j++ )
524  X->p[j / ciL] |= ((t_uint) buf[i - 1]) << ((j % ciL) << 3);
525 
526 cleanup:
527 
528  return( ret );
529 }
530 
531 /*
532  * Export X into unsigned binary data, big endian
533  */
534 int mpi_write_binary( const mpi *X, unsigned char *buf, size_t buflen )
535 {
536  size_t i, j, n;
537 
538  n = mpi_size( X );
539 
540  if( buflen < n )
542 
543  memset( buf, 0, buflen );
544 
545  for( i = buflen - 1, j = 0; n > 0; i--, j++, n-- )
546  buf[i] = (unsigned char)( X->p[j / ciL] >> ((j % ciL) << 3) );
547 
548  return( 0 );
549 }
550 
551 /*
552  * Left-shift: X <<= count
553  */
554 int mpi_shift_l( mpi *X, size_t count )
555 {
556  int ret;
557  size_t i, v0, t1;
558  t_uint r0 = 0, r1;
559 
560  v0 = count / (biL );
561  t1 = count & (biL - 1);
562 
563  i = mpi_msb( X ) + count;
564 
565  if( X->n * biL < i )
566  MPI_CHK( mpi_grow( X, BITS_TO_LIMBS( i ) ) );
567 
568  ret = 0;
569 
570  /*
571  * shift by count / limb_size
572  */
573  if( v0 > 0 )
574  {
575  for( i = X->n; i > v0; i-- )
576  X->p[i - 1] = X->p[i - v0 - 1];
577 
578  for( ; i > 0; i-- )
579  X->p[i - 1] = 0;
580  }
581 
582  /*
583  * shift by count % limb_size
584  */
585  if( t1 > 0 )
586  {
587  for( i = v0; i < X->n; i++ )
588  {
589  r1 = X->p[i] >> (biL - t1);
590  X->p[i] <<= t1;
591  X->p[i] |= r0;
592  r0 = r1;
593  }
594  }
595 
596 cleanup:
597 
598  return( ret );
599 }
600 
601 /*
602  * Right-shift: X >>= count
603  */
604 int mpi_shift_r( mpi *X, size_t count )
605 {
606  size_t i, v0, v1;
607  t_uint r0 = 0, r1;
608 
609  v0 = count / biL;
610  v1 = count & (biL - 1);
611 
612  /*
613  * shift by count / limb_size
614  */
615  if( v0 > 0 )
616  {
617  for( i = 0; i < X->n - v0; i++ )
618  X->p[i] = X->p[i + v0];
619 
620  for( ; i < X->n; i++ )
621  X->p[i] = 0;
622  }
623 
624  /*
625  * shift by count % limb_size
626  */
627  if( v1 > 0 )
628  {
629  for( i = X->n; i > 0; i-- )
630  {
631  r1 = X->p[i - 1] << (biL - v1);
632  X->p[i - 1] >>= v1;
633  X->p[i - 1] |= r0;
634  r0 = r1;
635  }
636  }
637 
638  return( 0 );
639 }
640 
641 /*
642  * Compare unsigned values
643  */
644 int mpi_cmp_abs( const mpi *X, const mpi *Y )
645 {
646  size_t i, j;
647 
648  for( i = X->n; i > 0; i-- )
649  if( X->p[i - 1] != 0 )
650  break;
651 
652  for( j = Y->n; j > 0; j-- )
653  if( Y->p[j - 1] != 0 )
654  break;
655 
656  if( i == 0 && j == 0 )
657  return( 0 );
658 
659  if( i > j ) return( 1 );
660  if( j > i ) return( -1 );
661 
662  for( ; i > 0; i-- )
663  {
664  if( X->p[i - 1] > Y->p[i - 1] ) return( 1 );
665  if( X->p[i - 1] < Y->p[i - 1] ) return( -1 );
666  }
667 
668  return( 0 );
669 }
670 
671 /*
672  * Compare signed values
673  */
674 int mpi_cmp_mpi( const mpi *X, const mpi *Y )
675 {
676  size_t i, j;
677 
678  for( i = X->n; i > 0; i-- )
679  if( X->p[i - 1] != 0 )
680  break;
681 
682  for( j = Y->n; j > 0; j-- )
683  if( Y->p[j - 1] != 0 )
684  break;
685 
686  if( i == 0 && j == 0 )
687  return( 0 );
688 
689  if( i > j ) return( X->s );
690  if( j > i ) return( -Y->s );
691 
692  if( X->s > 0 && Y->s < 0 ) return( 1 );
693  if( Y->s > 0 && X->s < 0 ) return( -1 );
694 
695  for( ; i > 0; i-- )
696  {
697  if( X->p[i - 1] > Y->p[i - 1] ) return( X->s );
698  if( X->p[i - 1] < Y->p[i - 1] ) return( -X->s );
699  }
700 
701  return( 0 );
702 }
703 
704 /*
705  * Compare signed values
706  */
707 int mpi_cmp_int( const mpi *X, t_sint z )
708 {
709  mpi Y;
710  t_uint p[1];
711 
712  *p = ( z < 0 ) ? -z : z;
713  Y.s = ( z < 0 ) ? -1 : 1;
714  Y.n = 1;
715  Y.p = p;
716 
717  return( mpi_cmp_mpi( X, &Y ) );
718 }
719 
720 /*
721  * Unsigned addition: X = |A| + |B| (HAC 14.7)
722  */
723 int mpi_add_abs( mpi *X, const mpi *A, const mpi *B )
724 {
725  int ret;
726  size_t i, j;
727  t_uint *o, *p, c;
728 
729  if( X == B )
730  {
731  const mpi *T = A; A = X; B = T;
732  }
733 
734  if( X != A )
735  MPI_CHK( mpi_copy( X, A ) );
736 
737  /*
738  * X should always be positive as a result of unsigned additions.
739  */
740  X->s = 1;
741 
742  for( j = B->n; j > 0; j-- )
743  if( B->p[j - 1] != 0 )
744  break;
745 
746  MPI_CHK( mpi_grow( X, j ) );
747 
748  o = B->p; p = X->p; c = 0;
749 
750  for( i = 0; i < j; i++, o++, p++ )
751  {
752  *p += c; c = ( *p < c );
753  *p += *o; c += ( *p < *o );
754  }
755 
756  while( c != 0 )
757  {
758  if( i >= X->n )
759  {
760  MPI_CHK( mpi_grow( X, i + 1 ) );
761  p = X->p + i;
762  }
763 
764  *p += c; c = ( *p < c ); i++;
765  }
766 
767 cleanup:
768 
769  return( ret );
770 }
771 
772 /*
773  * Helper for mpi substraction
774  */
775 static void mpi_sub_hlp( size_t n, t_uint *s, t_uint *d )
776 {
777  size_t i;
778  t_uint c, z;
779 
780  for( i = c = 0; i < n; i++, s++, d++ )
781  {
782  z = ( *d < c ); *d -= c;
783  c = ( *d < *s ) + z; *d -= *s;
784  }
785 
786  while( c != 0 )
787  {
788  z = ( *d < c ); *d -= c;
789  c = z; i++; d++;
790  }
791 }
792 
793 /*
794  * Unsigned substraction: X = |A| - |B| (HAC 14.9)
795  */
796 int mpi_sub_abs( mpi *X, const mpi *A, const mpi *B )
797 {
798  mpi TB;
799  int ret;
800  size_t n;
801 
802  if( mpi_cmp_abs( A, B ) < 0 )
804 
805  mpi_init( &TB );
806 
807  if( X == B )
808  {
809  MPI_CHK( mpi_copy( &TB, B ) );
810  B = &TB;
811  }
812 
813  if( X != A )
814  MPI_CHK( mpi_copy( X, A ) );
815 
816  /*
817  * X should always be positive as a result of unsigned substractions.
818  */
819  X->s = 1;
820 
821  ret = 0;
822 
823  for( n = B->n; n > 0; n-- )
824  if( B->p[n - 1] != 0 )
825  break;
826 
827  mpi_sub_hlp( n, B->p, X->p );
828 
829 cleanup:
830 
831  mpi_free( &TB );
832 
833  return( ret );
834 }
835 
836 /*
837  * Signed addition: X = A + B
838  */
839 int mpi_add_mpi( mpi *X, const mpi *A, const mpi *B )
840 {
841  int ret, s = A->s;
842 
843  if( A->s * B->s < 0 )
844  {
845  if( mpi_cmp_abs( A, B ) >= 0 )
846  {
847  MPI_CHK( mpi_sub_abs( X, A, B ) );
848  X->s = s;
849  }
850  else
851  {
852  MPI_CHK( mpi_sub_abs( X, B, A ) );
853  X->s = -s;
854  }
855  }
856  else
857  {
858  MPI_CHK( mpi_add_abs( X, A, B ) );
859  X->s = s;
860  }
861 
862 cleanup:
863 
864  return( ret );
865 }
866 
867 /*
868  * Signed substraction: X = A - B
869  */
870 int mpi_sub_mpi( mpi *X, const mpi *A, const mpi *B )
871 {
872  int ret, s = A->s;
873 
874  if( A->s * B->s > 0 )
875  {
876  if( mpi_cmp_abs( A, B ) >= 0 )
877  {
878  MPI_CHK( mpi_sub_abs( X, A, B ) );
879  X->s = s;
880  }
881  else
882  {
883  MPI_CHK( mpi_sub_abs( X, B, A ) );
884  X->s = -s;
885  }
886  }
887  else
888  {
889  MPI_CHK( mpi_add_abs( X, A, B ) );
890  X->s = s;
891  }
892 
893 cleanup:
894 
895  return( ret );
896 }
897 
898 /*
899  * Signed addition: X = A + b
900  */
901 int mpi_add_int( mpi *X, const mpi *A, t_sint b )
902 {
903  mpi _B;
904  t_uint p[1];
905 
906  p[0] = ( b < 0 ) ? -b : b;
907  _B.s = ( b < 0 ) ? -1 : 1;
908  _B.n = 1;
909  _B.p = p;
910 
911  return( mpi_add_mpi( X, A, &_B ) );
912 }
913 
914 /*
915  * Signed substraction: X = A - b
916  */
917 int mpi_sub_int( mpi *X, const mpi *A, t_sint b )
918 {
919  mpi _B;
920  t_uint p[1];
921 
922  p[0] = ( b < 0 ) ? -b : b;
923  _B.s = ( b < 0 ) ? -1 : 1;
924  _B.n = 1;
925  _B.p = p;
926 
927  return( mpi_sub_mpi( X, A, &_B ) );
928 }
929 
930 /*
931  * Helper for mpi multiplication
932  */
933 static void mpi_mul_hlp( size_t i, t_uint *s, t_uint *d, t_uint b )
934 {
935  t_uint c = 0, t = 0;
936 
937 #if defined(MULADDC_HUIT)
938  for( ; i >= 8; i -= 8 )
939  {
941  MULADDC_HUIT
943  }
944 
945  for( ; i > 0; i-- )
946  {
950  }
951 #else
952  for( ; i >= 16; i -= 16 )
953  {
959 
965  }
966 
967  for( ; i >= 8; i -= 8 )
968  {
972 
976  }
977 
978  for( ; i > 0; i-- )
979  {
983  }
984 #endif
985 
986  t++;
987 
988  do {
989  *d += c; c = ( *d < c ); d++;
990  }
991  while( c != 0 );
992 }
993 
994 /*
995  * Baseline multiplication: X = A * B (HAC 14.12)
996  */
997 int mpi_mul_mpi( mpi *X, const mpi *A, const mpi *B )
998 {
999  int ret;
1000  size_t i, j;
1001  mpi TA, TB;
1002 
1003  mpi_init( &TA ); mpi_init( &TB );
1004 
1005  if( X == A ) { MPI_CHK( mpi_copy( &TA, A ) ); A = &TA; }
1006  if( X == B ) { MPI_CHK( mpi_copy( &TB, B ) ); B = &TB; }
1007 
1008  for( i = A->n; i > 0; i-- )
1009  if( A->p[i - 1] != 0 )
1010  break;
1011 
1012  for( j = B->n; j > 0; j-- )
1013  if( B->p[j - 1] != 0 )
1014  break;
1015 
1016  MPI_CHK( mpi_grow( X, i + j ) );
1017  MPI_CHK( mpi_lset( X, 0 ) );
1018 
1019  for( i++; j > 0; j-- )
1020  mpi_mul_hlp( i - 1, A->p, X->p + j - 1, B->p[j - 1] );
1021 
1022  X->s = A->s * B->s;
1023 
1024 cleanup:
1025 
1026  mpi_free( &TB ); mpi_free( &TA );
1027 
1028  return( ret );
1029 }
1030 
1031 /*
1032  * Baseline multiplication: X = A * b
1033  */
1034 int mpi_mul_int( mpi *X, const mpi *A, t_sint b )
1035 {
1036  mpi _B;
1037  t_uint p[1];
1038 
1039  _B.s = 1;
1040  _B.n = 1;
1041  _B.p = p;
1042  p[0] = b;
1043 
1044  return( mpi_mul_mpi( X, A, &_B ) );
1045 }
1046 
1047 /*
1048  * Division by mpi: A = Q * B + R (HAC 14.20)
1049  */
1050 int mpi_div_mpi( mpi *Q, mpi *R, const mpi *A, const mpi *B )
1051 {
1052  int ret;
1053  size_t i, n, t, k;
1054  mpi X, Y, Z, T1, T2;
1055 
1056  if( mpi_cmp_int( B, 0 ) == 0 )
1058 
1059  mpi_init( &X ); mpi_init( &Y ); mpi_init( &Z );
1060  mpi_init( &T1 ); mpi_init( &T2 );
1061 
1062  if( mpi_cmp_abs( A, B ) < 0 )
1063  {
1064  if( Q != NULL ) MPI_CHK( mpi_lset( Q, 0 ) );
1065  if( R != NULL ) MPI_CHK( mpi_copy( R, A ) );
1066  return( 0 );
1067  }
1068 
1069  MPI_CHK( mpi_copy( &X, A ) );
1070  MPI_CHK( mpi_copy( &Y, B ) );
1071  X.s = Y.s = 1;
1072 
1073  MPI_CHK( mpi_grow( &Z, A->n + 2 ) );
1074  MPI_CHK( mpi_lset( &Z, 0 ) );
1075  MPI_CHK( mpi_grow( &T1, 2 ) );
1076  MPI_CHK( mpi_grow( &T2, 3 ) );
1077 
1078  k = mpi_msb( &Y ) % biL;
1079  if( k < biL - 1 )
1080  {
1081  k = biL - 1 - k;
1082  MPI_CHK( mpi_shift_l( &X, k ) );
1083  MPI_CHK( mpi_shift_l( &Y, k ) );
1084  }
1085  else k = 0;
1086 
1087  n = X.n - 1;
1088  t = Y.n - 1;
1089  mpi_shift_l( &Y, biL * (n - t) );
1090 
1091  while( mpi_cmp_mpi( &X, &Y ) >= 0 )
1092  {
1093  Z.p[n - t]++;
1094  mpi_sub_mpi( &X, &X, &Y );
1095  }
1096  mpi_shift_r( &Y, biL * (n - t) );
1097 
1098  for( i = n; i > t ; i-- )
1099  {
1100  if( X.p[i] >= Y.p[t] )
1101  Z.p[i - t - 1] = ~0;
1102  else
1103  {
1104 #if defined(POLARSSL_HAVE_LONGLONG)
1105  t_udbl r;
1106 
1107  r = (t_udbl) X.p[i] << biL;
1108  r |= (t_udbl) X.p[i - 1];
1109  r /= Y.p[t];
1110  if( r > ((t_udbl) 1 << biL) - 1)
1111  r = ((t_udbl) 1 << biL) - 1;
1112 
1113  Z.p[i - t - 1] = (t_uint) r;
1114 #else
1115  /*
1116  * __udiv_qrnnd_c, from gmp/longlong.h
1117  */
1118  t_uint q0, q1, r0, r1;
1119  t_uint d0, d1, d, m;
1120 
1121  d = Y.p[t];
1122  d0 = ( d << biH ) >> biH;
1123  d1 = ( d >> biH );
1124 
1125  q1 = X.p[i] / d1;
1126  r1 = X.p[i] - d1 * q1;
1127  r1 <<= biH;
1128  r1 |= ( X.p[i - 1] >> biH );
1129 
1130  m = q1 * d0;
1131  if( r1 < m )
1132  {
1133  q1--, r1 += d;
1134  while( r1 >= d && r1 < m )
1135  q1--, r1 += d;
1136  }
1137  r1 -= m;
1138 
1139  q0 = r1 / d1;
1140  r0 = r1 - d1 * q0;
1141  r0 <<= biH;
1142  r0 |= ( X.p[i - 1] << biH ) >> biH;
1143 
1144  m = q0 * d0;
1145  if( r0 < m )
1146  {
1147  q0--, r0 += d;
1148  while( r0 >= d && r0 < m )
1149  q0--, r0 += d;
1150  }
1151  r0 -= m;
1152 
1153  Z.p[i - t - 1] = ( q1 << biH ) | q0;
1154 #endif
1155  }
1156 
1157  Z.p[i - t - 1]++;
1158  do
1159  {
1160  Z.p[i - t - 1]--;
1161 
1162  MPI_CHK( mpi_lset( &T1, 0 ) );
1163  T1.p[0] = (t < 1) ? 0 : Y.p[t - 1];
1164  T1.p[1] = Y.p[t];
1165  MPI_CHK( mpi_mul_int( &T1, &T1, Z.p[i - t - 1] ) );
1166 
1167  MPI_CHK( mpi_lset( &T2, 0 ) );
1168  T2.p[0] = (i < 2) ? 0 : X.p[i - 2];
1169  T2.p[1] = (i < 1) ? 0 : X.p[i - 1];
1170  T2.p[2] = X.p[i];
1171  }
1172  while( mpi_cmp_mpi( &T1, &T2 ) > 0 );
1173 
1174  MPI_CHK( mpi_mul_int( &T1, &Y, Z.p[i - t - 1] ) );
1175  MPI_CHK( mpi_shift_l( &T1, biL * (i - t - 1) ) );
1176  MPI_CHK( mpi_sub_mpi( &X, &X, &T1 ) );
1177 
1178  if( mpi_cmp_int( &X, 0 ) < 0 )
1179  {
1180  MPI_CHK( mpi_copy( &T1, &Y ) );
1181  MPI_CHK( mpi_shift_l( &T1, biL * (i - t - 1) ) );
1182  MPI_CHK( mpi_add_mpi( &X, &X, &T1 ) );
1183  Z.p[i - t - 1]--;
1184  }
1185  }
1186 
1187  if( Q != NULL )
1188  {
1189  mpi_copy( Q, &Z );
1190  Q->s = A->s * B->s;
1191  }
1192 
1193  if( R != NULL )
1194  {
1195  mpi_shift_r( &X, k );
1196  mpi_copy( R, &X );
1197 
1198  R->s = A->s;
1199  if( mpi_cmp_int( R, 0 ) == 0 )
1200  R->s = 1;
1201  }
1202 
1203 cleanup:
1204 
1205  mpi_free( &X ); mpi_free( &Y ); mpi_free( &Z );
1206  mpi_free( &T1 ); mpi_free( &T2 );
1207 
1208  return( ret );
1209 }
1210 
1211 /*
1212  * Division by int: A = Q * b + R
1213  *
1214  * Returns 0 if successful
1215  * 1 if memory allocation failed
1216  * POLARSSL_ERR_MPI_DIVISION_BY_ZERO if b == 0
1217  */
1218 int mpi_div_int( mpi *Q, mpi *R, const mpi *A, t_sint b )
1219 {
1220  mpi _B;
1221  t_uint p[1];
1222 
1223  p[0] = ( b < 0 ) ? -b : b;
1224  _B.s = ( b < 0 ) ? -1 : 1;
1225  _B.n = 1;
1226  _B.p = p;
1227 
1228  return( mpi_div_mpi( Q, R, A, &_B ) );
1229 }
1230 
1231 /*
1232  * Modulo: R = A mod B
1233  */
1234 int mpi_mod_mpi( mpi *R, const mpi *A, const mpi *B )
1235 {
1236  int ret;
1237 
1238  if( mpi_cmp_int( B, 0 ) < 0 )
1240 
1241  MPI_CHK( mpi_div_mpi( NULL, R, A, B ) );
1242 
1243  while( mpi_cmp_int( R, 0 ) < 0 )
1244  MPI_CHK( mpi_add_mpi( R, R, B ) );
1245 
1246  while( mpi_cmp_mpi( R, B ) >= 0 )
1247  MPI_CHK( mpi_sub_mpi( R, R, B ) );
1248 
1249 cleanup:
1250 
1251  return( ret );
1252 }
1253 
1254 /*
1255  * Modulo: r = A mod b
1256  */
1257 int mpi_mod_int( t_uint *r, const mpi *A, t_sint b )
1258 {
1259  size_t i;
1260  t_uint x, y, z;
1261 
1262  if( b == 0 )
1264 
1265  if( b < 0 )
1267 
1268  /*
1269  * handle trivial cases
1270  */
1271  if( b == 1 )
1272  {
1273  *r = 0;
1274  return( 0 );
1275  }
1276 
1277  if( b == 2 )
1278  {
1279  *r = A->p[0] & 1;
1280  return( 0 );
1281  }
1282 
1283  /*
1284  * general case
1285  */
1286  for( i = A->n, y = 0; i > 0; i-- )
1287  {
1288  x = A->p[i - 1];
1289  y = ( y << biH ) | ( x >> biH );
1290  z = y / b;
1291  y -= z * b;
1292 
1293  x <<= biH;
1294  y = ( y << biH ) | ( x >> biH );
1295  z = y / b;
1296  y -= z * b;
1297  }
1298 
1299  /*
1300  * If A is negative, then the current y represents a negative value.
1301  * Flipping it to the positive side.
1302  */
1303  if( A->s < 0 && y != 0 )
1304  y = b - y;
1305 
1306  *r = y;
1307 
1308  return( 0 );
1309 }
1310 
1311 /*
1312  * Fast Montgomery initialization (thanks to Tom St Denis)
1313  */
1314 static void mpi_montg_init( t_uint *mm, const mpi *N )
1315 {
1316  t_uint x, m0 = N->p[0];
1317 
1318  x = m0;
1319  x += ( ( m0 + 2 ) & 4 ) << 1;
1320  x *= ( 2 - ( m0 * x ) );
1321 
1322  if( biL >= 16 ) x *= ( 2 - ( m0 * x ) );
1323  if( biL >= 32 ) x *= ( 2 - ( m0 * x ) );
1324  if( biL >= 64 ) x *= ( 2 - ( m0 * x ) );
1325 
1326  *mm = ~x + 1;
1327 }
1328 
1329 /*
1330  * Montgomery multiplication: A = A * B * R^-1 mod N (HAC 14.36)
1331  */
1332 static void mpi_montmul( mpi *A, const mpi *B, const mpi *N, t_uint mm, const mpi *T )
1333 {
1334  size_t i, n, m;
1335  t_uint u0, u1, *d;
1336 
1337  memset( T->p, 0, T->n * ciL );
1338 
1339  d = T->p;
1340  n = N->n;
1341  m = ( B->n < n ) ? B->n : n;
1342 
1343  for( i = 0; i < n; i++ )
1344  {
1345  /*
1346  * T = (T + u0*B + u1*N) / 2^biL
1347  */
1348  u0 = A->p[i];
1349  u1 = ( d[0] + u0 * B->p[0] ) * mm;
1350 
1351  mpi_mul_hlp( m, B->p, d, u0 );
1352  mpi_mul_hlp( n, N->p, d, u1 );
1353 
1354  *d++ = u0; d[n + 1] = 0;
1355  }
1356 
1357  memcpy( A->p, d, (n + 1) * ciL );
1358 
1359  if( mpi_cmp_abs( A, N ) >= 0 )
1360  mpi_sub_hlp( n, N->p, A->p );
1361  else
1362  /* prevent timing attacks */
1363  mpi_sub_hlp( n, A->p, T->p );
1364 }
1365 
1366 /*
1367  * Montgomery reduction: A = A * R^-1 mod N
1368  */
1369 static void mpi_montred( mpi *A, const mpi *N, t_uint mm, const mpi *T )
1370 {
1371  t_uint z = 1;
1372  mpi U;
1373 
1374  U.n = U.s = z;
1375  U.p = &z;
1376 
1377  mpi_montmul( A, &U, N, mm, T );
1378 }
1379 
1380 /*
1381  * Sliding-window exponentiation: X = A^E mod N (HAC 14.85)
1382  */
1383 int mpi_exp_mod( mpi *X, const mpi *A, const mpi *E, const mpi *N, mpi *_RR )
1384 {
1385  int ret;
1386  size_t wbits, wsize, one = 1;
1387  size_t i, j, nblimbs;
1388  size_t bufsize, nbits;
1389  t_uint ei, mm, state;
1390  mpi RR, T, W[ 2 << POLARSSL_MPI_WINDOW_SIZE ];
1391 
1392  if( mpi_cmp_int( N, 0 ) < 0 || ( N->p[0] & 1 ) == 0 )
1394 
1395  /*
1396  * Init temps and window size
1397  */
1398  mpi_montg_init( &mm, N );
1399  mpi_init( &RR ); mpi_init( &T );
1400  memset( W, 0, sizeof( W ) );
1401 
1402  i = mpi_msb( E );
1403 
1404  wsize = ( i > 671 ) ? 6 : ( i > 239 ) ? 5 :
1405  ( i > 79 ) ? 4 : ( i > 23 ) ? 3 : 1;
1406 
1407  if( wsize > POLARSSL_MPI_WINDOW_SIZE )
1408  wsize = POLARSSL_MPI_WINDOW_SIZE;
1409 
1410  j = N->n + 1;
1411  MPI_CHK( mpi_grow( X, j ) );
1412  MPI_CHK( mpi_grow( &W[1], j ) );
1413  MPI_CHK( mpi_grow( &T, j * 2 ) );
1414 
1415  /*
1416  * If 1st call, pre-compute R^2 mod N
1417  */
1418  if( _RR == NULL || _RR->p == NULL )
1419  {
1420  MPI_CHK( mpi_lset( &RR, 1 ) );
1421  MPI_CHK( mpi_shift_l( &RR, N->n * 2 * biL ) );
1422  MPI_CHK( mpi_mod_mpi( &RR, &RR, N ) );
1423 
1424  if( _RR != NULL )
1425  memcpy( _RR, &RR, sizeof( mpi ) );
1426  }
1427  else
1428  memcpy( &RR, _RR, sizeof( mpi ) );
1429 
1430  /*
1431  * W[1] = A * R^2 * R^-1 mod N = A * R mod N
1432  */
1433  if( mpi_cmp_mpi( A, N ) >= 0 )
1434  mpi_mod_mpi( &W[1], A, N );
1435  else mpi_copy( &W[1], A );
1436 
1437  mpi_montmul( &W[1], &RR, N, mm, &T );
1438 
1439  /*
1440  * X = R^2 * R^-1 mod N = R mod N
1441  */
1442  MPI_CHK( mpi_copy( X, &RR ) );
1443  mpi_montred( X, N, mm, &T );
1444 
1445  if( wsize > 1 )
1446  {
1447  /*
1448  * W[1 << (wsize - 1)] = W[1] ^ (wsize - 1)
1449  */
1450  j = one << (wsize - 1);
1451 
1452  MPI_CHK( mpi_grow( &W[j], N->n + 1 ) );
1453  MPI_CHK( mpi_copy( &W[j], &W[1] ) );
1454 
1455  for( i = 0; i < wsize - 1; i++ )
1456  mpi_montmul( &W[j], &W[j], N, mm, &T );
1457 
1458  /*
1459  * W[i] = W[i - 1] * W[1]
1460  */
1461  for( i = j + 1; i < (one << wsize); i++ )
1462  {
1463  MPI_CHK( mpi_grow( &W[i], N->n + 1 ) );
1464  MPI_CHK( mpi_copy( &W[i], &W[i - 1] ) );
1465 
1466  mpi_montmul( &W[i], &W[1], N, mm, &T );
1467  }
1468  }
1469 
1470  nblimbs = E->n;
1471  bufsize = 0;
1472  nbits = 0;
1473  wbits = 0;
1474  state = 0;
1475 
1476  while( 1 )
1477  {
1478  if( bufsize == 0 )
1479  {
1480  if( nblimbs-- == 0 )
1481  break;
1482 
1483  bufsize = sizeof( t_uint ) << 3;
1484  }
1485 
1486  bufsize--;
1487 
1488  ei = (E->p[nblimbs] >> bufsize) & 1;
1489 
1490  /*
1491  * skip leading 0s
1492  */
1493  if( ei == 0 && state == 0 )
1494  continue;
1495 
1496  if( ei == 0 && state == 1 )
1497  {
1498  /*
1499  * out of window, square X
1500  */
1501  mpi_montmul( X, X, N, mm, &T );
1502  continue;
1503  }
1504 
1505  /*
1506  * add ei to current window
1507  */
1508  state = 2;
1509 
1510  nbits++;
1511  wbits |= (ei << (wsize - nbits));
1512 
1513  if( nbits == wsize )
1514  {
1515  /*
1516  * X = X^wsize R^-1 mod N
1517  */
1518  for( i = 0; i < wsize; i++ )
1519  mpi_montmul( X, X, N, mm, &T );
1520 
1521  /*
1522  * X = X * W[wbits] R^-1 mod N
1523  */
1524  mpi_montmul( X, &W[wbits], N, mm, &T );
1525 
1526  state--;
1527  nbits = 0;
1528  wbits = 0;
1529  }
1530  }
1531 
1532  /*
1533  * process the remaining bits
1534  */
1535  for( i = 0; i < nbits; i++ )
1536  {
1537  mpi_montmul( X, X, N, mm, &T );
1538 
1539  wbits <<= 1;
1540 
1541  if( (wbits & (one << wsize)) != 0 )
1542  mpi_montmul( X, &W[1], N, mm, &T );
1543  }
1544 
1545  /*
1546  * X = A^E * R * R^-1 mod N = A^E mod N
1547  */
1548  mpi_montred( X, N, mm, &T );
1549 
1550 cleanup:
1551 
1552  for( i = (one << (wsize - 1)); i < (one << wsize); i++ )
1553  mpi_free( &W[i] );
1554 
1555  mpi_free( &W[1] ); mpi_free( &T );
1556 
1557  if( _RR == NULL )
1558  mpi_free( &RR );
1559 
1560  return( ret );
1561 }
1562 
1563 /*
1564  * Greatest common divisor: G = gcd(A, B) (HAC 14.54)
1565  */
1566 int mpi_gcd( mpi *G, const mpi *A, const mpi *B )
1567 {
1568  int ret;
1569  size_t lz, lzt;
1570  mpi TG, TA, TB;
1571 
1572  mpi_init( &TG ); mpi_init( &TA ); mpi_init( &TB );
1573 
1574  MPI_CHK( mpi_copy( &TA, A ) );
1575  MPI_CHK( mpi_copy( &TB, B ) );
1576 
1577  lz = mpi_lsb( &TA );
1578  lzt = mpi_lsb( &TB );
1579 
1580  if ( lzt < lz )
1581  lz = lzt;
1582 
1583  MPI_CHK( mpi_shift_r( &TA, lz ) );
1584  MPI_CHK( mpi_shift_r( &TB, lz ) );
1585 
1586  TA.s = TB.s = 1;
1587 
1588  while( mpi_cmp_int( &TA, 0 ) != 0 )
1589  {
1590  MPI_CHK( mpi_shift_r( &TA, mpi_lsb( &TA ) ) );
1591  MPI_CHK( mpi_shift_r( &TB, mpi_lsb( &TB ) ) );
1592 
1593  if( mpi_cmp_mpi( &TA, &TB ) >= 0 )
1594  {
1595  MPI_CHK( mpi_sub_abs( &TA, &TA, &TB ) );
1596  MPI_CHK( mpi_shift_r( &TA, 1 ) );
1597  }
1598  else
1599  {
1600  MPI_CHK( mpi_sub_abs( &TB, &TB, &TA ) );
1601  MPI_CHK( mpi_shift_r( &TB, 1 ) );
1602  }
1603  }
1604 
1605  MPI_CHK( mpi_shift_l( &TB, lz ) );
1606  MPI_CHK( mpi_copy( G, &TB ) );
1607 
1608 cleanup:
1609 
1610  mpi_free( &TG ); mpi_free( &TA ); mpi_free( &TB );
1611 
1612  return( ret );
1613 }
1614 
1615 int mpi_fill_random( mpi *X, size_t size,
1616  int (*f_rng)(void *, unsigned char *, size_t),
1617  void *p_rng )
1618 {
1619  int ret;
1620 
1621  MPI_CHK( mpi_grow( X, CHARS_TO_LIMBS( size ) ) );
1622  MPI_CHK( mpi_lset( X, 0 ) );
1623 
1624  MPI_CHK( f_rng( p_rng, (unsigned char *) X->p, size ) );
1625 
1626 cleanup:
1627  return( ret );
1628 }
1629 
1630 #if defined(POLARSSL_GENPRIME)
1631 
1632 /*
1633  * Modular inverse: X = A^-1 mod N (HAC 14.61 / 14.64)
1634  */
1635 int mpi_inv_mod( mpi *X, const mpi *A, const mpi *N )
1636 {
1637  int ret;
1638  mpi G, TA, TU, U1, U2, TB, TV, V1, V2;
1639 
1640  if( mpi_cmp_int( N, 0 ) <= 0 )
1642 
1643  mpi_init( &TA ); mpi_init( &TU ); mpi_init( &U1 ); mpi_init( &U2 );
1644  mpi_init( &G ); mpi_init( &TB ); mpi_init( &TV );
1645  mpi_init( &V1 ); mpi_init( &V2 );
1646 
1647  MPI_CHK( mpi_gcd( &G, A, N ) );
1648 
1649  if( mpi_cmp_int( &G, 1 ) != 0 )
1650  {
1652  goto cleanup;
1653  }
1654 
1655  MPI_CHK( mpi_mod_mpi( &TA, A, N ) );
1656  MPI_CHK( mpi_copy( &TU, &TA ) );
1657  MPI_CHK( mpi_copy( &TB, N ) );
1658  MPI_CHK( mpi_copy( &TV, N ) );
1659 
1660  MPI_CHK( mpi_lset( &U1, 1 ) );
1661  MPI_CHK( mpi_lset( &U2, 0 ) );
1662  MPI_CHK( mpi_lset( &V1, 0 ) );
1663  MPI_CHK( mpi_lset( &V2, 1 ) );
1664 
1665  do
1666  {
1667  while( ( TU.p[0] & 1 ) == 0 )
1668  {
1669  MPI_CHK( mpi_shift_r( &TU, 1 ) );
1670 
1671  if( ( U1.p[0] & 1 ) != 0 || ( U2.p[0] & 1 ) != 0 )
1672  {
1673  MPI_CHK( mpi_add_mpi( &U1, &U1, &TB ) );
1674  MPI_CHK( mpi_sub_mpi( &U2, &U2, &TA ) );
1675  }
1676 
1677  MPI_CHK( mpi_shift_r( &U1, 1 ) );
1678  MPI_CHK( mpi_shift_r( &U2, 1 ) );
1679  }
1680 
1681  while( ( TV.p[0] & 1 ) == 0 )
1682  {
1683  MPI_CHK( mpi_shift_r( &TV, 1 ) );
1684 
1685  if( ( V1.p[0] & 1 ) != 0 || ( V2.p[0] & 1 ) != 0 )
1686  {
1687  MPI_CHK( mpi_add_mpi( &V1, &V1, &TB ) );
1688  MPI_CHK( mpi_sub_mpi( &V2, &V2, &TA ) );
1689  }
1690 
1691  MPI_CHK( mpi_shift_r( &V1, 1 ) );
1692  MPI_CHK( mpi_shift_r( &V2, 1 ) );
1693  }
1694 
1695  if( mpi_cmp_mpi( &TU, &TV ) >= 0 )
1696  {
1697  MPI_CHK( mpi_sub_mpi( &TU, &TU, &TV ) );
1698  MPI_CHK( mpi_sub_mpi( &U1, &U1, &V1 ) );
1699  MPI_CHK( mpi_sub_mpi( &U2, &U2, &V2 ) );
1700  }
1701  else
1702  {
1703  MPI_CHK( mpi_sub_mpi( &TV, &TV, &TU ) );
1704  MPI_CHK( mpi_sub_mpi( &V1, &V1, &U1 ) );
1705  MPI_CHK( mpi_sub_mpi( &V2, &V2, &U2 ) );
1706  }
1707  }
1708  while( mpi_cmp_int( &TU, 0 ) != 0 );
1709 
1710  while( mpi_cmp_int( &V1, 0 ) < 0 )
1711  MPI_CHK( mpi_add_mpi( &V1, &V1, N ) );
1712 
1713  while( mpi_cmp_mpi( &V1, N ) >= 0 )
1714  MPI_CHK( mpi_sub_mpi( &V1, &V1, N ) );
1715 
1716  MPI_CHK( mpi_copy( X, &V1 ) );
1717 
1718 cleanup:
1719 
1720  mpi_free( &TA ); mpi_free( &TU ); mpi_free( &U1 ); mpi_free( &U2 );
1721  mpi_free( &G ); mpi_free( &TB ); mpi_free( &TV );
1722  mpi_free( &V1 ); mpi_free( &V2 );
1723 
1724  return( ret );
1725 }
1726 
1727 static const int small_prime[] =
1728 {
1729  3, 5, 7, 11, 13, 17, 19, 23,
1730  29, 31, 37, 41, 43, 47, 53, 59,
1731  61, 67, 71, 73, 79, 83, 89, 97,
1732  101, 103, 107, 109, 113, 127, 131, 137,
1733  139, 149, 151, 157, 163, 167, 173, 179,
1734  181, 191, 193, 197, 199, 211, 223, 227,
1735  229, 233, 239, 241, 251, 257, 263, 269,
1736  271, 277, 281, 283, 293, 307, 311, 313,
1737  317, 331, 337, 347, 349, 353, 359, 367,
1738  373, 379, 383, 389, 397, 401, 409, 419,
1739  421, 431, 433, 439, 443, 449, 457, 461,
1740  463, 467, 479, 487, 491, 499, 503, 509,
1741  521, 523, 541, 547, 557, 563, 569, 571,
1742  577, 587, 593, 599, 601, 607, 613, 617,
1743  619, 631, 641, 643, 647, 653, 659, 661,
1744  673, 677, 683, 691, 701, 709, 719, 727,
1745  733, 739, 743, 751, 757, 761, 769, 773,
1746  787, 797, 809, 811, 821, 823, 827, 829,
1747  839, 853, 857, 859, 863, 877, 881, 883,
1748  887, 907, 911, 919, 929, 937, 941, 947,
1749  953, 967, 971, 977, 983, 991, 997, -103
1750 };
1751 
1752 /*
1753  * Miller-Rabin primality test (HAC 4.24)
1754  */
1755 int mpi_is_prime( mpi *X,
1756  int (*f_rng)(void *, unsigned char *, size_t),
1757  void *p_rng )
1758 {
1759  int ret, xs;
1760  size_t i, j, n, s;
1761  mpi W, R, T, A, RR;
1762 
1763  if( mpi_cmp_int( X, 0 ) == 0 ||
1764  mpi_cmp_int( X, 1 ) == 0 )
1766 
1767  if( mpi_cmp_int( X, 2 ) == 0 )
1768  return( 0 );
1769 
1770  mpi_init( &W ); mpi_init( &R ); mpi_init( &T ); mpi_init( &A );
1771  mpi_init( &RR );
1772 
1773  xs = X->s; X->s = 1;
1774 
1775  /*
1776  * test trivial factors first
1777  */
1778  if( ( X->p[0] & 1 ) == 0 )
1780 
1781  for( i = 0; small_prime[i] > 0; i++ )
1782  {
1783  t_uint r;
1784 
1785  if( mpi_cmp_int( X, small_prime[i] ) <= 0 )
1786  return( 0 );
1787 
1788  MPI_CHK( mpi_mod_int( &r, X, small_prime[i] ) );
1789 
1790  if( r == 0 )
1792  }
1793 
1794  /*
1795  * W = |X| - 1
1796  * R = W >> lsb( W )
1797  */
1798  MPI_CHK( mpi_sub_int( &W, X, 1 ) );
1799  s = mpi_lsb( &W );
1800  MPI_CHK( mpi_copy( &R, &W ) );
1801  MPI_CHK( mpi_shift_r( &R, s ) );
1802 
1803  i = mpi_msb( X );
1804  /*
1805  * HAC, table 4.4
1806  */
1807  n = ( ( i >= 1300 ) ? 2 : ( i >= 850 ) ? 3 :
1808  ( i >= 650 ) ? 4 : ( i >= 350 ) ? 8 :
1809  ( i >= 250 ) ? 12 : ( i >= 150 ) ? 18 : 27 );
1810 
1811  for( i = 0; i < n; i++ )
1812  {
1813  /*
1814  * pick a random A, 1 < A < |X| - 1
1815  */
1816  MPI_CHK( mpi_fill_random( &A, X->n * ciL, f_rng, p_rng ) );
1817 
1818  if( mpi_cmp_mpi( &A, &W ) >= 0 )
1819  {
1820  j = mpi_msb( &A ) - mpi_msb( &W );
1821  MPI_CHK( mpi_shift_r( &A, j + 1 ) );
1822  }
1823  A.p[0] |= 3;
1824 
1825  /*
1826  * A = A^R mod |X|
1827  */
1828  MPI_CHK( mpi_exp_mod( &A, &A, &R, X, &RR ) );
1829 
1830  if( mpi_cmp_mpi( &A, &W ) == 0 ||
1831  mpi_cmp_int( &A, 1 ) == 0 )
1832  continue;
1833 
1834  j = 1;
1835  while( j < s && mpi_cmp_mpi( &A, &W ) != 0 )
1836  {
1837  /*
1838  * A = A * A mod |X|
1839  */
1840  MPI_CHK( mpi_mul_mpi( &T, &A, &A ) );
1841  MPI_CHK( mpi_mod_mpi( &A, &T, X ) );
1842 
1843  if( mpi_cmp_int( &A, 1 ) == 0 )
1844  break;
1845 
1846  j++;
1847  }
1848 
1849  /*
1850  * not prime if A != |X| - 1 or A == 1
1851  */
1852  if( mpi_cmp_mpi( &A, &W ) != 0 ||
1853  mpi_cmp_int( &A, 1 ) == 0 )
1854  {
1856  break;
1857  }
1858  }
1859 
1860 cleanup:
1861 
1862  X->s = xs;
1863 
1864  mpi_free( &W ); mpi_free( &R ); mpi_free( &T ); mpi_free( &A );
1865  mpi_free( &RR );
1866 
1867  return( ret );
1868 }
1869 
1870 /*
1871  * Prime number generation
1872  */
1873 int mpi_gen_prime( mpi *X, size_t nbits, int dh_flag,
1874  int (*f_rng)(void *, unsigned char *, size_t),
1875  void *p_rng )
1876 {
1877  int ret;
1878  size_t k, n;
1879  mpi Y;
1880 
1881  if( nbits < 3 || nbits > POLARSSL_MPI_MAX_BITS )
1883 
1884  mpi_init( &Y );
1885 
1886  n = BITS_TO_LIMBS( nbits );
1887 
1888  MPI_CHK( mpi_fill_random( X, n * ciL, f_rng, p_rng ) );
1889 
1890  k = mpi_msb( X );
1891  if( k < nbits ) MPI_CHK( mpi_shift_l( X, nbits - k ) );
1892  if( k > nbits ) MPI_CHK( mpi_shift_r( X, k - nbits ) );
1893 
1894  X->p[0] |= 3;
1895 
1896  if( dh_flag == 0 )
1897  {
1898  while( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) != 0 )
1899  {
1900  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1901  goto cleanup;
1902 
1903  MPI_CHK( mpi_add_int( X, X, 2 ) );
1904  }
1905  }
1906  else
1907  {
1908  MPI_CHK( mpi_sub_int( &Y, X, 1 ) );
1909  MPI_CHK( mpi_shift_r( &Y, 1 ) );
1910 
1911  while( 1 )
1912  {
1913  if( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) == 0 )
1914  {
1915  if( ( ret = mpi_is_prime( &Y, f_rng, p_rng ) ) == 0 )
1916  break;
1917 
1918  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1919  goto cleanup;
1920  }
1921 
1922  if( ret != POLARSSL_ERR_MPI_NOT_ACCEPTABLE )
1923  goto cleanup;
1924 
1925  MPI_CHK( mpi_add_int( &Y, X, 1 ) );
1926  MPI_CHK( mpi_add_int( X, X, 2 ) );
1927  MPI_CHK( mpi_shift_r( &Y, 1 ) );
1928  }
1929  }
1930 
1931 cleanup:
1932 
1933  mpi_free( &Y );
1934 
1935  return( ret );
1936 }
1937 
1938 #endif
1939 
1940 #if defined(POLARSSL_SELF_TEST)
1941 
1942 #define GCD_PAIR_COUNT 3
1943 
1944 static const int gcd_pairs[GCD_PAIR_COUNT][3] =
1945 {
1946  { 693, 609, 21 },
1947  { 1764, 868, 28 },
1948  { 768454923, 542167814, 1 }
1949 };
1950 
1951 /*
1952  * Checkup routine
1953  */
1954 int mpi_self_test( int verbose )
1955 {
1956  int ret, i;
1957  mpi A, E, N, X, Y, U, V;
1958 
1959  mpi_init( &A ); mpi_init( &E ); mpi_init( &N ); mpi_init( &X );
1960  mpi_init( &Y ); mpi_init( &U ); mpi_init( &V );
1961 
1962  MPI_CHK( mpi_read_string( &A, 16,
1963  "EFE021C2645FD1DC586E69184AF4A31E" \
1964  "D5F53E93B5F123FA41680867BA110131" \
1965  "944FE7952E2517337780CB0DB80E61AA" \
1966  "E7C8DDC6C5C6AADEB34EB38A2F40D5E6" ) );
1967 
1968  MPI_CHK( mpi_read_string( &E, 16,
1969  "B2E7EFD37075B9F03FF989C7C5051C20" \
1970  "34D2A323810251127E7BF8625A4F49A5" \
1971  "F3E27F4DA8BD59C47D6DAABA4C8127BD" \
1972  "5B5C25763222FEFCCFC38B832366C29E" ) );
1973 
1974  MPI_CHK( mpi_read_string( &N, 16,
1975  "0066A198186C18C10B2F5ED9B522752A" \
1976  "9830B69916E535C8F047518A889A43A5" \
1977  "94B6BED27A168D31D4A52F88925AA8F5" ) );
1978 
1979  MPI_CHK( mpi_mul_mpi( &X, &A, &N ) );
1980 
1981  MPI_CHK( mpi_read_string( &U, 16,
1982  "602AB7ECA597A3D6B56FF9829A5E8B85" \
1983  "9E857EA95A03512E2BAE7391688D264A" \
1984  "A5663B0341DB9CCFD2C4C5F421FEC814" \
1985  "8001B72E848A38CAE1C65F78E56ABDEF" \
1986  "E12D3C039B8A02D6BE593F0BBBDA56F1" \
1987  "ECF677152EF804370C1A305CAF3B5BF1" \
1988  "30879B56C61DE584A0F53A2447A51E" ) );
1989 
1990  if( verbose != 0 )
1991  printf( " MPI test #1 (mul_mpi): " );
1992 
1993  if( mpi_cmp_mpi( &X, &U ) != 0 )
1994  {
1995  if( verbose != 0 )
1996  printf( "failed\n" );
1997 
1998  return( 1 );
1999  }
2000 
2001  if( verbose != 0 )
2002  printf( "passed\n" );
2003 
2004  MPI_CHK( mpi_div_mpi( &X, &Y, &A, &N ) );
2005 
2006  MPI_CHK( mpi_read_string( &U, 16,
2007  "256567336059E52CAE22925474705F39A94" ) );
2008 
2009  MPI_CHK( mpi_read_string( &V, 16,
2010  "6613F26162223DF488E9CD48CC132C7A" \
2011  "0AC93C701B001B092E4E5B9F73BCD27B" \
2012  "9EE50D0657C77F374E903CDFA4C642" ) );
2013 
2014  if( verbose != 0 )
2015  printf( " MPI test #2 (div_mpi): " );
2016 
2017  if( mpi_cmp_mpi( &X, &U ) != 0 ||
2018  mpi_cmp_mpi( &Y, &V ) != 0 )
2019  {
2020  if( verbose != 0 )
2021  printf( "failed\n" );
2022 
2023  return( 1 );
2024  }
2025 
2026  if( verbose != 0 )
2027  printf( "passed\n" );
2028 
2029  MPI_CHK( mpi_exp_mod( &X, &A, &E, &N, NULL ) );
2030 
2031  MPI_CHK( mpi_read_string( &U, 16,
2032  "36E139AEA55215609D2816998ED020BB" \
2033  "BD96C37890F65171D948E9BC7CBAA4D9" \
2034  "325D24D6A3C12710F10A09FA08AB87" ) );
2035 
2036  if( verbose != 0 )
2037  printf( " MPI test #3 (exp_mod): " );
2038 
2039  if( mpi_cmp_mpi( &X, &U ) != 0 )
2040  {
2041  if( verbose != 0 )
2042  printf( "failed\n" );
2043 
2044  return( 1 );
2045  }
2046 
2047  if( verbose != 0 )
2048  printf( "passed\n" );
2049 
2050 #if defined(POLARSSL_GENPRIME)
2051  MPI_CHK( mpi_inv_mod( &X, &A, &N ) );
2052 
2053  MPI_CHK( mpi_read_string( &U, 16,
2054  "003A0AAEDD7E784FC07D8F9EC6E3BFD5" \
2055  "C3DBA76456363A10869622EAC2DD84EC" \
2056  "C5B8A74DAC4D09E03B5E0BE779F2DF61" ) );
2057 
2058  if( verbose != 0 )
2059  printf( " MPI test #4 (inv_mod): " );
2060 
2061  if( mpi_cmp_mpi( &X, &U ) != 0 )
2062  {
2063  if( verbose != 0 )
2064  printf( "failed\n" );
2065 
2066  return( 1 );
2067  }
2068 
2069  if( verbose != 0 )
2070  printf( "passed\n" );
2071 #endif
2072 
2073  if( verbose != 0 )
2074  printf( " MPI test #5 (simple gcd): " );
2075 
2076  for ( i = 0; i < GCD_PAIR_COUNT; i++)
2077  {
2078  MPI_CHK( mpi_lset( &X, gcd_pairs[i][0] ) );
2079  MPI_CHK( mpi_lset( &Y, gcd_pairs[i][1] ) );
2080 
2081  MPI_CHK( mpi_gcd( &A, &X, &Y ) );
2082 
2083  if( mpi_cmp_int( &A, gcd_pairs[i][2] ) != 0 )
2084  {
2085  if( verbose != 0 )
2086  printf( "failed at %d\n", i );
2087 
2088  return( 1 );
2089  }
2090  }
2091 
2092  if( verbose != 0 )
2093  printf( "passed\n" );
2094 
2095 cleanup:
2096 
2097  if( ret != 0 && verbose != 0 )
2098  printf( "Unexpected error, return code = %08X\n", ret );
2099 
2100  mpi_free( &A ); mpi_free( &E ); mpi_free( &N ); mpi_free( &X );
2101  mpi_free( &Y ); mpi_free( &U ); mpi_free( &V );
2102 
2103  if( verbose != 0 )
2104  printf( "\n" );
2105 
2106  return( ret );
2107 }
2108 
2109 #endif
2110 
2111 #endif