Some relecture

This commit is contained in:
Yann Esposito 2012-06-14 16:30:46 +02:00
parent 505a5dd323
commit 8527baa019
8 changed files with 132 additions and 54 deletions

View file

@ -39,7 +39,9 @@ and something nice to see in 3D.
Here is what you'll end with:
blogimage("GoldenMandelbulb.png","A golden mandelbulb")
blogfigure("GoldenMandelbulb.png","The entire Mandelbulb")
blogfigure("3DMandelbulbDetail.png","A Mandelbulb detail")
blogfigure("3DMandelbulbDetail2.png","Another detail of the Mandelbulb")
And here are the intermediate steps:

View file

@ -140,7 +140,7 @@ Given two coordinates in pixels, it returns some integer value:
> f (complex r i) 0 64
It uses the main Mandelbrot function for each complex \\(c\\).
The Mandelbrot set is the set of complex number c such that the following sequence does not escape to infinity.
The Mandelbrot set is the set of complex number \\(c\\) such that the following sequence does not escape to infinity.
Let us define \\(f_c: \mathbb{C} \to \mathbb{C}\\)

View file

@ -55,7 +55,7 @@ The method I use is a rough approximation.
I consider the Mandelbrot set to be almost convex.
The result will be good enough.
We change slightly the drawMandelbrot function.
We change slightly the `drawMandelbrot` function.
We replace the `Points` by `LineLoop`
> drawMandelbrot =
@ -77,12 +77,12 @@ we will choose only point on the surface.
> map (\(x,y,c) -> (x,-y,c)) (reverse positivePoints)
We only need to compute the positive point.
The Mandelbrot set is symmetric on the abscisse axis.
The Mandelbrot set is symmetric relatively to the abscisse axis.
> positivePoints :: [(GLfloat,GLfloat,Color3 GLfloat)]
> positivePoints = do
> x <- [-width..width]
> let y = findMaxOrdFor (mandel x) 0 height (log2 height)
> let y = maxZeroIndex (mandel x) 0 height (log2 height)
> if y < 1 -- We don't draw point in the absciss
> then []
> else return (x/width,y/height,colorFromValue $ mandel x y)
@ -92,22 +92,30 @@ The Mandelbrot set is symmetric on the abscisse axis.
This function is interesting.
For those not used to the list monad here is a natural language version of this function:
~~~
<code class="no-highlight">
positivePoints =
for all x in the range [-width..width]
let y be smallest number s.t. mandel x y > 0
if y is on 0 then don't return a point
else return the value corresonding to (x,y,color for (x+iy))
~~~
</code>
In fact using the list monad you write like if you consider only one element at a time and the computation is done non deterministically.
To find the smallest number such that `mandel x y > 0` we use a simple dichotomy:
> findMaxOrdFor func minval maxval 0 = (minval+maxval)/2
> findMaxOrdFor func minval maxval n =
> -- given f min max nbtest,
> -- considering
> -- - f is an increasing function
> -- - f(min)=0
> -- - f(max)≠0
> -- then maxZeroIndex f min max nbtest returns x such that
> -- f(x - ε)=0 and f(x + ε)≠0
> -- where ε=(max-min)/2^(nbtest+1)
> maxZeroIndex func minval maxval 0 = (minval+maxval)/2
> maxZeroIndex func minval maxval n =
> if (func medpoint) /= 0
> then findMaxOrdFor func minval medpoint (n-1)
> else findMaxOrdFor func medpoint maxval (n-1)
> then maxZeroIndex func minval medpoint (n-1)
> else maxZeroIndex func medpoint maxval (n-1)
> where medpoint = (minval+maxval)/2
No rocket science here. See the result now:

View file

