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