```MODULE CryptoPrimes; (* 2002.9.11 g.f. based on 'bn_prime.c' Copyright of the original: <-- middle click here *) IMPORT B := CryptoBigNumbers, Out := KernelLog; CONST N = 2048; TYPE BigNumber = B.BigNumber; VAR one, two: BigNumber; primes: ARRAY N OF LONGINT; PROCEDURE NewPrime*( bits: INTEGER; safe: BOOLEAN ): BigNumber; VAR checks, i: INTEGER; t, p: BigNumber; BEGIN checks := Checks( bits ); i := 1; LOOP p := ProbablePrime( bits ); Out.Char( '.' ); IF IsPrime( p, checks, FALSE ) THEN Out.Char( 'P' ); IF ~safe THEN Out.Ln; RETURN p ELSE Out.Int( i, 0 ); INC( i ); (* for "safe prime" generation, check that (p-1)/2 is prime. *) B.Copy( p, t ); B.Shift( t, -1 ); IF IsPrime( t, checks, FALSE ) THEN Out.String( " (safe prime)" ); Out.Ln; RETURN p END; Out.String( " (not safe)" ); Out.Ln; END END; END END NewPrime; PROCEDURE NewDHPrime*( bits: INTEGER; safe: BOOLEAN; add, rem: BigNumber ): BigNumber; VAR checks: INTEGER; t, p: BigNumber; BEGIN checks := Checks( bits ); LOOP IF safe THEN p := ProbableDHPrimeSafe( bits, add, rem ) ELSE p := ProbableDHPrime( bits, add, rem ) END; IF IsPrime( p, checks, FALSE ) THEN IF ~safe THEN RETURN p ELSE (* for "safe prime" generation, check that (p-1)/2 is prime. *) B.Copy( p, t ); B.Shift( t, -1 ); IF IsPrime( t, checks, FALSE ) THEN RETURN p END END END END END NewDHPrime; PROCEDURE Checks( b: LONGINT ): INTEGER; (* number of Miller-Rabin iterations for an error rate of less than 2^-80 for random 'b'-bit input, b >= 100 (taken from table 4.4 in the Handbook of Applied Cryptography [Menezes, van Oorschot, Vanstone; CRC Press 1996] *) VAR t: INTEGER; BEGIN ASSERT( b >= 100 ); IF b >= 1300 THEN t := 2 ELSIF b >= 850 THEN t := 3 ELSIF b >= 650 THEN t := 4 ELSIF b >= 550 THEN t := 5 ELSIF b >= 450 THEN t := 6 ELSIF b >= 400 THEN t := 7 ELSIF b >= 350 THEN t := 8 ELSIF b >= 300 THEN t := 9 ELSIF b >= 250 THEN t := 12 ELSIF b >= 200 THEN t := 15 ELSIF b >= 150 THEN t := 18 ELSE t := 27 END; RETURN t END Checks; PROCEDURE ProbablePrime( bits: INTEGER ): BigNumber; VAR t: BigNumber; delta, i: LONGINT; p: BigNumber; mods: ARRAY N OF LONGINT; BEGIN LOOP p := B.NewRand( bits, 1, 1 ); FOR i := 0 TO N - 1 DO mods[i] := B.ModWord( p, primes[i] ) END; (* check that p is not a prime and also that gcd( p-1, primes) = 1 (except for 2) *) i := 0; delta := 0; LOOP INC( i ); IF i >= N THEN B.AssignInt( t, delta ); p := B.Add( p, t ); RETURN p END; IF (mods[i] + delta) MOD primes[i] <= 1 THEN INC( delta, 2 ); i := 0; IF delta < 0 THEN (* overfow! try new random *) EXIT END; END; END END END ProbablePrime; PROCEDURE ProbableDHPrime( bits: INTEGER; add, rem: BigNumber ): BigNumber; VAR d, r, p: BigNumber; i: LONGINT; BEGIN p := B.NewRand( bits, 0, 1 ); (* we need (p - rem) mod add = 0 *) r := B.Mod( p, add ); p := B.Sub( p, r ); IF B.Zero( rem ) THEN B.Inc( p ) ELSE p := B.Add( p, rem ) END; (* we now have a random number 'p' to test. *) i := 0; LOOP INC( i ); IF i >= N THEN RETURN p END; B.AssignInt( d, primes[i] ); r := B.Mod( p, d ); IF B.Zero( r ) THEN p := B.Add( p, add ); i := 0 END; END END ProbableDHPrime; PROCEDURE ProbableDHPrimeSafe( bits: INTEGER; padd, rem: BigNumber ): BigNumber; VAR d, q, r, qr, qadd, p: BigNumber; i: LONGINT; BEGIN B.Copy( padd, qadd ); B.Shift( qadd, -1 ); q := B.NewRand( bits, 0, 1 ); r := B.Mod( q, qadd ); q := B.Sub( q, r ); IF B.Zero( rem ) THEN B.Inc( q ) ELSE B.Copy( rem, r ); B.Shift( r, -1 ); q := B.Add( q, r ) END; (* we now have a random number to test. *) B.Copy( q, p ); B.Shift( p, 1 ); B.Inc( p ); i := 0; LOOP INC( i ); IF i >= N THEN RETURN p END; B.AssignInt( d, primes[i] ); r := B.Mod( p, d ); qr := B.Mod( q, d ); IF B.Zero( r ) OR B.Zero( qr ) THEN p := B.Add( p, padd ); q := B.Add( q, qadd ); i := 0 END; END END ProbableDHPrimeSafe; PROCEDURE IsPrime*( a: BigNumber; checks: INTEGER; doTrialDiv: BOOLEAN ): BOOLEAN; VAR i, k: INTEGER; A, A1, A1odd, check: BigNumber; BEGIN IF checks = 0 THEN checks := Checks( B.BitSize( a ) ) END; IF ~ODD( a.d[0] ) THEN RETURN FALSE END; IF doTrialDiv THEN FOR i := 1 TO N - 1 DO IF B.ModWord( a, primes[i] ) = 0 THEN RETURN FALSE END END END; B.Copy( a, A ); IF A.neg THEN B.Neg( A ) END; B.Copy( A, A1 ); B.Dec( A1 ); IF B.Zero( A1 ) THEN RETURN FALSE END; (* write A1 as A1odd * 2^k *) k := 1; WHILE ~B.BitSet( A1, k ) DO INC( k ) END; B.Copy( A1, A1odd ); B.Shift( A1odd, -k ); FOR i := 1 TO checks DO check := B.NewRand( B.BitSize( A1 ), 0, 0 ); IF B.Cmp( check, A1 ) >= 0 THEN check := B.Sub( check, A1 ) END; B.Inc( check ); (* now 1 <= check < A *) IF ~witness( check, A, A1, A1odd, k ) THEN RETURN FALSE END; END; RETURN TRUE END IsPrime; PROCEDURE witness( W, a, a1, a1odd: BigNumber; k: INTEGER): BOOLEAN; VAR w: BigNumber; BEGIN w := B.ModExp( W, a1odd, a ); IF B.Cmp( w, one ) = 0 THEN (* probably prime *) RETURN TRUE END; IF B.Cmp( w, a1 ) = 0 THEN RETURN TRUE END; (* w = -1 (mod a), a is probably prime *) WHILE k > 0 DO w := B.ModMul( w, w, a ); (* w = w^2 mod a *) IF B.Cmp( w, one ) = 0 THEN RETURN FALSE END; (* a is composite, otherwise a previous w would have been = -1 (mod a) *) IF B.Cmp( w, a1 ) = 0 THEN RETURN TRUE END; (* w = -1 (mod a), a is probably prime *) DEC( k ) END; (* If we get here, w is the (a-1)/2-th power of the original w, *) (* and it is neither -1 nor +1 -- so a cannot be prime *) RETURN FALSE END witness; PROCEDURE Init; VAR sieve: ARRAY N OF SET; i, j, p, n: LONGINT; BEGIN (* compute the first N small primes *) FOR i := 0 TO N - 1 DO sieve[i] := {0..31} END; primes[0] := 2; n := 1; i := 1; WHILE n < N DO IF i MOD 32 IN sieve[i DIV 32] THEN p := 2*i + 1; primes[n] := p; INC( n ); j := i; WHILE j DIV 32 < N DO EXCL( sieve[j DIV 32], j MOD 32 ); INC( j, p ) END; END; INC( i ) END; END Init; BEGIN Init; B.AssignInt( one, 1 ); B.AssignInt( two, 2 ); END CryptoPrimes. ```