module Main  where{
import Control.Exception;
import Ratio;
import IO;
import System.IO.Unsafe;
import List;
import Monad;
import Array;
import System;

show_list::(Show a)=>[](a)->String;
show_list l = (unlines (map show l));
quiet::Bool;
quiet = True;
cerr::[](String)->a->a;
cerr message x = (case quiet of {
  (True)->x;
  (_)->(unsafePerformIO (do{
  (hPutStrLn stderr (concat message));
  (return x);
}))
});
cerr_x::[](String)->a->a;
cerr_x _ x = x;
peek::(Show a)=>a->a;
peek x = (seq x (cerr [(show x)] x));
peek_more::(a->[](String))->a->a;
peek_more f x = (case quiet of {
  (True)->x;
  (_)->(seq x (cerr (f x) x))
});
peek_more_x::(a->[](String))->a->a;
peek_more_x _ x = x;
zip_map::(a->b)->[](a)->[]((a,b));
zip_map f l = (zip l (map f l));
show_and::(Show a,Show b)=>(a->b)->[](a)->IO(());
show_and f l = (putStr(show_list((zip_map f l))));
ntakes::Int->[](a)->[]([](a));
ntakes n l = ((:) (take n l) (ntakes n (tail l)));
zero_to::Int->[](Int);
zero_to n = (take n (enumFrom 0));
enum_from_by::(Num a)=>a->a->[](a);
enum_from_by x delta = ((:) x (enum_from_by ((+) x delta) delta));
rcs_code::String;
rcs_code = "$Id: primes-2.ll,v 1.312 2004/10/13 07:53:14 kenta Exp kenta $";
n_log_n::Int->Int;
n_log_n n = (floor ((*) (fromIntegral n) (log (fromIntegral n))));
is_prime::Integer->Bool;
is_prime x = (case (compare x 2) of {
  (LT)->False;
  (EQ)->True;
  (GT)->(case (mod x 2) of {
  (0)->False;
  (_)->(miller_rabin x miller_rabin_bases)
})
});
miller_rabin_bases::[](Integer);
miller_rabin_bases = (take 20 primes);
miller_rabin::Integer->[](Integer)->Bool;
miller_rabin n bases = (and (map (miller_rabin_1 n) bases));
miller_rabin_1::Integer->Integer->Bool;
miller_rabin_1 n a = ((assert ((<) 2 n))((assert ((==) 1 (mod n 2)))((case (compare a ((-) n 1)) of {
  (LT)->(case (pow_2_special a ((-) n 1) n) of {
  (Probable_prime)->True;
  (Composite)->False;
  (N(1))->True;
  (_)->False
});
  (_)->True
}))));
pow_2_special::Integer->Integer->Integer->Miller_rabin_result;
pow_2_special a n modulo_p = (let {
end_show::Miller_rabin_result->[](String);
end_show x = ["pow-2-special ",(show a)," ",(show n)," ",(show x)]
}
 in (peek_more end_show (case (mod n 2) of {
  (0)->(let {
result::Miller_rabin_result;
result = (modular_square_miller_rabin (pow_2_special a (div n 2) modulo_p) modulo_p)
}
 in (case result of {
  (N(b))->(case (classify b modulo_p) of {
  (Negative_unity)->Probable_prime;
  (Unity)->(cerr ["nontrivial sqrt of unity ",(show (div n 2))] Composite);
  (Number)->result
});
  (_)->result
}));
  (_)->(let {
result::Integer;
result = (modular_exponentiation a n modulo_p)
}
 in (case (classify result modulo_p) of {
  (Negative_unity)->Probable_prime;
  (Unity)->Probable_prime;
  (Number)->(N result)
}))
})));
data Miller_rabin_result = Probable_prime | Composite | N(Integer) deriving (Show);
data Modular_classify = Negative_unity | Unity | Number;
classify::Integer->Integer->Modular_classify;
classify n modulo_p = (case ((==) n ((-) modulo_p 1)) of {
  (True)->Negative_unity;
  (_)->(case ((==) n 1) of {
  (True)->Unity;
  (_)->Number
})
});
modular_exponentiation::Integer->Integer->Integer->Integer;
modular_exponentiation a n modulo_p = (let {
end_show::Integer->[](String);
end_show x = ["modular-exponentiation ",(show a)," ",(show n)," = ",(show x)]
}
 in (peek_more end_show (case (compare n 1) of {
  (EQ)->a;
  (GT)->(case (mod n 2) of {
  (0)->(let {
result::Integer;
result = (modular_exponentiation a (div n 2) modulo_p)
}
 in (modular_square result modulo_p));
  (_)->(mod ((*) a (modular_exponentiation a ((-) n 1) modulo_p)) modulo_p)
})
})));
modular_square::Integer->Integer->Integer;
modular_square result modulo_p = (mod ((*) result result) modulo_p);
modular_square_miller_rabin::Miller_rabin_result->Integer->Miller_rabin_result;
modular_square_miller_rabin x modulo_p = (case x of {
  (N(n))->(N (modular_square n modulo_p));
  (_)->x
});
primes::[](Integer);
primes = (2:eratosthenes_sieve);
eratosthenes_sieve::[](Integer);
eratosthenes_sieve = (filter by_eratosthenes (enumFrom 3));
by_eratosthenes::Integer->Bool;
by_eratosthenes x = (and (map (indivisible x) (takeWhile (lt_sqrt x) primes)));
lt_sqrt::Integer->Integer->Bool;
lt_sqrt x p = ((<=) ((*) p p) x);
indivisible::Integer->Integer->Bool;
indivisible big small = ((/=) 0 (mod big small));
exhaustive_miller_rabin::Integer->Integer;
exhaustive_miller_rabin x = (genericLength (filter id (map (miller_rabin_1 x) (enumFromTo 2 ((-) x 2)))));
is_interesting::(Integer,Integer)->Bool;
is_interesting exhaustive_miller_rabin_result = (case exhaustive_miller_rabin_result of {
  ((a),(b))->((&&) ((/=) b 1) ((/=) a (b+2 )))
});
ip_fast::Integer->Bool;
ip_fast x = (miller_rabin_1 x 2);
good::Integer->Bool;
good x = (and [(ip_fast ((+) 1 ((*) 2 x))),(ip_fast ((+) 1 ((*) 4 x)))]);
make_good::Integer->Integer;
make_good x = ((*) ((+) 1 ((*) 2 x)) ((+) 1 ((*) 4 x)));
high_witnesses::[](Integer)->IO(());
high_witnesses l = (putStr (show_list (filter is_interesting (zip_map exhaustive_miller_rabin l))));
double_is_prime::(Integer,Integer)->Bool;
double_is_prime x = ((&&) (miller_rabin_1 (fst x) 2) (miller_rabin_1 (snd x) 2));
diff_sub::[](Integer)->[](Integer)->[](Integer);
diff_sub many some = (case (many,some) of {
  (((:)(h_many) (t_many)),((:)(h_some) (t_some)))->(case (compare h_many h_some) of {
  (LT)->((:) h_many (diff_sub t_many some));
  (EQ)->(diff_sub t_many t_some);
  (GT)->(diff_sub many t_some)
})
});
in_desired_range::Integer->Integer->Bool;
in_desired_range p num_wit = ((&&) ((<) ((*) 2 num_wit) p) ((<) p ((*) 5 num_wit)));
mersenne::Integer->Integer;
mersenne x = ((-) ((^) 2 x) 1);
sub_fib::[](Integer);
sub_fib = ((:) 1 ((:) 1 ((:) 1 (zipWith (+) sub_fib (tail (tail sub_fib))))));
sub_fib_numbered::[]((Integer,Int));
sub_fib_numbered = (zip sub_fib (enumFrom 0));
two_powers::Int->Int->IO(());
two_powers start next = (mapM_ single_two_power (enumFromThen start next));
two_powers_sieve::Int->Int->IO(());
two_powers_sieve start next = (mapM_ sievers_two_power (enumFromThen start next));
run_sequence::(Int->IO(()))->Int->Int->IO(());
run_sequence fun start by = (run_sequence_on fun (enumFromThen start ((+) start by)));
run_sequence_on::(Int->IO(()))->[](Int)->IO(());
run_sequence_on fun sequence = (mapM_ fun sequence);
logn_sequence_from_by::Int->Int->[](Int);
logn_sequence_from_by start by = ((dropWhile ((>) start))((map n_log_n (enumFromThen 1 ((+) 1 by)))));
is_prime_offset::Integer->Integer->Bool;
is_prime_offset base offset = (is_prime ((+) base offset));
prime_filter_offset::Integer->[](Integer)->[](Integer);
prime_filter_offset base l = (filter (is_prime_offset base) l);
single_two_power::Int->IO(());
single_two_power s = (let {
a::[](Integer);
a = (take 10 (prime_filter_offset (two_pow s) (enumFrom 0)));
b::[](Integer);
b = (take 10 (prime_filter_offset (two_pow s) (map negate (enumFrom 1))))
}
 in (do{
  (putStr (show s));
  (putStr " ");
  (putStrLn((show )(((++) (reverse b) a))));
}));
num_neighbors::Int;
num_neighbors = 3;
sievers_two_power::Int->IO(());
sievers_two_power exponent = (let {
base_point::Integer;
base_point = ((^) 2 (fromIntegral exponent));
a::[](Integer);
a = ((take num_neighbors)((prime_filter_offset base_point)((map fst)((filter snd)((sieve_forward base_point))))));
b::[](Integer);
b = ((take num_neighbors)((prime_filter_offset base_point)((map fst)((filter snd)((sieve_backward base_point))))))
}
 in (do{
  (putStr (show exponent));
  (putStr " ");
  (putStrLn((show )(((++) (reverse b) a))));
}));
test_print::Int->IO(());
test_print exponent = (let {
base_point::Integer;
base_point = ((^) 2 (fromIntegral exponent));
a::[](Integer);
a = ((twin_process base_point)((sieve_forward base_point)))
}
 in (do{
  (putStr (show exponent));
  (putStr " ");
  (putStrLn (show a));
}));
powers_of_two::[](Int);
powers_of_two = ((:) 4 (map ((*) 2) powers_of_two));
sophie_germain_sieve::(Integer->[]((Integer,Bool)))->Int->IO(());
sophie_germain_sieve sieve_direction exponent = (let {
base_point::Integer;
base_point = (two_pow exponent)
}
 in (do{
  (putStr "2^");
  (putStr (show exponent));
  (putStr " + ");
  (putStrLn(show((map fst)((take 1)((filter snd)((zip_map (sg_delta base_point))((map fst)((filter snd)((sieve_direction base_point))))))))));
}));
two_pow::Int->Integer;
two_pow exponent = ((^) 2 (fromIntegral exponent));
sophie_germain::[](Integer)->Int->IO(());
sophie_germain direction exponent = (do{
  (putStr "2^");
  (putStr (show exponent));
  (putStr " + ");
  (putStrLn(show((map fst)((take 2)((filter snd)((zip_map (sg_delta (two_pow exponent)))(direction)))))));
});
sg_delta::Integer->Integer->Bool;
sg_delta base_point delta = (is_sophie_germaine_prime_minor ((+) delta base_point));
twin_primes::Int->IO(());
twin_primes exponent = (let {
base_point::Integer;
base_point = ((^) 2 (fromIntegral exponent));
a::[](Integer);
a = ((twin_process base_point)((sieve_forward base_point)));
b::[](Integer);
b = ((twin_process base_point)((sieve_backward base_point)))
}
 in (do{
  (putStr (show exponent));
  (putStr " ");
  (putStrLn(show(((++) (reverse b) a))));
}));
twin_process_old::Integer->[]((Integer,Bool))->[](Integer);
twin_process_old base_point candidates = ((take 10)((map ((flip (!!)) 1))((filter (is_really_twin_prime base_point))((map (map fst))((filter twin_positive)((ntakes 3)(candidates)))))));
twin_process_2::Integer->[]((Int,Bool))->[](Int);
twin_process_2 base_point candidates = ((take 10)((map fst)((filter snd)((map (twin_real_test_2 base_point))(tails(candidates))))));
twin_process::Integer->[]((Integer,Bool))->[](Integer);
twin_process base_point candidates = ((take 10)((map fst)((filter snd)((map (twin_real_test base_point))(tails(candidates))))));
twin_real_test::Integer->[]((Integer,Bool))->(Integer,Bool);
twin_real_test base_point candidates_tail = (case candidates_tail of {
  ((:)((a),(True)) ((:)((middle),(_)) ((:)((b),(True)) (_))))->(middle,((&&) (is_prime_offset base_point a) (is_prime_offset base_point b)));
  (_)->(0,False)
});
twin_real_test_2::Integer->[]((Int,Bool))->(Int,Bool);
twin_real_test_2 base_point candidates_tail = (case candidates_tail of {
  ((:)((a),(True)) ((:)((middle),(_)) ((:)((b),(True)) (_))))->(middle,((&&) ((<) 10000000 a) ((<) 10000000 b)));
  (_)->(0,False)
});
sieve_backward::Integer->[]((Integer,Bool));
sieve_backward base_point = ((zip (map negate (enumFrom 1)))((map and)((many_sievers_reverse base_point num_primes_to_sieve))));
twin_positive::[]((a,Bool))->Bool;
twin_positive x = (case x of {
  [((_),(True)),(_),((_),(True))]->True;
  (_)->False
});
sieve_forward::Integer->[]((Integer,Bool));
sieve_forward base_point = ((zip (enumFrom 0))((map and)((many_sievers base_point num_primes_to_sieve))));
num_primes_to_sieve::Int;
num_primes_to_sieve = 5;
is_really_twin_prime::Integer->[](Integer)->Bool;
is_really_twin_prime base l = ((&&) (is_prime_offset base ((!!) l 0)) (is_prime_offset base ((!!) l 2)));
siever::Integer->Int->[](Bool);
siever start modulo = ((take modulo)((map (((/=) 0) . ((flip mod) modulo)))((enumFrom (fromInteger (mod start (fromIntegral modulo)))))));
siever_by::Integer->Int->Int->[](Bool);
siever_by start by_count modulo = ((take modulo)((map (((/=0) ) . ((flip mod) modulo)) (enum_from_by (fromInteger (mod start (fromIntegral modulo))) by_count))));
many_sievers::Integer->Int->[]([](Bool));
many_sievers start num = (let {
first_num_primes::[](Int);
first_num_primes = (map fromInteger (take num primes))
}
 in (transpose((map (cycle . (siever start)))(first_num_primes))));
many_sievers_reverse::Integer->Int->[]([](Bool));
many_sievers_reverse start num = (let {
first_num_primes::[](Int);
first_num_primes = (map fromInteger (take num primes))
}
 in (transpose((map (cycle . reverse . (siever start)))(first_num_primes))));
is_sophie_germaine_prime_minor::Integer->Bool;
is_sophie_germaine_prime_minor p = ((&&) (is_prime (div ((-) p 1) 2)) (is_prime p));
is_sophie_germaine_prime::Integer->Bool;
is_sophie_germaine_prime p = ((&&) (is_prime p) (is_prime ((+) 1 ((*) 2 p))));
negatives::[](Integer);
negatives = (enumFromThen (negate 1) (negate 2));
positives::[](Integer);
positives = (enumFrom 1);
halfways::[](Int)->[](Int);
halfways l = (case l of {
  ((:)(a) rest@((:)(b) (_)))->((:) a ((:) (div ((+) a b) 2) (halfways rest)))
});
hh_powers_of_two::[](Int);
hh_powers_of_two = (halfways (halfways powers_of_two));
sophie_germain_caveat::IO(());
sophie_germain_caveat = (do{
  (putStrLn "Here are some primes p such that p and (p-1)/2 are both prime. The corresponding");
  (putStrLn "Sophie Germain prime is (p-1)/2.");
});
hh_powers_of_two_from::Int->[](Int);
hh_powers_of_two_from minimum = (dropWhile ((>=) minimum) hh_powers_of_two);
smaller_factors::Int->[](Int);
smaller_factors n = (let {
sqrt_smaller::Int->Bool;
sqrt_smaller p = ((<=) ((*) p p) n)
}
 in (takeWhile sqrt_smaller (map fromInteger primes)));
factor n = (let {
sndpositive x = ((<) 0 (snd x))
}
 in ((filter sndpositive)((factor_l (smaller_factors n) n))));
factor_l list_factors victim = (case victim of {
  (1)->[];
  (_)->(case list_factors of {
  ([])->[(victim,1)];
  (_)->(let {
p = (divide_through victim (head list_factors))
}
 in ((:) ((head list_factors),(fst p)) (factor_l (tail list_factors) (snd p))))
})
});
divide_through victim primefactor = (case (mod victim primefactor) of {
  (0)->(let {
final = (divide_through (div victim primefactor) primefactor)
}
 in (((+) 1 (fst final)),(snd final)));
  (_)->(0,victim)
});
big_product_prime::Integer;
big_product_prime = (product (take 1000 primes));
gcd_test::Integer->Integer;
gcd_test n = (gcd n big_product_prime);
mersenne_gcd_test::Integer->Integer;
mersenne_gcd_test n = (gcd_test (mod (pred ((+) big_product_prime (modular_exponentiation 2 n big_product_prime))) big_product_prime));
main::IO(());
main = (do{
  (hSetBuffering stdout LineBuffering);
  (hPutStrLn stdout rcs_code);
  args<-getArgs;
  (hPutStrLn stdout (show args));
  (putStrLn ("% num-primes-to-sieve " ++ (show num_primes_to_sieve)));
  (case args of {
  ["primes",(num)]->(putStr(show_list((take (read num) primes))));
  ["sub-fib"]->(putStr(show_list((filter (is_prime . fst))(sub_fib_numbered))));
  ["sub-fib-squares"]->(error "no");
  ["two-powers",(start),(next)]->(two_powers (read start) (read next));
  ["two-powers-sieve",(start),(next)]->(two_powers_sieve (read start) (read next));
  ["single-two-power",(x)]->(single_two_power (read x));
  ["many-sievers",(x)]->(putStr(show_list((take 10)((filter (is_prime . ((+) ((^) 2 1000)) . fst))((filter snd)((zip (enumFrom 0))((map and)((many_sievers ((^) 2 1000) (read x))))))))));
  ["sievers-two-power",(s)]->(sievers_two_power (read s));
  ["twin-primes",(start),(by)]->(run_sequence twin_primes (read start) (read by));
  ["twin-primes-l",(start),(by)]->(run_sequence_on twin_primes (logn_sequence_from_by (read start) (read by)));
  ["sophie-germain",(start),(by)]->(run_sequence (sophie_germain negatives) (read start) (read by));
  ["sophie-germain-powers-of-two-minus",(minimum)]->(do{
  sophie_germain_caveat;
  (mapM_ (sophie_germain_sieve sieve_backward) (hh_powers_of_two_from (read minimum)));
});
  ["test-sophie",(value)]->(sophie_germain_sieve sieve_backward (read value));
  ["sophie-germain-powers-of-two-plus",(minimum)]->(do{
  sophie_germain_caveat;
  (mapM_ (sophie_germain_sieve sieve_forward) (hh_powers_of_two_from (read minimum)));
});
  ["primes-hh",(minimum)]->(do{
  (run_sequence_on sievers_two_power (hh_powers_of_two_from (read minimum)));
});
  ["test",(x)]->(twin_primes (read x));
  ["test-is-prime",(base),(exponent),(delta)]->(putStrLn(show(is_prime(((+) ((^) (read base) (read exponent)) (read delta))))));
  ["test-print",(start),(by)]->(run_sequence test_print (read start) (read by));
  ["many-sievers-test"]->(print(length((filter id)((map and)((many_sievers 0 num_primes_to_sieve))))));
  ["time-modexp",(bits)]->(let {
mers::Integer;
mers = ((-) ((^) 2 (read bits)) 1)
}
 in (putStrLn(show((modular_exponentiation 2 ((-) mers 1) mers)))));
  ["sierp",(test_val),(k),(n)]->(let {
p::Integer;
p = ((+) 1 ((*) m v));
m::Integer;
m = (read k);
v::Integer;
v = ((^) 2 (read n))
}
 in (do{
  (putStrLn(show(((==) 1)((modular_exponentiation (read test_val) v p)))));
  (putStrLn(show(((==) 1)((modular_exponentiation (read test_val) (div ((-) p 1) 2) p)))));
  (putStrLn "==");
  (putStrLn(show((miller_rabin_1 p 2))));
}))
});
})
}