Concurrency demos/Zeta

A simple example of parallelism in Haskell

This little piece of code computes an approximation of Riemann's zeta function, balancing the work to be done between N threads.

import Control.Concurrent
import Control.Concurrent.MVar
import Control.Monad
import Data.Complex
import System.Environment

-- Return the list of the terms of the zeta function for the given range.
-- We don't sum the terms here but let the main thread sum the lists returned
-- by all the other threads so as to avoid accumulating rounding imprecisions.
zetaRange :: (Floating a, Integral b) => a -> (b, b) -> [a]
zetaRange s (x,y) = [ (fromIntegral n) ** (-s) | n <- [x..y] ]

cut :: (Integral a) => (a, a) -> a -> [(a, a)]
cut (x,y) n = (x, x + mine - 1) : cut' (x + mine) size (y - mine)
 where
  (size, modulo)   = y `divMod` n
  mine             = size + modulo

  cut' _ _ 0       = []
  cut' x' size' n' = (x', x' + size' - 1) : cut' (x' + size') size' (n' - size') 

getParams :: IO (Int, Int, Complex Double)
getParams = do
  argv <- getArgs
  case argv of
    (t:n:s:[]) -> return (read t, read n, read s)
    _          -> error "usage: zeta <nthreads> <boundary> <s>"

main :: IO ()
main = do
  (t, n, s) <- getParams
  childs    <- mapM (thread s) (cut (1, n) t)
  results   <- mapM takeMVar childs
  print (sum (concat results))
 where
  thread s range = do
    putStrLn ("Starting thread for range " ++ show range)
    mvar <- newEmptyMVar
    forkIO (putMVar mvar (zetaRange s range))
    return mvar

Or using Strategies

Replace the Control.Concurrent... imports by

import Control.Parallel.Strategies

and replace main by

main :: IO ()
main = do
  (t, n, s) <- getParams
  let ranges    = cut (1, n) t
      results   = map (zetaRange s) ranges `using` parList rnf
  putStr $ unlines [ "Starting thread for range " ++ show r | r <- ranges ]
  print (sum (concat results))

Using a Chan instead of MVars

This isn't working yet...

import Control.Concurrent
import Control.Monad
import Data.Complex
import System.Environment

zetaRange :: (Floating a, Integral b) => a -> (b, b) -> [a]
zetaRange s (x,y) = [ (fromIntegral n) ** (-s) | n <- [x..y] ]

cut :: (Integral a) => (a, a) -> a -> [(a, a)]
cut (x,y) n = (x, x + mine - 1) : cut' (x + mine) size (y - mine)
 where
  (size, modulo)   = y `divMod` n
  mine             = size + modulo

  cut' _ _ 0       = []
  cut' x' size' n' = (x', x' + size' - 1) : cut' (x' + size') size' (n' - size')

getParams :: IO (Int, Int, Complex Double)
getParams = do
  argv <- getArgs
  case argv of
    (t:n:s:[]) -> return (read t, read n, read s)
    _          -> error "usage: zeta <nthreads> <boundary> <s>"

main :: IO ()
main = do
  (t, n, s) <- getParams
  chan      <- newChan
  terms     <- getChanContents chan

  forM_ (cut (1,n) t) $ thread chan s

  let wait xs i result
        | i >= n    = print result  -- Done.
        | otherwise = case xs of
           Nothing : rest -> wait rest (i + 1) result
           Just x : rest  -> wait rest i (result + x)

  wait terms 0 0
 where
  thread chan s range = do
    putStrLn ("Starting thread for range " ++ show range)
    forkIO $ do
      writeList2Chan chan (map Just (zetaRange s range))
      writeChan chan Nothing

Benchmarks

Insert benchmarks here! :-)