@ -10,7 +10,7 @@ But it will be enough for us to create something that look nice.
This section is quite long, but don't be afraid,
most of the code is some OpenGL boilerplate.
For those you want to skim,
If you just want to skim this section,
here is a high level representation:
> - OpenGL Boilerplate
@ -263,7 +263,7 @@ depthPoints = do
x <- [-width..width]
y <- [-height..height]
let
depthOf x' y' = findMaxOrdFor (mandel x' y') 0 deep logdeep
depthOf x' y' = maxZeroIndex (mandel x' y') 0 deep logdeep
logdeep = floor ((log deep) / log 2)
z1 = depthOf x y
z2 = depthOf (x+1) y
@ -292,7 +292,7 @@ Here is a harder to read but shorter and more generic rewritten function:
> y <- [-height..height]
> let
> neighbors = [(x,y),(x+1,y),(x+1,y+1),(x,y+1)]
> depthOf (u,v) = findMaxOrdFor (mandel u v) 0 deep logdeep
> depthOf (u,v) = maxZeroIndex (mandel u v) 0 deep logdeep
> logdeep = floor ((log deep) / log 2)
> -- zs are 3D points with found depth
> zs = map (\(u,v) -> (u,v,depthOf (u,v))) neighbors
@ -324,11 +324,21 @@ The rest of the program is very close to the preceding one.
<div style="display:none">
> findMaxOrdFor func minval maxval 0 = (minval+maxval)/2
> findMaxOrdFor func minval maxval n =
> -- given f min max nbtest,
> -- considering
> -- - f is an increasing function
> -- - f(min)=0
> -- - f(max)≠0
> -- then maxZeroIndex f min max nbtest returns x such that
> -- f(x - ε)=0 and f(x + ε)≠0
> -- where ε=(max-min)/2^(nbtest+1)
> maxZeroIndex :: (Fractional a,Num a,Num b,Eq b) =>
> (a -> b) -> a -> a -> Int -> a
> maxZeroIndex func minval maxval 0 = (minval+maxval)/2
> maxZeroIndex func minval maxval n =
> if (func medpoint) /= 0
> then findMaxOrdFor func minval medpoint (n-1)
> else findMaxOrdFor func medpoint maxval (n-1)
> then maxZeroIndex func minval medpoint (n-1)
> else maxZeroIndex func medpoint maxval (n-1)
> where medpoint = (minval+maxval)/2
I made the color slightly brighter

View file

@ -42,7 +42,7 @@ This is similar to the preceding section.
> y <- [0..height]
> let
> neighbors = [(x,y),(x+1,y),(x+1,y+1),(x,y+1)]
> depthOf (u,v) = findMaxOrdFor (ymandel u v) 0 deep 7
> depthOf (u,v) = maxZeroIndex (ymandel u v) 0 deep 7
> -- zs are 3D points with found depth
> zs = map (\(u,v) -> (u,v,depthOf (u,v))) neighbors
> -- ts are 3D pixels + mandel value
@ -56,11 +56,20 @@ This is similar to the preceding section.
> -- Draw two triangles
> else [ps!!0,ps!!1,ps!!2,ps!!0,ps!!2,ps!!3]
>
> findMaxOrdFor func minval maxval 0 = (minval+maxval)/2
> findMaxOrdFor func minval maxval n =
>
> -- given f min max nbtest,
> -- considering
> -- - f is an increasing function
> -- - f(min)=0
> -- - f(max)≠0
> -- then maxZeroIndex f min max nbtest returns x such that
> -- f(x - ε)=0 and f(x + ε)≠0
> -- where ε=(max-min)/2^(nbtest+1)
> maxZeroIndex func minval maxval 0 = (minval+maxval)/2
> maxZeroIndex func minval maxval n =
> if (func medpoint) /= 0
> then findMaxOrdFor func minval medpoint (n-1)
> else findMaxOrdFor func medpoint maxval (n-1)
> then maxZeroIndex func minval medpoint (n-1)
> else maxZeroIndex func medpoint maxval (n-1)
> where medpoint = (minval+maxval)/2
>
> colorFromValue n =

View file

