[Shootout-list] Improved GHC Mandelbrot

[Alson Kemp]

  Some small mods to the GHC Mandelbrot code sped it up quite a bit.
Changes
   * Added signatures
   * handcoded 'magnitude z' because magnitude is slow
   * calculate z' before evaluating |z|>limit to match C version

  'fractal 400' now runs in about 1.2s instead of 40s

	- Alson

--- CUT HERE ---
-- Fractal.hs mandlebrot example
-- dumps a *.ppm graphics file to stdout
-- compile:  ghc -O2 -o fractal fractal.hs
-- run: fractal 600 >mandel.pbm
-- by Greg Buchholz
-- modified by Alson Kemp

import Complex
import System(getArgs)
import Char(chr)
import System.IO

limit :: Double
limit  = 4.0

iter :: Int
iter   = 50

main :: IO ()
main = do   [arg] <- getArgs
            let width = read arg
            --AK:optional;prevent newline mangle on PC
            hSetBinaryMode stdout True	 
            putStr $ "P4\n" ++ (show width) ++ " " ++ (show width) ++ "\n"
            putStr $ concatMap (makePBM 0 0) (chunk width (fractal (points
width width)))

points :: Int -> Int -> [Complex Double]
points width height = [(2.0*x/w - 1.5) :+ (2.0*y/h - 1) | y<-[0..h-1],
x<-[0..w-1]]
                    where w = fromIntegral width
                          h = fromIntegral height

fractal :: [Complex Double] -> [Int]
fractal pts = map (fractal' (0.0 :+ 0.0) iter) pts

--magnitude is sloooooowwwwww, so hand code abs^2
fractal' :: Complex Double -> Int -> Complex Double -> Int
fractal' z i c | (((realPart z')*(realPart z') + (imagPart z')*(imagPart
z')) > limit ) = 0
               | (i == 0) = iter
               | otherwise = fractal' z' (i-1) c
			     where z' = z*z+c
 
makePBM :: Int -> Int -> [Int] -> [Char]
makePBM i acc []     = (chr (acc * 2^(8-i))) : []
makePBM i acc (x:xs) | i==8      = (chr acc) : makePBM 0 0 (x:xs)
                     | otherwise = makePBM (i+1) n xs
                                    where
                                      n = if x==iter then (acc*2+1) else
(acc*2)

chunk :: Int -> [Int] -> [[Int]]
chunk width [] = []
chunk width c  = (fst split) : chunk width (snd split)
                 where split = splitAt width c