start restructuring project
This commit is contained in:
43
src/Factorization/PollardPminus1.hs
Normal file
43
src/Factorization/PollardPminus1.hs
Normal file
@@ -0,0 +1,43 @@
|
||||
import Text.Read (readMaybe)
|
||||
import System.Exit (exitSuccess)
|
||||
import Utils (askNumber)
|
||||
import ModularArithmeticUtils (modExp)
|
||||
|
||||
main :: IO ()
|
||||
main = do
|
||||
n <- askNumber "Enter an integer n>1 to factor:"
|
||||
b <- askNumber "Enter the bound B:"
|
||||
putStrLn ("n = " ++ show n)
|
||||
|
||||
case pollardP1 n b of
|
||||
Just factor -> do
|
||||
putStrLn ("Found factor: " ++ show factor)
|
||||
let otherFactor = n `div` factor
|
||||
putStrLn ("Other factor: " ++ show otherFactor)
|
||||
Nothing -> do
|
||||
putStrLn "Failed to find a factor using Pollard p-1 method. n could be prime or different parameters needed."
|
||||
|
||||
-- Pollard p-1 factorization algorithm
|
||||
pollardP1 :: Integer -> Integer -> Maybe Integer
|
||||
pollardP1 n b = tryBases [2..5]
|
||||
where
|
||||
tryBases [] = Nothing
|
||||
tryBases (a:as) =
|
||||
case pollardP1WithParams n a b of
|
||||
Just factor -> Just factor
|
||||
Nothing -> tryBases as
|
||||
|
||||
-- Pollard p-1 with configurable parameters
|
||||
pollardP1WithParams :: Integer -> Integer -> Integer -> Maybe Integer
|
||||
pollardP1WithParams n a b = go a 2
|
||||
where
|
||||
go currentA k
|
||||
| k > b = Nothing -- Bound exceeded, failed to find factor
|
||||
| otherwise =
|
||||
let newA = modExp currentA k n
|
||||
d = gcd (newA - 1) n
|
||||
in if d > 1 && d < n
|
||||
then Just d -- Found a non-trivial factor
|
||||
else if d == n
|
||||
then Nothing -- Found n itself, try with different parameters
|
||||
else go newA (k + 1)
|
||||
24
src/Factorization/fermat-factorization.hs
Normal file
24
src/Factorization/fermat-factorization.hs
Normal file
@@ -0,0 +1,24 @@
|
||||
import Text.Read (readMaybe)
|
||||
import System.Exit (exitSuccess)
|
||||
import Utils (askNumber)
|
||||
|
||||
main :: IO ()
|
||||
main = do
|
||||
n <- askNumber "Enter an integer:"
|
||||
let (d, root) = findIntegerSqrt n
|
||||
q = (root - d)
|
||||
p = (root + d)
|
||||
putStrLn ("n = " ++ show n)
|
||||
putStrLn ("Found d = " ++ show d ++ ", sqrt(n + d^2) = " ++ show root)
|
||||
putStrLn ("q = " ++ show q ++ ", p = " ++ show p)
|
||||
|
||||
-- Find the smallest integer d >= 0 such that sqrt(n + d^2) is an integer
|
||||
findIntegerSqrt :: Integer -> (Integer, Integer)
|
||||
findIntegerSqrt n = go 0
|
||||
where
|
||||
go d =
|
||||
let val = n + d*d
|
||||
root = floor (sqrt (fromIntegral val))
|
||||
in if root * root == val
|
||||
then (d, root)
|
||||
else go (d + 1)
|
||||
60
src/ModularArithmeticUtils.hs
Normal file
60
src/ModularArithmeticUtils.hs
Normal file
@@ -0,0 +1,60 @@
|
||||
module ModularArithmeticUtils (modExp, modMul, legendre, factorOutTwos, jacobi) where
|
||||
|
||||
import Data.Bits (testBit, shiftR)
|
||||
|
||||
-- Modular exponentiation: compute a^b mod m
|
||||
modExp :: Integer -> Integer -> Integer -> Integer
|
||||
modExp b 0 m = 1
|
||||
modExp b e m = t * modExp ((b * b) `mod` m) (shiftR e 1) m `mod` m
|
||||
where t = if testBit e 0 then b `mod` m else 1
|
||||
|
||||
-- Compute (a * b) mod m
|
||||
modMul :: Integer -> Integer -> Integer -> Integer
|
||||
modMul a b m = (a * b) `mod` m
|
||||
|
||||
-- Compute the Legendre symbol (a/p)
|
||||
-- (a/p) = a^((p-1)/2) mod p
|
||||
-- (a/p) in {-1, 0, 1}
|
||||
legendre :: Integer -> Integer -> Integer
|
||||
legendre a p = modExp a ((p - 1) `div` 2) p
|
||||
|
||||
-- Factor p-1 as q * 2^s, where q is odd
|
||||
factorOutTwos :: Integer -> (Integer, Integer)
|
||||
factorOutTwos n = go n 0
|
||||
where
|
||||
go x s
|
||||
| even x = go (x `div` 2) (s + 1)
|
||||
| otherwise = (x, s)
|
||||
|
||||
-- jacobi symbol (a/n)
|
||||
-- (a/n) in {-1, 0, 1}
|
||||
jacobi :: Integer -> Integer -> Integer
|
||||
jacobi a n
|
||||
| n <= 0 = error "jacobi: n must be positive"
|
||||
| even n = error "jacobi: n must be odd"
|
||||
| a `mod` n == 0 = 0
|
||||
| otherwise = go (a `mod` n) n 1
|
||||
where
|
||||
go :: Integer -> Integer -> Integer -> Integer
|
||||
go 0 _ _ = 0
|
||||
go 1 _ s = s
|
||||
go x m s =
|
||||
let (xOdd, e) = factorOutTwos x
|
||||
s' = if even e then s else s * jacobi2 m
|
||||
in
|
||||
if xOdd == 1
|
||||
then s'
|
||||
else
|
||||
let s'' = if (xOdd `mod` 4 == 3) && (m `mod` 4 == 3) then -s' else s'
|
||||
xNext = m `mod` xOdd
|
||||
in go xNext xOdd s''
|
||||
|
||||
-- Jacobi(2, m), helper function
|
||||
jacobi2 :: Integer -> Integer
|
||||
jacobi2 m =
|
||||
case m `mod` 8 of
|
||||
1 -> 1
|
||||
7 -> 1
|
||||
3 -> -1
|
||||
5 -> -1
|
||||
_ -> error "jacobi2: unexpected (m mod 8) for odd m"
|
||||
17
src/Primes/FermatPrimeTest.hs
Normal file
17
src/Primes/FermatPrimeTest.hs
Normal file
@@ -0,0 +1,17 @@
|
||||
module Primes.FermatPrimeTest (fermatPrimeTest) where
|
||||
|
||||
import System.Random (randomRIO)
|
||||
import ModularArithmeticUtils (modExp)
|
||||
|
||||
fermatPrimeTest :: Integer -> Integer -> IO Bool
|
||||
fermatPrimeTest n k
|
||||
| n <= 3 = return (n == 2 || n == 3)
|
||||
| even n = return False
|
||||
| otherwise = go k
|
||||
where
|
||||
go 0 = return True
|
||||
go i = do
|
||||
a <- randomRIO (2, n - 2)
|
||||
if modExp a (n - 1) n /= 1
|
||||
then return False
|
||||
else go (i - 1)
|
||||
30
src/Primes/MillerRabin.hs
Normal file
30
src/Primes/MillerRabin.hs
Normal file
@@ -0,0 +1,30 @@
|
||||
module Primes.MillerRabin (millerRabin) where
|
||||
|
||||
import System.Random (randomRIO)
|
||||
import ModularArithmeticUtils (modExp, factorOutTwos)
|
||||
|
||||
millerRabinWitness :: Integer -> Integer -> Integer -> Integer -> Bool
|
||||
millerRabinWitness n a d s =
|
||||
let x0 = modExp a d n
|
||||
in x0 == 1 || x0 == n - 1 || loop x0 (s - 1)
|
||||
where
|
||||
loop _ 0 = False
|
||||
loop x r =
|
||||
let x' = (x * x) `mod` n
|
||||
in x' == n - 1 || loop x' (r - 1)
|
||||
|
||||
millerRabin :: Integer -> Integer -> IO Bool
|
||||
millerRabin n k
|
||||
| n < 2 = return False
|
||||
| n == 2 = return True
|
||||
| even n = return False
|
||||
| otherwise = do
|
||||
let (d, s) = factorOutTwos (n - 1)
|
||||
go k d s
|
||||
where
|
||||
go 0 _ _ = return True
|
||||
go i d s = do
|
||||
a <- randomRIO (2, n - 2)
|
||||
if millerRabinWitness n a d s
|
||||
then go (i - 1) d s
|
||||
else return False
|
||||
28
src/Primes/SoloveyStrassen.hs
Normal file
28
src/Primes/SoloveyStrassen.hs
Normal file
@@ -0,0 +1,28 @@
|
||||
module Primes.SoloveyStrassen (soloveyStrassen) where
|
||||
|
||||
import ModularArithmeticUtils (modExp, jacobi)
|
||||
import System.Random (randomRIO)
|
||||
|
||||
|
||||
soloveyStrassen :: Integer -> Integer -> IO Bool
|
||||
soloveyStrassen n k
|
||||
| n < 2 = return False
|
||||
| n == 2 = return True
|
||||
| even n = return False
|
||||
| otherwise = go k
|
||||
where
|
||||
expHalf = (n - 1) `div` 2
|
||||
|
||||
go 0 = return True
|
||||
go i = do
|
||||
a <- randomRIO (2, n - 2)
|
||||
if gcd a n /= 1
|
||||
then return False
|
||||
else do
|
||||
let x = modExp a expHalf n
|
||||
j = jacobi a n `mod` n
|
||||
|
||||
if x /= j
|
||||
then return False
|
||||
else go (i - 1)
|
||||
|
||||
43
src/TonelliShanks.hs
Normal file
43
src/TonelliShanks.hs
Normal file
@@ -0,0 +1,43 @@
|
||||
module TonelliShanks (tonelliShanks) where
|
||||
|
||||
import ModularArithmeticUtils (modExp, modMul, legendre, factorOutTwos)
|
||||
|
||||
-- Tonelli-Shanks algorithm: find x such that x^2 = n (mod p)
|
||||
-- Returns Just x if it exists, Nothing otherwise.
|
||||
tonelliShanks :: Integer -> Integer -> Maybe Integer
|
||||
tonelliShanks n p
|
||||
| legendre n p /= 1 = Nothing -- no square root
|
||||
| p == 2 = Just n -- special case
|
||||
| otherwise = Just x
|
||||
where
|
||||
(q, s) = factorOutTwos (p - 1)
|
||||
|
||||
-- find z which is a quadratic non-residue mod p
|
||||
z = head [z' | z' <- [2..p-1], legendre z' p == p - 1]
|
||||
|
||||
m0 = s
|
||||
c0 = modExp z q p
|
||||
t0 = modExp n q p
|
||||
r0 = modExp n ((q + 1) `div` 2) p
|
||||
|
||||
(x, _, _) = loop m0 c0 t0 r0
|
||||
|
||||
loop :: Integer -> Integer -> Integer -> Integer -> (Integer, Integer, Integer)
|
||||
loop m c t r
|
||||
| t == 0 = (0, c, t)
|
||||
| t == 1 = (r, c, t)
|
||||
| otherwise =
|
||||
let
|
||||
i = smallestI 0 t m
|
||||
b = modExp c (2 ^ (m - i - 1)) p
|
||||
m' = i
|
||||
c' = modMul b b p
|
||||
t' = modMul t c' p
|
||||
r' = modMul r b p
|
||||
in loop m' c' t' r'
|
||||
|
||||
-- find smallest i (0 <= i < m) such that t^(2^i) = 1 mod p
|
||||
smallestI i t m
|
||||
| i >= m = error "no valid i found"
|
||||
| modExp t (2 ^ i) p == 1 = i
|
||||
| otherwise = smallestI (i + 1) t m
|
||||
15
src/Utils.hs
Normal file
15
src/Utils.hs
Normal file
@@ -0,0 +1,15 @@
|
||||
module Utils (askNumber) where
|
||||
|
||||
import Text.Read (readMaybe)
|
||||
import System.Exit (exitSuccess)
|
||||
|
||||
-- Ask user for an integer > 1, or exit on invalid input
|
||||
askNumber :: String -> IO Integer
|
||||
askNumber s = do
|
||||
putStrLn s
|
||||
input <- getLine
|
||||
case readMaybe input of
|
||||
Just n | n > 1 -> return n
|
||||
_ -> do
|
||||
putStrLn "Not a valid integer"
|
||||
exitSuccess
|
||||
Reference in New Issue
Block a user