@ -2,21 +2,23 @@
Some points:
1. OpenGL and GLUT are linked to the C library.
1. OpenGL and GLUT is done in C.
In particular the `mainLoop` function is a direct link to the C library (FFI).
This function if so far from the pure spirit of functional languages.
This function is clearly far from the functional paradigm.
Could we make this better?
We will have two choices, or create our own `mainLoop` function to make it more functional.
Or deal with the imperative nature of the GLUT `mainLoop` function.
As a goal of this article is to understand how to deal with existing library and particularly the one coming from imperative language we will continue to use the `mainLoop` function.
2. Or main problem come from user interaction.
If you ask the Internet, about how to deal with user interaction with a functional paradigm, the main answer is to use _functional reactive programming_ (FRP).
I read very few about FRP, and I might be completely wrong when I say that it is about creating a DSL where atoms are time functions.
While I'm writing these lines, I don't know if I'll do something looking close to that.
For now I'll simply try to resolve the first problem.
As one of the goal of this article is to understand how to deal with existing libraries and particularly the one coming from imperative languages, we will continue to use the `mainLoop` function.
2. Our main problem come from user interaction.
If you ask "the Internet",
about how to deal with user interaction with a functional paradigm,
the main answer is to use _functional reactive programming_ (FRP).
I won't use FRP in this article.
Instead, I'll use a simpler while less effective way to deal with user interaction.
But The method I'll use will be as pure and functional as possible.
Then here is how I imagine things should go.
First, what the main loop should look like:
Here is how I imagine things should go.
First, what the main loop should look like if we could make our own:
<code class="no-highlight">
functionalMainLoop =
@ -89,6 +91,14 @@ We simply have to provide the definition of some functions.
> camPos = position w,
> camDir = angle w,
> camZoom = scale w }
> -- objects for world w
> -- is the list of one unique element
> -- The element is an YObject
> -- more precisely the XYFunc Function3D Box3D
> -- where the Function3D is the type
> -- Point -> Point -> Maybe (Point,Color)
> -- and its value here is ((shape w) res)
> -- and the Box3D value is defbox
> objects w = [XYFunc ((shape w) res) defbox]
> where
> res = resolution $ box w
@ -110,8 +120,8 @@ These function are used to update the world state.
> zdir :: Point3D
> zdir = makePoint3D (0,0,1)
Note `(-*<)` is scalar product.
Also note we could add Point3D as numbers.
Note `(-*<)` is the scalar product (`α -*< (x,y,z) = (αx,αy,αz)`).
Also note we could add two Point3D.
> rotate :: Point3D -> Scalar -> World -> World
> rotate dir angleValue world =
@ -188,9 +198,9 @@ Because we consider partial functions
> shapeFunc :: Scalar -> Function3D
> shapeFunc res x y =
> let
> z = findMaxOrdFor (ymandel x y) 0 1 20
> z = maxZeroIndex (ymandel x y) 0 1 20
> in
> if and [ findMaxOrdFor (ymandel (x+xeps) (y+yeps)) 0 1 20 < 0.000001 |
> if and [ maxZeroIndex (ymandel (x+xeps) (y+yeps)) 0 1 20 < 0.000001 |
> val <- [res], xeps <- [-val,val], yeps<-[-val,val]]
> then Nothing
> else Just (z,colorFromValue ((ymandel x y z) * 64))
@ -207,13 +217,21 @@ With the color function.
The rest is similar to the preceding sections.
> findMaxOrdFor :: (Fractional a,Num a,Num b,Eq b) =>
> -- given f min max nbtest,
> -- considering
> -- - f is an increasing function
> -- - f(min)=0
> -- - f(max)≠0
> -- then maxZeroIndex f min max nbtest returns x such that
> -- f(x - ε)=0 and f(x + ε)≠0
> -- where ε=(max-min)/2^(nbtest+1)
> maxZeroIndex :: (Fractional a,Num a,Num b,Eq b) =>
> (a -> b) -> a -> a -> Int -> a
> findMaxOrdFor _ minval maxval 0 = (minval+maxval)/2
> findMaxOrdFor func minval maxval n =
> if func medpoint /= 0
> then findMaxOrdFor func minval medpoint (n-1)
> else findMaxOrdFor func medpoint maxval (n-1)
> maxZeroIndex func minval maxval 0 = (minval+maxval)/2
> maxZeroIndex func minval maxval n =
> if (func medpoint) /= 0
> then maxZeroIndex func minval medpoint (n-1)
> else maxZeroIndex func medpoint maxval (n-1)
> where medpoint = (minval+maxval)/2
>
> ymandel :: Point -> Point -> Point -> Point

View file

@ -1,7 +1,8 @@
## Optimization
From the architecture stand point all is clear.
If you read the code, you'll see I didn't made everything perfect, for example, I didn't coded nicely the lights.
If you read the code of `YGL.hs`, you'll see I didn't made everything perfect.
For example, I didn't finished the code of the lights.
But I believe it is a good first step and it will be easy to go further.
The separation between rendering and world behavior is clear.
Unfortunately the program of the preceding session is extremely slow.
@ -124,13 +125,14 @@ Our initial world state is slightly changed:
> , anglePerSec = 5.0
> , position = makePoint3D (0,0,0)
> , scale = 1.0
> , box = Box3D { minPoint = makePoint3D (-2,-2,-2)
> , maxPoint = makePoint3D (2,2,2)
> , resolution = 0.03 }
> , box = Box3D { minPoint = makePoint3D (0-eps, 0-eps, 0-eps)
> , maxPoint = makePoint3D (0+eps, 0+eps, 0+eps)
> , resolution = 0.02 }
> , told = 0
> -- We declare cache directly this time
> , cache = objectFunctionFromWorld initialWorld
> }
> where eps=2
We use the `YGL.getObject3DFromShapeFunction` function directly.
This way instead of providing `XYFunc`, we provide directly a list of Atoms.
@ -183,9 +185,9 @@ All the rest is exactly the same.
> shapeFunc :: Scalar -> Function3D
> shapeFunc res x y =
> let
> z = findMaxOrdFor (ymandel x y) 0 1 20
> z = maxZeroIndex (ymandel x y) 0 1 20
> in
> if and [ findMaxOrdFor (ymandel (x+xeps) (y+yeps)) 0 1 20 < 0.000001 |
> if and [ maxZeroIndex (ymandel (x+xeps) (y+yeps)) 0 1 20 < 0.000001 |
> val <- [res], xeps <- [-val,val], yeps<-[-val,val]]
> then Nothing
> else Just (z,colorFromValue 0)
@ -198,13 +200,21 @@ All the rest is exactly the same.
> in
> makeColor (t n) (t (n+5)) (t (n+10))
>
> findMaxOrdFor :: (Fractional a,Num a,Num b,Eq b) =>
> -- given f min max nbtest,
> -- considering
> -- - f is an increasing function
> -- - f(min)=0
> -- - f(max)≠0
> -- then maxZeroIndex f min max nbtest returns x such that
> -- f(x - ε)=0 and f(x + ε)≠0
> -- where ε=(max-min)/2^(nbtest+1)
> maxZeroIndex :: (Fractional a,Num a,Num b,Eq b) =>
> (a -> b) -> a -> a -> Int -> a
> findMaxOrdFor _ minval maxval 0 = (minval+maxval)/2
> findMaxOrdFor func minval maxval n =
> maxZeroIndex _ minval maxval 0 = (minval+maxval)/2
> maxZeroIndex func minval maxval n =
> if func medpoint /= 0
> then findMaxOrdFor func minval medpoint (n-1)
> else findMaxOrdFor func medpoint maxval (n-1)
> then maxZeroIndex func minval medpoint (n-1)
> else maxZeroIndex func medpoint maxval (n-1)
> where medpoint = (minval+maxval)/2
>
> ymandel :: Point -> Point -> Point -> Point
@ -212,6 +222,8 @@ All the rest is exactly the same.
</div>
And you can also consider small changes in other source files.
- [`YGL.hs`](code/06_Mandelbulb/YGL.hs), the 3D rendering framework
- [`Mandel`](code/06_Mandelbulb/Mandel.hs), the mandel function
- [`ExtComplex`](code/06_Mandelbulb/ExtComplex.hs), the extended complexes

View file

@ -0,0 +1,19 @@
## Conclusion
In conclusion, if you add the courage to read this long article,
you should see how to organize your code in a functional way.
As we can use imperative style in a functional language,
know you can use functional style in imperative languages.
The usage of Haskell made it very simple to use our own data type.
The gigantic advantage is to code while being perfectly pure.
In the two last sections, the code is completely pure and functional.
Furthermore I don't use `GLfloat`, `Color3` or any other OpenGL type.
Therefore, if tomorrow I want to use another library I could do this simply by updating the `YGL` module.
The `YGL` module is a "wrapper".
It is a clean separator between the imperative paradigm and functional paradigm.
Furthermore I demonstrated you can make user interaction and display 3D objects in real time using only a pure functional paradigm.