{-# LANGUAGE TemplateHaskell #-}
module Main where
import Data.Bits
import Data.IORef
import Data.Word
$( let wbytes' = sizeOf (undefined :: Word)
wbits' = 8*wbytes'
wbits1' = wbits' - 1
wlog' = popCount wbits1'
in [d| wbytes = wbytes'
wbits = wbits'
wbits1 = wbits1'
wlog = wlog'
|])
{-# INLINE soeIO #-}
soeIO :: Int -> (Int -> IO ()) -> IO ()
soeIO n usePrime =
allocaArray (maxW + 1) $ \arr -> do
let readB ix = do
w <- peekElemOff arr (shiftR ix wlog)
return (w .&. bit (ix .&. wbits1) /= 0)
clearB ix = do
let y = shiftR ix wlog
w <- peekElemOff arr y
pokeElemOff arr y (w .&. complement (bit (ix .&. wbits1)))
forM_ [0..maxW] (flip (pokeElemOff arr) (-1 :: Word))
forM_ [2..maxN] $ \i -> do
b <- readB i
when b $ do
usePrime i
let j0 = i*i
j0 `seq` mapM_ clearB [j0, j0 + i .. maxN]
where
maxN = n - 1
maxW = shiftR (n + wbits1) wlog - 1
main :: IO ()
main = do
var <- newIORef (0 :: Int)
soeIO 50000 (\_ -> modifyIORef' var (+ 1))
ey0jIExBTkdVQUdFIFRlbXBsYXRlSGFza2VsbCAjLX0KCm1vZHVsZSBNYWluIHdoZXJlCgppbXBvcnQgQ29udHJvbC5Nb25hZAppbXBvcnQgRGF0YS5CaXRzCmltcG9ydCBEYXRhLklPUmVmCmltcG9ydCBEYXRhLldvcmQKaW1wb3J0IEZvcmVpZ24uTWFyc2hhbAppbXBvcnQgRm9yZWlnbi5TdG9yYWJsZQoKCiQoIGxldCB3Ynl0ZXMnID0gc2l6ZU9mICh1bmRlZmluZWQgOjogV29yZCkKICAgICAgIHdiaXRzJyAgPSA4KndieXRlcycKICAgICAgIHdiaXRzMScgPSB3Yml0cycgLSAxCiAgICAgICB3bG9nJyAgID0gcG9wQ291bnQgd2JpdHMxJwoKICAgaW4gW2R8IHdieXRlcyA9IHdieXRlcycKICAgICAgICAgIHdiaXRzICA9IHdiaXRzJwogICAgICAgICAgd2JpdHMxID0gd2JpdHMxJwogICAgICAgICAgd2xvZyAgID0gd2xvZycKICAgICAgIHxdKQoKCnstIyBJTkxJTkUgc29lSU8gIy19Cgpzb2VJTyA6OiBJbnQgLT4gKEludCAtPiBJTyAoKSkgLT4gSU8gKCkKc29lSU8gbiB1c2VQcmltZSA9CiAgICBhbGxvY2FBcnJheSAobWF4VyArIDEpICQgXGFyciAtPiBkbwogICAgICAgIGxldCByZWFkQiBpeCA9IGRvCiAgICAgICAgICAgICAgICB3IDwtIHBlZWtFbGVtT2ZmIGFyciAoc2hpZnRSIGl4IHdsb2cpCiAgICAgICAgICAgICAgICByZXR1cm4gKHcgLiYuIGJpdCAoaXggLiYuIHdiaXRzMSkgLz0gMCkKCiAgICAgICAgICAgIGNsZWFyQiBpeCA9IGRvCiAgICAgICAgICAgICAgICBsZXQgeSA9IHNoaWZ0UiBpeCB3bG9nCiAgICAgICAgICAgICAgICB3IDwtIHBlZWtFbGVtT2ZmIGFyciB5CiAgICAgICAgICAgICAgICBwb2tlRWxlbU9mZiBhcnIgeSAodyAuJi4gY29tcGxlbWVudCAoYml0IChpeCAuJi4gd2JpdHMxKSkpCgogICAgICAgIGZvck1fIFswLi5tYXhXXSAoZmxpcCAocG9rZUVsZW1PZmYgYXJyKSAoLTEgOjogV29yZCkpCgogICAgICAgIGZvck1fIFsyLi5tYXhOXSAkIFxpIC0+IGRvCiAgICAgICAgICAgIGIgPC0gcmVhZEIgaQogICAgICAgICAgICB3aGVuIGIgJCBkbwogICAgICAgICAgICAgICAgdXNlUHJpbWUgaQogICAgICAgICAgICAgICAgbGV0IGowID0gaSppCiAgICAgICAgICAgICAgICBqMCBgc2VxYCBtYXBNXyBjbGVhckIgW2owLCBqMCArIGkgLi4gbWF4Tl0KCiAgICB3aGVyZQogICAgbWF4TiA9IG4gLSAxCiAgICBtYXhXID0gc2hpZnRSIChuICsgd2JpdHMxKSB3bG9nIC0gMQoKCm1haW4gOjogSU8gKCkKbWFpbiA9IGRvCiAgICB2YXIgPC0gbmV3SU9SZWYgKDAgOjogSW50KQogICAgc29lSU8gNTAwMDAgKFxfIC0+IG1vZGlmeUlPUmVmJyB2YXIgKCsgMSkpCiAgICByZWFkSU9SZWYgdmFyID4+PSBwcmludA==