2012-06-01 13:24:58 +00:00
|
|
|
## First version
|
|
|
|
|
|
|
|
We can consider two parts.
|
2012-06-12 14:17:25 +00:00
|
|
|
The first being mostly some boilerplate[^011].
|
|
|
|
And the second part more focused on OpenGL and content.
|
2012-06-01 13:24:58 +00:00
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
[^011]: Generally in Haskell you need to declare a lot of import lines.
|
2012-06-01 13:24:58 +00:00
|
|
|
This is something I find annoying.
|
|
|
|
In particular, it should be possible to create a special file, Import.hs
|
|
|
|
which make all the necessary import for you, as you generally need them all.
|
|
|
|
I understand why this is cleaner to force the programmer not to do so,
|
|
|
|
but, each time I do a copy/paste, I feel something is wrong.
|
|
|
|
I believe this concern can be generalized to the lack of namespace in Haskell.
|
|
|
|
|
|
|
|
### Let's play the song of our people
|
|
|
|
|
|
|
|
> import Graphics.Rendering.OpenGL
|
|
|
|
> import Graphics.UI.GLUT
|
|
|
|
> import Data.IORef
|
|
|
|
|
|
|
|
For efficiency reason, I won't use the default Haskell `Complex` data type.
|
|
|
|
|
|
|
|
> newtype Complex = C (Float,Float) deriving (Show,Eq)
|
|
|
|
|
|
|
|
|
|
|
|
> instance Num Complex where
|
|
|
|
> fromInteger n = C (fromIntegral n,0.0)
|
|
|
|
> C (x,y) * C (z,t) = C (z*x - y*t, y*z + x*t)
|
|
|
|
> C (x,y) + C (z,t) = C (x+z, y+t)
|
|
|
|
> abs (C (x,y)) = C (sqrt (x*x + y*y),0.0)
|
|
|
|
> signum (C (x,y)) = C (signum x , 0.0)
|
|
|
|
|
|
|
|
We declare some useful functions for manipulating complex numbers:
|
|
|
|
|
|
|
|
> complex :: Float -> Float -> Complex
|
|
|
|
> complex x y = C (x,y)
|
|
|
|
>
|
|
|
|
> real :: Complex -> Float
|
|
|
|
> real (C (x,y)) = x
|
|
|
|
>
|
|
|
|
> im :: Complex -> Float
|
|
|
|
> im (C (x,y)) = y
|
|
|
|
>
|
|
|
|
> magnitude :: Complex -> Float
|
|
|
|
> magnitude = real.abs
|
|
|
|
|
|
|
|
|
|
|
|
### Let us start
|
|
|
|
|
|
|
|
We start by giving the main architecture of our program:
|
|
|
|
|
|
|
|
> main :: IO ()
|
|
|
|
> main = do
|
|
|
|
> -- GLUT need to be initialized
|
|
|
|
> (progname,_) <- getArgsAndInitialize
|
|
|
|
> -- We will use the double buffered mode (GL constraint)
|
|
|
|
> initialDisplayMode $= [DoubleBuffered]
|
|
|
|
> -- We create a window with some title
|
|
|
|
> createWindow "Mandelbrot Set with Haskell and OpenGL"
|
|
|
|
> -- Each time we will need to update the display
|
|
|
|
> -- we will call the function 'display'
|
|
|
|
> displayCallback $= display
|
|
|
|
> -- We enter the main loop
|
|
|
|
> mainLoop
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
Mainly, we initialize our OpenGL application.
|
|
|
|
We declared that the function `display` will be used to render the graphics:
|
2012-06-01 13:24:58 +00:00
|
|
|
|
|
|
|
> display = do
|
|
|
|
> clear [ColorBuffer] -- make the window black
|
|
|
|
> loadIdentity -- reset any transformation
|
|
|
|
> preservingMatrix drawMandelbrot
|
|
|
|
> swapBuffers -- refresh screen
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
Also here, there is only one interesting line;
|
|
|
|
the draw will occur in the function `drawMandelbrot`.
|
2012-06-01 13:24:58 +00:00
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
This function will provide a list of draw actions.
|
|
|
|
Remember that OpenGL is imperative by design.
|
|
|
|
Then, one of the consequence is you must write the actions in the right order.
|
2012-06-01 13:24:58 +00:00
|
|
|
No easy parallel drawing here.
|
|
|
|
Here is the function which will render something on the screen:
|
|
|
|
|
|
|
|
> drawMandelbrot =
|
|
|
|
> -- We will print Points (not triangles for example)
|
|
|
|
> renderPrimitive Points $ do
|
|
|
|
> mapM_ drawColoredPoint allPoints
|
|
|
|
> where
|
|
|
|
> drawColoredPoint (x,y,c) = do
|
|
|
|
> color c -- set the current color to c
|
|
|
|
> -- then draw the point at position (x,y,0)
|
|
|
|
> -- remember we're in 3D
|
|
|
|
> vertex $ Vertex3 x y 0
|
|
|
|
|
|
|
|
The `mapM_` function is mainly the same as map but inside a monadic context.
|
|
|
|
More precisely, this can be transformed as a list of actions where the order is important:
|
|
|
|
|
|
|
|
~~~
|
|
|
|
drawMandelbrot =
|
|
|
|
renderPrimitive Points $ do
|
|
|
|
color color1
|
|
|
|
vertex $ Vertex3 x1 y1 0
|
|
|
|
...
|
|
|
|
color colorN
|
|
|
|
vertex $ Vertex3 xN yN 0
|
|
|
|
~~~
|
|
|
|
|
|
|
|
We also need some kind of global variables.
|
2012-06-12 14:17:25 +00:00
|
|
|
In fact, global variable are a proof of a design problem.
|
|
|
|
We will get rid of them later.
|
2012-06-01 13:24:58 +00:00
|
|
|
|
|
|
|
> width = 320 :: GLfloat
|
|
|
|
> height = 320 :: GLfloat
|
|
|
|
|
|
|
|
And of course our list of colored points.
|
|
|
|
In OpenGL the default coordinate are from -1 to 1.
|
|
|
|
|
|
|
|
> allPoints :: [(GLfloat,GLfloat,Color3 GLfloat)]
|
|
|
|
> allPoints = [ (x/width,y/height,colorFromValue $ mandel x y) |
|
|
|
|
> x <- [-width..width],
|
|
|
|
> y <- [-height..height]]
|
|
|
|
>
|
|
|
|
|
|
|
|
We need a function which transform an integer value to some color:
|
|
|
|
|
|
|
|
> colorFromValue n =
|
|
|
|
> let
|
|
|
|
> t :: Int -> GLfloat
|
|
|
|
> t i = 0.5 + 0.5*cos( fromIntegral i / 10 )
|
|
|
|
> in
|
|
|
|
> Color3 (t n) (t (n+5)) (t (n+10))
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
And now the `mandel` function.
|
2012-06-01 13:24:58 +00:00
|
|
|
Given two coordinates in pixels, it returns some integer value:
|
|
|
|
|
|
|
|
> mandel x y =
|
|
|
|
> let r = 2.0 * x / width
|
|
|
|
> i = 2.0 * y / height
|
|
|
|
> in
|
|
|
|
> f (complex r i) 0 64
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
It uses the main Mandelbrot function for each complex \\(c\\).
|
2012-06-14 16:04:16 +00:00
|
|
|
The Mandelbrot set is the set of complex number \\(c\\) such that the following sequence does not escape to infinity.
|
2012-06-01 13:24:58 +00:00
|
|
|
|
|
|
|
Let us define \\(f_c: \mathbb{C} \to \mathbb{C}\\)
|
|
|
|
|
|
|
|
$$ f_c(z) = z^2 + c $$
|
|
|
|
|
|
|
|
The sequence is:
|
|
|
|
|
|
|
|
$$ 0 \rightarrow f_c(0) \rightarrow f_c(f_c(0)) \rightarrow \cdots \rightarrow f^n_c(0) \rightarrow \cdots $$
|
|
|
|
|
|
|
|
Of course, instead of trying to test the real limit, we just make a test after a finite number of occurrences.
|
|
|
|
|
|
|
|
> f :: Complex -> Complex -> Int -> Int
|
|
|
|
> f c z 0 = 0
|
|
|
|
> f c z n = if (magnitude z > 2 )
|
|
|
|
> then n
|
|
|
|
> else f c ((z*z)+c) (n-1)
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
Well, if you download this file (look at the bottom of this section), compile it and run it this is the result:
|
2012-06-01 13:24:58 +00:00
|
|
|
|
|
|
|
blogimage("hglmandel_v01.png","The mandelbrot set version 1")
|
|
|
|
|
|
|
|
A first very interesting property of this program is that the computation for all the points is done only once.
|
2012-06-12 14:17:25 +00:00
|
|
|
It is a bit long before the first image appears, but if you resize the window, it updates instantaneously.
|
2012-06-01 13:24:58 +00:00
|
|
|
This property is a direct consequence of purity.
|
|
|
|
If you look closely, you see that `allPoints` is a pure list.
|
2012-06-12 14:17:25 +00:00
|
|
|
Therefore, calling `allPoints` will always render the same result and Haskell is clever enough to use this property.
|
2012-06-01 13:24:58 +00:00
|
|
|
While Haskell doesn't garbage collect `allPoints` the result is reused for free.
|
|
|
|
We didn't specified this value should be saved for later use.
|
|
|
|
It is saved for us.
|
|
|
|
|
|
|
|
See what occurs if we make the window bigger:
|
|
|
|
|
|
|
|
blogimage("hglmandel_v01_too_wide.png","The mandelbrot too wide, black lines and columns")
|
|
|
|
|
2012-06-12 14:17:25 +00:00
|
|
|
We see some black lines because we drawn less point than there is on the surface.
|
2012-06-01 13:24:58 +00:00
|
|
|
We can repair this by drawing little squares instead of just points.
|
|
|
|
But, instead we will do something a bit different and unusual.
|