re-add the primitives to generate primes and test for prime.

This commit is contained in:
Vincent Hanquez 2015-03-29 10:55:46 +01:00
parent d9b16a529e
commit c3d9570881

View File

@ -1,3 +1,10 @@
-- |
-- Module : Crypto.Number.Prime
-- License : BSD-style
-- Maintainer : Vincent Hanquez <vincent@snarc.org>
-- Stability : experimental
-- Portability : Good
{-# LANGUAGE CPP #-} {-# LANGUAGE CPP #-}
{-# LANGUAGE BangPatterns #-} {-# LANGUAGE BangPatterns #-}
#ifndef MIN_VERSION_integer_gmp #ifndef MIN_VERSION_integer_gmp
@ -6,31 +13,25 @@
#if MIN_VERSION_integer_gmp(0,5,1) #if MIN_VERSION_integer_gmp(0,5,1)
{-# LANGUAGE MagicHash #-} {-# LANGUAGE MagicHash #-}
#endif #endif
-- |
-- Module : Crypto.Number.Prime
-- License : BSD-style
-- Maintainer : Vincent Hanquez <vincent@snarc.org>
-- Stability : experimental
-- Portability : Good
module Crypto.Number.Prime module Crypto.Number.Prime
( (
{- generatePrime
generatePrime
, generateSafePrime , generateSafePrime
, isProbablyPrime , isProbablyPrime
, findPrimeFrom , findPrimeFrom
, findPrimeFromWith , findPrimeFromWith
, primalityTestMillerRabin , primalityTestMillerRabin
-} , primalityTestNaive
primalityTestNaive
, primalityTestFermat , primalityTestFermat
, isCoprime , isCoprime
) where ) where
import Control.Applicative
import Crypto.Number.Generate import Crypto.Number.Generate
import Crypto.Number.Basic (sqrti, gcde_binary) import Crypto.Number.Basic (sqrti, gcde_binary)
import Crypto.Number.ModArithmetic (exponantiation) import Crypto.Number.ModArithmetic (exponantiation)
import Crypto.Random.Types
#if MIN_VERSION_integer_gmp(0,5,1) #if MIN_VERSION_integer_gmp(0,5,1)
import GHC.Integer.GMP.Internals import GHC.Integer.GMP.Internals
@ -39,101 +40,103 @@ import GHC.Base
import Data.Bits import Data.Bits
#endif #endif
{-
-- | returns if the number is probably prime. -- | returns if the number is probably prime.
-- first a list of small primes are implicitely tested for divisibility, -- first a list of small primes are implicitely tested for divisibility,
-- then a fermat primality test is used with arbitrary numbers and -- then a fermat primality test is used with arbitrary numbers and
-- then the Miller Rabin algorithm is used with an accuracy of 30 recursions -- then the Miller Rabin algorithm is used with an accuracy of 30 recursions
isProbablyPrime :: CPRG g => g -> Integer -> (Bool, g) isProbablyPrime :: MonadRandom m => Integer -> m Bool
isProbablyPrime rng !n isProbablyPrime !n
| any (\p -> p `divides` n) (filter (< n) firstPrimes) = (False, rng) | any (\p -> p `divides` n) (filter (< n) firstPrimes) = return False
| primalityTestFermat 50 (n`div`2) n = primalityTestMillerRabin rng 30 n | primalityTestFermat 50 (n`div`2) n = primalityTestMillerRabin 30 n
| otherwise = (False, rng) | otherwise = return False
-- | generate a prime number of the required bitsize -- | generate a prime number of the required bitsize
generatePrime :: CPRG g => g -> Int -> (Integer, g) generatePrime :: MonadRandom m => Int -> m Integer
generatePrime rng bits = generatePrime bits = do
let (sp, rng') = generateOfSize rng bits sp <- generateOfSize bits
in findPrimeFrom rng' sp findPrimeFrom sp
-- | generate a prime number of the form 2p+1 where p is also prime. -- | generate a prime number of the form 2p+1 where p is also prime.
-- it is also knowed as a Sophie Germaine prime or safe prime. -- it is also knowed as a Sophie Germaine prime or safe prime.
-- --
-- The number of safe prime is significantly smaller to the number of prime, -- The number of safe prime is significantly smaller to the number of prime,
-- as such it shouldn't be used if this number is supposed to be kept safe. -- as such it shouldn't be used if this number is supposed to be kept safe.
generateSafePrime :: CPRG g => g -> Int -> (Integer, g) generateSafePrime :: MonadRandom m => Int -> m Integer
generateSafePrime rng bits = generateSafePrime bits = do
let (sp, rng') = generateOfSize rng bits sp <- generateOfSize bits
(p, rng'') = findPrimeFromWith rng' (\g i -> isProbablyPrime g (2*i+1)) (sp `div` 2) p <- findPrimeFromWith (\i -> isProbablyPrime (2*i+1)) (sp `div` 2)
in (2*p+1, rng'') return (2*p+1)
-- | find a prime from a starting point where the property hold. -- | find a prime from a starting point where the property hold.
findPrimeFromWith :: CPRG g => g -> (g -> Integer -> (Bool,g)) -> Integer -> (Integer, g) findPrimeFromWith :: MonadRandom m => (Integer -> m Bool) -> Integer -> m Integer
findPrimeFromWith rng prop !n findPrimeFromWith prop !n
| even n = findPrimeFromWith rng prop (n+1) | even n = findPrimeFromWith prop (n+1)
| otherwise = case isProbablyPrime rng n of | otherwise = do
(False, rng') -> findPrimeFromWith rng' prop (n+2) primed <- isProbablyPrime n
(True, rng') -> if not primed
case prop rng' n of then findPrimeFromWith prop (n+2)
(False, rng'') -> findPrimeFromWith rng'' prop (n+2) else do
(True, rng'') -> (n, rng'') validate <- prop n
if validate
then return n
else findPrimeFromWith prop (n+2)
-- | find a prime from a starting point with no specific property. -- | find a prime from a starting point with no specific property.
findPrimeFrom :: CPRG g => g -> Integer -> (Integer, g) findPrimeFrom :: MonadRandom m => Integer -> m Integer
findPrimeFrom rng n = findPrimeFrom n =
#if MIN_VERSION_integer_gmp(0,5,1) #if MIN_VERSION_integer_gmp(0,5,1)
(nextPrimeInteger n, rng) return $ nextPrimeInteger n
#else #else
findPrimeFromWith rng (\g _ -> (True, g)) n findPrimeFromWith (\_ -> return True) n
#endif #endif
-- | Miller Rabin algorithm return if the number is probably prime or composite. -- | Miller Rabin algorithm return if the number is probably prime or composite.
-- the tries parameter is the number of recursion, that determines the accuracy of the test. -- the tries parameter is the number of recursion, that determines the accuracy of the test.
primalityTestMillerRabin :: CPRG g => g -> Int -> Integer -> (Bool, g) primalityTestMillerRabin :: MonadRandom m => Int -> Integer -> m Bool
#if MIN_VERSION_integer_gmp(0,5,1) #if MIN_VERSION_integer_gmp(0,5,1)
primalityTestMillerRabin rng (I# tries) !n = primalityTestMillerRabin (I# tries) !n =
case testPrimeInteger n tries of case testPrimeInteger n tries of
0# -> (False, rng) 0# -> return False
_ -> (True, rng) _ -> return True
#else #else
primalityTestMillerRabin rng tries !n primalityTestMillerRabin tries !n
| n <= 3 = error "Miller-Rabin requires tested value to be > 3" | n <= 3 = error "Miller-Rabin requires tested value to be > 3"
| even n = (False, rng) | even n = return False
| tries <= 0 = error "Miller-Rabin tries need to be > 0" | tries <= 0 = error "Miller-Rabin tries need to be > 0"
| otherwise = let (witnesses, rng') = generateTries tries rng | otherwise = loop <$> generateTries tries
in (loop witnesses, rng') where
where !nm1 = n-1 !nm1 = n-1
!nm2 = n-2 !nm2 = n-2
(!s,!d) = (factorise 0 nm1) (!s,!d) = (factorise 0 nm1)
generateTries 0 g = ([], g) generateTries 0 = return []
generateTries t g = let (v,g') = generateBetween g 2 nm2 generateTries t = do
(vs,g'') = generateTries (t-1) g' v <- generateBetween 2 nm2
in (v:vs, g'') vs <- generateTries (t-1)
return (v:vs)
-- factorise n-1 into the form 2^s*d -- factorise n-1 into the form 2^s*d
factorise :: Integer -> Integer -> (Integer, Integer) factorise :: Integer -> Integer -> (Integer, Integer)
factorise !si !vi factorise !si !vi
| vi `testBit` 0 = (si, vi) | vi `testBit` 0 = (si, vi)
| otherwise = factorise (si+1) (vi `shiftR` 1) -- probably faster to not shift v continously, but just once. | otherwise = factorise (si+1) (vi `shiftR` 1) -- probably faster to not shift v continously, but just once.
expmod = exponantiation expmod = exponantiation
-- when iteration reach zero, we have a probable prime -- when iteration reach zero, we have a probable prime
loop [] = True loop [] = True
loop (w:ws) = let x = expmod w d n loop (w:ws) = let x = expmod w d n
in if x == (1 :: Integer) || x == nm1 in if x == (1 :: Integer) || x == nm1
then loop ws then loop ws
else loop' ws ((x*x) `mod` n) 1 else loop' ws ((x*x) `mod` n) 1
-- loop from 1 to s-1. if we reach the end then it's composite -- loop from 1 to s-1. if we reach the end then it's composite
loop' ws !x2 !r loop' ws !x2 !r
| r == s = False | r == s = False
| x2 == 1 = False | x2 == 1 = False
| x2 /= nm1 = loop' ws ((x2*x2) `mod` n) (r+1) | x2 /= nm1 = loop' ws ((x2*x2) `mod` n) (r+1)
| otherwise = loop ws | otherwise = loop ws
#endif #endif
-}
{- {-
n < z -> witness to test n < z -> witness to test