This was fun for a while, but eventually I got bored of DOS programming and decided to try some bare metal stuff – and somehow, perhaps vaguely influenced by demoscene<\/a> style demos, I ended up setting myself the challenge of writing a bare-metal Mandelbrot generator within a 512-byte bootsector.<\/p>\n In fact, if we plot each iteration, we can see that it’s going around in a strange, almost spirograph-like shape: To take another example, let’s try it with the point (0.5,0.5). The iterations go:<\/p>\n So, this looks like it’s shooting off to infinity… and fast.<\/p>\n But this raises an interesting question: at what point do we decide whether the sequence is or isn’t bounded? In general we’d have to continue iterating to infinity to know for sure. Perhaps we should set a threshold, beyond which we declare it to be unbounded? But where should it be? And how do we know that sequences which exceed that threshold aren’t coming back?<\/p>\n Well, as it turns out, it can be proven<\/a> (see A6d) that if any point leaves a circle of radius 2, (ie |zn<\/sub>|>2), it will never return. So a reasonable approach is to compute a fixed number of iterations, and then check whether the result is less than 2 in absolute value.<\/p>\n Here are two simple implementations I wrote just to demonstrate the idea. First, this is for the C++-minded people:<\/p>\n This goes through each pixel one-by-one and for each one runs 64 iterations. The output is a PBM file to stdout, so you’ll have to pipe it and open it in a compatible viewer (GIMP works) or convert it to a more common format.<\/p>\n Or, if you’re more of a MATLAB fan:<\/p>\n The same idea as the C++ example, but using vectorisation rather than iterating through pixels, and using MATLAB’s internal imshow command to display the image.<\/p>\n With all that setting the scene done, let’s get down to how I wrote an 88-byte 8086 bare-metal mandelbrot generator.<\/p>\n The general structure of the C++ example can be carried over, but the arithmetic gets difficult. Firstly, the complex number datatype in C++ is an abstraction. The processor doesn’t understand complex numbers (not an 8086 at least). So we’ll have to treat them as two reals. Addition and subtraction are trivial. For multiplication we just have to remember that: This is conceptually not too difficult, but it does mean that complex operations take up a lot of program space.<\/p>\n The second thing to consider is that floating point operations don’t exist in 16-bit intel. We’re going to have to work in integer instructions instead. With a little bit of trickery though, we can create a very simple fixed-point system. Essentially, we just represent every number as a fixed multiple of its true value. Say, this could be 10. Addition (and by extension, subtraction) is consistent: Multiplication works too, with the caveat that we have to divide by the scale factor after: It’s for this reason that it’s best for that multiple to be a power of two: this reduces the division to a right shift. Note the use of sar rather than shr to preserve the sign.<\/p>\n So, without further ado, let’s start deconstructing my code, a section at a time:<\/p>\n First a bit of initialisation – instructions to the assembler that the code should be 16-bit, and will be loaded into memory at location 0x7c00. This video mode will fit a nice big 512×256 canvas on the screen.<\/p>\n Total assembled size so far: 5 bytes<\/p>\n So, here we initialise the two variables for our c value. DX will contain the y value and DI the x value. I’m starting at y=-1 and x=-2, using a scaling factor of 1\/128. My canvas is to be y from -1 to 1 and x from -2 to 2. The z variable is stored with its real part in SI and its imaginary part in BP. We skip out the zeroth iteration here by initialising z=c.<\/p>\n Later, the interrupt to write a pixel to the screen will read the pixel value from AL. This wants to be initialised to zero (black). If the pixel is in the set, this will be changed to one (white) later. Otherwise it will be unchanged. This block: 12 bytes, So far: 17 bytes.<\/p>\n At each iteration step we need to calculate As the comment suggests, here we’re calculating Now, back to the jump. This is how we catch the sequences which go off to infinity – we’ll see where it goes later. In the C++ example, we used abs(z) to check whether we were inside an r=2 circle. But calculating an abs value is expensive. It requires us to square two numbers and sum them. But consider this: any point which leaves the r=2 circle will eventually go off to infinity. If I wanted to I could use an r=3 circle, or an r=1,000,000 circle – it would just potentially take more iterations. Any region which completely contains the r=2 circle will do the job. So how about a rectangle?<\/p>\n My jump is conditional on So, provided that |Im(c)|<6 - which it always will be since we're only working with y in the range of -1 to 1 - this will contain the whole of the r=2 circle. Perfect.\n\nThis block: 10 bytes, So far: 27 bytes\n\n[code];calculate Re(z^2)=x^2-y^2\nmov bx,si ;we will work in the BX register for now\nadd bx,bp ;bx contains x+y\nsub si,bp ;si contains x-y\nimul si,bx\njo overflow\nsar si,7 ;si contains Re(z^2)[\/code]\n\n After we’ve calculated Again, we check for overflow and SAR 7 to take out the factor of 128. Here, the overflow checks if The code is compiled with nasm:<\/p>\n And run in qemu:<\/p>\n Here’s a screenshot and a video of it running:The Mandelbrot Set<\/h2>\n
![](\"https:\/\/upload.wikimedia.org\/wikipedia\/commons\/thumb\/2\/21\/Mandel_zoom_00_mandelbrot_set.jpg\/320px-Mandel_zoom_00_mandelbrot_set.jpg\")
\nThe Mandelbrot Set<\/a> is the set of points, c, on the complex plane such that the sequence zn+1<\/sub>=zn<\/sub>2<\/sup>+c (where z0<\/sub>=0) remains bounded. To clarify, here’s an example:<\/p>\n\n
\n![](\"http:\/\/blog.nothinguntoward.eu\/wp-content\/uploads\/2015\/01\/mandelbrotiter-300x225.png\" "\\"z1")
<\/a>
\nBut overall, the general trend seems to be going around in a closed circle, rather than shooting off to infinity, so it would probably be fairly safe to say here that the sequence is bounded. If you’re not convinced, here’s a plot of the magnitude at each iteration:
\n![](\"http:\/\/blog.nothinguntoward.eu\/wp-content\/uploads\/2015\/01\/mandelbrotiter-mag-300x225.png\" "\\"all")
<\/a><\/p>\n\n
#include <iostream>\r\n#include <complex>\r\nusing namespace std;\r\nint main(){\r\n\tcout << "P1 800 800 ";\r\n\tfor(int i=0;i<800;i++){\r\n\t\tfor(int j=0;j<800;j++){\r\n\t\t\tcomplex<float> c((j-400)\/200.0,(i-400)\/200.0); \/\/scale the indicies to floats in the range -2..2\r\n\t\t\tcomplex<float> z(0,0); \/\/ initialise z to 0\r\n\t\t\tfor (int iter=0;iter<64;iter++) \/\/iterate 64 times\r\n\t\t\t\tz=z*z+c;\r\n\t\t\tcout << (abs(z)<2)?"1 ":"0 "; \/\/output a black pixel if in the set, white otherwise\r\n\t\t\t}\r\n\t\t}\r\n\tcout<<endl;\r\n\treturn 0;\r\n\t}<\/pre>\nc=ones(800,1)*linspace(-2,2,800);\r\nc=complex(c,c');\r\nz=zeros(800,800);\r\nfor i=1:50\r\n\tz=z.^2+c;\r\nend\r\nimshow(abs(z)<2)<\/pre>\n
The Bare-Metal Version<\/h2>\n
A note on arithmetic<\/h3>\n
\n
\nand specifically:
\n
\n
\n
\nFor example, here is how we might multiply two numbers with a scaling of 1\/32:<\/p>\nmov al,73 ; move the number 2.28 (73\/32) into al\r\nmov bl,-193 ;move -6.03 (-193\/32) into bl\r\nimul bl ; multiply bl (implicitly, with al). stores result in ax.\r\nsar ax,5 ; shift right by 5 places to divide by 2^5=32<\/pre>\n
The code<\/h3>\n
[bits 16]\r\n[org 7c00h]\r\nmov ax,0fh\r\nint 10h<\/pre>\n
\nThen, we call BIOS interrupt 10h<\/a> with AH=0 (change video mode) and AL=0fh (640×350 mono<\/a> resolution).<\/p>\n\r\nmov dx,-128 ;start on the first line\r\nyloop:\r\nmov di,-256 ;go to the beginning of the line\r\nxloop:\r\n;set z=c initially - start on iteration 1\r\nmov si,di\r\nmov bp,dx\r\nxor ax,ax ; set colour to black (al=0) and iteration counter to zero (ah=0)<\/pre>\n
\nAt the same time, I want to set an iteration counter to zero. Sneakily, I’ve chosen to use the other byte of AX for this, so that I can set both the pixel value and the iteration counter to zero in a single instruction. This single xor saves two bytes vs two 8-bit literal moves.<\/p>\niterate:\r\n;calculate Im(z^2)=2xy\r\nmov cx,bp\r\nimul cx,si\r\njo overflow\r\nsar cx,6 ;maybe mov 6 into a register first?\r\n;cx contains Im(z^2)<\/pre>\n
![\"areas\"](\"http:\/\/blog.nothinguntoward.eu\/wp-content\/uploads\/2015\/01\/areas-180x300.png\")
<\/a>And now to calculate 
The result<\/h3>\n
$ nasm -o mandelbrot mandelbrot-bootsector.asm<\/pre>\n
$ qemu-system-i386 -video cirrus -fda mandelbrot -net none<\/pre>\n
\n![\"mandelbrotscreenshot\"](\"http:\/\/blog.nothinguntoward.eu\/wp-content\/uploads\/2015\/01\/mandelbrotscreenshot-300x176.png\")
<\/a>
\n
