/home/dko/projects/mobilec/trunk/src/security/xyssl-0.7/library/bignum.c

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

Generated on Fri May 16 14:49:55 2008 for Mobile-C by  doxygen 1.5.4