Hill Climbing in Haskell with Accelerate

climb.hs

-- Parallel Bit Climbing
-- climb.hs
-- Mike Vollmer
--
--   This program uses Accelerate to attempt to find the minimum
--   point of the sinbowl function. It uses parallel bit climbing,
--   a form of local search.
--
--   It takes three command-line parameters:
--
--   ./climb <num> <size> <iters>
--
--   Where <num> is the number of simultaneous solutions to search,
--   <size> is the number of bits in each potential solution, and 
--   <iters> is the number of iterations to search before returning
--   the best answer.
--
--   With reasonable input (like 50 10 20 ,for example) the
--   program should return something around -0.84.

import Prelude                    as P
import Data.Array.Accelerate      as A
import Data.Array.Accelerate.CUDA as R
import System.Random
import System.Environment   
import Control.Monad (replicateM)

-- I got tired of typing these so I made shortcuts
type TwoD     = Array DIM2 Int
type OneD     = Array DIM1 Int
type OneDF    = Array DIM1 Float
type EvalFunc = Acc OneD -> Acc OneDF

-- Convert arrays of 1s and 0s into integers.
translateNums :: Acc TwoD -> Acc OneD
translateNums = process
  where process     = powers >-> sum
        f strs ix   = let i = indexHead ix
                          n = strs ! ix
                      in shiftL n i
        powers strs = generate (shape strs) (f strs)
        sum         = A.fold1 (+)

-- Generate a new list of neighboring numbers
makeNeighbors :: Acc (Scalar Int) -> Acc TwoD -> Acc TwoD
makeNeighbors place strs = flip
  where place' = the place
        len    = indexHead (shape strs)
        f ix   = let Z :. j :. i = unlift ix
                     modEq       = (i `mod` len) ==* place'
                     value       = strs ! index2 j i
                     result      = complementBit value 0
                 in modEq ? (result, value)
        flip   = generate (shape strs) f                   

-- Pick new neighbors that are better (minimizing eval score).
compareNeighbors :: EvalFunc -> Acc TwoD -> Acc TwoD -> Acc TwoD
compareNeighbors evalFunc str1 str2 = finalStr
  where process  = translateNums >-> evalFunc
        str1e    = process str1
        str2e    = process str2
        f ix     = let Z :. j :. i = unlift ix
                       e1          = str1e ! index1 j
                       e2          = str2e ! index1 j
                       v1          = str1 ! index2 j i
                       v2          = str2 ! index2 j i
                   in e1 <* e2 ? (v1, v2)
        finalStr = generate (shape str1) f                     

-- Pass numbers to sinbowl function.
evaluateNum :: Acc (Scalar Float) -> Acc (Scalar Float) ->
               Acc (Scalar Float) -> Acc (Scalar Float) ->
               Acc OneD -> Acc OneDF
evaluateNum min max a b = process
  where process = A.map scale >-> A.map f
        min'    = the min
        max'    = the max
        a'      = the a
        b'      = the b
        scale x = let x' = A.fromIntegral x
                  in (((b' - a') * (x' - min')) / (max' - min')) + a'
        f x     = abs x * 0.1 - sin x

-- Random numbers are awkward in Haskell
makeRandomString :: Int -> IO [Int]
makeRandomString size = do
  g <- newStdGen
  return $ P.take size (randomRs (0,1) g :: [Int])

-- Use makeRandomString to generate all the random strings.
makeAllStrings :: Int -> Int -> IO [[Int]]
makeAllStrings num size = replicateM num $ makeRandomString size

-- Given a list of random numbers, an evaluation function, and an
-- initial population, run the search using foldr. The number of
-- iterations will be the length of rn.
runSearch :: [Int] -> EvalFunc -> Acc TwoD -> Acc TwoD
runSearch rn f strs = foldr compute strs rn
  where compute r strs = let flipBit   = unit $ constant r
                             neighbors = makeNeighbors flipBit strs
                         in compareNeighbors f strs neighbors

-- This function launches the search.
startSearch :: Int -> Int -> Int -> IO OneDF
startSearch num size iters = do
  g <- newStdGen
  strs <- makeAllStrings num size
  let arrs     = fromList (Z:.num:.size) (concat strs)
      useArrs  = use arrs
      passNum  = unit . constant
      evalFunc = evaluateNum -- pass number parameters as unit arrays
                 (passNum 0) (passNum (2^size - 1))
                 (passNum (-20)) (passNum 20) -- scale number to [-20,20]
      rns      = P.take iters $ randomRs (0, size-1) g
      result   = runSearch rns evalFunc useArrs
  -- compute the f(x) values for each result and return them
  return $ run $ evalFunc $ translateNums result

main :: IO ()
main = do
  args <- getArgs
  let [a,b,c] = args
  values <- startSearch (read a :: Int) (read b :: Int) (read c :: Int)
  print $ P.minimum $ A.toList values