{-# LANGUAGE BangPatterns, CPP #-}
{-# OPTIONS -O2 -funbox-strict-fields -fvia-C -optc-O2 -fexcess-precision #-}
--
-- The Computer Language Benchmarks Game
-- http://shootout.alioth.debian.org/
--
-- Contributed by Spencer Janssen, Trevor McCort, Christophe Poucet and Don Stewart
-- Parallelised by Tim Newsham
--
-- This version: no threads, Eden parMap, by Jost Berthold
--
-- Must be compiled with the -fexcess-precision flag as a pragma. GHC
-- currently doesn't recognise the -fexcess-precision flag on the command
-- line (!).
--
-- The following flags are suggested when compiling:
--
--      ghc -optc-march=pentium4 -optc-mfpmath=sse -optc-msse2 -parmpi --make
--
--      $ time ./A 6400 +RTS -N6
--

import System.Environment
import System.IO
import Foreign
import Foreign.Marshal.Array
import Control.Monad

import Control.Parallel.Eden

import ParMap

main = do
    args <- getArgs
    -- width in pixels, parmap version
    let (w,v) = case args of
                  []   -> (1000,3)
                  [s]  -> (trunc8 (read s), 3)
                  (s:k:_) -> (trunc8 (read s), read k)
        trunc8 n = 8 * (n `div` 8)
        -- parmap version:
        myParMap = (!! v) 
           (cycle [map, parMap, parMap2]) --, parMapWP ])
        -- width in bytes
        n      = w `div` 8
        -- width of a pixel in the complex plane
        m    = 2 / fromIntegral w
        rows = [T 1 0 y (fromIntegral y * m - 1) | y <- [0..w-1]]

    putStrLn ("P4\n"++show w++" "++show w)

    let results = myParMap (computeRow w m n) rows
                  --------
    mapM_ (\r -> writeBytes r n) results

--------------------------------

-- computes one line of the image, starting at the given pixel
-- w - width in pixels
-- m - width of a pixel in the complex plane
-- n - width in bytes
-- start 
computeRow :: Int -> Double -> Int -> T -> [Word8]
computeRow w m n (T _  x y ci) -- coord 
    = map (compute8 w m n)
         [ T bx' x' y ci | bx' <- [ 1..n ], let x' = x+8*(bx'-1) ]

----------------------------------------

-- T bs x y ci
--    bx - x position in bytes
--    x  - x position in pixels
--    y  - y position in pixels
--    ci - y position in complex plane
data T = T !Int !Int !Int !Double

instance NFData T
instance Trans T
instance Trans Word8

-- w - image width in pixels
-- iw - pixel width in the complex plane
-- bw - image width in bytes
compute8 !w !iw !bw (T bx x y ci) = (loop_x w x 8 iw ci 0)

-- w - image width in pixels
-- x - current x coordinate in pixels
-- n - bit positition from 8 to 0
-- iw - pixel width in the complex plane
-- ci - current y coordinate in complex plane
-- b - accumulated bit value.
loop_x !w !x !n !iw !ci !b
    | x < w = if n == 0
                    then b
                    else loop_x w (x+1) (n-1) iw ci (b+b+v)
    | otherwise = b `shiftL` n
  where
    v = fractal 0 0 (fromIntegral x * iw - 1.5) ci 50

-- julia function (r :+ i) (cr :+ ci) with max iterations k.
fractal :: Double -> Double -> Double -> Double -> Int -> Word8
fractal !r !i !cr !ci !k
    | r2 + i2 > 4 = 0
    | k == 0      = 1
    | otherwise   = fractal (r2-i2+cr) ((r+r)*i+ci) cr ci (k-1)
  where
    (!r2,!i2) = (r*r,i*i)


writeBytes :: [Word8] -> Int -> IO ()
writeBytes list size = do a <- newArray list
                          hPutBuf stdout a size
          
