%%HP: T(3)A(R)F(.); DIR @ FRACTALS, by Dan Ciarniello @ @ These programs will generate the Mandelbrot set. To speed up calculations, a @ couple of properties of the Mandelbrot set are used. @ @ Since the Mandelbrot set is symmetric about the x-axis, if the x-axis is @ present, the program will generate the image on only one side and then map @ the image to the other side of the axis. @ @ The Mandelbrot set is a closed connected set. This means that if all points @ on the border are in the Mandelbrot set, then all points interior to the @ rectangle are also in the Mandelbrot set. Ditto if all points are not in the @ Mandelbrot set. Thus the program calculates all points around the border of @ a rectangle. If all points on the border are the same (in or out), then fill @ the rectangle the appropriate colour and go to the next rectangle. @ Otherwise, cut the rectangle in half and check the border of the new @ rectangle. This is called the Mariani technique. Note that loop counters @ run through pixel coordinates. These pixel coordinates are converted to user @ coordinates using PX\->C which uses PPAR to determine conversion factors. CST { MANDEL SAVE @ Custom menu. LOAD { DITH @ DITH turns dithering on \<< 3 SF @ NODITH turns dithering off \>> } { NODITH @ \<< 3 CF @ \>> } } @ LOAD @ LOAD restores a previously saved \<< 2 MENU @ image and its associated PPAR. "Enter Filename" "" @ INPUT OBJ\-> OBJ\-> @ DROP PICT STO @ 'PPAR' STO 0 MENU @ GRAPH @ \>> @ SAVE @ SAVE an image along with its PPAR. \<< PPAR PICT RCL @ PPAR is necessary so that the correct 2 \->LIST @ user units are restored when the "Enter Filename" { @ image is restored with LOAD. "" \Ga } INPUT OBJ\-> @ STO @ \>> @ MANDEL @ The main program \<< ERASE { # 0d @ Erase PICT and view it as program # 0d } PVIEW RCLF @ runs; save flag status; start TICKS \-> c1 c2 f t @ timing; save initial window corner \<< 4 CF 5 CF @ values. IF 'nITTR' @ Check whether nITTR exists. If VTYPE -1 == @ not, create it with a default value THEN 100 @ of 100. 'nITTR' STO @ END @ IF c2 IM c1 @ Correct the aspect ratio of the IM - c2 RE c1 RE - @ image. 64 * 131 / > @ THEN c2 IM @ If y2-y1 > (64/131)(x2-x1) then set c1 IM - 131 * 64 / @ x-value of upper right corner (x2) to c1 RE + c2 IM @ x1+(y2-y1)*131/64 ELSE c2 RE @ DUP c1 RE - 64 * @ If y2-y1 < (64/131)(x2-x1) then set 131 / c1 IM + @ y2 to y1+(x2-x1)*64/131 END R\->C @ 'c2' STO @ IF 'PPAR' @ VTYPE -1 == @ Check whether PPAR exists. If not, THEN { @ create it. (-2,-1) (2.09375,1) @ X 0 (100,100) TRUTH @ Y } 'PPAR' STO @ END PPAR 2 @ c2 PUT 'PPAR' STO @ Store the window corners in PPAR to PPAR 1 c1 PUT @ set the user coordinates. 'PPAR' STO @ IF c1 IM c2 @ IM * 0 < @ Determine whether the x-axis is THEN c1 RE @ present. If so determine the pixel 0 R\->C C\->PX c2 RE 0 @ coordinates of the left and right R\->C C\->PX \-> lside @ end. rside @ \<< @ IF c2 @ Determine whether the x-axis is in IM c1 IM ABS > @ the upper or lower half of the THEN 0 @ screen. If its in the lower half, 130 @ then generate the image from top of FOR x @ the screen down to the x-axis. x R\->B # 0d 2 \->LIST @ ITERATE around the first rectangle PX\->C ITERATE @ starting at the top of the screen. NEXT @ 1 rside 2 GET B\->R @ FOR y @ ITERATE down right side to x-axis. # 130d y R\->B 2 @ \->LIST PX\->C ITERATE @ NEXT @ 129 0 @ FOR x @ ITERATE along the x-axis from right x R\->B lside 2 GET 2 @ to left. \->LIST PX\->C ITERATE @ -1 @ STEP @ lside 2 GET B\->R 1 @ FOR y @ ITERATE up the left side of the # 0d y R\->B 2 \->LIST @ screen from the x-axis to the top PX\->C ITERATE -1 @ of the screen. STEP @ { # 0d # 0d } rside @ MARIANI 63 lside 2 @ Pass the upper left and lower right GET B\->R @ pixel coordinates of the rectangle FOR y @ to MARIANI. # 0d y R\->B 2 \->LIST @ PICT OVER PX\->C CONJ @ Map the area above the x-axis to C\->PX DUP 1 # 130d @ that below the x-axis line by line. PUT SUB PICT 3 @ ROLLD REPL -1 @ STEP @ ELSE 0 @ 130 @ If the x-axis is in the upper half FOR x @ of the screen, generate the image x R\->B lside 2 GET 2 @ from the x-axis to the bottom of \->LIST PX\->C ITERATE @ the screen. NEXT @ Start with the x-axis from left to lside 2 GET 1 + B\->R @ right. 63 @ FOR y @ # 130d y R\->B 2 @ Down left side to bottom of screen. \->LIST PX\->C ITERATE @ NEXT @ 129 0 @ FOR x @ x R\->B # 63d 2 \->LIST @ Along bottom of screen right to PX\->C ITERATE -1 @ left. STEP @ 62 lside 2 GET 1 + @ B\->R @ FOR y @ # 0d y R\->B 2 \->LIST @ Up right side to x-axis. PX\->C ITERATE -1 @ STEP @ lside { # 130d @ Pass upper left and lower right # 63d } MARIANI 0 @ corners of rectangle to MARIANI. lside 2 GET B\->R @ FOR y @ Map area below the x-axis to that # 0d y R\->B 2 \->LIST @ above the x-axis line by line. PICT OVER PX\->C CONJ @ C\->PX DUP 1 # 130d @ PUT SUB PICT 3 @ ROLLD REPL @ NEXT @ END @ \>> @ ELSE 1 130 @ FOR x x @ If the x-axis is not present, just R\->B # 0d 2 \->LIST @ ITERATE around the edges of the PX\->C ITERATE @ screen to start. NEXT 1 63 @ Along top. FOR y @ # 130d y R\->B 2 @ \->LIST PX\->C ITERATE @ Down right side NEXT 129 @ 0 @ FOR x x @ R\->B # 63d 2 \->LIST @ PX\->C ITERATE -1 @ Along bottom STEP 62 1 @ FOR y @ # 0d y R\->B 2 \->LIST @ PX\->C ITERATE -1 @ Up left side STEP { @ # 0d # 0d } { @ Pass upper left and lower right # 130d # 63d } @ coordinates of the screen to MARIANI @ MARIANI END TICKS t @ - B\->R 29491200 / f @ Done! Calculate execution time. STOF @ Restore flags. \>> @ \>> @ MARIANI @ @ MARIANI is a recursive routine. @ It takes a rectangle (defined by its upper left corner (ulc) and lower right @ corner (rlc), divides it in two and iterates along the dividing line. It @ then checks the border of each new rectangle in turn to determine whether all @ points on the border have the same status (either on or off). If so, fill @ the rectangle. If not, MARIANI calls itself with the coordinates of the @ smaller rectangle. This contues until a rectangle is filled or until the @ current rectangle is less than 6 pixels in area in which case the state of @ each pixel is determined individually and no further division occurs. \<< \-> ulc lrc @ \<< lrc 1 GET @ ulc 1 GET - lrc 2 @ GET ulc 2 GET - * @ IF # 6d < @ THEN ulc 1 @ If area less than 6 pixels, GET 1 + B\->R lrc 1 @ determine state of each pixel in GET 1 - B\->R @ rectangle individually. FOR x ulc @ 2 GET 1 + B\->R lrc 2 @ GET 1 - B\->R @ FOR y x @ R\->B y R\->B 2 \->LIST @ PX\->C ITERATE @ NEXT @ NEXT @ ELSE lrc 1 @ If area greater than 6, divide GET ulc 1 GET - lrc @ rectangle in two. 2 GET ulc 2 GET - @ IF < @ THEN lrc @ If the rectangle is taller than it 2 GET ulc 2 GET + 2 @ is wide, determine line which / DUP ulc 2 ROT PUT @ divides it into upper and lower SWAP lrc 2 ROT PUT @ halves. Set flag 4 to indicate 4 SF @ this. ELSE lrc @ 1 GET ulc 1 GET + 2 @ If the rectangle is wider than it / DUP ulc 1 ROT PUT @ is tall, determine the line which SWAP lrc 1 ROT PUT @ divides it into left and right END \-> @ halves. lines linet @ \<< @ lines: start coordinate of line IF 4 @ linet: end coordinate of line FS?C @ THEN @ lines 1 GET 1 + B\->R @ linet 1 GET B\->R 1 - @ If flag 4 set, draw the horizontal FOR x @ dividing line. lines 1 x R\->B PUT @ PX\->C ITERATE @ NEXT @ ELSE @ lines 2 GET 1 + B\->R @ linet 2 GET B\->R 1 - @ FOR y @ Otherwise draw the vertical one. lines 2 y R\->B PUT @ PX\->C ITERATE @ NEXT @ END ulc @ linet CHECKBORDER @ Check the border of the left (or IF 1 @ upper) rectangle (defined by FS? 2 FS? OR @ corners ulc and linet). THEN @ If either flag 1 or flag 2 is set, PICT ulc 1 GET 1 + @ then all pixels around the border ulc 2 GET 1 + 2 @ of the rectangle have the same \->LIST linet 1 GET @ status. ulc 1 GET - 1 - @ linet 2 GET ulc 2 @ GET - 1 - BLANK @ IF 2 @ FS? @ Fill black if flag 2 is set. THEN @ Fill white if flag 1 is set. NEG @ END @ REPL 1 CF 2 CF @ ELSE @ ulc linet MARIANI @ If border pixels are not all one END @ state, call MARIANI with corners lines lrc @ ulc and linet. CHECKBORDER @ IF 1 @ Check border of right (or lower) FS? 2 FS? OR @ rectangle (defined by lines and THEN @ lrc). PICT lines 1 GET 1 @ + lines 2 GET 1 + 2 @ \->LIST lrc 1 GET @ lines 1 GET - 1 - @ Fill according to flag status as lrc 2 GET lines 2 @ above. GET - 1 - BLANK @ IF 2 @ FS? @ THEN @ NEG @ END @ REPL 1 CF 2 CF @ ELSE @ lines lrc MARIANI @ Or call MARIANI with corners lines END @ and lrc. \>> @ END @ \>> @ \>> @ CHECKBORDER @ \<< \-> ulc lrc @ Determine status of rectangle \<< ulc PIX? \-> @ border. state @ First determine state of ulc pixel. \<< ulc 1 GET @ 1 + B\->R lrc 1 GET @ B\->R @ FOR x ulc @ 1 x R\->B PUT PIX? @ Now check each pixel around border IF @ starting with top line and compare state \=/ @ to state of ulc. If they are the THEN 5 @ same continue around border SF lrc 1 GET B\->R @ otherwise break out of the loop, 'x' STO @ there is no need to check further. END @ Set flag 5 to stop further NEXT @ checking. IF 5 FC? @ THEN ulc @ 2 GET 1 + B\->R lrc 2 @ If flag 5 still clear check down GET B\->R @ right side. FOR y @ lrc 2 y R\->B PUT @ PIX? @ IF @ state \=/ @ THEN @ 5 SF lrc 2 GET B\->R @ 'y' STO @ END @ NEXT @ END @ IF 5 FC? @ THEN ulc @ If flag 5 still clear, check left 2 GET 1 + B\->R lrc 2 @ side. GET B\->R @ FOR y @ ulc 2 y R\->B PUT @ PIX? @ IF @ state \=/ @ THEN @ 5 SF lrc 2 GET B\->R @ 'y' STO @ END @ NEXT @ END @ IF 5 FC? @ THEN ulc @ 1 GET 1 + B\->R lrc 1 @ If flag 5 still clear, check GET 1 - B\->R @ bottom. FOR x @ lrc 1 x R\->B PUT @ PIX? @ IF @ state \=/ @ THEN @ 5 SF lrc 1 GET B\->R @ 'x' STO @ END @ NEXT @ END @ IF 5 FC?C @ THEN @ If flag 5 clear then border all one state 1 + SF @ state. Set flag 1 for all pixels END @ off. Set flag 2 for all pixels \>> @ on. \>> @ \>> @ ITERATE @ ITERATE determines whether a \<< DUP DUP PIXON @ coordinate is in the Mandelbrot set PIXOFF C\->R @ or not. Toggle the current pixel 0 \-> x y i @ as a visual aid in following \<< 0 0 0 0 @ progress of program. WHILE SQ @ Initialize Z(0) as (0,0) and SWAP SQ SWAP DUP2 + @ duplicate. Real and imaginary 4 < 'i' INCR nITTR @ parts are treated separately on the < AND @ stack. The loop calculates REPEAT - x @ Z(n+1)=SQ(Z(n))+c where c=x+iy. + 3 ROLLD * DUP + y @ It repeats until SQ(Z)>4 or nITTR + DUP2 @ iterations are completed. END 4 DROPN @ x y R\->C i @ IF nITTR \>= @ THEN PIXON @ If nITTR iterations are completed, ELSE @ the coordinate is in the Mandelbrot IF i 20 < @ set. Turn the pixel on. 3 FS? AND @ THEN i 2 @ If dithering is on, check to see MOD @ whether the point diverged in less IF @ than 20 iterations. If so and it THEN @ is odd, turn the pixel on. PIXON @ ELSE @ PIXOFF @ END @ ELSE @ PIXOFF @ END @ END @ \>> @ \>> @ END @