Skip to content

Mandelbrot Sets (2)

Vector 7:4, page 110

The text of this article is not yet online. No complete scan of this issue has been found. The only copy found, on vector.org.uk, is damaged: its page images are lost, but the text machine-read from them survives, shown below (help find a scan).

Unedited OCR text, machine-read from a scan now lost: it contains errors, and its APL is wrong. Shown for searching.

VECTOR Vol.7 No.4

Mandelbrot Sets

by Ray Cannon

Let me say first of all, that this article has very little to do with APL. I did write some Mandelbrot code in APL for an IBM PC compatible, but it took too long to calculate even the simplest areas.

The reason it takes so long is image 480 x 480 pixels in size, requires 1,200,000,000 odd floating points instructions.

Question: What is the Mandelbrot set? Answer: It is the set of all the connected Julian sets. (Ask a silly question and...)

OK, think of a number, square it, add the number you originally thought of (2,4,6), square it, add first number again (36,38), square it , add the first number (1444,1446)... Nothing much strange there, the result rapidly rises. Try -1 (- 1,1,0,0,-1,1,0...), or zero (0,0,0,0...). All numbers seem to fall into these 2 groups. Hither the number shoots off towards infinity, or the result stabilises cycle (which may only have one element).

Now try it with complex numbers, plotting each starting point on the complex plane (x/y graph) as either black (cycles) or white (shoots off). The result is a Mandelbrot diagram. (To speed things up, please note that ALL points outside the circle of radius 2 starting at the origin 0,0 are white, so only points INSIDE this circle need be considered.)

Thus the Mandelbrot set can be described via the transformation of z goes to z squared plus c, where ¢ is a

To produce a picture of part of the Mandelbrot set, follow the following steps:

1) Decide on the number of (pixels) points in your picture (say 480x480) and map these onto a the complex plane” within the region -2.0 to +0.5 on the real (x) axis, and -1.25 to +1.25 on the imaginary (y) axis.

2) Each pixel represents a point within a small section on the complex plane. Loop for each point in turn, taking its co-ordinates as the value for c , witha starting value for z

3) Por each ce , calculate the value of z1 by squaring z and adding c.

4) Reset z to the value of z1 and repeat step 3).

110

VECTOR Vol.7 No.4

5) Continue doing steps 3) and 4) until either the absolute value of z is greater than 2, or the loop count is greater than an arbitrary value (say 256).

6) colour of the pixel dependent on the loop count. Plot the point in black if it repeated, or reached the maximum value chosen.

7) Go on to the next point. See

the function called MANDOT for a specification in APL of steps 3, 4 and 5. (I have written MANDO? in a very simple manner which may become apparent later.)

One method of speeding things up is to utilise the fact areas in any one picture contain the same value, and as the MANDELBROT set is connected there are no local maxima or minima”. (Please don’t ask me to prove that statement.)

From this it became apparent to me, that if all the points on the boundary of an area have the same value, then all the points within the boundary have the SAME value as the boundary. This can be used in a kind of binary search” type algorithm.

So, define a the area to be mapped, (say 128 x 128 pixels). Calculate the values on the boundary, and if they are all the same, fill the complete area with that value.

If boundary has not the same value, divide the square into 4 smaller squares (64 x 64) and repeat the process on each of the smaller squares. Although the process can be repeated until the square is 2x2, due to the overheads involved with this process, squares of about 8 x 8 may not be worth sub-dividing.

This method reduces the number of calculations required in the BLACK area from an X squared to an X (plus) where X is the size of the square filled . (I know the square has 4 sides, not 2, but boundaries between adjacent squares are common and only need to be calculated once.)

