Haskell Fractal Plotter

articles ✒ haskell-fractals

In this article, I'll show you how to plot the Mandelbrot set with a short Haskell program.

For speed you'll probably want to use ghc to compile the program rather than hugs or ghci which will be very slow!

module Main where
import Complex

type Image = [Colour]
type Colour = [Integer]

-- Create a string suitable for saving as a PPM file to represent a given image
createppm :: Integer -> Integer -> Image -> String
createppm width height im = "P3\n" ++ show width ++ " " ++ show height ++ "\n255" ++ concat (map (('\n':) . spacejoin . map show) im) ++ "\n"

spacejoin (x:[]) = x
spacejoin (x:xs) = x ++ " " ++ spacejoin xs

We use our familiar PPM bitmap file creation function again. See my raytracer article for an explanation.

To find the colour of a co-ordinate in the complex plane we find the number of iterations for the recurrence $ z := z^2 + c $ to produce a z with modulus > 2.

mandelbrot :: Double -> Double -> Integer
mandelbrot x y = mandelbrot' (x :+ y) (0 :+ 0) 0

mandelbrot' :: Complex Double -> Complex Double -> Integer -> Integer
mandelbrot' c x iter | magnitude x > 2 = iter
                     | iter >= 255 = 255
                     | otherwise = mandelbrot' c (c + (x * x)) (iter + 1)

As you can see above, we give up after 255 iterations (otherwise for some points we'd never finish!).

And we just apply this calculation to lots of points to create an image.

evaluate2D :: (Double -> Double -> Integer) -> Double -> Double -> Double -> Double -> Integer -> Integer -> [[Integer]]
evaluate2D f x1 x2 y1 y2 width height = [ [ f x y | x <- map fromPixelX [0..width - 1] ] | y <- map fromPixelY [0..height - 1] ]
    where fromPixelX x = x1 + (fromInteger x * xscale)
          fromPixelY y = y1 + (fromInteger y * yscale)
          xscale = (x2 - x1) / fromInteger width
          yscale = (y2 - y1) / fromInteger height

To make the image we change the integers produced into colours with an arbitrary formula. Make your own if you wish.

tocol x = [x,x,x]

Then we bring it all together to make our PPM file.

plot = createppm 400 400 (map tocol (concat (evaluate2D mandelbrot (-2.0) 2.0 (-2.0) 2.0 400 400)))

main = writeFile "mandelbrot.ppm" plot

Here's the kind of output you'll get.

Mandelbrot fractal

Note: The image above wasn't actually produced by the haskell program. ghc is very slow compared to C++ code compiled with g++, so I rewrote above code in C++, rendered it at 3200x3200, then downscaled it with the GIMP to improve its appearance (making it smoother).

See the C++ version of the graph plotter.

comment on this page