マルチスレッドでマルチコア

Haskellではマルチスレッドにするだけでマルチコアで動く、との情報を得た。
さっそくやってみる。
お題は平方剰余の解を求めるプログラム。

-- module Repunit06 where

import Control.Concurrent
import Data.Time


r19 = 1111111111111111111
ei = r19`div`8
a = 10000000000

lucas :: MVar Int -> Integer -> Integer -> Integer -> Integer -> IO ()
lucas ref _ x y _ | x == y = do{ threadDelay 1000000; print 0; putMVar ref 1; }
lucas ref a x y g = do
  if (x`mod`100000000==0)
    then getZonedTime >>= \k -> print k
    else return()
  if (z^2)`mod`r19 == 10+a
    then do{ print z; putMVar ref (-1); }
    else lucas ref a (x+1) y g'
      where
        z = g
        g' = ((69069*g +1)`mod`(2^64))`mod`(y-x+1)+x

main :: IO ()
main = do
  t1 <- getZonedTime
  print t1
  ref <- newMVar 0
  s1 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092553
  s2 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092554
  s3 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092555
  s4 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092556
  s5 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092557
  s6 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092558
  s7 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092559
  s8 <- forkIO $ lucas ref a 1054092553 (r19-1) 1054092510
{--
  s1 <- forkIO $ lucas ref a 0 10000
  s2 <- forkIO $ lucas ref a 10001 20000
  s3 <- forkIO $ lucas ref a 20001 30000
  s4 <- forkIO $ lucas ref a 30001 40000
  s5 <- forkIO $ lucas ref a 40001 50000
  s6 <- forkIO $ lucas ref a 50001 60000
  s7 <- forkIO $ lucas ref a 60001 70000
  s8 <- forkIO $ lucas ref a 70001 80000
--}
  go ref t1 s1 s2 s3 s4 s5 s6 s7 s8 0
  where
    go ref t1 s1 s2 s3 s4 s5 s6 s7 s8 eflg = do
      tf <- takeMVar ref
      case tf of
        -1 -> do
          killThread s1
          killThread s2
          killThread s3
          killThread s4
          killThread s5
          killThread s6
          killThread s7
          killThread s8
          t2 <- getCurrentTime
          print $ diffUTCTime t2 $ zonedTimeToUTC t1
          putStrLn "found."
        1 -> do
          if eflg==7 then do
            t2 <- getCurrentTime
            print $ diffUTCTime t2 $ zonedTimeToUTC t1
            putStrLn "not found."
          else go ref t1 s1 s2 s3 s4 s5 s6 s7 s8 (eflg+1)
        _ -> go ref t1 s1 s2 s3 s4 s5 s6 s7 s8 eflg

実行してみる。

C:\me\Haskell\repunit\03>ghc -threaded --make Repunit06.hs -O

C:\me\Haskell\repunit\03>Repunit06.exe +RTS -N
2016-04-29 14:13:30.5929137 東京 (標準時)
2016-04-29 14:14:10.4418537 東京 (標準時)
2016-04-29 14:14:10.4428544 東京 (標準時)
2016-04-29 14:14:10.45086 東京 (標準時)
2016-04-29 14:14:10.4558639 東京 (標準時)
2016-04-29 14:14:10.4598664 東京 (標準時)
2016-04-29 14:14:10.4738763 東京 (標準時)
2016-04-29 14:14:10.4858848 東京 (標準時)
2016-04-29 14:14:10.5109025 東京 (標準時)

タスクマネージャーを見る。

うおお、マルチスレッドにするだけでCPU使用率100%になったぞー。

• • • • -

罠が1つあって、それは乱数です。
スレッド内で乱数を使用すると、CPU使用率が30%ぐらいまで落ちます。
それで線形合同法もどきを使っている訳ですが。

• • • • -

オチとしては、WolframAlphaで

PowerMod[10000000010, 1/2, 1111111111111111111]

とか打って、 decimalを押すと、一瞬で平方剰余の解を計算してくれることでした。

• • • • -

追記
こちらに平方剰余の解の求め方がありました。
Tech Tips: 平方剰余の根