My Weblog

SPOJ 6560. N DIV PHI_N (Hard)

SPOJ 6560. N DIV PHI_N (Hard) also accepted using the same idea accept using a bit memoization. Precalculation of product of primes in prime table $<=$ 10^25000 and linear search is enough to get accepted.Its first 6000 primes whose products are sufficient . Accepted Haskell code. I also tried a slightly improved Rabin Miller primality testing with idea give topcoder forum for faster modular exponentiation for large numbers and code with isProbable function got accepted with 3.95 and time limit exceed other time.

```import Data.List
import Data.Maybe
import Data.Bits
import qualified Data.ByteString.Char8 as BS

powM::Integer->Integer->Integer->Integer
powM a d n
| d == 0 = 1
| d == 1 = mod a n
| otherwise = mod q n  where
p = powM  ( mod ( a^2 ) n ) ( shiftR d 1 ) n
q = if (.&.) d 1 == 1 then mod ( a * p ) n else p

powL :: Integer -> Integer -> Integer -> Integer
powL a d n
| d <= 10^10 = powM a d n
| otherwise =  mod ( powM ( powL a ( div d 1000000 ) n ) 1000000 n   *  powM a ( mod d 1000000 ) n ) n

isProbable::Integer-> Bool
isProbable n |n<=1 = False
|n==2 = True
|even n = False
|otherwise	= rabinMiller [2,3] n s d where
s = until ( \x -> testBit ( n - 1) ( fromIntegral x ) )  ( +1 ) 0
d = div ( n - 1 ) (  shiftL 1 ( fromIntegral s )  )
--d=until (\x->mod x 2/=0) (`div` 2) (n-1) --(n-1=2^s*d)
--s=until (\x->d*2^x==pred n) (+1) 0  --counts the power

rabinMiller::[Integer]->Integer->Integer->Integer-> Bool
rabinMiller [] _ _ _ = True
rabinMiller (a:xs) n s d =
if a == n then True
else  case x==1 of
True -> rabinMiller xs n s d
_ -> if any ( == pred n) . take ( fromIntegral s ) . iterate ( \e -> mod (e^2) n ) \$  x then rabinMiller xs n s d
else  False
where x=powL a d n

prime :: [Integer]
prime = 2:3:filter isPrime [2*k+1 | k<-[2..]]

isPrime :: Integer -> Bool
isPrime n = all ((/=0). mod n  ).takeWhile ( <= ( truncate.sqrt.fromIntegral \$ n) ) \$ prime

anstable ::[Integer]
anstable =  scanl1 (*)   prime

helpSolve :: [Integer] -> Integer -> Integer -> Integer
helpSolve (x:xs) n ret
| x > n = ret
| otherwise = helpSolve xs n x

solve :: Integer -> BS.ByteString
solve 1 = BS.pack "1"
solve n = BS.pack.show \$ a  where
a = helpSolve ( tail anstable ) n ( head anstable )

main = BS.interact \$ BS.unlines.map ( solve .readInt ) .  tail . BS.lines
```