Mobile-C: /home/dko/projects/mobilec/trunk/src/security/xyssl-0.9/library/bignum.c Source File

00001 /*
00002  *  Multi-precision integer library
00003  *
00004  *  Copyright (C) 2006-2007  Christophe Devine
00005  *
00006  *  This program is free software; you can redistribute it and/or modify
00007  *  it under the terms of the GNU General Public License as published by
00008  *  the Free Software Foundation; either version 2 of the License, or
00009  *  (at your option) any later version.
00010  *
00011  *  This program is distributed in the hope that it will be useful,
00012  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00014  *  GNU General Public License for more details.
00015  *
00016  *  You should have received a copy of the GNU General Public License along
00017  *  with this program; if not, write to the Free Software Foundation, Inc.,
00018  *  51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
00019  */
00020 /*
00021  *  This MPI implementation is based on:
00022  *
00023  *  http://www.cacr.math.uwaterloo.ca/hac/about/chap14.pdf
00024  *  http://www.stillhq.com/extracted/gnupg-api/mpi/
00025  *  http://math.libtomcrypt.com/files/tommath.pdf
00026  */
00027 
00028 #include "xyssl/config.h"
00029 
00030 #if defined(XYSSL_BIGNUM_C)
00031 
00032 #include "xyssl/bignum.h"
00033 #include "xyssl/bn_mul.h"
00034 
00035 #include <string.h>
00036 #include <stdlib.h>
00037 #include <stdarg.h>
00038 
00039 #define ciL    ((int) sizeof(t_int))    /* chars in limb  */
00040 #define biL    (ciL << 3)               /* bits  in limb  */
00041 #define biH    (ciL << 2)               /* half limb size */
00042 
00043 /*
00044  * Convert between bits/chars and number of limbs
00045  */
00046 #define BITS_TO_LIMBS(i)  (((i) + biL - 1) / biL)
00047 #define CHARS_TO_LIMBS(i) (((i) + ciL - 1) / ciL)
00048 
00049 /*
00050  * Initialize one or more mpi
00051  */
00052 void mpi_init( mpi *X, ... )
00053 {
00054     va_list args;
00055 
00056     va_start( args, X );
00057 
00058     while( X != NULL )
00059     {
00060         X->s = 1;
00061         X->n = 0;
00062         X->p = NULL;
00063 
00064         X = va_arg( args, mpi* );
00065     }
00066 
00067     va_end( args );
00068 }
00069 
00070 /*
00071  * Unallocate one or more mpi
00072  */
00073 void mpi_free( mpi *X, ... )
00074 {
00075     va_list args;
00076 
00077     va_start( args, X );
00078 
00079     while( X != NULL )
00080     {
00081         if( X->p != NULL )
00082         {
00083             memset( X->p, 0, X->n * ciL );
00084             free( X->p );
00085         }
00086 
00087         X->s = 1;
00088         X->n = 0;
00089         X->p = NULL;
00090 
00091         X = va_arg( args, mpi* );
00092     }
00093 
00094     va_end( args );
00095 }
00096 
00097 /*
00098  * Enlarge to the specified number of limbs
00099  */
00100 int mpi_grow( mpi *X, int nblimbs )
00101 {
00102     t_int *p;
00103 
00104     if( X->n < nblimbs )
00105     {
00106         if( ( p = (t_int *) malloc( nblimbs * ciL ) ) == NULL )
00107             return( 1 );
00108 
00109         memset( p, 0, nblimbs * ciL );
00110 
00111         if( X->p != NULL )
00112         {
00113             memcpy( p, X->p, X->n * ciL );
00114             memset( X->p, 0, X->n * ciL );
00115             free( X->p );
00116         }
00117 
00118         X->n = nblimbs;
00119         X->p = p;
00120     }
00121 
00122     return( 0 );
00123 }
00124 
00125 /*
00126  * Copy the contents of Y into X
00127  */
00128 int mpi_copy( mpi *X, mpi *Y )
00129 {
00130     int ret, i;
00131 
00132     if( X == Y )
00133         return( 0 );
00134 
00135     for( i = Y->n - 1; i > 0; i-- )
00136         if( Y->p[i] != 0 )
00137             break;
00138     i++;
00139 
00140     X->s = Y->s;
00141 
00142     MPI_CHK( mpi_grow( X, i ) );
00143 
00144     memset( X->p, 0, X->n * ciL );
00145     memcpy( X->p, Y->p, i * ciL );
00146 
00147 cleanup:
00148 
00149     return( ret );
00150 }
00151 
00152 /*
00153  * Swap the contents of X and Y
00154  */
00155 void mpi_swap( mpi *X, mpi *Y )
00156 {
00157     mpi T;
00158 
00159     memcpy( &T,  X, sizeof( mpi ) );
00160     memcpy(  X,  Y, sizeof( mpi ) );
00161     memcpy(  Y, &T, sizeof( mpi ) );
00162 }
00163 
00164 /*
00165  * Set value from integer
00166  */
00167 int mpi_lset( mpi *X, int z )
00168 {
00169     int ret;
00170 
00171     MPI_CHK( mpi_grow( X, 1 ) );
00172     memset( X->p, 0, X->n * ciL );
00173 
00174     X->p[0] = ( z < 0 ) ? -z : z;
00175     X->s    = ( z < 0 ) ? -1 : 1;
00176 
00177 cleanup:
00178 
00179     return( ret );
00180 }
00181 
00182 /*
00183  * Return the number of least significant bits
00184  */
00185 int mpi_lsb( mpi *X )
00186 {
00187     int i, j, count = 0;
00188 
00189     for( i = 0; i < X->n; i++ )
00190         for( j = 0; j < (int) biL; j++, count++ )
00191             if( ( ( X->p[i] >> j ) & 1 ) != 0 )
00192                 return( count );
00193 
00194     return( 0 );
00195 }
00196 
00197 /*
00198  * Return the number of most significant bits
00199  */
00200 int mpi_msb( mpi *X )
00201 {
00202     int i, j;
00203 
00204     for( i = X->n - 1; i > 0; i-- )
00205         if( X->p[i] != 0 )
00206             break;
00207 
00208     for( j = biL - 1; j >= 0; j-- )
00209         if( ( ( X->p[i] >> j ) & 1 ) != 0 )
00210             break;
00211 
00212     return( ( i * biL ) + j + 1 );
00213 }
00214 
00215 /*
00216  * Return the total size in bytes
00217  */
00218 int mpi_size( mpi *X )
00219 {
00220     return( ( mpi_msb( X ) + 7 ) >> 3 );
00221 }
00222 
00223 /*
00224  * Convert an ASCII character to digit value
00225  */
00226 static int mpi_get_digit( t_int *d, int radix, char c )
00227 {
00228     *d = 255;
00229 
00230     if( c >= 0x30 && c <= 0x39 ) *d = c - 0x30;
00231     if( c >= 0x41 && c <= 0x46 ) *d = c - 0x37;
00232     if( c >= 0x61 && c <= 0x66 ) *d = c - 0x57;
00233 
00234     if( *d >= (t_int) radix )
00235         return( XYSSL_ERR_MPI_INVALID_CHARACTER );
00236 
00237     return( 0 );
00238 }
00239 
00240 /*
00241  * Import from an ASCII string
00242  */
00243 int mpi_read_string( mpi *X, int radix, char *s )
00244 {
00245     int ret, i, j, n;
00246     t_int d;
00247     mpi T;
00248 
00249     if( radix < 2 || radix > 16 )
00250         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
00251 
00252     mpi_init( &T, NULL );
00253 
00254     if( radix == 16 )
00255     {
00256         n = BITS_TO_LIMBS( strlen( s ) << 2 );
00257 
00258         MPI_CHK( mpi_grow( X, n ) );
00259         MPI_CHK( mpi_lset( X, 0 ) );
00260 
00261         for( i = strlen( s ) - 1, j = 0; i >= 0; i--, j++ )
00262         {
00263             if( i == 0 && s[i] == '-' )
00264             {
00265                 X->s = -1;
00266                 break;
00267             }
00268 
00269             MPI_CHK( mpi_get_digit( &d, radix, s[i] ) );
00270             X->p[j / (2 * ciL)] |= d << ( (j % (2 * ciL)) << 2 );
00271         }
00272     }
00273     else
00274     {
00275         MPI_CHK( mpi_lset( X, 0 ) );
00276 
00277         for( i = 0; i < (int) strlen( s ); i++ )
00278         {
00279             if( i == 0 && s[i] == '-' )
00280             {
00281                 X->s = -1;
00282                 continue;
00283             }
00284 
00285             MPI_CHK( mpi_get_digit( &d, radix, s[i] ) );
00286             MPI_CHK( mpi_mul_int( &T, X, radix ) );
00287             MPI_CHK( mpi_add_int( X, &T, d ) );
00288         }
00289     }
00290 
00291 cleanup:
00292 
00293     mpi_free( &T, NULL );
00294 
00295     return( ret );
00296 }
00297 
00298 /*
00299  * Helper to write the digits high-order first
00300  */
00301 static int mpi_write_hlp( mpi *X, int radix, char **p )
00302 {
00303     int ret;
00304     t_int r;
00305 
00306     if( radix < 2 || radix > 16 )
00307         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
00308 
00309     MPI_CHK( mpi_mod_int( &r, X, radix ) );
00310     MPI_CHK( mpi_div_int( X, NULL, X, radix ) );
00311 
00312     if( mpi_cmp_int( X, 0 ) != 0 )
00313         MPI_CHK( mpi_write_hlp( X, radix, p ) );
00314 
00315     if( r < 10 )
00316         *(*p)++ = (char)( r + 0x30 );
00317     else
00318         *(*p)++ = (char)( r + 0x37 );
00319 
00320 cleanup:
00321 
00322     return( ret );
00323 }
00324 
00325 /*
00326  * Export into an ASCII string
00327  */
00328 int mpi_write_string( mpi *X, int radix, char *s, int *slen )
00329 {
00330     int ret = 0, n;
00331     char *p;
00332     mpi T;
00333 
00334     if( radix < 2 || radix > 16 )
00335         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
00336 
00337     n = mpi_msb( X );
00338     if( radix >=  4 ) n >>= 1;
00339     if( radix >= 16 ) n >>= 1;
00340     n += 3;
00341 
00342     if( *slen < n )
00343     {
00344         *slen = n;
00345         return( XYSSL_ERR_MPI_BUFFER_TOO_SMALL );
00346     }
00347 
00348     p = s;
00349     mpi_init( &T, NULL );
00350 
00351     if( X->s == -1 )
00352         *p++ = '-';
00353 
00354     if( radix == 16 )
00355     {
00356         int c, i, j, k;
00357 
00358         for( i = X->n - 1, k = 0; i >= 0; i-- )
00359         {
00360             for( j = ciL - 1; j >= 0; j-- )
00361             {
00362                 c = ( X->p[i] >> (j << 3) ) & 0xFF;
00363 
00364                 if( c == 0 && k == 0 && (i + j) != 0 )
00365                     continue;
00366 
00367                 p += sprintf( p, "%02X", c );
00368                 k = 1;
00369             }
00370         }
00371     }
00372     else
00373     {
00374         MPI_CHK( mpi_copy( &T, X ) );
00375         MPI_CHK( mpi_write_hlp( &T, radix, &p ) );
00376     }
00377 
00378     *p++ = '\0';
00379     *slen = p - s;
00380 
00381 cleanup:
00382 
00383     mpi_free( &T, NULL );
00384 
00385     return( ret );
00386 }
00387 
00388 /*
00389  * Read X from an opened file
00390  */
00391 int mpi_read_file( mpi *X, int radix, FILE *fin )
00392 {
00393     t_int d;
00394     int slen;
00395     char *p;
00396     char s[1024];
00397 
00398     memset( s, 0, sizeof( s ) );
00399     if( fgets( s, sizeof( s ) - 1, fin ) == NULL )
00400         return( XYSSL_ERR_MPI_FILE_IO_ERROR );
00401 
00402     slen = strlen( s );
00403     if( s[slen - 1] == '\n' ) { slen--; s[slen] = '\0'; }
00404     if( s[slen - 1] == '\r' ) { slen--; s[slen] = '\0'; }
00405 
00406     p = s + slen;
00407     while( --p >= s )
00408         if( mpi_get_digit( &d, radix, *p ) != 0 )
00409             break;
00410 
00411     return( mpi_read_string( X, radix, p + 1 ) );
00412 }
00413 
00414 /*
00415  * Write X into an opened file (or stdout if fout == NULL)
00416  */
00417 int mpi_write_file( char *p, mpi *X, int radix, FILE *fout )
00418 {
00419     int n, ret;
00420     size_t slen;
00421     size_t plen;
00422     char s[1024];
00423 
00424     n = sizeof( s );
00425     memset( s, 0, n );
00426     n -= 2;
00427 
00428     MPI_CHK( mpi_write_string( X, radix, s, (int *) &n ) );
00429 
00430     if( p == NULL ) p = "";
00431 
00432     plen = strlen( p );
00433     slen = strlen( s );
00434     s[slen++] = '\r';
00435     s[slen++] = '\n';
00436 
00437     if( fout != NULL )
00438     {
00439         if( fwrite( p, 1, plen, fout ) != plen ||
00440             fwrite( s, 1, slen, fout ) != slen )
00441             return( XYSSL_ERR_MPI_FILE_IO_ERROR );
00442     }
00443     else
00444         printf( "%s%s", p, s );
00445 
00446 cleanup:
00447 
00448     return( ret );
00449 }
00450 
00451 /*
00452  * Import X from unsigned binary data, big endian
00453  */
00454 int mpi_read_binary( mpi *X, unsigned char *buf, int buflen )
00455 {
00456     int ret, i, j, n;
00457 
00458     for( n = 0; n < buflen; n++ )
00459         if( buf[n] != 0 )
00460             break;
00461 
00462     MPI_CHK( mpi_grow( X, CHARS_TO_LIMBS( buflen - n ) ) );
00463     MPI_CHK( mpi_lset( X, 0 ) );
00464 
00465     for( i = buflen - 1, j = 0; i >= n; i--, j++ )
00466         X->p[j / ciL] |= ((t_int) buf[i]) << ((j % ciL) << 3);
00467 
00468 cleanup:
00469 
00470     return( ret );
00471 }
00472 
00473 /*
00474  * Export X into unsigned binary data, big endian
00475  */
00476 int mpi_write_binary( mpi *X, unsigned char *buf, int buflen )
00477 {
00478     int i, j, n;
00479 
00480     n = mpi_size( X );
00481 
00482     if( buflen < n )
00483         return( XYSSL_ERR_MPI_BUFFER_TOO_SMALL );
00484 
00485     memset( buf, 0, buflen );
00486 
00487     for( i = buflen - 1, j = 0; n > 0; i--, j++, n-- )
00488         buf[i] = (unsigned char)( X->p[j / ciL] >> ((j % ciL) << 3) );
00489 
00490     return( 0 );
00491 }
00492 
00493 /*
00494  * Left-shift: X <<= count
00495  */
00496 int mpi_shift_l( mpi *X, int count )
00497 {
00498     int ret, i, v0, t1;
00499     t_int r0 = 0, r1;
00500 
00501     v0 = count / (biL    );
00502     t1 = count & (biL - 1);
00503 
00504     i = mpi_msb( X ) + count;
00505 
00506     if( X->n * (int) biL < i )
00507         MPI_CHK( mpi_grow( X, BITS_TO_LIMBS( i ) ) );
00508 
00509     ret = 0;
00510 
00511     /*
00512      * shift by count / limb_size
00513      */
00514     if( v0 > 0 )
00515     {
00516         for( i = X->n - 1; i >= v0; i-- )
00517             X->p[i] = X->p[i - v0];
00518 
00519         for( ; i >= 0; i-- )
00520             X->p[i] = 0;
00521     }
00522 
00523     /*
00524      * shift by count % limb_size
00525      */
00526     if( t1 > 0 )
00527     {
00528         for( i = v0; i < X->n; i++ )
00529         {
00530             r1 = X->p[i] >> (biL - t1);
00531             X->p[i] <<= t1;
00532             X->p[i] |= r0;
00533             r0 = r1;
00534         }
00535     }
00536 
00537 cleanup:
00538 
00539     return( ret );
00540 }
00541 
00542 /*
00543  * Right-shift: X >>= count
00544  */
00545 int mpi_shift_r( mpi *X, int count )
00546 {
00547     int i, v0, v1;
00548     t_int r0 = 0, r1;
00549 
00550     v0 = count /  biL;
00551     v1 = count & (biL - 1);
00552 
00553     /*
00554      * shift by count / limb_size
00555      */
00556     if( v0 > 0 )
00557     {
00558         for( i = 0; i < X->n - v0; i++ )
00559             X->p[i] = X->p[i + v0];
00560 
00561         for( ; i < X->n; i++ )
00562             X->p[i] = 0;
00563     }
00564 
00565     /*
00566      * shift by count % limb_size
00567      */
00568     if( v1 > 0 )
00569     {
00570         for( i = X->n - 1; i >= 0; i-- )
00571         {
00572             r1 = X->p[i] << (biL - v1);
00573             X->p[i] >>= v1;
00574             X->p[i] |= r0;
00575             r0 = r1;
00576         }
00577     }
00578 
00579     return( 0 );
00580 }
00581 
00582 /*
00583  * Compare unsigned values
00584  */
00585 int mpi_cmp_abs( mpi *X, mpi *Y )
00586 {
00587     int i, j;
00588 
00589     for( i = X->n - 1; i >= 0; i-- )
00590         if( X->p[i] != 0 )
00591             break;
00592 
00593     for( j = Y->n - 1; j >= 0; j-- )
00594         if( Y->p[j] != 0 )
00595             break;
00596 
00597     if( i < 0 && j < 0 )
00598         return( 0 );
00599 
00600     if( i > j ) return(  1 );
00601     if( j > i ) return( -1 );
00602 
00603     for( ; i >= 0; i-- )
00604     {
00605         if( X->p[i] > Y->p[i] ) return(  1 );
00606         if( X->p[i] < Y->p[i] ) return( -1 );
00607     }
00608 
00609     return( 0 );
00610 }
00611 
00612 /*
00613  * Compare signed values
00614  */
00615 int mpi_cmp_mpi( mpi *X, mpi *Y )
00616 {
00617     int i, j;
00618 
00619     for( i = X->n - 1; i >= 0; i-- )
00620         if( X->p[i] != 0 )
00621             break;
00622 
00623     for( j = Y->n - 1; j >= 0; j-- )
00624         if( Y->p[j] != 0 )
00625             break;
00626 
00627     if( i < 0 && j < 0 )
00628         return( 0 );
00629 
00630     if( i > j ) return(  X->s );
00631     if( j > i ) return( -X->s );
00632 
00633     if( X->s > 0 && Y->s < 0 ) return(  1 );
00634     if( Y->s > 0 && X->s < 0 ) return( -1 );
00635 
00636     for( ; i >= 0; i-- )
00637     {
00638         if( X->p[i] > Y->p[i] ) return(  X->s );
00639         if( X->p[i] < Y->p[i] ) return( -X->s );
00640     }
00641 
00642     return( 0 );
00643 }
00644 
00645 /*
00646  * Compare signed values
00647  */
00648 int mpi_cmp_int( mpi *X, int z )
00649 {
00650     mpi Y;
00651     t_int p[1];
00652 
00653     *p  = ( z < 0 ) ? -z : z;
00654     Y.s = ( z < 0 ) ? -1 : 1;
00655     Y.n = 1;
00656     Y.p = p;
00657 
00658     return( mpi_cmp_mpi( X, &Y ) );
00659 }
00660 
00661 /*
00662  * Unsigned addition: X = |A| + |B|  (HAC 14.7)
00663  */
00664 int mpi_add_abs( mpi *X, mpi *A, mpi *B )
00665 {
00666     int ret, i, j;
00667     t_int *o, *p, c;
00668 
00669     if( X == B )
00670     {
00671         mpi *T = A; A = X; B = T;
00672     }
00673 
00674     if( X != A )
00675         MPI_CHK( mpi_copy( X, A ) );
00676 
00677     for( j = B->n - 1; j >= 0; j-- )
00678         if( B->p[j] != 0 )
00679             break;
00680 
00681     MPI_CHK( mpi_grow( X, j + 1 ) );
00682 
00683     o = B->p; p = X->p; c = 0;
00684 
00685     for( i = 0; i <= j; i++, o++, p++ )
00686     {
00687         *p +=  c; c  = ( *p <  c );
00688         *p += *o; c += ( *p < *o );
00689     }
00690 
00691     while( c != 0 )
00692     {
00693         if( i >= X->n )
00694         {
00695             MPI_CHK( mpi_grow( X, i + 1 ) );
00696             p = X->p + i;
00697         }
00698 
00699         *p += c; c = ( *p < c ); i++;
00700     }
00701 
00702 cleanup:
00703 
00704     return( ret );
00705 }
00706 
00707 /*
00708  * Helper for mpi substraction
00709  */
00710 static void mpi_sub_hlp( int n, t_int *s, t_int *d )
00711 {
00712     int i;
00713     t_int c, z;
00714 
00715     for( i = c = 0; i < n; i++, s++, d++ )
00716     {
00717         z = ( *d <  c );     *d -=  c;
00718         c = ( *d < *s ) + z; *d -= *s;
00719     }
00720 
00721     while( c != 0 )
00722     {
00723         z = ( *d < c ); *d -= c;
00724         c = z; i++; d++;
00725     }
00726 }
00727 
00728 /*
00729  * Unsigned substraction: X = |A| - |B|  (HAC 14.9)
00730  */
00731 int mpi_sub_abs( mpi *X, mpi *A, mpi *B )
00732 {
00733     mpi TB;
00734     int ret, n;
00735 
00736     if( mpi_cmp_abs( A, B ) < 0 )
00737         return( XYSSL_ERR_MPI_NEGATIVE_VALUE );
00738 
00739     mpi_init( &TB, NULL );
00740 
00741     if( X == B )
00742     {
00743         MPI_CHK( mpi_copy( &TB, B ) );
00744         B = &TB;
00745     }
00746 
00747     if( X != A )
00748         MPI_CHK( mpi_copy( X, A ) );
00749 
00750     ret = 0;
00751 
00752     for( n = B->n - 1; n >= 0; n-- )
00753         if( B->p[n] != 0 )
00754             break;
00755 
00756     mpi_sub_hlp( n + 1, B->p, X->p );
00757 
00758 cleanup:
00759 
00760     mpi_free( &TB, NULL );
00761 
00762     return( ret );
00763 }
00764 
00765 /*
00766  * Signed addition: X = A + B
00767  */
00768 int mpi_add_mpi( mpi *X, mpi *A, mpi *B )
00769 {
00770     int ret, s = A->s;
00771 
00772     if( A->s * B->s < 0 )
00773     {
00774         if( mpi_cmp_abs( A, B ) >= 0 )
00775         {
00776             MPI_CHK( mpi_sub_abs( X, A, B ) );
00777             X->s =  s;
00778         }
00779         else
00780         {
00781             MPI_CHK( mpi_sub_abs( X, B, A ) );
00782             X->s = -s;
00783         }
00784     }
00785     else
00786     {
00787         MPI_CHK( mpi_add_abs( X, A, B ) );
00788         X->s = s;
00789     }
00790 
00791 cleanup:
00792 
00793     return( ret );
00794 }
00795 
00796 /*
00797  * Signed substraction: X = A - B
00798  */
00799 int mpi_sub_mpi( mpi *X, mpi *A, mpi *B )
00800 {
00801     int ret, s = A->s;
00802 
00803     if( A->s * B->s > 0 )
00804     {
00805         if( mpi_cmp_abs( A, B ) >= 0 )
00806         {
00807             MPI_CHK( mpi_sub_abs( X, A, B ) );
00808             X->s =  s;
00809         }
00810         else
00811         {
00812             MPI_CHK( mpi_sub_abs( X, B, A ) );
00813             X->s = -s;
00814         }
00815     }
00816     else
00817     {
00818         MPI_CHK( mpi_add_abs( X, A, B ) );
00819         X->s = s;
00820     }
00821 
00822 cleanup:
00823 
00824     return( ret );
00825 }
00826 
00827 /*
00828  * Signed addition: X = A + b
00829  */
00830 int mpi_add_int( mpi *X, mpi *A, int b )
00831 {
00832     mpi _B;
00833     t_int p[1];
00834 
00835     p[0] = ( b < 0 ) ? -b : b;
00836     _B.s = ( b < 0 ) ? -1 : 1;
00837     _B.n = 1;
00838     _B.p = p;
00839 
00840     return( mpi_add_mpi( X, A, &_B ) );
00841 }
00842 
00843 /*
00844  * Signed substraction: X = A - b
00845  */
00846 int mpi_sub_int( mpi *X, mpi *A, int b )
00847 {
00848     mpi _B;
00849     t_int p[1];
00850 
00851     p[0] = ( b < 0 ) ? -b : b;
00852     _B.s = ( b < 0 ) ? -1 : 1;
00853     _B.n = 1;
00854     _B.p = p;
00855 
00856     return( mpi_sub_mpi( X, A, &_B ) );
00857 }
00858 
00859 /*
00860  * Helper for mpi multiplication
00861  */ 
00862 static void mpi_mul_hlp( int i, t_int *s, t_int *d, t_int b )
00863 {
00864     t_int c = 0, t = 0;
00865 
00866 #if defined(MULADDC_HUIT)
00867     for( ; i >= 8; i -= 8 )
00868     {
00869         MULADDC_INIT
00870         MULADDC_HUIT
00871         MULADDC_STOP
00872     }
00873 
00874     for( ; i > 0; i-- )
00875     {
00876         MULADDC_INIT
00877         MULADDC_CORE
00878         MULADDC_STOP
00879     }
00880 #else
00881     for( ; i >= 16; i -= 16 )
00882     {
00883         MULADDC_INIT
00884         MULADDC_CORE   MULADDC_CORE
00885         MULADDC_CORE   MULADDC_CORE
00886         MULADDC_CORE   MULADDC_CORE
00887         MULADDC_CORE   MULADDC_CORE
00888 
00889         MULADDC_CORE   MULADDC_CORE
00890         MULADDC_CORE   MULADDC_CORE
00891         MULADDC_CORE   MULADDC_CORE
00892         MULADDC_CORE   MULADDC_CORE
00893         MULADDC_STOP
00894     }
00895 
00896     for( ; i >= 8; i -= 8 )
00897     {
00898         MULADDC_INIT
00899         MULADDC_CORE   MULADDC_CORE
00900         MULADDC_CORE   MULADDC_CORE
00901 
00902         MULADDC_CORE   MULADDC_CORE
00903         MULADDC_CORE   MULADDC_CORE
00904         MULADDC_STOP
00905     }
00906 
00907     for( ; i > 0; i-- )
00908     {
00909         MULADDC_INIT
00910         MULADDC_CORE
00911         MULADDC_STOP
00912     }
00913 #endif
00914 
00915     t++;
00916 
00917     do {
00918         *d += c; c = ( *d < c ); d++;
00919     }
00920     while( c != 0 );
00921 }
00922 
00923 /*
00924  * Baseline multiplication: X = A * B  (HAC 14.12)
00925  */
00926 int mpi_mul_mpi( mpi *X, mpi *A, mpi *B )
00927 {
00928     int ret, i, j;
00929     mpi TA, TB;
00930 
00931     mpi_init( &TA, &TB, NULL );
00932 
00933     if( X == A ) { MPI_CHK( mpi_copy( &TA, A ) ); A = &TA; }
00934     if( X == B ) { MPI_CHK( mpi_copy( &TB, B ) ); B = &TB; }
00935 
00936     for( i = A->n - 1; i >= 0; i-- )
00937         if( A->p[i] != 0 )
00938             break;
00939 
00940     for( j = B->n - 1; j >= 0; j-- )
00941         if( B->p[j] != 0 )
00942             break;
00943 
00944     MPI_CHK( mpi_grow( X, i + j + 2 ) );
00945     MPI_CHK( mpi_lset( X, 0 ) );
00946 
00947     for( i++; j >= 0; j-- )
00948         mpi_mul_hlp( i, A->p, X->p + j, B->p[j] );
00949 
00950     X->s = A->s * B->s;
00951 
00952 cleanup:
00953 
00954     mpi_free( &TB, &TA, NULL );
00955 
00956     return( ret );
00957 }
00958 
00959 /*
00960  * Baseline multiplication: X = A * b
00961  */
00962 int mpi_mul_int( mpi *X, mpi *A, t_int b )
00963 {
00964     mpi _B;
00965     t_int p[1];
00966 
00967     _B.s = 1;
00968     _B.n = 1;
00969     _B.p = p;
00970     p[0] = b;
00971 
00972     return( mpi_mul_mpi( X, A, &_B ) );
00973 }
00974 
00975 /*
00976  * Division by mpi: A = Q * B + R  (HAC 14.20)
00977  */
00978 int mpi_div_mpi( mpi *Q, mpi *R, mpi *A, mpi *B )
00979 {
00980     int ret, i, n, t, k;
00981     mpi X, Y, Z, T1, T2;
00982 
00983     if( mpi_cmp_int( B, 0 ) == 0 )
00984         return( XYSSL_ERR_MPI_DIVISION_BY_ZERO );
00985 
00986     mpi_init( &X, &Y, &Z, &T1, &T2, NULL );
00987 
00988     if( mpi_cmp_abs( A, B ) < 0 )
00989     {
00990         if( Q != NULL ) MPI_CHK( mpi_lset( Q, 0 ) );
00991         if( R != NULL ) MPI_CHK( mpi_copy( R, A ) );
00992         return( 0 );
00993     }
00994 
00995     MPI_CHK( mpi_copy( &X, A ) );
00996     MPI_CHK( mpi_copy( &Y, B ) );
00997     X.s = Y.s = 1;
00998 
00999     MPI_CHK( mpi_grow( &Z, A->n + 2 ) );
01000     MPI_CHK( mpi_lset( &Z,  0 ) );
01001     MPI_CHK( mpi_grow( &T1, 2 ) );
01002     MPI_CHK( mpi_grow( &T2, 3 ) );
01003 
01004     k = mpi_msb( &Y ) % biL;
01005     if( k < (int) biL - 1 )
01006     {
01007         k = biL - 1 - k;
01008         MPI_CHK( mpi_shift_l( &X, k ) );
01009         MPI_CHK( mpi_shift_l( &Y, k ) );
01010     }
01011     else k = 0;
01012 
01013     n = X.n - 1;
01014     t = Y.n - 1;
01015     mpi_shift_l( &Y, biL * (n - t) );
01016 
01017     while( mpi_cmp_mpi( &X, &Y ) >= 0 )
01018     {
01019         Z.p[n - t]++;
01020         mpi_sub_mpi( &X, &X, &Y );
01021     }
01022     mpi_shift_r( &Y, biL * (n - t) );
01023 
01024     for( i = n; i > t ; i-- )
01025     {
01026         if( X.p[i] >= Y.p[t] )
01027             Z.p[i - t - 1] = ~0;
01028         else
01029         {
01030 #if defined(XYSSL_HAVE_LONGLONG)
01031             t_dbl r;
01032 
01033             r  = (t_dbl) X.p[i] << biL;
01034             r |= (t_dbl) X.p[i - 1];
01035             r /= Y.p[t];
01036             if( r > ((t_dbl) 1 << biL) - 1)
01037                 r = ((t_dbl) 1 << biL) - 1;
01038 
01039             Z.p[i - t - 1] = (t_int) r;
01040 #else
01041             /*
01042              * __udiv_qrnnd_c, from gmp/longlong.h
01043              */
01044             t_int q0, q1, r0, r1;
01045             t_int d0, d1, d, m;
01046 
01047             d  = Y.p[t];
01048             d0 = ( d << biH ) >> biH;
01049             d1 = ( d >> biH );
01050 
01051             q1 = X.p[i] / d1;
01052             r1 = X.p[i] - d1 * q1;
01053             r1 <<= biH;
01054             r1 |= ( X.p[i - 1] >> biH );
01055 
01056             m = q1 * d0;
01057             if( r1 < m )
01058             {
01059                 q1--, r1 += d;
01060                 while( r1 >= d && r1 < m )
01061                     q1--, r1 += d;
01062             }
01063             r1 -= m;
01064 
01065             q0 = r1 / d1;
01066             r0 = r1 - d1 * q0;
01067             r0 <<= biH;
01068             r0 |= ( X.p[i - 1] << biH ) >> biH;
01069 
01070             m = q0 * d0;
01071             if( r0 < m )
01072             {
01073                 q0--, r0 += d;
01074                 while( r0 >= d && r0 < m )
01075                     q0--, r0 += d;
01076             }
01077             r0 -= m;
01078 
01079             Z.p[i - t - 1] = ( q1 << biH ) | q0;
01080 #endif
01081         }
01082 
01083         Z.p[i - t - 1]++;
01084         do
01085         {
01086             Z.p[i - t - 1]--;
01087 
01088             MPI_CHK( mpi_lset( &T1, 0 ) );
01089             T1.p[0] = (t < 1) ? 0 : Y.p[t - 1];
01090             T1.p[1] = Y.p[t];
01091             MPI_CHK( mpi_mul_int( &T1, &T1, Z.p[i - t - 1] ) );
01092 
01093             MPI_CHK( mpi_lset( &T2, 0 ) );
01094             T2.p[0] = (i < 2) ? 0 : X.p[i - 2];
01095             T2.p[1] = (i < 1) ? 0 : X.p[i - 1];
01096             T2.p[2] = X.p[i];
01097         }
01098         while( mpi_cmp_mpi( &T1, &T2 ) > 0 );
01099 
01100         MPI_CHK( mpi_mul_int( &T1, &Y, Z.p[i - t - 1] ) );
01101         MPI_CHK( mpi_shift_l( &T1,  biL * (i - t - 1) ) );
01102         MPI_CHK( mpi_sub_mpi( &X, &X, &T1 ) );
01103 
01104         if( mpi_cmp_int( &X, 0 ) < 0 )
01105         {
01106             MPI_CHK( mpi_copy( &T1, &Y ) );
01107             MPI_CHK( mpi_shift_l( &T1, biL * (i - t - 1) ) );
01108             MPI_CHK( mpi_add_mpi( &X, &X, &T1 ) );
01109             Z.p[i - t - 1]--;
01110         }
01111     }
01112 
01113     if( Q != NULL )
01114     {
01115         mpi_copy( Q, &Z );
01116         Q->s = A->s * B->s;
01117     }
01118 
01119     if( R != NULL )
01120     {
01121         mpi_shift_r( &X, k );
01122         mpi_copy( R, &X );
01123 
01124         R->s = A->s;
01125         if( mpi_cmp_int( R, 0 ) == 0 )
01126             R->s = 1;
01127     }
01128 
01129 cleanup:
01130 
01131     mpi_free( &X, &Y, &Z, &T1, &T2, NULL );
01132 
01133     return( ret );
01134 }
01135 
01136 /*
01137  * Division by int: A = Q * b + R
01138  *
01139  * Returns 0 if successful
01140  *         1 if memory allocation failed
01141  *         XYSSL_ERR_MPI_DIVISION_BY_ZERO if b == 0
01142  */
01143 int mpi_div_int( mpi *Q, mpi *R, mpi *A, int b )
01144 {
01145     mpi _B;
01146     t_int p[1];
01147 
01148     p[0] = ( b < 0 ) ? -b : b;
01149     _B.s = ( b < 0 ) ? -1 : 1;
01150     _B.n = 1;
01151     _B.p = p;
01152 
01153     return( mpi_div_mpi( Q, R, A, &_B ) );
01154 }
01155 
01156 /*
01157  * Modulo: R = A mod B
01158  */
01159 int mpi_mod_mpi( mpi *R, mpi *A, mpi *B )
01160 {
01161     int ret;
01162 
01163     MPI_CHK( mpi_div_mpi( NULL, R, A, B ) );
01164 
01165     while( mpi_cmp_int( R, 0 ) < 0 )
01166       MPI_CHK( mpi_add_mpi( R, R, B ) );
01167 
01168     while( mpi_cmp_mpi( R, B ) >= 0 )
01169       MPI_CHK( mpi_sub_mpi( R, R, B ) );
01170 
01171 cleanup:
01172 
01173     return( ret );
01174 }
01175 
01176 /*
01177  * Modulo: r = A mod b
01178  */
01179 int mpi_mod_int( t_int *r, mpi *A, int b )
01180 {
01181     int i;
01182     t_int x, y, z;
01183 
01184     if( b == 0 )
01185         return( XYSSL_ERR_MPI_DIVISION_BY_ZERO );
01186 
01187     if( b < 0 )
01188         b = -b;
01189 
01190     /*
01191      * handle trivial cases
01192      */
01193     if( b == 1 )
01194     {
01195         *r = 0;
01196         return( 0 );
01197     }
01198 
01199     if( b == 2 )
01200     {
01201         *r = A->p[0] & 1;
01202         return( 0 );
01203     }
01204 
01205     /*
01206      * general case
01207      */
01208     for( i = A->n - 1, y = 0; i >= 0; i-- )
01209     {
01210         x  = A->p[i];
01211         y  = ( y << biH ) | ( x >> biH );
01212         z  = y / b;
01213         y -= z * b;
01214 
01215         x <<= biH;
01216         y  = ( y << biH ) | ( x >> biH );
01217         z  = y / b;
01218         y -= z * b;
01219     }
01220 
01221     *r = y;
01222 
01223     return( 0 );
01224 }
01225 
01226 /*
01227  * Fast Montgomery initialization (thanks to Tom St Denis)
01228  */
01229 static void mpi_montg_init( t_int *mm, mpi *N )
01230 {
01231     t_int x, m0 = N->p[0];
01232 
01233     x  = m0;
01234     x += ( ( m0 + 2 ) & 4 ) << 1;
01235     x *= ( 2 - ( m0 * x ) );
01236 
01237     if( biL >= 16 ) x *= ( 2 - ( m0 * x ) );
01238     if( biL >= 32 ) x *= ( 2 - ( m0 * x ) );
01239     if( biL >= 64 ) x *= ( 2 - ( m0 * x ) );
01240 
01241     *mm = ~x + 1;
01242 }
01243 
01244 /*
01245  * Montgomery multiplication: A = A * B * R^-1 mod N  (HAC 14.36)
01246  */
01247 static void mpi_montmul( mpi *A, mpi *B, mpi *N, t_int mm, mpi *T )
01248 {
01249     int i, n, m;
01250     t_int u0, u1, *d;
01251 
01252     memset( T->p, 0, T->n * ciL );
01253 
01254     d = T->p;
01255     n = N->n;
01256     m = ( B->n < n ) ? B->n : n;
01257 
01258     for( i = 0; i < n; i++ )
01259     {
01260         /*
01261          * T = (T + u0*B + u1*N) / 2^biL
01262          */
01263         u0 = A->p[i];
01264         u1 = ( d[0] + u0 * B->p[0] ) * mm;
01265 
01266         mpi_mul_hlp( m, B->p, d, u0 );
01267         mpi_mul_hlp( n, N->p, d, u1 );
01268 
01269         *d++ = u0; d[n + 1] = 0;
01270     }
01271 
01272     memcpy( A->p, d, (n + 1) * ciL );
01273 
01274     if( mpi_cmp_abs( A, N ) >= 0 )
01275         mpi_sub_hlp( n, N->p, A->p );
01276     else
01277         /* prevent timing attacks */
01278         mpi_sub_hlp( n, A->p, T->p );
01279 }
01280 
01281 /*
01282  * Montgomery reduction: A = A * R^-1 mod N
01283  */
01284 static void mpi_montred( mpi *A, mpi *N, t_int mm, mpi *T )
01285 {
01286     t_int z = 1;
01287     mpi U;
01288 
01289     U.n = U.s = z;
01290     U.p = &z;
01291 
01292     mpi_montmul( A, &U, N, mm, T );
01293 }
01294 
01295 /*
01296  * Sliding-window exponentiation: X = A^E mod N  (HAC 14.85)
01297  */
01298 int mpi_exp_mod( mpi *X, mpi *A, mpi *E, mpi *N, mpi *_RR )
01299 {
01300     int ret, i, j, wsize, wbits;
01301     int bufsize, nblimbs, nbits;
01302     t_int ei, mm, state;
01303     mpi RR, T, W[64];
01304 
01305     if( mpi_cmp_int( N, 0 ) < 0 || ( N->p[0] & 1 ) == 0 )
01306         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
01307 
01308     /*
01309      * Init temps and window size
01310      */
01311     mpi_montg_init( &mm, N );
01312     mpi_init( &RR, &T, NULL );
01313     memset( W, 0, sizeof( W ) );
01314 
01315     i = mpi_msb( E );
01316 
01317     wsize = ( i > 671 ) ? 6 : ( i > 239 ) ? 5 :
01318             ( i >  79 ) ? 4 : ( i >  23 ) ? 3 : 1;
01319 
01320     j = N->n + 1;
01321     MPI_CHK( mpi_grow( X, j ) );
01322     MPI_CHK( mpi_grow( &W[1],  j ) );
01323     MPI_CHK( mpi_grow( &T, j * 2 ) );
01324 
01325     /*
01326      * If 1st call, pre-compute R^2 mod N
01327      */
01328     if( _RR == NULL || _RR->p == NULL )
01329     {
01330         MPI_CHK( mpi_lset( &RR, 1 ) );
01331         MPI_CHK( mpi_shift_l( &RR, N->n * 2 * biL ) );
01332         MPI_CHK( mpi_mod_mpi( &RR, &RR, N ) );
01333 
01334         if( _RR != NULL )
01335             memcpy( _RR, &RR, sizeof( mpi ) );
01336     }
01337     else
01338         memcpy( &RR, _RR, sizeof( mpi ) );
01339 
01340     /*
01341      * W[1] = A * R^2 * R^-1 mod N = A * R mod N
01342      */
01343     if( mpi_cmp_mpi( A, N ) >= 0 )
01344         mpi_mod_mpi( &W[1], A, N );
01345     else   mpi_copy( &W[1], A );
01346 
01347     mpi_montmul( &W[1], &RR, N, mm, &T );
01348 
01349     /*
01350      * X = R^2 * R^-1 mod N = R mod N
01351      */
01352     MPI_CHK( mpi_copy( X, &RR ) );
01353     mpi_montred( X, N, mm, &T );
01354 
01355     if( wsize > 1 )
01356     {
01357         /*
01358          * W[1 << (wsize - 1)] = W[1] ^ (wsize - 1)
01359          */
01360         j =  1 << (wsize - 1);
01361 
01362         MPI_CHK( mpi_grow( &W[j], N->n + 1 ) );
01363         MPI_CHK( mpi_copy( &W[j], &W[1]    ) );
01364 
01365         for( i = 0; i < wsize - 1; i++ )
01366             mpi_montmul( &W[j], &W[j], N, mm, &T );
01367     
01368         /*
01369          * W[i] = W[i - 1] * W[1]
01370          */
01371         for( i = j + 1; i < (1 << wsize); i++ )
01372         {
01373             MPI_CHK( mpi_grow( &W[i], N->n + 1 ) );
01374             MPI_CHK( mpi_copy( &W[i], &W[i - 1] ) );
01375 
01376             mpi_montmul( &W[i], &W[1], N, mm, &T );
01377         }
01378     }
01379 
01380     nblimbs = E->n;
01381     bufsize = 0;
01382     nbits   = 0;
01383     wbits   = 0;
01384     state   = 0;
01385 
01386     while( 1 )
01387     {
01388         if( bufsize == 0 )
01389         {
01390             if( nblimbs-- == 0 )
01391                 break;
01392 
01393             bufsize = sizeof( t_int ) << 3;
01394         }
01395 
01396         bufsize--;
01397 
01398         ei = (E->p[nblimbs] >> bufsize) & 1;
01399 
01400         /*
01401          * skip leading 0s
01402          */
01403         if( ei == 0 && state == 0 )
01404             continue;
01405 
01406         if( ei == 0 && state == 1 )
01407         {
01408             /*
01409              * out of window, square X
01410              */
01411             mpi_montmul( X, X, N, mm, &T );
01412             continue;
01413         }
01414 
01415         /*
01416          * add ei to current window
01417          */
01418         state = 2;
01419 
01420         nbits++;
01421         wbits |= (ei << (wsize - nbits));
01422 
01423         if( nbits == wsize )
01424         {
01425             /*
01426              * X = X^wsize R^-1 mod N
01427              */
01428             for( i = 0; i < wsize; i++ )
01429                 mpi_montmul( X, X, N, mm, &T );
01430 
01431             /*
01432              * X = X * W[wbits] R^-1 mod N
01433              */
01434             mpi_montmul( X, &W[wbits], N, mm, &T );
01435 
01436             state--;
01437             nbits = 0;
01438             wbits = 0;
01439         }
01440     }
01441 
01442     /*
01443      * process the remaining bits
01444      */
01445     for( i = 0; i < nbits; i++ )
01446     {
01447         mpi_montmul( X, X, N, mm, &T );
01448 
01449         wbits <<= 1;
01450 
01451         if( (wbits & (1 << wsize)) != 0 )
01452             mpi_montmul( X, &W[1], N, mm, &T );
01453     }
01454 
01455     /*
01456      * X = A^E * R * R^-1 mod N = A^E mod N
01457      */
01458     mpi_montred( X, N, mm, &T );
01459 
01460 cleanup:
01461 
01462     for( i = (1 << (wsize - 1)); i < (1 << wsize); i++ )
01463         mpi_free( &W[i], NULL );
01464 
01465     if( _RR != NULL )
01466          mpi_free( &W[1], &T, NULL );
01467     else mpi_free( &W[1], &T, &RR, NULL );
01468 
01469     return( ret );
01470 }
01471 
01472 #if defined(XYSSL_GENPRIME)
01473 
01474 /*
01475  * Greatest common divisor: G = gcd(A, B)  (HAC 14.54)
01476  */
01477 int mpi_gcd( mpi *G, mpi *A, mpi *B )
01478 {
01479     int ret;
01480     mpi TG, TA, TB;
01481 
01482     mpi_init( &TG, &TA, &TB, NULL );
01483 
01484     MPI_CHK( mpi_lset( &TG, 1 ) );
01485     MPI_CHK( mpi_copy( &TA, A ) );
01486     MPI_CHK( mpi_copy( &TB, B ) );
01487 
01488     TA.s = TB.s = 1;
01489 
01490     while( mpi_cmp_int( &TA, 0 ) != 0 )
01491     {
01492         while( ( TA.p[0] & 1 ) == 0 ) MPI_CHK( mpi_shift_r( &TA, 1 ) );
01493         while( ( TB.p[0] & 1 ) == 0 ) MPI_CHK( mpi_shift_r( &TB, 1 ) );
01494 
01495         if( mpi_cmp_mpi( &TA, &TB ) >= 0 )
01496         {
01497             MPI_CHK( mpi_sub_abs( &TA, &TA, &TB ) );
01498             MPI_CHK( mpi_shift_r( &TA, 1 ) );
01499         }
01500         else
01501         {
01502             MPI_CHK( mpi_sub_abs( &TB, &TB, &TA ) );
01503             MPI_CHK( mpi_shift_r( &TB, 1 ) );
01504         }
01505     }
01506 
01507     MPI_CHK( mpi_mul_mpi( G, &TG, &TB ) );
01508 
01509 cleanup:
01510 
01511     mpi_free( &TB, &TA, &TG, NULL );
01512 
01513     return( ret );
01514 }
01515 
01516 /*
01517  * Modular inverse: X = A^-1 mod N  (HAC 14.61 / 14.64)
01518  */
01519 int mpi_inv_mod( mpi *X, mpi *A, mpi *N )
01520 {
01521     int ret;
01522     mpi G, TA, TU, U1, U2, TB, TV, V1, V2;
01523 
01524     if( mpi_cmp_int( N, 0 ) <= 0 )
01525         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
01526 
01527     mpi_init( &TA, &TU, &U1, &U2, &G,
01528               &TB, &TV, &V1, &V2, NULL );
01529 
01530     MPI_CHK( mpi_gcd( &G, A, N ) );
01531 
01532     if( mpi_cmp_int( &G, 1 ) != 0 )
01533     {
01534         ret = XYSSL_ERR_MPI_NOT_ACCEPTABLE;
01535         goto cleanup;
01536     }
01537 
01538     MPI_CHK( mpi_mod_mpi( &TA, A, N ) );
01539     MPI_CHK( mpi_copy( &TU, &TA ) );
01540     MPI_CHK( mpi_copy( &TB, N ) );
01541     MPI_CHK( mpi_copy( &TV, N ) );
01542 
01543     MPI_CHK( mpi_lset( &U1, 1 ) );
01544     MPI_CHK( mpi_lset( &U2, 0 ) );
01545     MPI_CHK( mpi_lset( &V1, 0 ) );
01546     MPI_CHK( mpi_lset( &V2, 1 ) );
01547 
01548     do
01549     {
01550         while( ( TU.p[0] & 1 ) == 0 )
01551         {
01552             MPI_CHK( mpi_shift_r( &TU, 1 ) );
01553 
01554             if( ( U1.p[0] & 1 ) != 0 || ( U2.p[0] & 1 ) != 0 )
01555             {
01556                 MPI_CHK( mpi_add_mpi( &U1, &U1, &TB ) );
01557                 MPI_CHK( mpi_sub_mpi( &U2, &U2, &TA ) );
01558             }
01559 
01560             MPI_CHK( mpi_shift_r( &U1, 1 ) );
01561             MPI_CHK( mpi_shift_r( &U2, 1 ) );
01562         }
01563 
01564         while( ( TV.p[0] & 1 ) == 0 )
01565         {
01566             MPI_CHK( mpi_shift_r( &TV, 1 ) );
01567 
01568             if( ( V1.p[0] & 1 ) != 0 || ( V2.p[0] & 1 ) != 0 )
01569             {
01570                 MPI_CHK( mpi_add_mpi( &V1, &V1, &TB ) );
01571                 MPI_CHK( mpi_sub_mpi( &V2, &V2, &TA ) );
01572             }
01573 
01574             MPI_CHK( mpi_shift_r( &V1, 1 ) );
01575             MPI_CHK( mpi_shift_r( &V2, 1 ) );
01576         }
01577 
01578         if( mpi_cmp_mpi( &TU, &TV ) >= 0 )
01579         {
01580             MPI_CHK( mpi_sub_mpi( &TU, &TU, &TV ) );
01581             MPI_CHK( mpi_sub_mpi( &U1, &U1, &V1 ) );
01582             MPI_CHK( mpi_sub_mpi( &U2, &U2, &V2 ) );
01583         }
01584         else
01585         {
01586             MPI_CHK( mpi_sub_mpi( &TV, &TV, &TU ) );
01587             MPI_CHK( mpi_sub_mpi( &V1, &V1, &U1 ) );
01588             MPI_CHK( mpi_sub_mpi( &V2, &V2, &U2 ) );
01589         }
01590     }
01591     while( mpi_cmp_int( &TU, 0 ) != 0 );
01592 
01593     while( mpi_cmp_int( &V1, 0 ) < 0 )
01594         MPI_CHK( mpi_add_mpi( &V1, &V1, N ) );
01595 
01596     while( mpi_cmp_mpi( &V1, N ) >= 0 )
01597         MPI_CHK( mpi_sub_mpi( &V1, &V1, N ) );
01598 
01599     MPI_CHK( mpi_copy( X, &V1 ) );
01600 
01601 cleanup:
01602 
01603     mpi_free( &V2, &V1, &TV, &TB, &G,
01604               &U2, &U1, &TU, &TA, NULL );
01605 
01606     return( ret );
01607 }
01608 
01609 static const int small_prime[] =
01610 {
01611         3,    5,    7,   11,   13,   17,   19,   23,
01612        29,   31,   37,   41,   43,   47,   53,   59,
01613        61,   67,   71,   73,   79,   83,   89,   97,
01614       101,  103,  107,  109,  113,  127,  131,  137,
01615       139,  149,  151,  157,  163,  167,  173,  179,
01616       181,  191,  193,  197,  199,  211,  223,  227,
01617       229,  233,  239,  241,  251,  257,  263,  269,
01618       271,  277,  281,  283,  293,  307,  311,  313,
01619       317,  331,  337,  347,  349,  353,  359,  367,
01620       373,  379,  383,  389,  397,  401,  409,  419,
01621       421,  431,  433,  439,  443,  449,  457,  461,
01622       463,  467,  479,  487,  491,  499,  503,  509,
01623       521,  523,  541,  547,  557,  563,  569,  571,
01624       577,  587,  593,  599,  601,  607,  613,  617,
01625       619,  631,  641,  643,  647,  653,  659,  661,
01626       673,  677,  683,  691,  701,  709,  719,  727,
01627       733,  739,  743,  751,  757,  761,  769,  773,
01628       787,  797,  809,  811,  821,  823,  827,  829,
01629       839,  853,  857,  859,  863,  877,  881,  883,
01630       887,  907,  911,  919,  929,  937,  941,  947,
01631       953,  967,  971,  977,  983,  991,  997, -103
01632 };
01633 
01634 /*
01635  * Miller-Rabin primality test  (HAC 4.24)
01636  */
01637 int mpi_is_prime( mpi *X, int (*f_rng)(void *), void *p_rng )
01638 {
01639     int ret, i, j, n, s, xs;
01640     mpi W, R, T, A, RR;
01641     unsigned char *p;
01642 
01643     if( mpi_cmp_int( X, 0 ) == 0 )
01644         return( 0 );
01645 
01646     mpi_init( &W, &R, &T, &A, &RR, NULL );
01647 
01648     xs = X->s; X->s = 1;
01649 
01650     /*
01651      * test trivial factors first
01652      */
01653     if( ( X->p[0] & 1 ) == 0 )
01654         return( XYSSL_ERR_MPI_NOT_ACCEPTABLE );
01655 
01656     for( i = 0; small_prime[i] > 0; i++ )
01657     {
01658         t_int r;
01659 
01660         if( mpi_cmp_int( X, small_prime[i] ) <= 0 )
01661             return( 0 );
01662 
01663         MPI_CHK( mpi_mod_int( &r, X, small_prime[i] ) );
01664 
01665         if( r == 0 )
01666             return( XYSSL_ERR_MPI_NOT_ACCEPTABLE );
01667     }
01668 
01669     /*
01670      * W = |X| - 1
01671      * R = W >> lsb( W )
01672      */
01673     s = mpi_lsb( &W );
01674     MPI_CHK( mpi_sub_int( &W, X, 1 ) );
01675     MPI_CHK( mpi_copy( &R, &W ) );
01676     MPI_CHK( mpi_shift_r( &R, s ) );
01677 
01678     i = mpi_msb( X );
01679     /*
01680      * HAC, table 4.4
01681      */
01682     n = ( ( i >= 1300 ) ?  2 : ( i >=  850 ) ?  3 :
01683           ( i >=  650 ) ?  4 : ( i >=  350 ) ?  8 :
01684           ( i >=  250 ) ? 12 : ( i >=  150 ) ? 18 : 27 );
01685 
01686     for( i = 0; i < n; i++ )
01687     {
01688         /*
01689          * pick a random A, 1 < A < |X| - 1
01690          */
01691         MPI_CHK( mpi_grow( &A, X->n ) );
01692 
01693         p = (unsigned char *) A.p;
01694         for( j = 0; j < A.n * ciL; j++ )
01695             *p++ = (unsigned char) f_rng( p_rng );
01696 
01697         j = mpi_msb( &A ) - mpi_msb( &W );
01698         MPI_CHK( mpi_shift_r( &A, j + 1 ) );
01699         A.p[0] |= 3;
01700 
01701         /*
01702          * A = A^R mod |X|
01703          */
01704         MPI_CHK( mpi_exp_mod( &A, &A, &R, X, &RR ) );
01705 
01706         if( mpi_cmp_mpi( &A, &W ) == 0 ||
01707             mpi_cmp_int( &A,  1 ) == 0 )
01708             continue;
01709 
01710         j = 1;
01711         while( j < s && mpi_cmp_mpi( &A, &W ) != 0 )
01712         {
01713             /*
01714              * A = A * A mod |X|
01715              */
01716             MPI_CHK( mpi_mul_mpi( &T, &A, &A ) );
01717             MPI_CHK( mpi_mod_mpi( &A, &T, X  ) );
01718 
01719             if( mpi_cmp_int( &A, 1 ) == 0 )
01720                 break;
01721 
01722             j++;
01723         }
01724 
01725         /*
01726          * not prime if A != |X| - 1 or A == 1
01727          */
01728         if( mpi_cmp_mpi( &A, &W ) != 0 ||
01729             mpi_cmp_int( &A,  1 ) == 0 )
01730         {
01731             ret = XYSSL_ERR_MPI_NOT_ACCEPTABLE;
01732             break;
01733         }
01734     }
01735 
01736 cleanup:
01737 
01738     X->s = xs;
01739 
01740     mpi_free( &RR, &A, &T, &R, &W, NULL );
01741 
01742     return( ret );
01743 }
01744 
01745 /*
01746  * Prime number generation
01747  */
01748 int mpi_gen_prime( mpi *X, int nbits, int dh_flag,
01749                    int (*f_rng)(void *), void *p_rng )
01750 {
01751     int ret, k, n;
01752     unsigned char *p;
01753     mpi Y;
01754 
01755     if( nbits < 3 )
01756         return( XYSSL_ERR_MPI_BAD_INPUT_DATA );
01757 
01758     mpi_init( &Y, NULL );
01759 
01760     n = BITS_TO_LIMBS( nbits );
01761 
01762     MPI_CHK( mpi_grow( X, n ) );
01763     MPI_CHK( mpi_lset( X, 0 ) );
01764 
01765     p = (unsigned char *) X->p;
01766     for( k = 0; k < X->n * ciL; k++ )
01767         *p++ = (unsigned char) f_rng( p_rng );
01768 
01769     k = mpi_msb( X );
01770     if( k < nbits ) MPI_CHK( mpi_shift_l( X, nbits - k ) );
01771     if( k > nbits ) MPI_CHK( mpi_shift_r( X, k - nbits ) );
01772 
01773     X->p[0] |= 3;
01774 
01775     if( dh_flag == 0 )
01776     {
01777         while( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) != 0 )
01778         {
01779             if( ret != XYSSL_ERR_MPI_NOT_ACCEPTABLE )
01780                 goto cleanup;
01781 
01782             MPI_CHK( mpi_add_int( X, X, 2 ) );
01783         }
01784     }
01785     else
01786     {
01787         MPI_CHK( mpi_sub_int( &Y, X, 1 ) );
01788         MPI_CHK( mpi_shift_r( &Y, 1 ) );
01789 
01790         while( 1 )
01791         {
01792             if( ( ret = mpi_is_prime( X, f_rng, p_rng ) ) == 0 )
01793             {
01794                 if( ( ret = mpi_is_prime( &Y, f_rng, p_rng ) ) == 0 )
01795                     break;
01796 
01797                 if( ret != XYSSL_ERR_MPI_NOT_ACCEPTABLE )
01798                     goto cleanup;
01799             }
01800 
01801             if( ret != XYSSL_ERR_MPI_NOT_ACCEPTABLE )
01802                 goto cleanup;
01803 
01804             MPI_CHK( mpi_add_int( &Y, X, 1 ) );
01805             MPI_CHK( mpi_add_int(  X, X, 2 ) );
01806             MPI_CHK( mpi_shift_r( &Y, 1 ) );
01807         }
01808     }
01809 
01810 cleanup:
01811 
01812     mpi_free( &Y, NULL );
01813 
01814     return( ret );
01815 }
01816 
01817 #endif
01818 
01819 #if defined(XYSSL_SELF_TEST)
01820 
01821 /*
01822  * Checkup routine
01823  */
01824 int mpi_self_test( int verbose )
01825 {
01826     int ret;
01827     mpi A, E, N, X, Y, U, V;
01828 
01829     mpi_init( &A, &E, &N, &X, &Y, &U, &V, NULL );
01830 
01831     MPI_CHK( mpi_read_string( &A, 16,
01832         "EFE021C2645FD1DC586E69184AF4A31E" \
01833         "D5F53E93B5F123FA41680867BA110131" \
01834         "944FE7952E2517337780CB0DB80E61AA" \
01835         "E7C8DDC6C5C6AADEB34EB38A2F40D5E6" ) );
01836 
01837     MPI_CHK( mpi_read_string( &E, 16,
01838         "B2E7EFD37075B9F03FF989C7C5051C20" \
01839         "34D2A323810251127E7BF8625A4F49A5" \
01840         "F3E27F4DA8BD59C47D6DAABA4C8127BD" \
01841         "5B5C25763222FEFCCFC38B832366C29E" ) );
01842 
01843     MPI_CHK( mpi_read_string( &N, 16,
01844         "0066A198186C18C10B2F5ED9B522752A" \
01845         "9830B69916E535C8F047518A889A43A5" \
01846         "94B6BED27A168D31D4A52F88925AA8F5" ) );
01847 
01848     MPI_CHK( mpi_mul_mpi( &X, &A, &N ) );
01849 
01850     MPI_CHK( mpi_read_string( &U, 16,
01851         "602AB7ECA597A3D6B56FF9829A5E8B85" \
01852         "9E857EA95A03512E2BAE7391688D264A" \
01853         "A5663B0341DB9CCFD2C4C5F421FEC814" \
01854         "8001B72E848A38CAE1C65F78E56ABDEF" \
01855         "E12D3C039B8A02D6BE593F0BBBDA56F1" \
01856         "ECF677152EF804370C1A305CAF3B5BF1" \
01857         "30879B56C61DE584A0F53A2447A51E" ) );
01858 
01859     if( verbose != 0 )
01860         printf( "  MPI test #1 (mul_mpi): " );
01861 
01862     if( mpi_cmp_mpi( &X, &U ) != 0 )
01863     {
01864         if( verbose != 0 )
01865             printf( "failed\n" );
01866 
01867         return( 1 );
01868     }
01869 
01870     if( verbose != 0 )
01871         printf( "passed\n" );
01872 
01873     MPI_CHK( mpi_div_mpi( &X, &Y, &A, &N ) );
01874 
01875     MPI_CHK( mpi_read_string( &U, 16,
01876         "256567336059E52CAE22925474705F39A94" ) );
01877 
01878     MPI_CHK( mpi_read_string( &V, 16,
01879         "6613F26162223DF488E9CD48CC132C7A" \
01880         "0AC93C701B001B092E4E5B9F73BCD27B" \
01881         "9EE50D0657C77F374E903CDFA4C642" ) );
01882 
01883     if( verbose != 0 )
01884         printf( "  MPI test #2 (div_mpi): " );
01885 
01886     if( mpi_cmp_mpi( &X, &U ) != 0 ||
01887         mpi_cmp_mpi( &Y, &V ) != 0 )
01888     {
01889         if( verbose != 0 )
01890             printf( "failed\n" );
01891 
01892         return( 1 );
01893     }
01894 
01895     if( verbose != 0 )
01896         printf( "passed\n" );
01897 
01898     MPI_CHK( mpi_exp_mod( &X, &A, &E, &N, NULL ) );
01899 
01900     MPI_CHK( mpi_read_string( &U, 16,
01901         "36E139AEA55215609D2816998ED020BB" \
01902         "BD96C37890F65171D948E9BC7CBAA4D9" \
01903         "325D24D6A3C12710F10A09FA08AB87" ) );
01904 
01905     if( verbose != 0 )
01906         printf( "  MPI test #3 (exp_mod): " );
01907 
01908     if( mpi_cmp_mpi( &X, &U ) != 0 )
01909     {
01910         if( verbose != 0 )
01911             printf( "failed\n" );
01912 
01913         return( 1 );
01914     }
01915 
01916     if( verbose != 0 )
01917         printf( "passed\n" );
01918 
01919     MPI_CHK( mpi_inv_mod( &X, &A, &N ) );
01920 
01921     MPI_CHK( mpi_read_string( &U, 16,
01922         "003A0AAEDD7E784FC07D8F9EC6E3BFD5" \
01923         "C3DBA76456363A10869622EAC2DD84EC" \
01924         "C5B8A74DAC4D09E03B5E0BE779F2DF61" ) );
01925 
01926     if( verbose != 0 )
01927         printf( "  MPI test #4 (inv_mod): " );
01928 
01929     if( mpi_cmp_mpi( &X, &U ) != 0 )
01930     {
01931         if( verbose != 0 )
01932             printf( "failed\n" );
01933 
01934         return( 1 );
01935     }
01936 
01937     if( verbose != 0 )
01938         printf( "passed\n" );
01939 
01940 cleanup:
01941 
01942     if( ret != 0 && verbose != 0 )
01943         printf( "Unexpected error, return code = %08X\n", ret );
01944 
01945     mpi_free( &V, &U, &Y, &X, &N, &E, &A, NULL );
01946 
01947     if( verbose != 0 )
01948         printf( "\n" );
01949 
01950     return( ret );
01951 }
01952 
01953 #endif
01954 
01955 #endif