As [have said, APL is to slow to do this work reasonable length of time. So I wrote my MANDELBROT calculator in C with the MANDOT function written in PC assembler using the floating point maths co-processor.

I did find it useful however, to model the 80387 floating point co-processor assembler code in APL. This let me write and test the basic assembler code within an APL environment. The top level function MAN can be compared with the resulting assembler listing.

111

VECTOR Vol.7 No.4

(NOTE on 80387 maths chip. The 80387 chip works on the principle of a floating

stack, and has 8 internal stack registers. Most instructions act between the top of

the stack and other stack registers. Programming it is a bit like programming a

pocket calculator in reverse polish notation” .)

The APL listing and the resulting assembler code are shown below.

Notes on Functions

Comment Functions

MANDOT Noxmal (but very simple) APL MAN Produces same result as MANDOT, calls FP fns FADD Floating point FADD machine code instruction simulator FLD Floating point FLD (load from stack) FLDM Floating point FLD (load from memory) etc.

assembler to do the MANDELBROT calculation. called from ¢ #C calling syntax man(x,y)" }xeturn mandot value for given U and V values @CODE SEGMENT CODE  MAN PROC near push bp nov bp, sp Eninit jelear stack Kor ax,ax solear AX yatack constant 4 fld q{val4] #Push 4 into stack 4 iget x,y u and v from parameter string fld q{bp+oc] : atack ye fld q[bp+04] : stack xyda fld it : stack vxyd fld 1 ? XOX Cx,cx :Clear registers xor dx,dx Initialise loop control mov bx,1 Main loop stack uyayd tpush y stack yuvay4 ¥2 ptack {y2) uvxy4 scale x2 stack x(y2) uvxy4 1x2 atack (x2) (y2) uvay4 vexchange stack xlyZ)uv(x2)y4 new ¥ value calculations

wnew X value calculations spush x2 stack feale x2-y2 sealc x2-y2+u stack

check answer and loop or exit fnstew ax

112 calc 2x stack (2x) (y2) av (xady4 teale 2xy stack (2xy) (y2)uv(x2)y4 scale ¥(y2)uv (x2) y4 jetore new Y atack (y2) uv (x2) ¥4

(x2) (y2)uv (xa) ¥4 stack (x2-y2) (y2)uv(x2)¥4 Xty2)uv(x2)¥4 store new K stack (x2) (y2) uvxya feompare x2+y2 > 4.0 reale (y2+x2) uvxY4 jcomp stack uvxYa

VECTOR Vol.7 No.4

and ax,041¢0 jz FOUND :Found to be 4 or more inc loop control cmp bx,dx est against loop max ext loop jot found go return 0

j7Found so return count

;eleax return value ignore top byte

APL Functions

Vv FADD ARGS Vv FADDP ARGS;RES (1] #30387 FADD [1] 930387 FADDP [2] ARGS+~2+S5P[01,STLARGS] [2] RES++/STCARGS] (3] SPACK[0]++/ARGS [3] STACKCARGS[0]]+RES (4] SP+s+STACK C4] STACK+14SPACK v (5) ST+8+STACK v v

FLD MEM V FLDM VAL [1] 8 30387 PLD {1] 930387 LOAD FROM MEMORY [2] STACK+ST(MEM] , STACK (2] SPACK+VAL, STACK [3] SP+atSPACK [3] ST+B4STACK v v

Vv FLDZ v FMUL ARGS [1] 30387 FLDZ [1] 930387 FMUL [2] STACK+0, STACK [2} ARGS+2+ST[01, STLARGS] [3] SP+e+STACK [3] SPACKLO]+*/ARGS v [4] ST+B+STACK v

v FNINIT v [1] 030387 FINIT [4] 030387 FSPP [2] aINIT FLOATING POINT (2] SPACKCARG]<STACK[0] [3] SPACK+.0 [3] SPACK+1+STACK [4] ST+8t+STACK [4] SP«B+STACK v v

113

VECTOR Vol.7 No.4

Cal [1] vy FSTPM NAME n30397 FSTP STR TO MEMORY 4&NAME,'«ST[O]' STACK+1+STACK ST+84+STACK

v a30387 FXCH {11

STACK (ARGS[0]3+SP[ARGS(1]] STACK LARGS(1]}+ST[ARGS(O]} ST+8+ STACK v V

R<U MAN V;MAT; COUNT; ST; STACK C13

cs] 6] {7J C8] C9] [10]

C15] (16)

£18] (19) £20] (21) £22] {23]

C25] (26) [27] [28]

[30] [34] ASIMULATE ASS MANDELBROT Reo MAT+*O COUNT+O FNINIT 4 FLIDZ FLD 0 FLDM U FLEDM ¥ LP: FLD 3 FMUL FLD 3 FMUL PXCH FADD FMUL FADD FSTP FLD 3 FSUB FADD 2 FXCH O 4 FADDP 1 0 FSUB 5 FSTPM' MAT! + (MAT20) /FOUND COUNT+COUNT+1 +(COUNT<256)/LP FOUND: R+COUNT v uadveoso nore v FSUB ARGS A30387 FSUB ARGS+” 24+ST[0],ST[LARGS] STACK(OJ+-/ARGS SI+84+ STACK v

Y RY MANDO? U;X;¥;X2;Y2;AL C17 (2) [3] £4]

c6] C7] CsJ

C10)

[13] C14] [15]

£17] 18] C19]

[22] [23]

144 ARETURN MBP VALUE FOR U V ALONG WINDED LISTING - RCOMPARE AGAINST MAN AL+t XeX24U Y¥Y2+V

EP:X2*+X2X Y2+Y2xy Y+yxX X+Y Y«¥+X YeY+V X+X2 X+X-Y¥2 X+X+0 X2*X2+Y2 X2+X2-4 +(X220)/END AL*AL+1 >+(AL>255)/END X2+X Y2ey +LP END: R+AL v

VECTOR Vol.7 No.4

TECHNICAL SECTION This section of VECTOR is aimed principally at those of our readers who already know APL. It will people with differing degrees of fluency in APL.

Contents NED: A nesting editor for APL*PLUS/PC Olle Evero Thoughts ond+f Maurice Jordan Full Screen Methods with APL2 Peter Branson A Note on the Match Function in APL Joseph L.F. De Kerf The Steam-Hammer and the Fly Gerard A Langlet 116 119 126 138

115