Difference between revisions of "Mandelbrot set"
(5 intermediate revisions by the same user not shown) | |||
Line 1: | Line 1: | ||
[[image:Mandelbrotset.png|right|thumb|The Mandelbrot set]] | [[image:Mandelbrotset.png|right|thumb|The Mandelbrot set]] | ||
+ | [[image:Mandelbrot.png|right|thumb|Detail of the Mandelbrot set, each pixel is coloured by the number of iterations it took to escape]] | ||
The formula for the mandelbrot set using complex numbers is: | The formula for the mandelbrot set using complex numbers is: | ||
Line 35: | Line 36: | ||
Zr = Cr | Zr = Cr | ||
Zi = Ci | Zi = Ci | ||
+ | |||
For i = 1 To 100 'Give up after 100 iterations to trade accuracy for speed | For i = 1 To 100 'Give up after 100 iterations to trade accuracy for speed | ||
Zrt = (Zr * Zr - Zi * Zi) + Cr | Zrt = (Zr * Zr - Zi * Zi) + Cr | ||
Line 52: | Line 54: | ||
− | 'Integer version with 24 bit fractional resolution | + | ' Integer version with 24 bit fractional resolution |
For ci = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.01 * (2 ^ 24) | For ci = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.01 * (2 ^ 24) | ||
For cr = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.008 * (2 ^ 24) | For cr = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.008 * (2 ^ 24) | ||
Line 72: | Line 74: | ||
If i = 101 Then | If i = 101 Then | ||
− | complexPlane. | + | complexPlane.Pset(cr, ci), RGB(0, 0, 0) |
Else | Else | ||
− | complexPlane. | + | complexPlane.Pset(cr, ci), RGB(255, 255, 255) |
End If | End If | ||
Line 80: | Line 82: | ||
Next ci | Next ci | ||
+ | // F# version using complex datatype | ||
+ | open System.Numerics | ||
+ | |||
+ | let maxIteration = 5000 | ||
+ | |||
+ | let mandelbrot c = | ||
+ | let rec loop z iteration = | ||
+ | if Complex.Abs(Complex.Pow(z,2.0)) >= 4.0 then iteration | ||
+ | elif iteration = maxIteration then 0 | ||
+ | else loop ((z * z) + c) (iteration + 1) | ||
+ | loop c 0 | ||
+ | |||
+ | for y in -1.2..0.15..1.2 do | ||
+ | for x in -2.0..0.07..0.9 do | ||
+ | printf "%c" " .:-;!+>|)&IH8%#".[mandelbrot(Complex(x, y)) &&& 15] | ||
+ | printfn "" | ||
[[Category:General information]] | [[Category:General information]] |
Latest revision as of 14:35, 13 November 2015
The formula for the mandelbrot set using complex numbers is:
- Z = Z2 + C
If Z escapes to infinity then we are outside the set, otherwise we are inside.
C is a point in the complex plane where i is on the vertical axis and real numbers are on the horizontal axis. A complex number has the form of x+iy where x and y are real numbers and i is the sqare root of -1. so the formula can also be written:
- (Zr+Zi) = ((Zr+Zi) * (Zr+Zi)) + (Cr+Ci)
Complex addition has the form:
- (a+bi)+(c+di) = (a+c)+i(b+d)
Complex multiplication has the form:
- (a+bi)*(c+di) = (ac-bd)+i(ad+bc)
so
- Z2 = ((Zr+Zi) * (Zr+Zi)) = (Zr*Zr - Zir*Zir) + i(Zr*Zir + Zir*Zr)
and
- Z2 + C = (Zr*Zr - Zir*Zir) + i(Zr*Zir + Zir*Zr) + (Cr+Ci) = (Zr*Zr - Zir*Zir) + Cr + i((Zr*Zir + Zir*Zr) + Cir)
that gives us
- Zr = (Zr*Zr - Zir*Zir) + Cr
- Zi = i((Zr*Zir + Zir*Zr) + Cir)
we don't need the imaginary part to determine if Z escapes to infinity so we can use
- Zir = (Zr*Zir + Zir*Zr) + Cir
we must remember to use temporary values to calculate Zr so it does not affect the result of Zi.
Z = C when we enter the loop. If Z escapes to infinity at any time during the iterations then C is outside the mandelbrot set. When C is confirmed to be outside or inside we plot a corresponding white or black pixel at coordinate C. After N iterations we give up and assume that Z will never escape to infinity to make sure the program finishes in a reasonable time.
' Floating point version For Ci = -2 To 2 Step 0.01 For Cr = -2 To 2 Step 0.008 Zr = Cr Zi = Ci For i = 1 To 100 'Give up after 100 iterations to trade accuracy for speed Zrt = (Zr * Zr - Zi * Zi) + Cr Zi = (Zr * Zi + Zi * Zr) + Ci Zr = Zrt If (Zr * Zr + Zi * Zi) > 4 Then Exit For Next i If i = 101 Then complexPlane.PSet (Cr, Ci), RGB(0, 0, 0) 'Looks like we are inside Else complexPlane.PSet (Cr, Ci), RGB(255, 255, 255) 'Infinity End If Next Cr Next Ci
' Integer version with 24 bit fractional resolution For ci = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.01 * (2 ^ 24) For cr = -2 * (2 ^ 24) To 2 * (2 ^ 24) Step 0.008 * (2 ^ 24) zr = cr zi = ci For i = 1 To 100 'Give up after 100 iterations to trade accuracy for speed zr2 = zr >> 12 zr2 *= zr2 zi2 = zi >> 12 zi2 *= zi2 zi = (((zr >> 12) * (zi >> 12)) * 2) + ci zr = (zr2 - zi2) + cr If (zr2 + zi2) > (4 << 24) Then Exit For Next i If i = 101 Then complexPlane.Pset(cr, ci), RGB(0, 0, 0) Else complexPlane.Pset(cr, ci), RGB(255, 255, 255) End If Next cr Next ci
// F# version using complex datatype open System.Numerics let maxIteration = 5000 let mandelbrot c = let rec loop z iteration = if Complex.Abs(Complex.Pow(z,2.0)) >= 4.0 then iteration elif iteration = maxIteration then 0 else loop ((z * z) + c) (iteration + 1) loop c 0 for y in -1.2..0.15..1.2 do for x in -2.0..0.07..0.9 do printf "%c" " .:-;!+>|)&IH8%#".[mandelbrot(Complex(x, y)) &&& 15] printfn ""