Mandelbrot set: Difference between revisions
Updated first D entry |
|||
Line 1,090:
=={{header|D}}==
===Textual
This uses <code>std.complex.Complex</code> because D built-in complex numbers are deprecated.
<lang d>import std.stdio, std.
void main() {
enum maxIter =
foreach (y; -39 .. 39) {
foreach (x; -39 .. 39) {
auto c = Complex!double(y/40.0 - 0.5
z = Complex!double(0.0
i = 0;
for (; i < maxIter && z.abs(
z = z ^^ 2 + c;
write(i == maxIter ? '#' : ' ');
Line 1,107 ⟶ 1,108:
}
}</lang>
{{out}}
<pre> #
#
#
#
#
###
#####
#####
###
#
#########
#############
###############
#####################
#####################
###################
###################
###################
###################
#######################
###################
###################
#####################
###################
###################
#################
###############
#############
#########
#
###############
#######################
# ######################### #
#############################
# ############################### #
#################################
###################################
#######################################
### ######################################### ###
#################################################
###############################################
#############################################
#############################################
###############################################
###############################################
###################################################
#################################################
#################################################
###################################################
###################################################
# ################################################### #
##### ################################################### #####
###### ################################################### ######
####### ################################################### #######
#######################################################################
######### ################################################### #########
###### ################################################### ######
##### ################################################### #####
###################################################
###################################################
###################################################
###################################################
#################################################
#################################################
###################################################
###############################################
###############################################
###########################################
#########################################
#############################################
#### ################## ################## ####
### ################ ################ ###
# ############## ############## #
########### ###########
## ##### ##### ##
# # # #
</pre>
===Graphical version===
{{libheader|QD}} {{libheader|SDL}} {{libheader|Phobos}}
|
Revision as of 23:04, 4 April 2012
This page uses content from Wikipedia. The original article was at Mandelbrot_set. The list of authors can be seen in the page history. As with Rosetta Code, the text of Wikipedia is available under the GNU FDL. (See links for details on variance) |
You are encouraged to solve this task according to the task description, using any language you may know.
Generate and draw the Mandelbrot set. Note that there are many algorithms to draw Mandelbrot set and there are many functions which generate it .
Ada
mandelbrot.adb: <lang Ada>with Lumen.Binary; package body Mandelbrot is
function Create_Image (Width, Height : Natural) return Lumen.Image.Descriptor is use type Lumen.Binary.Byte; Result : Lumen.Image.Descriptor; X0, Y0 : Float; X, Y, Xtemp : Float; Iteration : Float; Max_Iteration : constant Float := 1000.0; Color : Lumen.Binary.Byte; begin Result.Width := Width; Result.Height := Height; Result.Complete := True; Result.Values := new Lumen.Image.Pixel_Matrix (1 .. Width, 1 .. Height); for Screen_X in 1 .. Width loop for Screen_Y in 1 .. Height loop X0 := -2.5 + (3.5 / Float (Width) * Float (Screen_X)); Y0 := -1.0 + (2.0 / Float (Height) * Float (Screen_Y)); X := 0.0; Y := 0.0; Iteration := 0.0; while X * X + Y * Y <= 4.0 and then Iteration < Max_Iteration loop Xtemp := X * X - Y * Y + X0; Y := 2.0 * X * Y + Y0; X := Xtemp; Iteration := Iteration + 1.0; end loop; if Iteration = Max_Iteration then Color := 255; else Color := 0; end if; Result.Values (Screen_X, Screen_Y) := (R => Color, G => Color, B => Color, A => 0); end loop; end loop; return Result; end Create_Image;
end Mandelbrot;</lang>
mandelbrot.ads: <lang Ada>with Lumen.Image;
package Mandelbrot is
function Create_Image (Width, Height : Natural) return Lumen.Image.Descriptor;
end Mandelbrot;</lang>
test_mandelbrot.adb: <lang Ada>with System.Address_To_Access_Conversions; with Lumen.Window; with Lumen.Image; with Lumen.Events; with GL; with Mandelbrot;
procedure Test_Mandelbrot is
Program_End : exception;
Win : Lumen.Window.Handle; Image : Lumen.Image.Descriptor; Tx_Name : aliased GL.GLuint; Wide, High : Natural := 400;
-- Create a texture and bind a 2D image to it procedure Create_Texture is use GL;
package GLB is new System.Address_To_Access_Conversions (GLubyte);
IP : GLpointer; begin -- Create_Texture -- Allocate a texture name glGenTextures (1, Tx_Name'Unchecked_Access);
-- Bind texture operations to the newly-created texture name glBindTexture (GL_TEXTURE_2D, Tx_Name);
-- Select modulate to mix texture with color for shading glTexEnvi (GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_MODULATE);
-- Wrap textures at both edges glTexParameteri (GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_REPEAT); glTexParameteri (GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_REPEAT);
-- How the texture behaves when minified and magnified glTexParameteri (GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST); glTexParameteri (GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
-- Create a pointer to the image. This sort of horror show is going to -- be disappearing once Lumen includes its own OpenGL bindings. IP := GLB.To_Pointer (Image.Values.all'Address).all'Unchecked_Access;
-- Build our texture from the image we loaded earlier glTexImage2D (GL_TEXTURE_2D, 0, GL_RGBA, GLsizei (Image.Width), GLsizei (Image.Height), 0, GL_RGBA, GL_UNSIGNED_BYTE, IP); end Create_Texture;
-- Set or reset the window view parameters procedure Set_View (W, H : in Natural) is use GL; begin -- Set_View GL.glEnable (GL.GL_TEXTURE_2D); glClearColor (0.8, 0.8, 0.8, 1.0);
glMatrixMode (GL_PROJECTION); glLoadIdentity; glViewport (0, 0, GLsizei (W), GLsizei (H)); glOrtho (0.0, GLdouble (W), GLdouble (H), 0.0, -1.0, 1.0);
glMatrixMode (GL_MODELVIEW); glLoadIdentity; end Set_View;
-- Draw our scene procedure Draw is use GL; begin -- Draw -- clear the screen glClear (GL_COLOR_BUFFER_BIT or GL_DEPTH_BUFFER_BIT); GL.glBindTexture (GL.GL_TEXTURE_2D, Tx_Name);
-- fill with a single textured quad glBegin (GL_QUADS); begin glTexCoord2f (1.0, 0.0); glVertex2i (GLint (Wide), 0);
glTexCoord2f (0.0, 0.0); glVertex2i (0, 0);
glTexCoord2f (0.0, 1.0); glVertex2i (0, GLint (High));
glTexCoord2f (1.0, 1.0); glVertex2i (GLint (Wide), GLint (High)); end; glEnd;
-- flush rendering pipeline glFlush;
-- Now show it Lumen.Window.Swap (Win); end Draw;
-- Simple event handler routine for keypresses and close-window events procedure Quit_Handler (Event : in Lumen.Events.Event_Data) is begin -- Quit_Handler raise Program_End; end Quit_Handler;
-- Simple event handler routine for Exposed events procedure Expose_Handler (Event : in Lumen.Events.Event_Data) is pragma Unreferenced (Event); begin -- Expose_Handler Draw; end Expose_Handler;
-- Simple event handler routine for Resized events procedure Resize_Handler (Event : in Lumen.Events.Event_Data) is begin -- Resize_Handler Wide := Event.Resize_Data.Width; High := Event.Resize_Data.Height; Set_View (Wide, High);
-- Image := Mandelbrot.Create_Image (Width => Wide, Height => High); -- Create_Texture;
Draw; end Resize_Handler;
begin
-- Create Lumen window, accepting most defaults; turn double buffering off -- for simplicity Lumen.Window.Create (Win => Win, Name => "Mandelbrot fractal", Width => Wide, Height => High, Events => (Lumen.Window.Want_Exposure => True, Lumen.Window.Want_Key_Press => True, others => False));
-- Set up the viewport and scene parameters Set_View (Wide, High);
-- Now create the texture and set up to use it Image := Mandelbrot.Create_Image (Width => Wide, Height => High); Create_Texture;
-- Enter the event loop declare use Lumen.Events; begin Select_Events (Win => Win, Calls => (Key_Press => Quit_Handler'Unrestricted_Access, Exposed => Expose_Handler'Unrestricted_Access, Resized => Resize_Handler'Unrestricted_Access, Close_Window => Quit_Handler'Unrestricted_Access, others => No_Callback)); end;
exception
when Program_End => null;
end Test_Mandelbrot;</lang>
Output:
ALGOL 68
Plot part of the Mandelbrot set as a pseudo-gif image.
<lang algol68> INT pix = 300, max iter = 256, REAL zoom = 0.33 / pix; [-pix : pix, -pix : pix] INT plane; COMPL ctr = 0.05 I 0.75 # center of set #;
- Compute the length of an orbit. #
PROC iterate = (COMPL z0) INT:
BEGIN COMPL z := 0, INT iter := 1; WHILE (iter +:= 1) < max iter # not converged # AND ABS z < 2 # not diverged # DO z := z * z + z0 OD; iter END;
- Compute set and find maximum orbit length. #
INT max col := 0; FOR x FROM -pix TO pix DO FOR y FROM -pix TO pix
DO COMPL z0 = ctr + (x * zoom) I (y * zoom); IF (plane [x, y] := iterate (z0)) < max iter THEN (plane [x, y] > max col | max col := plane [x, y]) FI OD
OD;
- Make a plot. #
FILE plot; INT num pix = 2 * pix + 1; make device (plot, "gif", whole (num pix, 0) + "x" + whole (num pix, 0)); open (plot, "mandelbrot.gif", stand draw channel); FOR x FROM -pix TO pix DO FOR y FROM -pix TO pix
DO INT col = (plane [x, y] > max col | max col | plane [x, y]); REAL c = sqrt (1- col / max col); # sqrt to enhance contrast # draw colour (plot, c, c, c); draw point (plot, (x + pix) / (num pix - 1), (y + pix) / (num pix - 1)) OD
OD; close (plot) </lang>
Applesoft BASIC
This version takes into account the Apple II's funky 280×192 6-color display, which has an effective resolution of only 140×192 in color.
<lang basic> 10 HGR2 20 XC = -0.5 : REM CENTER COORD X 30 YC = 0 : REM " " Y 40 S = 2 : REM SCALE 45 IT = 20 : REM ITERATIONS 50 XR = S * 4 / 3 : REM TOTAL RANGE OF X 60 YR = S : REM " " " Y 70 X0 = XC - (XR/2) : REM MIN VALUE OF X 80 X1 = XC + (XR/2) : REM MAX " " X 90 Y0 = YC - (YR/2) : REM MIN " " Y 100 Y1 = YC - (YR/2) : REM MAX " " Y 110 XM = XR / 279 : REM SCALING FACTOR FOR X 120 YM = YR / 191 : REM " " " Y 130 FOR YI = 0 TO 3 : REM INTERLEAVE 140 FOR YS = 0+YI TO 188+YI STEP 4 : REM Y SCREEN COORDINATE 145 HCOLOR=3 : HPLOT 0,YS TO 279,YS 150 FOR XS = 0 TO 278 STEP 2 : REM X SCREEN COORDINATE 170 X = XS * XM + X0 : REM TRANSL SCREEN TO TRUE X 180 Y = YS * YM + Y0 : REM TRANSL SCREEN TO TRUE Y 190 ZX = 0 200 ZY = 0 210 XX = 0 220 YY = 0 230 FOR I = 0 TO IT 240 ZY = 2 * ZX * ZY + Y 250 ZX = XX - YY + X 260 XX = ZX * ZX 270 YY = ZY * ZY 280 C = IT-I 290 IF XX+YY >= 4 GOTO 301 300 NEXT I 301 IF C >= 8 THEN C = C - 8 : GOTO 301 310 HCOLOR = C : HPLOT XS, YS TO XS+1, YS 320 NEXT XS 330 NEXT YS 340 NEXT YI </lang>
By making the following modifications, the same code will render the Mandelbrot set in monochrome at full 280×192 resolution.
<lang basic> 150 FOR XS = 0 TO 279 301 C = (C - INT(C/2)*2)*3 310 HCOLOR = C: HPLOT XS, YS </lang>
AutoHotkey
<lang autohotkey>Max_Iteration := 256 Width := Height := 400
File := "MandelBrot." Width ".bmp" Progress, b2 w400 fs9, Creating Colours ... Gosub, CreateColours Gosub, CreateBitmap Progress, Off Gui, -Caption Gui, Margin, 0, 0 Gui, Add, Picture,, %File% Gui, Show,, MandelBrot Return
GuiClose: GuiEscape: ExitApp
- ---------------------------------------------------------------------------
CreateBitmap: ; create and save a 32bit bitmap file
- ---------------------------------------------------------------------------
; define header details HeaderBMP := 14 HeaderDIB := 40 DataOffset := HeaderBMP + HeaderDIB ImageSize := Width * Height * 4 ; 32bit FileSize := DataOffset + ImageSize Resolution := 3780 ; from mspaint
; create bitmap header VarSetCapacity(IMAGE, FileSize, 0) NumPut(Asc("B") , IMAGE, 0x00, "Char") NumPut(Asc("M") , IMAGE, 0x01, "Char") NumPut(FileSize , IMAGE, 0x02, "UInt") NumPut(DataOffset , IMAGE, 0x0A, "UInt") NumPut(HeaderDIB , IMAGE, 0x0E, "UInt") NumPut(Width , IMAGE, 0x12, "UInt") NumPut(Height , IMAGE, 0x16, "UInt") NumPut(1 , IMAGE, 0x1A, "Short") ; Planes NumPut(32 , IMAGE, 0x1C, "Short") ; Bits per Pixel NumPut(ImageSize , IMAGE, 0x22, "UInt") NumPut(Resolution , IMAGE, 0x26, "UInt") NumPut(Resolution , IMAGE, 0x2A, "UInt")
; fill in Data Gosub, CreatePixels
; save Bitmap to file FileDelete, %File% Handle := DllCall("CreateFile", "Str", File, "UInt", 0x40000000 , "UInt", 0, "UInt", 0, "UInt", 2, "UInt", 0, "UInt", 0) DllCall("WriteFile", "UInt", Handle, "UInt", &IMAGE, "UInt" , FileSize, "UInt *", Bytes, "UInt", 0) DllCall("CloseHandle", "UInt", Handle)
Return
- ---------------------------------------------------------------------------
CreatePixels: ; create pixels for [-2 < x < 1] [-1.5 < y < 1.5]
- ---------------------------------------------------------------------------
Loop, % Height // 2 + 1 { yi := A_Index - 1 y0 := -1.5 + yi / Height * 3 ; range -1.5 .. +1.5 Progress, % 200*yi // Height, % "Current line: " 2*yi " / " Height Loop, %Width% { xi := A_Index - 1 x0 := -2 + xi / Width * 3 ; range -2 .. +1 Gosub, Mandelbrot p1 := DataOffset + 4 * (Width * yi + xi) NumPut(Colour, IMAGE, p1, "UInt") p2 := DataOffset + 4 * (Width * (Height-yi) + xi) NumPut(Colour, IMAGE, p2, "UInt") } }
Return
- ---------------------------------------------------------------------------
Mandelbrot: ; calculate a colour for each pixel
- ---------------------------------------------------------------------------
x := y := Iteration := 0 While, (x*x + y*y <= 4) And (Iteration < Max_Iteration) { xtemp := x*x - y*y + x0 y := 2*x*y + y0 x := xtemp Iteration++ } Colour := Iteration = Max_Iteration ? 0 : Colour_%Iteration%
Return
- ---------------------------------------------------------------------------
CreateColours: ; borrowed from PureBasic example
- ---------------------------------------------------------------------------
Loop, 64 { i4 := (i3 := (i2 := (i1 := A_Index - 1) + 64) + 64) + 64 Colour_%i1% := RGB(4*i1 + 128, 4*i1, 0) Colour_%i2% := RGB(64, 255, 4*i1) Colour_%i3% := RGB(64, 255 - 4*i1, 255) Colour_%i4% := RGB(64, 0, 255 - 4*i1) }
Return
- ---------------------------------------------------------------------------
RGB(r, g, b) { ; return 24bit color value
- ---------------------------------------------------------------------------
Return, (r&0xFF)<<16 | g<<8 | b
}</lang>
AWK
<lang AWK> BEGIN {
XSize=59; YSize=21; MinIm=-1.0; MaxIm=1.0;MinRe=-2.0; MaxRe=1.0; StepX=(MaxRe-MinRe)/XSize; StepY=(MaxIm-MinIm)/YSize; for(y=0;y<YSize;y++) { Im=MinIm+StepY*y; for(x=0;x<XSize;x++) { Re=MinRe+StepX*x; Zr=Re; Zi=Im; for(n=0;n<30;n++) { a=Zr*Zr; b=Zi*Zi; if(a+b>4.0) break; Zi=2*Zr*Zi+Im; Zr=a-b+Re; } printf "%c",62-n; } print ""; } exit;
}
>>>>>>=====<<<<<<<<<<<<<<<;;;;;;:::96032:;;;;<<<<========== >>>>>===<<<<<<<<<<<<<<<<;;;;;;;:::873*079::;;;;<<<<<======= >>>>===<<<<<<<<<<<<<<<;;;;;;;::9974 (.9::::;;<<<<<====== >>>==<<<<<<<<<<<<<<<;;;;;;:98888764 5789999:;;<<<<<==== >>==<<<<<<<<<<<<<;;;;::::996. &2 45335:;<<<<<<=== >>=<<<<<<<<<<<;;;::::::999752 *79:;<<<<<<== >=<<<<<<<<;;;:599999999886 %78:;;<<<<<<= ><<<<;;;;;:::972456-567763 +9;;<<<<<<< ><;;;;;;::::9875& .3 *9;;;<<<<<< >;;;;;;::997564' ' 8:;;;<<<<<< >::988897735/ &89:;;;<<<<<< >::988897735/ &89:;;;<<<<<< >;;;;;;::997564' ' 8:;;;<<<<<< ><;;;;;;::::9875& .3 *9;;;<<<<<< ><<<<;;;;;:::972456-567763 +9;;<<<<<<< >=<<<<<<<<;;;:599999999886 %78:;;<<<<<<= >>=<<<<<<<<<<<;;;::::::999752 *79:;<<<<<<== >>==<<<<<<<<<<<<<;;;;::::996. &2 45335:;<<<<<<=== >>>==<<<<<<<<<<<<<<<;;;;;;:98888764 5789999:;;<<<<<==== >>>>===<<<<<<<<<<<<<<<;;;;;;;::9974 (.9::::;;<<<<<====== >>>>>===<<<<<<<<<<<<<<<<;;;;;;;:::873*079::;;;;<<<<<======= </lang>
BASIC
This is almost exactly the same as the pseudocode from the Wikipedia entry's "For programmers" section (which it's closely based on, of course). The image generated is very blocky ("low-res") due to the selected video mode, but it's fairly accurate.
<lang qbasic>SCREEN 13 WINDOW (-2, 1.5)-(2, -1.5) FOR x0 = -2 TO 2 STEP .01
FOR y0 = -1.5 TO 1.5 STEP .01 x = 0 y = 0
iteration = 0 maxIteration = 223
WHILE (x * x + y * y <= (2 * 2) AND iteration < maxIteration) xtemp = x * x - y * y + x0 y = 2 * x * y + y0
x = xtemp
iteration = iteration + 1 WEND
IF iteration <> maxIteration THEN c = iteration ELSE c = 0 END IF
PSET (x0, y0), c + 32 NEXT
NEXT</lang>
Brace
This is a simple Mandelbrot plotter. A longer version based on this smooths colors, and avoids calculating the time-consuming black pixels: http://sam.ai.ki/brace/examples/mandelbrot.d/1 <lang brace>#!/usr/bin/env bx use b
Main(): num outside = 16, ox = -0.5, oy = 0, r = 1.5 long i, max_i = 100, rb_i = 30 space() uint32_t *px = pixel() num d = 2*r/h, x0 = ox-d*w_2, y0 = oy+d*h_2 for(y, 0, h): cmplx c = x0 + (y0-d*y)*I repeat(w): cmplx w = 0 for i=0; i < max_i && cabs(w) < outside; ++i w = w*w + c *px++ = i < max_i ? rainbow(i*359 / rb_i % 360) : black c += d</lang>
An example plot from the longer version:
C
PPM non interactive
Here is one file program. It directly creates ppm file. <lang C> /*
c program: -------------------------------- 1. draws Mandelbrot set for Fc(z)=z*z +c using Mandelbrot algorithm ( boolean escape time ) ------------------------------- 2. technique of creating ppm file is based on the code of Claudio Rocchini http://en.wikipedia.org/wiki/Image:Color_complex_plot.jpg create 24 bit color graphic file , portable pixmap file = PPM see http://en.wikipedia.org/wiki/Portable_pixmap to see the file use external application ( graphic viewer) */ #include <stdio.h> int main() { /* screen ( integer) coordinate */ int iX,iY; const int iXmax = 800; const int iYmax = 800; /* world ( double) coordinate = parameter plane*/ double Cx,Cy; const double CxMin=-2.5; const double CxMax=1.5; const double CyMin=-2.0; const double CyMax=2.0; /* */ double PixelWidth=(CxMax-CxMin)/iXmax; double PixelHeight=(CyMax-CyMin)/iYmax; /* color component ( R or G or B) is coded from 0 to 255 */ /* it is 24 bit color RGB file */ const int MaxColorComponentValue=255; FILE * fp; char *filename="new1.ppm"; char *comment="# ";/* comment should start with # */ static unsigned char color[3]; /* Z=Zx+Zy*i ; Z0 = 0 */ double Zx, Zy; double Zx2, Zy2; /* Zx2=Zx*Zx; Zy2=Zy*Zy */ /* */ int Iteration; const int IterationMax=200; /* bail-out value , radius of circle ; */ const double EscapeRadius=2; double ER2=EscapeRadius*EscapeRadius; /*create new file,give it a name and open it in binary mode */ fp= fopen(filename,"wb"); /* b - binary mode */ /*write ASCII header to the file*/ fprintf(fp,"P6\n %s\n %d\n %d\n %d\n",comment,iXmax,iYmax,MaxColorComponentValue); /* compute and write image data bytes to the file*/ for(iY=0;iY<iYmax;iY++) { Cy=CyMin + iY*PixelHeight; if (fabs(Cy)< PixelHeight/2) Cy=0.0; /* Main antenna */ for(iX=0;iX<iXmax;iX++) { Cx=CxMin + iX*PixelWidth; /* initial value of orbit = critical point Z= 0 */ Zx=0.0; Zy=0.0; Zx2=Zx*Zx; Zy2=Zy*Zy; /* */ for (Iteration=0;Iteration<IterationMax && ((Zx2+Zy2)<ER2);Iteration++) { Zy=2*Zx*Zy + Cy; Zx=Zx2-Zy2 +Cx; Zx2=Zx*Zx; Zy2=Zy*Zy; }; /* compute pixel color (24 bit = 3 bytes) */ if (Iteration==IterationMax) { /* interior of Mandelbrot set = black */ color[0]=0; color[1]=0; color[2]=0; } else { /* exterior of Mandelbrot set = white */ color[0]=255; /* Red*/ color[1]=255; /* Green */ color[2]=255;/* Blue */ }; /*write color to the file*/ fwrite(color,1,3,fp); } } fclose(fp); return 0; }</lang >
PPM Interactive
Infinitely zoomable OpenGL program. Adjustable colors, max iteration, black and white, screen dump, etc. Compile with gcc mandelbrot.c -lglut -lGLU -lGL -lm
- OpenBSD users, install freeglut package, and compile with
make mandelbrot CPPFLAGS='-I/usr/local/include `pkg-config glu --cflags`' LDLIBS='-L/usr/local/lib -lglut `pkg-config glu --libs` -lm'
<lang c>#include <stdio.h>
- include <stdlib.h>
- include <math.h>
- include <GL/glut.h>
- include <GL/gl.h>
- include <GL/glu.h>
void set_texture();
typedef struct {unsigned char r, g, b;} rgb_t; rgb_t **tex = 0; int gwin; GLuint texture; int width, height; int tex_w, tex_h; double scale = 1./256; double cx = -.6, cy = 0; int color_rotate = 0; int saturation = 1; int invert = 0; int max_iter = 256;
void render() { double x = (double)width /tex_w, y = (double)height/tex_h;
glClear(GL_COLOR_BUFFER_BIT); glTexEnvi(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_REPLACE);
glBindTexture(GL_TEXTURE_2D, texture);
glBegin(GL_QUADS);
glTexCoord2f(0, 0); glVertex2i(0, 0); glTexCoord2f(x, 0); glVertex2i(width, 0); glTexCoord2f(x, y); glVertex2i(width, height); glTexCoord2f(0, y); glVertex2i(0, height);
glEnd();
glFlush(); glFinish(); }
int dump = 1; void screen_dump() { char fn[100]; int i; sprintf(fn, "screen%03d.ppm", dump++); FILE *fp = fopen(fn, "w"); fprintf(fp, "P6\n%d %d\n255\n", width, height); for (i = height - 1; i >= 0; i--) fwrite(tex[i], 1, width * 3, fp); fclose(fp); printf("%s written\n", fn); }
void keypress(unsigned char key, int x, int y) { switch(key) { case 'q': glFinish(); glutDestroyWindow(gwin); return; case 27: scale = 1./256; cx = -.6; cy = 0; break;
case 'r': color_rotate = (color_rotate + 1) % 6; break;
case '>': case '.': max_iter += 128; if (max_iter > 1 << 15) max_iter = 1 << 15; printf("max iter: %d\n", max_iter); break;
case '<': case ',': max_iter -= 128; if (max_iter < 128) max_iter = 128; printf("max iter: %d\n", max_iter); break;
case 'c': saturation = 1 - saturation; break;
case 's': screen_dump(); return; case 'z': max_iter = 4096; break; case 'x': max_iter = 128; break; case ' ': invert = !invert; } set_texture(); }
void hsv_to_rgb(int hue, int min, int max, rgb_t *p) { if (min == max) max = min + 1; if (invert) hue = max - (hue - min); if (!saturation) { p->r = p->g = p->b = 255 * (max - hue) / (max - min); return; } double h = fmod(color_rotate + 1e-4 + 4.0 * (hue - min) / (max - min), 6);
- define VAL 255
double c = VAL * saturation; double X = c * (1 - fabs(fmod(h, 2) - 1));
p->r = p->g = p->b = 0;
switch((int)h) { case 0: p->r = c; p->g = X; return; case 1: p->r = X; p->g = c; return; case 2: p->g = c; p->b = X; return; case 3: p->g = X; p->b = c; return; case 4: p->r = X; p->b = c; return; default:p->r = c; p->b = X; } }
void calc_mandel() { int i, j, iter, min, max; rgb_t *px; double x, y, zx, zy, zx2, zy2; min = max_iter; max = 0; for (i = 0; i < height; i++) { px = tex[i]; y = (i - height/2) * scale + cy; for (j = 0; j < width; j++, px++) { x = (j - width/2) * scale + cx; iter = 0;
zx = hypot(x - .25, y); if (x < zx - 2 * zx * zx + .25) iter = max_iter; if ((x + 1)*(x + 1) + y * y < 1/16) iter = max_iter;
zx = zy = zx2 = zy2 = 0; for (; iter < max_iter && zx2 + zy2 < 4; iter++) { zy = 2 * zx * zy + y; zx = zx2 - zy2 + x; zx2 = zx * zx; zy2 = zy * zy; } if (iter < min) min = iter; if (iter > max) max = iter; *(unsigned short *)px = iter; } }
for (i = 0; i < height; i++) for (j = 0, px = tex[i]; j < width; j++, px++) hsv_to_rgb(*(unsigned short*)px, min, max, px); }
void alloc_tex() { int i, ow = tex_w, oh = tex_h;
for (tex_w = 1; tex_w < width; tex_w <<= 1); for (tex_h = 1; tex_h < height; tex_h <<= 1);
if (tex_h != oh || tex_w != ow) tex = realloc(tex, tex_h * tex_w * 3 + tex_h * sizeof(rgb_t*));
for (tex[0] = (rgb_t *)(tex + tex_h), i = 1; i < tex_h; i++) tex[i] = tex[i - 1] + tex_w; }
void set_texture() { alloc_tex(); calc_mandel();
glEnable(GL_TEXTURE_2D); glBindTexture(GL_TEXTURE_2D, texture); glTexImage2D(GL_TEXTURE_2D, 0, 3, tex_w, tex_h, 0, GL_RGB, GL_UNSIGNED_BYTE, tex[0]);
glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST); glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST); render(); }
void mouseclick(int button, int state, int x, int y) { if (state != GLUT_UP) return;
cx += (x - width / 2) * scale; cy -= (y - height/ 2) * scale;
switch(button) { case GLUT_LEFT_BUTTON: /* zoom in */ if (scale > fabs(x) * 1e-16 && scale > fabs(y) * 1e-16) scale /= 2; break; case GLUT_RIGHT_BUTTON: /* zoom out */ scale *= 2; break; /* any other button recenters */ } set_texture(); }
void resize(int w, int h)
{
printf("resize %d %d\n", w, h);
width = w;
height = h;
glViewport(0, 0, w, h); glOrtho(0, w, 0, h, -1, 1);
set_texture(); }
void init_gfx(int *c, char **v) { glutInit(c, v); glutInitDisplayMode(GLUT_RGB); glutInitWindowSize(640, 480); glutDisplayFunc(render);
gwin = glutCreateWindow("Mandelbrot");
glutKeyboardFunc(keypress); glutMouseFunc(mouseclick); glutReshapeFunc(resize); glGenTextures(1, &texture); set_texture(); }
int main(int c, char **v) { init_gfx(&c, v); printf("keys:\n\tr: color rotation\n\tc: monochrome\n\ts: screen dump\n\t" "<, >: decrease/increase max iteration\n\tq: quit\n\tmouse buttons to zoom\n");
glutMainLoop(); return 0; }</lang>
C++
This generic function assumes that the image can be accessed like a two-dimensional array of colors. It may be passed a true array (in which case the Mandelbrot set will simply be drawn into that array, which then might be saved as image file), or a class which maps the subscript operator to the pixel drawing routine of some graphics library. In the latter case, there must be functions get_first_dimension and get_second_dimension defined for that type, to be found by argument dependent lookup. The code provides those functions for built-in arrays. <lang cpp>#include <cstdlib>
- include <complex>
// get dimensions for arrays template<typename ElementType, std::size_t dim1, std::size_t dim2>
std::size_t get_first_dimension(ElementType (&a)[dim1][dim2])
{
return dim1;
}
template<typename ElementType, std::size_t dim1, std::size_t dim2>
std::size_t get_second_dimension(ElementType (&a)[dim1][dim2])
{
return dim2;
}
template<typename ColorType, typename ImageType>
void draw_Mandelbrot(ImageType& image, //where to draw the image ColorType set_color, ColorType non_set_color, //which colors to use for set/non-set points double cxmin, double cxmax, double cymin, double cymax,//the rect to draw in the complex plane unsigned int max_iterations) //the maximum number of iterations
{
std::size_t const ixsize = get_first_dimension(ImageType); std::size_t const iysize = get_first_dimension(ImageType); for (std::size_t ix = 0; ix < ixsize; ++ix) for (std::size_t iy = 0; iy < iysize; ++iy) { std::complex<double> c(cxmin + ix/(ixsize-1.0)*(cxmax-cxmin), cymin + iy/(iysize-1.0)*(cymax-cymin)); std::complex<double> z = 0; unsigned int iterations;
for (iterations = 0; iterations < max_iterations && std::abs(z) < 2.0; ++iterations) z = z*z + c;
image[ix][iy] = (iterations == max_iterations) ? set_color : non_set_color;
}
}</lang>
C#
<lang csharp>using System; using System.Drawing; using System.Drawing.Imaging; using System.Threading; using System.Windows.Forms;
/// <summary> /// Generates bitmap of Mandelbrot Set and display it on the form. /// </summary> public class MandelbrotSetForm : Form {
const double MaxValueExtent = 2.0; Thread thread;
static double CalcMandelbrotSetColor(ComplexNumber c) { // from http://en.wikipedia.org/w/index.php?title=Mandelbrot_set const int MaxIterations = 1000; const double MaxNorm = MaxValueExtent * MaxValueExtent;
int iteration = 0; ComplexNumber z = new ComplexNumber(); do { z = z * z + c; iteration++; } while (z.Norm() < MaxNorm && iteration < MaxIterations); if (iteration < MaxIterations) return (double)iteration / MaxIterations; else return 0; // black }
static void GenerateBitmap(Bitmap bitmap) { double scale = 2 * MaxValueExtent / Math.Min(bitmap.Width, bitmap.Height); for (int i = 0; i < bitmap.Height; i++) { double y = (bitmap.Height / 2 - i) * scale; for (int j = 0; j < bitmap.Width; j++) { double x = (j - bitmap.Width / 2) * scale; double color = CalcMandelbrotSetColor(new ComplexNumber(x, y)); bitmap.SetPixel(j, i, GetColor(color)); } } }
static Color GetColor(double value) { const double MaxColor = 256; const double ContrastValue = 0.2; return Color.FromArgb(0, 0, (int)(MaxColor * Math.Pow(value, ContrastValue))); } public MandelbrotSetForm() { // form creation this.Text = "Mandelbrot Set Drawing"; this.BackColor = System.Drawing.Color.Black; this.BackgroundImageLayout = System.Windows.Forms.ImageLayout.Stretch; this.MaximizeBox = false; this.StartPosition = FormStartPosition.CenterScreen; this.FormBorderStyle = FormBorderStyle.FixedDialog; this.ClientSize = new Size(640, 640); this.Load += new System.EventHandler(this.MainForm_Load); }
void MainForm_Load(object sender, EventArgs e) { thread = new Thread(thread_Proc); thread.IsBackground = true; thread.Start(this.ClientSize); }
void thread_Proc(object args) { // start from small image to provide instant display for user Size size = (Size)args; int width = 16; while (width * 2 < size.Width) { int height = width * size.Height / size.Width; Bitmap bitmap = new Bitmap(width, height, PixelFormat.Format24bppRgb); GenerateBitmap(bitmap); this.BeginInvoke(new SetNewBitmapDelegate(SetNewBitmap), bitmap); width *= 2; Thread.Sleep(200); } // then generate final image Bitmap finalBitmap = new Bitmap(size.Width, size.Height, PixelFormat.Format24bppRgb); GenerateBitmap(finalBitmap); this.BeginInvoke(new SetNewBitmapDelegate(SetNewBitmap), finalBitmap); }
void SetNewBitmap(Bitmap image) { if (this.BackgroundImage != null) this.BackgroundImage.Dispose(); this.BackgroundImage = image; }
delegate void SetNewBitmapDelegate(Bitmap image);
static void Main() { Application.Run(new MandelbrotSetForm()); }
}
struct ComplexNumber {
public double Re; public double Im;
public ComplexNumber(double re, double im) { this.Re = re; this.Im = im; }
public static ComplexNumber operator +(ComplexNumber x, ComplexNumber y) { return new ComplexNumber(x.Re + y.Re, x.Im + y.Im); }
public static ComplexNumber operator *(ComplexNumber x, ComplexNumber y) { return new ComplexNumber(x.Re * y.Re - x.Im * y.Im, x.Re * y.Im + x.Im * y.Re); }
public double Norm() { return Re * Re + Im * Im; }
}</lang>
Clojure
<lang lisp>(ns mandelbrot
(:refer-clojure :exclude [+ * <]) (:use (clojure.contrib complex-numbers) (clojure.contrib.generic [arithmetic :only [+ *]] [comparison :only [<]] [math-functions :only [abs]])))
(defn mandelbrot? [z]
(loop [c 1 m (iterate #(+ z (* % %)) 0)] (if (and (> 20 c) (< (abs (first m)) 2) ) (recur (inc c) (rest m)) (if (= 20 c) true false))))
(defn mandelbrot []
(for [y (range 1 -1 -0.05)
x (range -2 0.5 0.0315)]
(if (mandelbrot? (complex x y)) "#" " ")))
(println (interpose \newline (map #(apply str %) (partition 80 (mandelbrot))))) </lang>
D
Textual Version
This uses std.complex.Complex
because D built-in complex numbers are deprecated.
<lang d>import std.stdio, std.complex;
void main() {
enum maxIter = 1_000; foreach (y; -39 .. 39) { foreach (x; -39 .. 39) { auto c = Complex!double(y/40.0 - 0.5, x/40.0), z = Complex!double(0.0, 0.0), i = 0; for (; i < maxIter && z.abs() < 4; i++) z = z ^^ 2 + c; write(i == maxIter ? '#' : ' '); } writeln(); }
}</lang>
- Output:
# # # # # ### ##### ##### ### # ######### ############# ############### ##################### ##################### ################### ################### ################### ################### ####################### ################### ################### ##################### ################### ################### ################# ############### ############# ######### # ############### ####################### # ######################### # ############################# # ############################### # ################################# ################################### ####################################### ### ######################################### ### ################################################# ############################################### ############################################# ############################################# ############################################### ############################################### ################################################### ################################################# ################################################# ################################################### ################################################### # ################################################### # ##### ################################################### ##### ###### ################################################### ###### ####### ################################################### ####### ####################################################################### ######### ################################################### ######### ###### ################################################### ###### ##### ################################################### ##### ################################################### ################################################### ################################################### ################################################### ################################################# ################################################# ################################################### ############################################### ############################################### ########################################### ######################################### ############################################# #### ################## ################## #### ### ################ ################ ### # ############## ############## # ########### ########### ## ##### ##### ## # # # #
Graphical version
<lang d>import qd;
double lensqr(cdouble c) { return c.re * c.re + c.im * c.im; }
const Limit = 150;
void main() {
screen(640, 480); for (int y = 0; y < screen.h; ++y) { flip; events; for (int x = 0; x < screen.w; ++x) { auto c_x = x * 1.0 / screen.w - 0.5, c_y = y * 1.0 / screen.h - 0.5, c = c_y * 2.0i + c_x * 3.0 - 1.0, z = 0.0i + 0.0, i = 0; for (; i < Limit; ++i) { z = z * z + c; if (lensqr(z) > 4) break; } auto value = cast(ubyte) (i * 255.0 / Limit); pset(x, y, rgb(value, value, value)); } } while (true) { flip; events; }
}</lang>
Dart
Implementation in Google Dart works on http://try.dartlang.org/ (as of 10/18/2011) since the language is very new, it may break in the future. The implementation uses a incomplete Complex class supporting operator overloading. <lang>class Complex {
double _r,_i;
Complex(this._r,this._i); double get r() => _r; double get i() => _i; String toString() => "($r,$i)";
Complex operator +(Complex other) => new Complex(r+other.r,i+other.i); Complex operator *(Complex other) => new Complex(r*other.r-i*other.i,r*other.i+other.r*i); double abs() => r*r+i*i;
}
void main() {
double start_x=-1.5; double start_y=-1.0; double step_x=0.03; double step_y=0.1;
for(int y=0;y<20;y++) { String line=""; for(int x=0;x<70;x++) { Complex c=new Complex(start_x+step_x*x,start_y+step_y*y); Complex z=new Complex(0,0); for(int i=0;i<100;i++) { z=z*(z)+c; if(z.abs()>2) { break; } } line+=z.abs()>2 ? " " : "*"; } print(line); }
}</lang>
DWScript
<lang delphi>const maxIter = 256;
var x, y, i : Integer; for y:=-39 to 39 do begin
for x:=-39 to 39 do begin var c := Complex(y/40-0.5, x/40); var z := Complex(0, 0); for i:=1 to maxIter do begin z := z*z + c; if Abs(z)>=4 then Break; end; if i>=maxIter then Print('#') else Print('.'); end; PrintLn();
end;</lang>
F#
<lang fsharp>open System.Drawing open System.Windows.Forms type Complex =
{ re : float; im : float }
let cplus (x:Complex) (y:Complex) : Complex =
{ re = x.re + y.re; im = x.im + y.im }
let cmult (x:Complex) (y:Complex) : Complex =
{ re = x.re * y.re - x.im * y.im; im = x.re * y.im + x.im * y.re; }
let norm (x:Complex) : float =
x.re*x.re + x.im*x.im
type Mandel = class
inherit Form static member xPixels = 500 static member yPixels = 500 val mutable bmp : Bitmap member x.mandelbrot xMin xMax yMin yMax maxIter = let rec mandelbrotIterator z c n = if (norm z) > 2.0 then false else match n with | 0 -> true | n -> let z' = cplus ( cmult z z ) c in mandelbrotIterator z' c (n-1) let dx = (xMax - xMin) / (float (Mandel.xPixels)) let dy = (yMax - yMin) / (float (Mandel.yPixels)) in for xi = 0 to Mandel.xPixels-1 do for yi = 0 to Mandel.yPixels-1 do let c = {re = xMin + (dx * float(xi) ) ; im = yMin + (dy * float(yi) )} in if (mandelbrotIterator {re=0.;im=0.;} c maxIter) then x.bmp.SetPixel(xi,yi,Color.Azure) else x.bmp.SetPixel(xi,yi,Color.Black) done done
member public x.generate () = x.mandelbrot (-1.5) 0.5 (-1.0) 1.0 200 ; x.Refresh()
new() as x = {bmp = new Bitmap(Mandel.xPixels , Mandel.yPixels)} then x.Text <- "Mandelbrot set" ; x.Width <- Mandel.xPixels ; x.Height <- Mandel.yPixels ; x.BackgroundImage <- x.bmp; x.generate(); x.Show();
end
let f = new Mandel() do Application.Run(f)</lang>
Forth
This uses grayscale image utilities. <lang Forth>500 value max-iter
- mandel ( gmp F: imin imax rmin rmax -- )
0e 0e { F: imin F: imax F: rmin F: rmax F: Zr F: Zi } dup bheight 0 do i s>f dup bheight s>f f/ imax imin f- f* imin f+ TO Zi dup bwidth 0 do i s>f dup bwidth s>f f/ rmax rmin f- f* rmin f+ TO Zr Zr Zi max-iter begin 1- dup while fover fdup f* fover fdup f* fover fover f+ 4e f< while f- Zr f+ frot frot f* 2e f* Zi f+ repeat fdrop fdrop drop 0 \ for a pretty grayscale image, replace with: 255 max-iter */ else drop 255 then fdrop fdrop over i j rot g! loop loop drop ;
80 24 graymap dup -1e 1e -2e 1e mandel dup gshow free bye</lang>
Fortran
<lang fortran>program mandelbrot
implicit none integer , parameter :: rk = selected_real_kind (9, 99) integer , parameter :: i_max = 800 integer , parameter :: j_max = 600 integer , parameter :: n_max = 100 real (rk), parameter :: x_centre = -0.5_rk real (rk), parameter :: y_centre = 0.0_rk real (rk), parameter :: width = 4.0_rk real (rk), parameter :: height = 3.0_rk real (rk), parameter :: dx_di = width / i_max real (rk), parameter :: dy_dj = -height / j_max real (rk), parameter :: x_offset = x_centre - 0.5_rk * (i_max + 1) * dx_di real (rk), parameter :: y_offset = y_centre - 0.5_rk * (j_max + 1) * dy_dj integer, dimension (i_max, j_max) :: image integer :: i integer :: j integer :: n real (rk) :: x real (rk) :: y real (rk) :: x_0 real (rk) :: y_0 real (rk) :: x_sqr real (rk) :: y_sqr
do j = 1, j_max y_0 = y_offset + dy_dj * j do i = 1, i_max x_0 = x_offset + dx_di * i x = 0.0_rk y = 0.0_rk n = 0 do x_sqr = x ** 2 y_sqr = y ** 2 if (x_sqr + y_sqr > 4.0_rk) then image (i, j) = 255 exit end if if (n == n_max) then image (i, j) = 0 exit end if y = y_0 + 2.0_rk * x * y x = x_0 + x_sqr - y_sqr n = n + 1 end do end do end do open (10, file = 'out.pgm') write (10, '(a/ i0, 1x, i0/ i0)') 'P2', i_max, j_max, 255 write (10, '(i0)') image close (10)
end program mandelbrot</lang>
gnuplot
The output from gnuplot is controlled by setting the appropriate values for the options terminal
and output
.
<lang gnuplot>set terminal png
set output 'mandelbrot.png'</lang>
The following script draws an image of the number of iterations it takes to escape the circle with radius rmax
with a maximum of nmax
.
<lang gnuplot>rmax = 2
nmax = 100
complex (x, y) = x * {1, 0} + y * {0, 1}
mandelbrot (z, z0, n) = n == nmax || abs (z) > rmax ? n : mandelbrot (z ** 2 + z0, z0, n + 1)
set samples 200
set isosamples 200
set pm3d map
set size square
splot [-2 : .8] [-1.4 : 1.4] mandelbrot (complex (0, 0), complex (x, y), 0) notitle</lang>Output:
Go
Text only, prints an 80-char by 41-line depiction. <lang go>package main
import "fmt" import "math/cmplx"
func mandelbrot(a complex128) (z complex128) {
for i := 0; i < 50; i++ { z = z*z + a } return
}
func main() {
for y := 1.0; y >= -1.0; y -= 0.05 { for x := -2.0; x <= 0.5; x += 0.0315 { if cmplx.Abs(mandelbrot(complex(x, y))) < 2 { fmt.Print("*") } else { fmt.Print(" ") } } fmt.Println("") }
}</lang>
Haskell
<lang haskell>import Data.Complex
mandelbrot a = iterate (\z -> z^2 + a) 0 !! 50
main = mapM_ putStrLn [[if magnitude (mandelbrot (x :+ y)) < 2 then '*' else ' '
| x <- [-2, -1.9685 .. 0.5]] | y <- [1, 0.95 .. -1]]</lang>
Save the code to file m.hs and run :
runhaskell m.hs
haXe
This version compiles for flash version 9 or greater. The compilation command is <lang haxe>haxe -swf9 mandelbrot.swf -main Mandelbrot</lang>
<lang haxe>class Mandelbrot extends flash.display.Sprite {
inline static var MAX_ITER = 255;
public static function main() { var w = flash.Lib.current.stage.stageWidth; var h = flash.Lib.current.stage.stageHeight; var mandelbrot = new Mandelbrot(w, h); flash.Lib.current.stage.addChild(mandelbrot); mandelbrot.drawMandelbrot(); }
var image:flash.display.BitmapData; public function new(width, height) { super(); var bitmap:flash.display.Bitmap; image = new flash.display.BitmapData(width, height, false); bitmap = new flash.display.Bitmap(image); this.addChild(bitmap); }
public function drawMandelbrot() { image.lock(); var step_x = 3.0 / (image.width-1); var step_y = 2.0 / (image.height-1); for (i in 0...image.height) { var ci = i * step_y - 1.0; for (j in 0...image.width) { var k = 0; var zr = 0.0; var zi = 0.0; var cr = j * step_x - 2.0; while (k <= MAX_ITER && (zr*zr + zi*zi) <= 4) { var temp = zr*zr - zi*zi + cr; zi = 2*zr*zi + ci; zr = temp; k ++; } paint(j, i, k); } } image.unlock(); }
inline function paint(x, y, iter) { var color = iter > MAX_ITER? 0 : iter * 0x100; image.setPixel(x, y, color); }
}</lang>
Icon and Unicon
<lang Icon>link graphics
procedure main()
width := 750 height := 600 limit := 100 WOpen("size="||width||","||height) every x:=1 to width & y:=1 to height do { z:=complex(0,0) c:=complex(2.5*x/width-2.0,(2.0*y/height-1.0)) j:=0 while j<limit & cAbs(z)<2.0 do { z := cAdd(cMul(z,z),c) j+:= 1 } Fg(mColor(j,limit)) DrawPoint(x,y) } WriteImage("./mandelbrot.gif") WDone()
end
procedure mColor(x,limit)
max_color := 2^16-1 color := integer(max_color*(real(x)/limit))
return(if x=limit then "black" else color||","||color||",0")
end
record complex(r,i)
procedure cAdd(x,y)
return complex(x.r+y.r,x.i+y.i)
end
procedure cMul(x,y)
return complex(x.r*y.r-x.i*y.i,x.r*y.i+x.i*y.r)
end
procedure cAbs(x)
return sqrt(x.r*x.r+x.i*x.i)
end</lang>
IDL
IDL - Interactive Data Language (free implementation: GDL - GNU Data Language http://gnudatalanguage.sourceforge.net) <lang IDL> PRO Mandelbrot,xRange,yRange,xPixels,yPixels,iterations
xPixelstartVec = Lindgen( xPixels) * Float(xRange[1]-xRange[0]) / $
xPixels + xRange[0]
yPixelstartVec = Lindgen( yPixels) * Float(YRANGE[1]-yrange[0])$
/ yPixels + yRange[0]
constArr = Complex( Rebin( xPixelstartVec, xPixels, yPixels),$
Rebin( Transpose(yPixelstartVec), xPixels, yPixels))
valArr = ComplexArr( xPixels, yPixels)
res = IntArr( xPixels, yPixels)
oriIndex = Lindgen( Long(xPixels) * yPixels)
FOR i = 0, iterations-1 DO BEGIN ; only one loop needed
; calculation for whole array at once valArr = valArr^2 - constArr
whereIn = Where( Abs( valArr) LE 4.0d, COMPLEMENT=whereOut)
IF whereIn[0] EQ -1 THEN BREAK
valArr = valArr[ whereIn]
constArr = constArr[ whereIn]
IF whereOut[0] NE -1 THEN BEGIN
res[ oriIndex[ whereOut]] = i+1
oriIndex = oriIndex[ whereIn] ENDIF
ENDFOR
tv,res ; open a window and show the result
END
Mandelbrot,[-1.,2.3],[-1.3,1.3],640,512,200
END
</lang> from the command line: <lang IDL> GDL>.run mandelbrot </lang> or <lang IDL> GDL> Mandelbrot,[-1.,2.3],[-1.3,1.3],640,512,200 </lang>
Inform 7
<lang inform7>"Mandelbrot"
The story headline is "A Non-Interactive Set".
Include Glimmr Drawing Commands by Erik Temple.
[Q20 fixed-point or floating-point: see definitions below] Use floating-point math.
Finished is a room.
The graphics-window is a graphics g-window spawned by the main-window. The position is g-placeabove.
When play begins: let f10 be 10 as float; now min re is ( -20 as float ) fdiv f10; now max re is ( 6 as float ) fdiv f10; now min im is ( -12 as float ) fdiv f10; now max im is ( 12 as float ) fdiv f10; now max iterations is 100; add color g-Black to the palette; add color g-Red to the palette; add hex "#FFA500" to the palette; add color g-Yellow to the palette; add color g-Green to the palette; add color g-Blue to the palette; add hex "#4B0082" to the palette; add hex "#EE82EE" to the palette; open up the graphics-window.
Min Re is a number that varies. Max Re is a number that varies. Min Im is a number that varies. Max Im is a number that varies.
Max Iterations is a number that varies.
Min X is a number that varies. Max X is a number that varies. Min Y is a number that varies. Max Y is a number that varies.
The palette is a list of numbers that varies.
[vertically mirrored version] Window-drawing rule for the graphics-window when max im is fneg min im: clear the graphics-window; let point be { 0, 0 }; now min X is 0 as float; now min Y is 0 as float; let mX be the width of the graphics-window minus 1; let mY be the height of the graphics-window minus 1; now max X is mX as float; now max Y is mY as float; let L be the column order with max mX; repeat with X running through L: now entry 1 in point is X; repeat with Y running from 0 to mY / 2: now entry 2 in point is Y; let the scaled point be the complex number corresponding to the point; let V be the Mandelbrot result for the scaled point; let C be the color corresponding to V; if C is 0, next; draw a rectangle (C) in the graphics-window at the point with size 1 by 1; now entry 2 in point is mY - Y; draw a rectangle (C) in the graphics-window at the point with size 1 by 1; yield to VM; rule succeeds.
[slower non-mirrored version] Window-drawing rule for the graphics-window: clear the graphics-window; let point be { 0, 0 }; now min X is 0 as float; now min Y is 0 as float; let mX be the width of the graphics-window minus 1; let mY be the height of the graphics-window minus 1; now max X is mX as float; now max Y is mY as float; let L be the column order with max mX; repeat with X running through L: now entry 1 in point is X; repeat with Y running from 0 to mY: now entry 2 in point is Y; let the scaled point be the complex number corresponding to the point; let V be the Mandelbrot result for the scaled point; let C be the color corresponding to V; if C is 0, next; draw a rectangle (C) in the graphics-window at the point with size 1 by 1; yield to VM; rule succeeds.
To decide which list of numbers is column order with max (N - number): let L be a list of numbers; let L2 be a list of numbers; let D be 64; let rev be false; while D > 0: let X be 0; truncate L2 to 0 entries; while X <= N: if D is 64 or X / D is odd, add X to L2; increase X by D; if rev is true: reverse L2; let rev be false; otherwise: let rev be true; add L2 to L; let D be D / 2; decide on L.
To decide which list of numbers is complex number corresponding to (P - list of numbers): let R be a list of numbers; extend R to 2 entries; let X be entry 1 in P as float; let X be (max re fsub min re) fmul (X fdiv max X); let X be X fadd min re; let Y be entry 2 in P as float; let Y be (max im fsub min im) fmul (Y fdiv max Y); let Y be Y fadd min im; now entry 1 in R is X; now entry 2 in R is Y; decide on R.
To decide which number is Mandelbrot result for (P - list of numbers): let c_re be entry 1 in P; let c_im be entry 2 in P; let z_re be 0 as float; let z_im be z_re; let threshold be 4 as float; let runs be 0; while 1 is 1: [ z = z * z ] let r2 be z_re fmul z_re; let i2 be z_im fmul z_im; let ri be z_re fmul z_im; let z_re be r2 fsub i2; let z_im be ri fadd ri; [ z = z + c ] let z_re be z_re fadd c_re; let z_im be z_im fadd c_im; let norm be (z_re fmul z_re) fadd (z_im fmul z_im); increase runs by 1; if norm is greater than threshold, decide on runs; if runs is max iterations, decide on 0.
To decide which number is color corresponding to (V - number): let L be the number of entries in the palette; let N be the remainder after dividing V by L; decide on entry (N + 1) in the palette.
Section - Fractional numbers (for Glulx only)
To decide which number is (N - number) as float: (- (numtof({N})) -). To decide which number is (N - number) fadd (M - number): (- (fadd({N}, {M})) -). To decide which number is (N - number) fsub (M - number): (- (fsub({N}, {M})) -). To decide which number is (N - number) fmul (M - number): (- (fmul({N}, {M})) -). To decide which number is (N - number) fdiv (M - number): (- (fdiv({N}, {M})) -). To decide which number is fneg (N - number): (- (fneg({N})) -). To yield to VM: (- glk_select_poll(gg_event); -).
Use Q20 fixed-point math translates as (- Constant Q20_MATH; -). Use floating-point math translates as (- Constant FLOAT_MATH; -).
Include (-
- ifdef Q20_MATH;
! Q11.20 format: 1 sign bit, 11 integer bits, 20 fraction bits [ numtof n r; @shiftl n 20 r; return r; ]; [ fadd n m; return n+m; ]; [ fsub n m; return n-m; ]; [ fmul n m; n = n + $$1000000000; @sshiftr n 10 n; m = m + $$1000000000; @sshiftr m 10 m; return n * m; ]; [ fdiv n m; @sshiftr m 20 m; return n / m; ]; [ fneg n; return -n; ];
- endif;
- ifdef FLOAT_MATH;
[ numtof f; @"S2:400" f f; return f; ]; [ fadd n m; @"S3:416" n m n; return n; ]; [ fsub n m; @"S3:417" n m n; return n; ]; [ fmul n m; @"S3:418" n m n; return n; ]; [ fdiv n m; @"S3:419" n m n; return n; ]; [ fneg n; @bitxor n $80000000 n; return n; ];
- endif;
-).</lang>
Newer Glulx interpreters provide 32-bit floating-point operations, but this solution also supports fixed-point math which is more widely supported and accurate enough for a zoomed-out view. Inform 6 inclusions are used for the low-level math functions in either case. The rendering process is extremely slow, since the graphics system is not optimized for pixel-by-pixel drawing, so this solution includes an optimization for vertical symmetry (as in the default view) and also includes extra logic to draw the lines in a more immediately useful order.
J
The characteristic function of the Mandelbrot can be defined as follows: <lang j>mcf=. (<: 2:)@|@(] ((*:@] + [)^:((<: 2:)@|@])^:1000) 0:) NB. 1000 iterations test</lang> The Mandelbrot set can be drawn as follows: <lang j>domain=. |.@|:@({.@[ + ] *~ j./&i.&>/@+.@(1j1 + ] %~ -~/@[))&>/
load'graph' viewmat mcf "0 @ domain (_2j_1 1j1) ; 0.01 NB. Complex interval and resolution</lang>
A smaller version, based on a black&white implementation of viewmat (and paraphrased, from html markup to wiki markup), is shown here:
<lang j> viewmat mcf "0 @ domain (_2j_1 1j1) ; 0.1 NB. Complex interval and resolution</lang>
The output is HTML-heavy and can be found here (split out to make editing this page easier).
Java
<lang java>import java.awt.Graphics; import java.awt.image.BufferedImage; import javax.swing.JFrame;
public class Mandelbrot extends JFrame {
private final int MAX_ITER = 570; private final double ZOOM = 150; private BufferedImage I; private double zx, zy, cX, cY, tmp;
public Mandelbrot() { super("Mandelbrot Set"); setBounds(100, 100, 800, 600); setResizable(false); setDefaultCloseOperation(EXIT_ON_CLOSE); I = new BufferedImage(getWidth(), getHeight(), BufferedImage.TYPE_INT_RGB); for (int y = 0; y < getHeight(); y++) { for (int x = 0; x < getWidth(); x++) { zx = zy = 0; cX = (x - 400) / ZOOM; cY = (y - 300) / ZOOM; int iter = MAX_ITER; while (zx * zx + zy * zy < 4 && iter > 0) { tmp = zx * zx - zy * zy + cX; zy = 2.0 * zx * zy + cY; zx = tmp; iter--; } I.setRGB(x, y, iter | (iter << 8)); } } }
@Override public void paint(Graphics g) { g.drawImage(I, 0, 0, this); }
public static void main(String[] args) { new Mandelbrot().setVisible(true); }
}</lang>
JavaScript
This needs the canvas tag of HTML 5 (it will not run on IE8 and lower or old browsers). <lang javascript> function Mandeliter(cx, cy, maxiter) {
var i; var x = 0.0; var y = 0.0; for (i = 0; i < maxiter && x*x + y*y <= 4; ++i) { var tmp = 2*x*y; x = x*x - y*y + cx; y = tmp + cy; } return i;
}
function Mandelbrot() {
var width = 900; var height = 600; var cd = document.getElementById('calcdata'); var xmin = parseFloat(cd.xmin.value); var xmax = parseFloat(cd.xmax.value); var ymin = parseFloat(cd.ymin.value); var ymax = parseFloat(cd.ymax.value); var iterations = parseInt(cd.iterations.value); var ctx = document.getElementById('mandelimage').getContext("2d"); var img = ctx.getImageData(0, 0, width, height); var pix = img.data; for (var ix = 0; ix < width; ++ix) for (var iy = 0; iy < height; ++iy) { var x = xmin + (xmax-xmin)*ix/(width-1); var y = ymin + (ymax-ymin)*iy/(height-1); var i = Mandeliter(x, y, iterations); var ppos = 4*(900*iy + ix); if (i == iterations) { pix[ppos] = 0; pix[ppos+1] = 0; pix[ppos+2] = 0; } else { var c = 3*Math.log(i)/Math.log(iterations - 1.0); if (c < 1) { pix[ppos] = 255*c; pix[ppos+1] = 0; pix[ppos+2] = 0; } else if (c < 2) { pix[ppos] = 255; pix[ppos+1] = 255*(c-1); pix[ppos+2] = 0; } else { pix[ppos] = 255; pix[ppos+1] = 255; pix[ppos+2] = 255*(c-2); } } pix[ppos+3] = 255; } ctx.putImageData(img,0,0);
}</lang> The corresponding HTML: <lang html> <!DOCTYPE html> <html>
<head> <title>Mandelbrot set</title> <script src="Mandelbrot.js" type="text/javascript"></script> </head>
<body onload="Mandelbrot()">
Mandelbrot set
<form id="calcdata" onsubmit="javascript:Mandelbrot(); return false;">
xmin = | <input name="xmin" type="text" size="10" value="-2"> | xmax = | <input name="xmax" type="text" size="10" value="1"> |
ymin = | <input name="ymin" type="text" size="10" value="-1"> | ymax = | <input name="ymax" type="text" size="10" value="1"> |
iterations = <input name="iterations" type="text" size="10" value="1000">
<input type="submit" value=" Calculate "> <input type="reset" value=" Reset form ">
</form>
<canvas id="mandelimage" width="900" height="600"> This page needs a browser with canvas support. </canvas> </body>
</html> </lang>
Output with default parameters:
LabVIEW
Liberty BASIC
Any words of description go outside of lang tags. <lang lb> nomainwin
WindowWidth =440 WindowHeight =460
open "Mandelbrot Set" for graphics_nsb_nf as #w
- w "trapclose [quit]"
- w "down"
for x0 = -2 to 1 step .0033
for y0 = -1.5 to 1.5 step .0075 x = 0 y = 0
iteration = 0 maxIteration = 255
while ( ( x *x +y *y) <=4) and ( iteration <maxIteration) xtemp =x *x -y *y +x0 y =2 *x *y +y0 x = xtemp iteration = iteration + 1 wend
if iteration <>maxIteration then c =iteration else c =0 end if
call pSet x0, y0, c scan next
next
- w "flush"
wait
sub pSet x, y, c
xScreen = 10 +( x +2) /3 *400 yScreen = 10 +( y +1.5) /3 *400 if c =0 then col$ ="red" else if c mod 2 =1 then col$ ="lightgray" else col$ ="white" end if #w "color "; col$ #w "set "; xScreen; " "; yScreen
end sub
[quit] close #w end </lang>
Logo
<lang logo>to mandelbrot :left :bottom :side :size
cs setpensize [1 1] make "inc :side/:size make "zr :left repeat :size [ make "zr :zr + :inc make "zi :bottom pu setxy repcount - :size/2 minus :size/2 pd repeat :size [ make "zi :zi + :inc setpencolor count.color calc :zr :zi fd 1 ] ]
end
to count.color :count
;op (list :count :count :count) if :count > 256 [op 0] ; black if :count > 128 [op 7] ; white if :count > 64 [op 5] ; magenta if :count > 32 [op 6] ; yellow if :count > 16 [op 4] ; red if :count > 8 [op 2] ; green if :count > 4 [op 1] ; blue op 3 ; cyan
end
to calc :zr :zi [:count 0] [:az 0] [:bz 0]
if :az*:az + :bz*:bz > 4 [op :count] if :count > 256 [op :count] op (calc :zr :zi (:count + 1) (:zr + :az*:az - :bz*:bz) (:zi + 2*:az*:bz))
end
mandelbrot -2 -1.25 2.5 400</lang>
Mathematica
The implementation could be better. But this is a start... <lang mathematica>eTime[z0_, maxIter_Integer: 100] := (Length@NestWhileList[(# + z0)^2 &, 0, (Abs@# <= 2) &, 1, maxIter]) - 1
DistributeDefinitions[eTime]; mesh = ParallelTable[eTime[(x + I*y), 1000], {y, 1.2, -1.2, -0.01}, {x, -1.72, 1, 0.01}]; ReliefPlot[mesh, Frame -> False]</lang>
Mathmap
filter mandelbrot (gradient coloration) c=ri:(xy/xy:[X,X]*1.5-xy:[0.5,0]); z=ri:[0,0]; # initial value z0 = 0 # iteration of z iter=0; while abs(z)<2 && iter<31 do z=z*z+c; # z(n+1) = fc(zn) iter=iter+1 end; coloration(iter/32) # color of pixel end
MATLAB
This solution uses the escape time algorithm to determine the coloring of the coordinates on the complex plane. The code can be reduced to a single line via vectorization after the Escape Time Algorithm function definition, but the code becomes unnecessarily obfuscated. Also, this code uses a lot of memory. You will need a computer with a lot of memory to compute the set with high resolution.
<lang MATLAB>function [theSet,realAxis,imaginaryAxis] = mandelbrotSet(start,gridSpacing,last,maxIteration)
%Define the escape time algorithm function escapeTime = escapeTimeAlgorithm(z0) escapeTime = 0; z = 0; while( (abs(z)<=2) && (escapeTime < maxIteration) ) z = (z + z0)^2; escapeTime = escapeTime + 1; end end %Define the imaginary axis imaginaryAxis = (imag(start):imag(gridSpacing):imag(last)); %Define the real axis realAxis = (real(start):real(gridSpacing):real(last)); %Construct the complex plane from the real and imaginary axes complexPlane = meshgrid(realAxis,imaginaryAxis) + meshgrid(imaginaryAxis(end:-1:1),realAxis)'.*i; %Apply the escape time algorithm to each point in the complex plane theSet = arrayfun(@escapeTimeAlgorithm, complexPlane);
%Draw the set pcolor(realAxis,imaginaryAxis,theSet); shading flat;
end</lang>
To use this function you must specify the:
- lower left hand corner of the complex plane from which to start the image,
- the grid spacing in both the imaginary and real directions,
- the upper right hand corner of the complex plane at which to end the image and
- the maximum iterations for the escape time algorithm.
For example:
- Lower Left Corner: -2.05-1.2i
- Grid Spacing: 0.004+0.0004i
- Upper Right Corner: 0.45+1.2i
- Maximum Iterations: 500
Sample usage: <lang MATLAB>mandelbrotSet(-2.05-1.2i,0.004+0.0004i,0.45+1.2i,500);</lang>
Modula-3
<lang modula3>MODULE Mandelbrot EXPORTS Main;
IMPORT Wr, Stdio, Fmt, Word;
CONST m = 50;
limit2 = 4.0;
TYPE UByte = BITS 8 FOR [0..16_FF];
VAR width := 200;
height := 200; bitnum: CARDINAL := 0; byteacc: UByte := 0; isOverLimit: BOOLEAN; Zr, Zi, Cr, Ci, Tr, Ti: REAL;
BEGIN
Wr.PutText(Stdio.stdout, "P4\n" & Fmt.Int(width) & " " & Fmt.Int(height) & "\n");
FOR y := 0 TO height - 1 DO FOR x := 0 TO width - 1 DO Zr := 0.0; Zi := 0.0; Cr := 2.0 * FLOAT(x) / FLOAT(width) - 1.5; Ci := 2.0 * FLOAT(y) / FLOAT(height) - 1.0; FOR i := 1 TO m + 1 DO Tr := Zr*Zr - Zi*Zi + Cr; Ti := 2.0*Zr*Zi + Ci; Zr := Tr; Zi := Ti; isOverLimit := Zr*Zr + Zi*Zi > limit2; IF isOverLimit THEN EXIT; END; END; IF isOverLimit THEN byteacc := Word.Xor(Word.LeftShift(byteacc, 1), 16_00); ELSE byteacc := Word.Xor(Word.LeftShift(byteacc, 1), 16_01); END;
INC(bitnum); IF bitnum = 8 THEN Wr.PutChar(Stdio.stdout, VAL(byteacc, CHAR)); byteacc := 0; bitnum := 0; ELSIF x = width - 1 THEN byteacc := Word.LeftShift(byteacc, 8 - (width MOD 8)); Wr.PutChar(Stdio.stdout, VAL(byteacc, CHAR)); byteacc := 0; bitnum := 0 END; Wr.Flush(Stdio.stdout); END; END;
END Mandelbrot.</lang>
OCaml
<lang ocaml>#load "graphics.cma";;
let mandelbrot xMin xMax yMin yMax xPixels yPixels maxIter =
let rec mandelbrotIterator z c n = if (Complex.norm z) > 2.0 then false else match n with | 0 -> true | n -> let z' = Complex.add (Complex.mul z z) c in mandelbrotIterator z' c (n-1) in Graphics.open_graph (" "^(string_of_int xPixels)^"x"^(string_of_int yPixels)); let dx = (xMax -. xMin) /. (float_of_int xPixels) and dy = (yMax -. yMin) /. (float_of_int yPixels) in for xi = 0 to xPixels - 1 do for yi = 0 to yPixels - 1 do let c = {Complex.re = xMin +. (dx *. float_of_int xi); Complex.im = yMin +. (dy *. float_of_int yi)} in if (mandelbrotIterator Complex.zero c maxIter) then (Graphics.set_color Graphics.white; Graphics.plot xi yi ) else (Graphics.set_color Graphics.black; Graphics.plot xi yi ) done done;;
mandelbrot (-1.5) 0.5 (-1.0) 1.0 500 500 200;;</lang>
Octave
This code runs rather slowly and produces coloured Mandelbrot set by accident (output image).
<lang octave>#! /usr/bin/octave -qf global width = 200; global height = 200; maxiter = 100;
z0 = 0; global cmax = 1 + i; global cmin = -2 - i;
function cs = pscale(c)
global cmax; global cmin; global width; global height; persistent px = (real(cmax-cmin))/width; persistent py = (imag(cmax-cmin))/height; cs = real(cmin) + px*real(c) + i*(imag(cmin) + py*imag(c));
endfunction
ms = zeros(width, height); for x = 0:width-1
for y = 0:height-1 z0 = 0; c = pscale(x+y*i); for ic = 1:maxiter z1 = z0^2 + c; if ( abs(z1) > 2 ) break; endif z0 = z1; endfor ms(x+1, y+1) = ic/maxiter; endfor
endfor
saveimage("mandel.ppm", round(ms .* 255).', "ppm");</lang>
Perl
translation / optimization of the ruby solution <lang perl>use Math::Complex;
sub mandelbrot {
my ($z, $c) = @_[0,0]; for (1 .. 20) { $z = $z * $z + $c; return $_ if abs $z > 2; }
}
for (my $y = 1; $y >= -1; $y -= 0.05) {
for (my $x = -2; $x <= 0.5; $x += 0.0315) {print mandelbrot($x + $y * i) ? ' ' : '#'} print "\n"
}</lang>
Perl 6
Variant of a Mandelbrot script from the Perl 6 ecosystem. Produces a Portable Pixel Map to STDOUT. Redirect into a file to save it. Converted to a .png file for display here.
<lang perl6>my $height = @*ARGS[0] // 31; $height = $height % 2 ?? +$height !! 1+$height; my $width = $height; my $max_iterations = 50;
my $upper-right = -2 + 5/4i; my $lower-left = 1/2 - 5/4i;
my @color_map = (
"0 0 0", "0 0 252", "64 0 252", "124 0 252", "188 0 252", "252 0 252", "252 0 188", "252 0 124", "252 0 64", "252 0 0", "252 64 0", "252 124 0", "252 188 0", "252 252 0", "188 252 0", "124 252 0", "64 252 0", "0 252 0", "0 252 64", "0 252 124", "0 252 188", "0 252 252", "0 188 252", "0 124 252", "0 64 252", "124 124 252", "156 124 252", "188 124 252", "220 124 252", "252 124 252", "252 124 220", "252 124 188", "252 124 156", "252 124 124", "252 156 124", "252 188 124", "252 220 124", "252 252 124", "220 252 124", "188 252 124", "156 252 124", "124 252 124", "124 252 156", "124 252 188", "124 252 220", "124 252 252", "124 220 252", "124 188 252", "124 156 252", "180 180 252", "196 180 252", "216 180 252", "232 180 252", "252 180 252", "252 180 232", "252 180 216", "252 180 196", "252 180 180", "252 196 180", "252 216 180", "252 232 180", "252 252 180", "232 252 180", "216 252 180", "196 252 180", "180 252 180", "180 252 196", "180 252 216", "180 252 232", "180 252 252", "180 232 252", "180 216 252", "180 196 252", "0 0 112", "28 0 112", "56 0 112", "84 0 112", "112 0 112", "112 0 84", "112 0 56", "112 0 28", "112 0 0", "112 28 0", "112 56 0", "112 84 0", "112 112 0", "84 112 0", "56 112 0", "28 112 0", "0 112 0", "0 112 28", "0 112 56", "0 112 84", "0 112 112", "0 84 112", "0 56 112", "0 28 112", "56 56 112", "68 56 112", "84 56 112", "96 56 112", "112 56 112", "112 56 96", "112 56 84", "112 56 68", "112 56 56", "112 68 56", "112 84 56", "112 96 56", "112 112 56", "96 112 56", "84 112 56", "68 112 56", "56 112 56", "56 112 68", "56 112 84", "56 112 96", "56 112 112", "56 96 112", "56 84 112", "56 68 112", "80 80 112", "88 80 112", "96 80 112", "104 80 112", "112 80 112", "112 80 104", "112 80 96", "112 80 88", "112 80 80", "112 88 80", "112 96 80", "112 104 80", "112 112 80", "104 112 80", "96 112 80", "88 112 80", "80 112 80", "80 112 88", "80 112 96", "80 112 104", "80 112 112", "80 104 112", "80 96 112", "80 88 112", "0 0 64", "16 0 64", "32 0 64", "48 0 64", "64 0 64", "64 0 48", "64 0 32", "64 0 16", "64 0 0", "64 16 0", "64 32 0", "64 48 0", "64 64 0", "48 64 0", "32 64 0", "16 64 0", "0 64 0", "0 64 16", "0 64 32", "0 64 48", "0 64 64", "0 48 64", "0 32 64", "0 16 64", "32 32 64", "40 32 64", "48 32 64", "56 32 64", "64 32 64", "64 32 56", "64 32 48", "64 32 40", "64 32 32", "64 40 32", "64 48 32", "64 56 32", "64 64 32", "56 64 32", "48 64 32", "40 64 32", "32 64 32", "32 64 40", "32 64 48", "32 64 56", "32 64 64", "32 56 64", "32 48 64", "32 40 64", "44 44 64", "48 44 64", "52 44 64", "60 44 64", "64 44 64", "64 44 60", "64 44 52", "64 44 48", "64 44 44", "64 48 44", "64 52 44", "64 60 44", "64 64 44", "60 64 44", "52 64 44", "48 64 44", "44 64 44", "44 64 48", "44 64 52", "44 64 60", "44 64 64", "44 60 64", "44 52 64", "44 48 64"
);
sub mandel (Complex $c) {
my $z = 0i; for ^$max_iterations -> $i { return $i + 1 if $z.abs > 2; $z = $z * $z + $c; } return 0;
}
sub subdivide($low, $high, $count) {
(^$count).map: { $low + ($_ / ($count - 1)) * ($high - $low) }
}
say "P3"; say "$width $height"; say "255";
for subdivide($upper-right.re, $lower-left.re, $height) -> $re {
my @line = subdivide($re + ($upper-right.im)i, $re + 0i, ($width + 1) / 2).map: { mandel($_) } my $middle = @line.pop; say ~(@line, $middle, @line.reverse).map: { @color_map[$_] }
}</lang>
PHP
<lang PHP>$min_x=-2; $max_x=1; $min_y=-1; $max_y=1;
$dim_x=400; $dim_y=300;
$im = @imagecreate($dim_x, $dim_y)
or die("Cannot Initialize new GD image stream");
header("Content-Type: image/png"); $black_color = imagecolorallocate($im, 0, 0, 0); $white_color = imagecolorallocate($im, 255, 255, 255);
for($y=0;$y<=$dim_y;$y++) {
for($x=0;$x<=$dim_x;$x++) { $c1=$min_x+($max_x-$min_x)/$dim_x*$x; $c2=$min_y+($max_y-$min_y)/$dim_y*$y;
$z1=0; $z2=0;
for($i=0;$i<100;$i++) { $new1=$z1*$z1-$z2*$z2+$c1; $new2=2*$z1*$z2+$c2; $z1=$new1; $z2=$new2; if($z1*$z1+$z2*$z2>=4) { break; } } if($i<100) { imagesetpixel ($im, $x, $y, $white_color); } }
}
imagepng($im); imagedestroy($im); </lang>
PicoLisp
<lang PicoLisp>(scl 6)
(let Ppm (make (do 300 (link (need 400))))
(for (Y . Row) Ppm (for (X . @) Row (let (ZX 0 ZY 0 CX (*/ (- X 250) 1.0 150) CY (*/ (- Y 150) 1.0 150) C 570) (while (and (> 4.0 (+ (*/ ZX ZX 1.0) (*/ ZY ZY 1.0))) (gt0 C)) (let Tmp (- (*/ ZX ZX 1.0) (*/ ZY ZY 1.0) (- CX)) (setq ZY (+ (*/ 2 ZX ZY 1.0) CY) ZX Tmp ) ) (dec 'C) ) (set (nth Ppm Y X) (list 0 C C)) ) ) ) (out "img.ppm" (prinl "P6") (prinl 400 " " 300) (prinl 255) (for Y Ppm (for X Y (apply wr X))) ) )</lang>
PostScript
<lang postscript>%!PS-Adobe-2.0 %%BoundingBox: 0 0 300 200 %%EndComments /origstate save def /ld {load def} bind def /m /moveto ld /g /setgray ld /dot { currentpoint 1 0 360 arc fill } bind def %%EndProlog % param /maxiter 200 def % complex manipulation /complex { 2 array astore } def /real { 0 get } def /imag { 1 get } def /cmul { /a exch def /b exch def
a real b real mul a imag b imag mul sub a real b imag mul a imag b real mul add 2 array astore
} def /cadd { aload pop 3 -1 roll aload pop
3 -1 roll add 3 1 roll add exch 2 array astore
} def /cconj { aload pop neg 2 array astore } def /cabs2 { dup cconj cmul 0 get} def % mandel 200 100 translate -200 1 100 { /x exch def
-100 1 100 { /y exch def /z0 0.0 0.0 complex def 0 1 maxiter { /iter exch def
x 100 div y 100 div complex z0 z0 cmul cadd dup /z0 exch def cabs2 4 gt {exit} if
} for iter maxiter div g x y m dot } for
} for % showpage origstate restore %%EOF</lang>
Prolog
SWI-Prolog has a graphic interface XPCE : <lang Prolog>:- use_module(library(pce)).
mandelbrot :-
new(D, window('Mandelbrot Set')), send(D, size, size(700, 650)), new(Img, image(@nil, width := 700, height := 650, kind := pixmap)),
forall(between(0,699, I), ( forall(between(0,649, J), ( get_RGB(I, J, R, G, B), R1 is (R * 256) mod 65536, G1 is (G * 256) mod 65536, B1 is (B * 256) mod 65536, send(Img, pixel(I, J, colour(@default, R1, G1, B1))))))), new(Bmp, bitmap(Img)), send(D, display, Bmp, point(0,0)), send(D, open).
get_RGB(X, Y, R, G, B) :-
CX is (X - 350) / 150, CY is (Y - 325) / 150, Iter = 570, compute_RGB(CX, CY, 0, 0, Iter, It), IterF is It \/ It << 15, R is IterF >> 16, Iter1 is IterF - R << 16, G is Iter1 >> 8, B is Iter1 - G << 8.
compute_RGB(CX, CY, ZX, ZY, Iter, IterF) :-
ZX * ZX + ZY * ZY < 4, Iter > 0, !, Tmp is ZX * ZX - ZY * ZY + CX, ZY1 is 2 * ZX * ZY + CY, Iter1 is Iter - 1, compute_RGB(CX, CY, Tmp, ZY1, Iter1, IterF).
compute_RGB(_CX, _CY, _ZX, _ZY, Iter, Iter).</lang>Example :
PureBasic
PureBasic forum: discussion <lang PureBasic>EnableExplicit
- Window1 = 0
- Image1 = 0
- ImgGadget = 0
- max_iteration = 64
- width = 800
- height = 600
Define.d x0 ,y0 ,xtemp ,cr, ci Define.i i, n, x, y ,Event ,color
Dim Color.l (255) For n = 0 To 63
Color( 0 + n ) = RGB( n*4+128, 4 * n, 0 ) Color( 64 + n ) = RGB( 64, 255, 4 * n ) Color( 128 + n ) = RGB( 64, 255 - 4 * n , 255 ) Color( 192 + n ) = RGB( 64, 0, 255 - 4 * n )
Next
If OpenWindow(#Window1, 0, 0, #width, #height, "'Mandelbrot set' PureBasic Example", #PB_Window_SystemMenu )
If CreateImage(#Image1, #width, #height) ImageGadget(#ImgGadget, 0, 0, #width, #height, ImageID(#Image1)) For y.i = 1 To #height -1 StartDrawing(ImageOutput(#Image1)) For x.i = 1 To #width -1 x0 = 0 y0 = 0; cr = (x / #width)*2.5 -2 ci = (y / #height)*2.5 -1.25 i = 0 While (x0*x0 + y0*y0 <= 4.0) And i < #max_iteration i +1 xtemp = x0*x0 - y0*y0 + cr y0 = 2*x0*y0 + ci x0 = xtemp Wend If i >= #max_iteration Plot(x, y, 0 ) Else Plot(x, y, Color(i & 255)) EndIf Next StopDrawing() SetGadgetState(#ImgGadget, ImageID(#Image1)) Repeat Event = WindowEvent() If Event = #PB_Event_CloseWindow End EndIf Until Event = 0 Next EndIf Repeat Event = WaitWindowEvent() Until Event = #PB_Event_CloseWindow EndIf</lang>Example:
Python
Translation of the ruby solution <lang python># Python 3.0+ and 2.5+ try:
from functools import reduce
except:
pass
def mandelbrot(a): return reduce(lambda z, _: z*z + a, range(50), 0)
def step(start, step, iterations): return (start + (i * step) for i in range(iterations))
rows = (('*' if abs(mandelbrot(complex(x, y))) < 2 else ' '
for x in step(-2.0, .0315, 80)) for y in step(1, -.05, 41))
print( '\n'.join(.join(row) for row in rows) ) </lang>
R
<lang R>createmandel <- function(zfrom, zto, inc=0.01, maxiter=100) {
x <- c() y <- c() cost <- zfrom while( Re(cost) <= Re(zto) ) { cost <- Re(cost) + Im(zfrom)*1i while ( Im(cost) <= Im(zto) ) { j <- 0 z <- 0 while( (abs(z) < 2) && (j < maxiter)) { z <- z^2 + cost j <- j + 1 } if ( j == maxiter ) { x <- c(x, Re(cost)) y <- c(y, Im(cost)) } cost <- cost + inc*1i } cost <- cost + inc } plot(x,y, pch=".")
}
createmandel(-2-1i, 1+1i)</lang>
Racket
<lang scheme>#lang racket
(require racket/draw)
(define (iterations a z i)
(define z′ (+ (* z z) a)) (if (or (= i 255) (> (magnitude z′) 2)) i (iterations a z′ (add1 i))))
(define (iter->color i)
(if (= i 255) (make-object color% "black") (make-object color% (* 5 (modulo i 15)) (* 32 (modulo i 7)) (* 8 (modulo i 31)))))
(define (mandelbrot width height)
(define target (make-bitmap width height)) (define dc (new bitmap-dc% [bitmap target])) (for* ([x width] [y height]) (define real-x (- (* 3.0 (/ x width)) 2.25)) (define real-y (- (* 2.5 (/ y height)) 1.25)) (send dc set-pen (iter->color (iterations (make-rectangular real-x real-y) 0 0)) 1 'solid) (send dc draw-point x y)) (send target save-file "mandelbrot.png" 'png))
(mandelbrot 300 200)</lang>
Ruby
Text only, prints an 80-char by 41-line depiction. Found here. <lang ruby>require 'complex'
def mandelbrot(a)
Array.new(50).inject(0) { |z,c| z*z + a }
end
(1.0).step(-1,-0.05) do |y|
(-2.0).step(0.5,0.0315) do |x| print mandelbrot(Complex(x,y)).abs < 2 ? '*' : ' ' end puts
end</lang>
Uses the text progress bar from Median filter#Ruby <lang ruby>class RGBColour
def self.mandel_colour(i) self.new( 16*(i % 15), 32*(i % 7), 8*(i % 31) ) end
end
class Pixmap
def self.mandelbrot(width, height) mandel = Pixmap.new(width,height) pb = ProgressBar.new(width) if $DEBUG width.times do |x| height.times do |y| x_ish = Float(x - width*11/15) / (width/3) y_ish = Float(y - height/2) / (height*3/10) mandel[x,y] = RGBColour.mandel_colour(mandel_iters(x_ish, y_ish)) end pb.update(x) if $DEBUG end pb.close if $DEBUG mandel end
def self.mandel_iters(cx,cy) x = y = 0.0 count = 0 while Math.hypot(x,y) < 2 and count < 255 x, y = (x**2 - y**2 + cx), (2*x*y + cy) count += 1 end count end
end
Pixmap.mandelbrot(300,300).save('mandel.ppm')</lang>
Scala
Uses RgbBitmap from Basic Bitmap Storage task and Complex number class from this programming task. <lang scala>import rosettacode.ArithmeticComplex._ import java.awt.Color
object Mandelbrot {
def generate(width:Int =600, height:Int =400)={ val bm=new RgbBitmap(width, height)
val maxIter=1000 val xMin = -2.0 val xMax = 1.0 val yMin = -1.0 val yMax = 1.0
val cx=(xMax-xMin)/width val cy=(yMax-yMin)/height
for(y <- 0 until bm.height; x <- 0 until bm.width){ val c=Complex(xMin+x*cx, yMin+y*cy) val iter=itMandel(c, maxIter, 4) bm.setPixel(x, y, getColor(iter, maxIter)) } bm }
def itMandel(c:Complex, imax:Int, bailout:Int):Int={ var z=Complex() for(i <- 0 until imax){ z=z*z+c; if(z.abs > bailout) return i } imax; }
def getColor(iter:Int, max:Int):Color={ if (iter==max) return Color.BLACK
var c=3*math.log(iter)/math.log(max-1.0) if(c<1) new Color((255*c).toInt, 0, 0) else if(c<2) new Color(255, (255*(c-1)).toInt, 0) else new Color(255, 255, (255*(c-2)).toInt) }
}</lang> <lang scala>val imgMandel=Mandelbrot.generate() val mainframe=new MainFrame(){title="Test"; visible=true
contents=new Label(){icon=new ImageIcon(imgMandel.image)}
}</lang>
Scheme
This implementation writes an image of the Mandelbrot set to a plain pgm file. The set itself is drawn in white, while the exterior is drawn in black. <lang scheme>(define x-centre -0.5) (define y-centre 0.0) (define width 4.0) (define i-max 800) (define j-max 600) (define n 100) (define r-max 2.0) (define file "out.pgm") (define colour-max 255) (define pixel-size (/ width i-max)) (define x-offset (- x-centre (* 0.5 pixel-size (+ i-max 1)))) (define y-offset (+ y-centre (* 0.5 pixel-size (+ j-max 1))))
(define (inside? z)
(define (*inside? z-0 z n) (and (< (magnitude z) r-max) (or (= n 0) (*inside? z-0 (+ (* z z) z-0) (- n 1))))) (*inside? z 0 n))
(define (boolean->integer b)
(if b colour-max 0))
(define (pixel i j)
(boolean->integer (inside? (make-rectangular (+ x-offset (* pixel-size i)) (- y-offset (* pixel-size j))))))
(define (plot)
(with-output-to-file file (lambda () (begin (display "P2") (newline) (display i-max) (newline) (display j-max) (newline) (display colour-max) (newline) (do ((j 1 (+ j 1))) ((> j j-max)) (do ((i 1 (+ i 1))) ((> i i-max)) (begin (display (pixel i j)) (newline))))))))
(plot)</lang>
Scratch
Seed7
$ include "seed7_05.s7i"; include "float.s7i"; include "complex.s7i"; include "draw.s7i"; include "keybd.s7i"; # Display the Mandelbrot set, that are points z[0] in the complex plane # for which the sequence z[n+1] := z[n] ** 2 + z[0] (n >= 0) is bounded. # Since this program is computing intensive it should be compiled with # hi comp -O2 mandelbr const integer: pix is 200; const integer: max_iter is 256; var array color: colorTable is max_iter times black; const func integer: iterate (in complex: z0) is func result var integer: iter is 1; local var complex: z is complex.value; begin z := z0; while sqrAbs(z) < 4.0 and # not diverged iter < max_iter do # not converged z *:= z; z +:= z0; incr(iter); end while; end func; const proc: displayMandelbrotSet (in complex: center, in float: zoom) is func local var integer: x is 0; var integer: y is 0; var complex: z0 is complex.value; begin for x range -pix to pix do for y range -pix to pix do z0 := center + complex(flt(x) * zoom, flt(y) * zoom); point(x + pix, y + pix, colorTable[iterate(z0)]); end for; end for; end func; const proc: main is func local const integer: num_pix is 2 * pix + 1; var integer: col is 0; begin screen(num_pix, num_pix); clear(curr_win, black); KEYBOARD := GRAPH_KEYBOARD; for col range 1 to pred(max_iter) do colorTable[col] := color(65535 - (col * 5003) mod 65535, (col * 257) mod 65535, (col * 2609) mod 65535); end for; displayMandelbrotSet(complex(-0.75, 0.0), 1.3 / flt(pix)); DRAW_FLUSH; readln(KEYBOARD); end func;
Original source: [1]
Tcl
This code makes extensive use of Tk's built-in photo image system, which provides a 32-bit RGBA plotting surface that can be then quickly drawn in any number of places in the application. It uses a computational color scheme that was easy to code... <lang tcl>package require Tk
proc mandelIters {cx cy} {
set x [set y 0.0] for {set count 0} {hypot($x,$y) < 2 && $count < 255} {incr count} { set x1 [expr {$x*$x - $y*$y + $cx}] set y1 [expr {2*$x*$y + $cy}] set x $x1; set y $y1 } return $count
} proc mandelColor {iter} {
set r [expr {16*($iter % 15)}] set g [expr {32*($iter % 7)}] set b [expr {8*($iter % 31)}] format "#%02x%02x%02x" $r $g $b
} image create photo mandel -width 300 -height 300
- Build picture in strips, updating as we go so we have "progress" monitoring
- Also set the cursor to tell the user to wait while we work.
pack [label .mandel -image mandel -cursor watch] update for {set x 0} {$x < 300} {incr x} {
for {set y 0} {$y < 300} {incr y} { set i [mandelIters [expr {($x-220)/100.}] [expr {($y-150)/90.}]] mandel put [mandelColor $i] -to $x $y } update
} .mandel configure -cursor {}</lang>
TI-83 BASIC
Based on the BASIC Version. Due to the TI-83's lack of power, it takes around 2 hours to complete at 16 iterations. <lang ti83b>PROGRAM:MANDELBR
- Input "ITER. ",D
- For(A,Xmin,Xmax,ΔX)
- For(B,Ymin,Ymax,ΔY)
- 0→X
- 0→Y
- 0→I
- D→M
- While X^2+Y^2≤4 and I<M
- X^2-Y^2+A→R
- 2XY+B→Y
- R→X
- I+1→I
- End
- If I≠M
- Then
- I→C
- Else
- 0→C
- End
- If C<1
- Pt-On(A,B)
- End
- End
- End
</lang>
XSLT
The fact that you can create an image of the Mandelbrot Set with XSLT is sometimes under-appreciated. However, it has been discussed extensively on the internet so is best reproduced here, and the code can be executed directly in your browser at that site.
<lang xml> <?xml version="1.0" encoding="UTF-8"?> <xsl:stylesheet version="1.0" xmlns:xsl="http://www.w3.org/1999/XSL/Transform">
<xsl:output method="html" indent="no"
doctype-public="-//W3C//DTD HTML 4.01//EN" doctype-system="http://www.w3.org/TR/REC-html40/strict.dtd" />
<xsl:template match="/fractal">
<html> <head> <title>XSLT fractal</title> <style type="text/css">
body { color:#55F; background:#000 } pre { font-family:monospace; font-size:7px } pre span { background:<xsl:value-of select="background" /> }
</style> </head> <body>
Copyright © 1992,2007 Joel Yliluoma (<a href="http://iki.fi/bisqwit/">http://iki.fi/bisqwit/</a>)
XSLT fractal
<xsl:call-template name="bisqwit-mandelbrot" />
</body> </html>
</xsl:template>
<xsl:template name="bisqwit-mandelbrot"
><xsl:call-template name="bisqwit-mandelbrot-line"> <xsl:with-param name="y" select="y/min"/> </xsl:call-template
></xsl:template>
<xsl:template name="bisqwit-mandelbrot-line"
><xsl:param name="y" /><xsl:call-template name="bisqwit-mandelbrot-column"> <xsl:with-param name="x" select="x/min"/> <xsl:with-param name="y" select="$y"/> </xsl:call-template ><xsl:if test="$y < y/max" >
<xsl:call-template name="bisqwit-mandelbrot-line"> <xsl:with-param name="y" select="$y + y/step"/> </xsl:call-template ></xsl:if
></xsl:template>
<xsl:template name="bisqwit-mandelbrot-column"
><xsl:param name="x" /><xsl:param name="y" /><xsl:call-template name="bisqwit-mandelbrot-slot"> <xsl:with-param name="x" select="$x" /> <xsl:with-param name="y" select="$y" /> <xsl:with-param name="zr" select="$x" /> <xsl:with-param name="zi" select="$y" /> </xsl:call-template ><xsl:if test="$x < x/max" ><xsl:call-template name="bisqwit-mandelbrot-column"> <xsl:with-param name="x" select="$x + x/step"/> <xsl:with-param name="y" select="$y" /> </xsl:call-template ></xsl:if
></xsl:template>
<xsl:template name="bisqwit-mandelbrot-slot" ><xsl:param name="x"
/><xsl:param name="y"
/><xsl:param name="zr"
/><xsl:param name="zi"
/><xsl:param name="iter" select="0"
/><xsl:variable name="zrsqr" select="($zr * $zr)"
/><xsl:variable name="zisqr" select="($zi * $zi)"
/><xsl:choose>
<xsl:when test="(4*scale*scale >= $zrsqr + $zisqr) and (maxiter > $iter+1)"
><xsl:call-template name="bisqwit-mandelbrot-slot">
<xsl:with-param name="x" select="$x" />
<xsl:with-param name="y" select="$y" />
<xsl:with-param name="zi" select="(2 * $zr * $zi) div scale + $y" />
<xsl:with-param name="zr" select="($zrsqr - $zisqr) div scale + $x" />
<xsl:with-param name="iter" select="$iter + 1" />
</xsl:call-template
></xsl:when>
<xsl:otherwise
><xsl:variable name="magnitude" select="magnitude[@value=$iter]"
/><xsl:value-of select="$magnitude/symbol"
/></xsl:otherwise>
</xsl:choose
></xsl:template>
</xsl:stylesheet> </lang>
- WikipediaSourced
- Raster graphics operations
- Programming Tasks
- Fractals
- Ada
- Lumen
- ALGOL 68
- Applesoft BASIC
- AutoHotkey
- AWK
- BASIC
- Brace
- C
- GLUT
- C++
- C sharp
- Clojure
- D
- QD
- SDL
- Phobos
- Dart
- DWScript
- F Sharp
- Forth
- Fortran
- Gnuplot
- Go
- Haskell
- HaXe
- Icon
- Unicon
- Icon Programming Library
- Unicon examples needing attention
- Examples needing attention
- IDL
- Inform 7
- Glimmr Drawing Commands by Erik Temple
- J
- Java
- Swing
- AWT
- JavaScript
- LabVIEW
- Liberty BASIC
- Logo
- Mathematica
- MATLAB
- Modula-3
- OCaml
- Octave
- Perl
- Perl 6
- PHP
- GD Graphics Library
- PicoLisp
- PostScript
- Prolog
- PureBasic
- Python
- R
- Racket
- Ruby
- Scala
- Scheme
- Scratch
- Seed7
- Tcl
- Tk
- TI-83 BASIC
- XSLT
- M4/Omit
- ML/I/Omit