📄 lib_kg.c
字号:
/****************************************************************************
* *
* cryptlib PKC Generation/Checking Routines *
* Copyright Peter Gutmann and Kevin Bluck 1997-2003 *
* *
****************************************************************************/
/* The Usenet Oracle has pondered your question deeply.
Your question was:
> O Oracle Most Wise,
>
> What is the largest prime number?
And in response, thus spake the Oracle:
} This is a question which has stumped some of the best minds in
} mathematics, but I will explain it so that even you can understand it.
} The first prime is 2, and the binary representation of 2 is 10.
} Consider the following series:
}
} Prime Decimal Representation Representation in its own base
} 1st 2 10
} 2nd 3 10
} 3rd 5 10
} 4th 7 10
} 5th 11 10
} 6th 13 10
} 7th 17 10
}
} From this demonstration you can see that there is only one prime, and
} it is ten. Therefore, the largest prime is ten.
-- The Usenet Oracle */
#include <stdlib.h>
#if defined( INC_ALL )
#include "crypt.h"
#include "context.h"
#elif defined( INC_CHILD )
#include "../crypt.h"
#include "../misc/context.h"
#else
#include "crypt.h"
#include "misc/context.h"
#endif /* Compiler-specific includes */
#if defined( INC_ALL )
#ifdef __TANDEMNSK__
#include "bnprime.h"
#else
#include "bn_prime.h"
#endif /* __TANDEMNSK__ */
#elif defined( INC_CHILD )
#include "../bn/bn_prime.h"
#else
#include "bn/bn_prime.h"
#endif /* Compiler-specific includes */
/****************************************************************************
* *
* Determine Discrete Log Exponent Bits *
* *
****************************************************************************/
/* The following function (provided by Colin Plumb) is used to calculate the
appropriate size exponent for a given prime size which is required to
provide equivalent security from small-exponent attacks
This is based on a paper by Michael Wiener on | The function defined
the difficulty of the two attacks, which has | below (not part of the
the following table: | original paper)
| produces the following
Table 1: Subgroup Sizes to Match Field Sizes | results:
|
Size of p Cost of each attack Size of q | Output Error
(bits) (instructions or (bits) | (+ is safe)
modular multiplies) |
|
512 9 x 10^17 119 | 137 +18
768 6 x 10^21 145 | 153 +8
1024 7 x 10^24 165 | 169 +4
1280 3 x 10^27 183 | 184 +1
1536 7 x 10^29 198 | 198 +0
1792 9 x 10^31 212 | 212 +0
2048 8 x 10^33 225 | 225 +0
2304 5 x 10^35 237 | 237 +0
2560 3 x 10^37 249 | 249 +0
2816 1 x 10^39 259 | 260 +1
3072 3 x 10^40 269 | 270 +1
3328 8 x 10^41 279 | 280 +1
3584 2 x 10^43 288 | 289 +1
3840 4 x 10^44 296 | 297 +1
4096 7 x 10^45 305 | 305 +0
4352 1 x 10^47 313 | 313 +0
4608 2 x 10^48 320 | 321 +1
4864 2 x 10^49 328 | 329 +1
5120 3 x 10^50 335 | 337 +2
This function fits a curve to this, which overestimates the size of the
exponent required, but by a very small amount in the important 1000-4000
bit range. It is a quadratic curve up to 3840 bits, and a linear curve
past that. They are designed to be C(1) (have the same value and the same
slope) at the point where they meet */
#define AN 1L /* a = -AN/AD/65536, the quadratic coefficient */
#define AD 3L
#define M 8L /* Slope = M/256, i.e. 1/32 where linear starts */
#define TX 3840L /* X value at the slope point, where linear starts */
#define TY 297L /* Y value at the slope point, where linear starts */
/* For a slope of M at the point (TX,TY), we only have one degree of freedom
left in a quadratic curve, so use the coefficient of x^2, namely a, as
that free parameter.
y = -AN/AD*((x-TX)/256)^2 + M*(x-TX)/256 + TY
= -AN*(x-TX)*(x-TX)/AD/256/256 + M*x/256 - M*TX/256 + TY
= -AN*x*x/AD/256/256 + 2*AN*x*TX/AD/256/256 - AN*TX*TX/AD/256/256 \
+ M*x/256 - M*TX/256 + TY
= -AN*(x/256)^2/AD + 2*AN*(TX/256)*(x/256)/AD + M*(x/256) \
- AN*(TX/256)^2/AD - M*(TX/256) + TY
= (AN*(2*TX/256 - x/256) + M*AD)*x/256/AD - (AN*(TX/256)/AD + M)*TX/256 \
+ TY
= (AN*(2*TX/256 - x/256) + M*AD)*x/256/AD \
- (AN*(TX/256) + M*AD)*TX/256/AD + TY
= ((M*AD + AN*(2*TX/256 - x/256))*x - (AN*(TX/256)+M*AD)*TX)/256/AD + TY
= ((M*AD + AN*(2*TX - x)/256)*x - (AN*(TX/256)+M*AD)*TX)/256/AD + TY
= ((M*AD + AN*(2*TX - x)/256)*x - (M*AD + AN*TX/256)*TX)/256/AD + TY
= (((256*M*AD+2*AN*TX-AN*x)/256)*x - (M*AD + AN*TX/256)*TX)/256/AD + TY
Since this is for the range 0...TX, in order to avoid having any
intermediate results less than 0, we need one final rearrangement, and a
compiler can easily take the constant-folding from there...
= TY + (((256*M*AD+2*AN*TX-AN*x)/256)*x - (M*AD + AN*TX/256)*TX)/256/AD
= TY - ((M*AD + AN*TX/256)*TX - ((256*M*AD+2*AN*TX-AN*x)/256)*x)/256/AD
*/
static int getDLPexpSize( const int primeBits )
{
long value; /* Necessary to avoid braindamage on 16-bit compilers */
/* If it's over TX bits, it's linear */
if( primeBits > TX )
value = M * primeBits / 256 - M * TX / 256 + TY;
else
/* It's quadratic */
value = TY - ( ( M * AD + AN * TX / 256 ) * TX - \
( ( 256 * M * AD + AN * 2 * TX - AN * primeBits ) / 256 ) * \
primeBits ) / ( AD * 256 );
/* Various standards require a minimum of 160 bits so we always return at
least that size even if it's not necessary */
return( value > 160 ? ( int ) value : 160 );
}
/****************************************************************************
* *
* Generate Random Bignum *
* *
****************************************************************************/
/* Generate a bignum of a specified length, with the given high and low 8
bits. 'high' is merged into the high 8 bits of the number (set it to 0x80
to ensure that the number is exactly 'bits' bits long, i.e. 2^(bits-1) <=
bn < 2^bits), 'low' is merged into the low 8 bits (set it to 1 to ensure
that the number is odd). In almost all cases used in cryptlib, 'high' is
set to 0xC0 and low is set to 0x01.
We don't need to pagelock the bignum buffer we're using because it's being
accessed continuously while there's data in it, so there's little chance
it'll be swapped unless the system is already thrashing */
int generateBignum( BIGNUM *bn, const int noBits, const BYTE high,
const BYTE low )
{
RESOURCE_DATA msgData;
BYTE buffer[ CRYPT_MAX_PKCSIZE ];
int noBytes = bitsToBytes( noBits ), status;
/* Clear the return value */
BN_zero( bn );
/* Load the random data into the bignum buffer */
setMessageData( &msgData, buffer, noBytes );
status = krnlSendMessage( SYSTEM_OBJECT_HANDLE, IMESSAGE_GETATTRIBUTE_S,
&msgData, CRYPT_IATTRIBUTE_RANDOM );
if( cryptStatusError( status ) )
{
zeroise( buffer, noBytes );
return( status );
}
/* Merge in the specified low bits, mask off any excess high bits, and
merge in the specified high bits. This is a bit more complex than
just masking in the byte values because the bignum may not be a
multiple of 8 bytes long */
buffer[ noBytes - 1 ] |= low;
buffer[ 0 ] &= 255 >> ( -noBits & 7 );
buffer[ 0 ] |= high >> ( -noBits & 7 );
if( noBytes > 1 && ( noBits & 7 ) )
buffer[ 1 ] |= high << ( noBits & 7 );
/* Turn the contents of the buffer into a bignum and zeroise the buffer */
status = ( BN_bin2bn( buffer, noBytes, bn ) == NULL ) ? \
CRYPT_ERROR_MEMORY : CRYPT_OK;
zeroise( buffer, noBytes );
return( status );
}
/****************************************************************************
* *
* Generate Non-specific Primes *
* *
****************************************************************************/
/* #include 4k of EAY copyright */
/* The following define is necessary in memory-starved environments. It
controls the size of the table used for the sieving */
#if defined( CONFIG_CONSERVE_MEMORY )
#define EIGHT_BIT
#endif /* CONFIG_CONSERVE_MEMORY */
/* The number of primes in the sieve (and their values) that result in a
given number of candidates remaining from 40,000. Even the first 100
primes weed out 91% of all the candidates, and after 500 you're only
removing a handful for each 100 extra primes.
Number Prime Candidates left
Values from 40,000
-------- --------- ---------------
0- 99 0- 541 3564
100-199 541-1223 3175
200-299 1223-1987 2969
300-399 1987-2741 2845
400-499 2741-3571 2755
500-599 3571-4409 2688
600-699 4409-5279 2629
700-799 5279-6133 2593
800-899 6133-6997 2555
900-999 6997-7919 2521
There is in fact an even faster prime tester due to Dan Piponi that uses
C++ templates as a universal computer and performs the primality test at
compile time, however this requires the use of a fairly advanced C++
compiler and isn't amenable to generating different primes */
/* The number of iterations of Miller-Rabin for an error probbility of
(1/2)^80, from HAC */
#define getNoPrimeChecks( noBits ) \
( ( noBits < 150 ) ? 18 : ( noBits < 200 ) ? 15 : \
( noBits < 250 ) ? 12 : ( noBits < 300 ) ? 9 : \
( noBits < 350 ) ? 8 : ( noBits < 400 ) ? 7 : \
( noBits < 500 ) ? 6 : ( noBits < 600 ) ? 5 : \
( noBits < 800 ) ? 4 : ( noBits < 1250 ) ? 3 : 2 )
/* Enable the following to cross-check the Miller-Rabin test using a Fermat
test and an alternative form of the Miller-Rabin test that merges the
test loop and the modexp at the start. Note that this displays
diagnostic timing output and expects to use Pentium performance counters
for timing, so it's only (optionally) enabled for Win32 debug */
#if defined( __WIN32__ ) && !defined( NDEBUG ) && 0
#define CHECK_PRIMETEST
#endif /* Win32 debug */
/* The size of the sieve array - 1 memory page (on most CPUs) = 4K candidate
values. When changing this value the LFSR parameters need to be adjusted
to match */
#define SIEVE_SIZE 4096
/* When we're doing a sieve of a singleton candidate, we don't run through
the whole range of sieve values since we run into the law of diminshing
returns after a certain point. The following value sieves with every
prime under 1000 */
#if NUMPRIMES < ( 21 * 8 )
#define FAST_SIEVE_NUMPRIMES NUMPRIMES
#else
#define FAST_SIEVE_NUMPRIMES ( 21 * 8 )
#endif /* Small prime table */
/* Set up the sieve array for the number. Every position that contains
a zero is non-divisible by all of the small primes */
static void initSieve( BOOLEAN *sieveArray, const BIGNUM *candidate )
{
int i;
memset( sieveArray, 0, SIEVE_SIZE * sizeof( BOOLEAN ) );
/* Walk down the list of primes marking the appropriate position in the
array as divisible by the prime. We start at index 1, since the
candidate will never be divisible by 2 (== primes[ 0 ]) */
for( i = 1; i < NUMPRIMES; i++ )
{
unsigned int step = primes[ i ];
int sieveIndex = ( int ) BN_mod_word( candidate, step );
/* Determine the correct start index for this value */
if( sieveIndex & 1 )
sieveIndex = ( step - sieveIndex ) / 2;
else
if( sieveIndex )
sieveIndex = ( ( step * 2 ) - sieveIndex ) / 2;
/* Mark each multiple of the divisor as being divisible */
while( sieveIndex < SIEVE_SIZE )
{
sieveArray[ sieveIndex ] = 1;
sieveIndex += step;
}
}
}
/* An LFSR to step through each entry in the sieve array. This isn't a true
pseudorandom selection since all it's really doing is going through the
numbers in a linear order with a different starting point, but it's good
enough as a randomiser */
#define LFSR_POLYNOMIAL 0x1053
#define LFSR_MASK 0x1000
static int nextEntry( int value )
{
assert( LFSR_MASK == SIEVE_SIZE );
/* Get the next value: Multiply by x and reduce by the polynomial */
value <<= 1;
if( value & LFSR_MASK )
value ^= LFSR_POLYNOMIAL;
return( value );
}
/* A one-off sieve check for when we're testing a singleton rather than
running over a range of values */
static BOOLEAN primeSieve( const BIGNUM *candidate )
{
int i;
for( i = 1; i < FAST_SIEVE_NUMPRIMES; i++ )
if( BN_mod_word( candidate, primes[ i ] ) == 0 )
return( FALSE );
return( TRUE );
}
#ifdef CHECK_PRIMETEST
/* Witness function, modified from original BN code. Found at a UFO crash
site. This looks nothing like a standard Miller-Rabin test because it
merges the modexp that usually needs to be performed as the first
portion of the test process and the remainder of the checking. Destroys
param6 + 7 */
static int witnessOld( PKC_INFO *pkcInfo, BIGNUM *a, BIGNUM *n, BIGNUM *n1,
BIGNUM *mont_n1, BIGNUM *mont_1,
BN_MONT_CTX *montCTX_n )
{
BIGNUM *y = &pkcInfo->param6;
BIGNUM *yPrime = &pkcInfo->param7; /* Safe to destroy */
BN_CTX *ctx = &pkcInfo->bnCTX;
BIGNUM *mont_a = &ctx->bn[ ctx->tos++ ];
const int k = BN_num_bits( n1 );
int i, bnStatus = BN_STATUS;
/* All values are manipulated in their Montgomery form, so before we
begin we have to convert a to this form as well */
if( !BN_to_montgomery( mont_a, a, montCTX_n, &pkcInfo->bnCTX ) )
{
ctx->tos--;
return( CRYPT_ERROR_FAILED );
}
CKPTR( BN_copy( y, mont_1 ) );
for ( i = k - 1; i >= 0; i-- )
{
/* Perform the y^2 mod n check. yPrime = y^2 mod n, if yPrime == 1
it's composite (this condition is virtually never met) */
CK( BN_mod_mul_montgomery( yPrime, y, y, montCTX_n,
&pkcInfo->bnCTX ) );
if( bnStatusError( bnStatus ) || \
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -