Generating Mandelbrot Set Fractals with Lisp.
At one of my previous jobs I was lucky to have a chance to do some work in Lisp. I tried then to use it here and there, but found that most of the time it's not really suitable to solve my tasks (apart from that most often we are not free to choose the tools to work with).Some Lisp publications claim that in many cases Lisp's performance can be comparable to that of modern languages, but unfortunately it was not my case - it seems you have to be a very experienced Lisper to achieve that (I'm not, so please don't be very critical of my Lisp code). The Mandelbrot set fractals generation program that you can download here generates an image in 800x600 resolution and is VERY hungry for CPU power. However, I hope it can be fun for someone to play with Lisp this way, and as a part of it, the program contains a set of easy-to-use functions to create image files in Microsoft BMP format (24-bit color, uncompressed).
I am using the very basic language elements so the code will be very easy to understand even for a beginner. The graphics functions use standard streams and should work on any Common Lisp platform. I ran it on CMU CL and CLISP. Compiled CMU CL is much faster; CLISP is slow but available for many OS.
To run the program, you can start an interactive session, then use (load "mandelbrot.lsp") to load the code, and use (gen-fractal) to run it. As it goes through the image rows, the current row number will be printed. Unless you change the file name in the code, the resulting image will be written in the /temp/testlisp.bmp. NOTE: There is no I/O error processing - please verify if the file name will be suitable for your environment and change it as needed.
When playing with the program, you will probably want to display different areas of the fractal or change the fractal size. As every calculation takes too long, you may want to temporarily change the maximum iteration number in the function num-iter, i.e., replace (cond ((> start 359) 359) with (cond ((> start 50) 359). This way you will get the draft result much faster; the image will contain almost no detail but you will be able to judge what can be the final calculation's result.
Any feedback or code improvements will be greatly appreciated!
The image fragment that corresponds to the code you can download (click to see it in full 800x600 resolution):
The code is below. You can copy it and save in a file named "mandelbrot.lsp":
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;; Functions for BMP file generation
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;; Write 16-bit number to a stream
(defun write-int16( number outstream )
(let ( (bbblo (byte 8 0))
(bbbhi (byte 8 8)) )
(write-byte (ldb bbblo number) outstream)
(write-byte (ldb bbbhi number) outstream)
)
)
;;; Write 32-bit number to a stream
(defun write-int32( number outstream )
(let ( (bbb0 (byte 8 0))
(bbb1 (byte 8 8))
(bbb2 (byte 8 16))
(bbb3 (byte 8 24)) )
(write-byte (ldb bbb0 number) outstream)
(write-byte (ldb bbb1 number) outstream)
(write-byte (ldb bbb2 number) outstream)
(write-byte (ldb bbb3 number) outstream)
)
)
;;; Create bitmap array
(defun make-bitmap( w h )
(make-array (* w h))
)
;;; Set pixel
(defun set-pixel1( bmp w h x y r g b )
(let ( (pix (+ r (* g 256) (* b 65536)))
)
(setf (aref bmp (+ (* y w) x)) pix)
)
)
;;; Get pixel
(defun get-pixel1( bmp w h x y )
(let ( (bbb0 (byte 8 0))
(bbb1 (byte 8 8))
(bbb2 (byte 8 16))
(pix (aref bmp (+ (* y w) x)))
)
(cond ((null pix) '(0 0 0))
(t (list (ldb bbb0 pix) (ldb bbb1 pix) (ldb bbb2 pix))))
)
)
;;; Write the bitmap array data to BMP file
(defun write-bitmap( bmp w h fname )
(let ( (bmp-file (open fname :direction :output :element-type 'unsigned-byte))
)
;; Create BMP file header according to MS specs
(write-byte #x42 bmp-file)
(write-byte #x4d bmp-file)
;; bfsize
(write-int32 (+ 54 (* 3 w h)) bmp-file)
;; 2x bfreserved
(write-int32 0 bmp-file)
;; offbits
(write-int32 54 bmp-file)
;; bmp info header
;; DWORD biSize;
(write-int32 40 bmp-file)
;; LONG biWidth;
(write-int32 w bmp-file)
;; LONG biHeight;
(write-int32 h bmp-file)
;; WORD biPlanes;
(write-int16 1 bmp-file)
;; WORD biBitCount
(write-int16 24 bmp-file)
;; DWORD biCompression;
(write-int32 0 bmp-file)
;; DWORD biSizeImage;
(write-int32 0 bmp-file)
;; LONG biXPelsPerMeter;
(write-int32 4000 bmp-file)
;; LONG biYPelsPerMeter;
(write-int32 4000 bmp-file)
;; DWORD biClrUsed;
(write-int32 0 bmp-file)
;; DWORD biClrImportant;
(write-int32 0 bmp-file)
;; write image bytes
(do ((y 0 (+ y 1)))
((= y h) t)
(do ((x 0 (+ x 1)))
((= x w) t)
(let ( (pix (get-pixel1 bmp w h x y)))
(write-byte (third pix) bmp-file)
(write-byte (second pix) bmp-file)
(write-byte (car pix) bmp-file)
)
)
)
;; close the file
(close bmp-file)
)
)
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;; Mandelbrot set generation
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;; Please refer to Internet for explanation on the
;;; mathematics of the Mandelbrot set and fractals.
;;; First, we need functions to convert a number from a range
;;; of [0;359) to a color to colorize the set image.
;;; HLS to RGB color conversion - helper function
(defun hls2rgb1( rn1 rn2 huei )
(let ( (hue huei)
)
(if (> hue 360.0) (setq hue (- hue 360.0)))
(if (< hue 0.0) (setq hue (+ hue 360.0)))
(cond ( (< hue 60.0) (+ rn1 (/ (* hue (- rn2 rn1) ) 60.0) ) )
( (< hue 180.0) rn2 )
( (< hue 240.0) (+ rn1 (/ (* (- 240.0 hue) (- rn2 rn1) ) 60.0) ) )
( t rn1 )
)
)
)
;;; HLS to RGB color conversion - helper function
(defun hls2rgb( hue light satur )
(let ( (rh 0.0) (rl 0.0) (rs 0.0) rm1 rm2 r g b)
(if (> hue 0.0) (setq rh hue))
(if (> rh 360.0) (setq rh 360.0))
(if (> light 0.0) (setq rl light))
(if (> rl 1.0) (setq rl 1.0))
(if (> satur 0.0) (setq rs satur))
(if (> rs 1.0) (setq rs 1.0))
(cond ((<= rl 0.5) (setq rm2 (* rl (+ 1.0 rs))) )
(t (setq rm2 (- (+ rl rs) (* rl rs) ) ) ) )
(setq rm1 (- (+ rl rl) rm2) )
(cond ((equal rs 0.0) (progn (setq r rl) (setq g rl) (setq b rl)))
(t (progn (setq r (hls2rgb1 rm1 rm2 (+ rh 120.0)))
(setq g (hls2rgb1 rm1 rm2 rh))
(setq b (hls2rgb1 rm1 rm2 (- rh 120.0))))))
(list r g b)
)
)
;;; Calculate iteration number for the point in the image.
;;; This will define the color of this point.
(defun num-iter( start x y xc yc)
(let* ( (x1 (+ (- (* x x) (* y y)) xc) )
(y1 (+ (* 2.0 x y) yc) )
(sq (+ (* x1 x1) (* y1 y1)))
)
(cond ((> start 359) 359)
((> sq 2.0) start)
(t (num-iter (+ start 1) x1 y1 xc yc))
)
)
)
(defun gen-fractal()
(let ( (bmp (make-bitmap 800 600)) )
(do ((y 0 (+ y 1)))
((= y 600) t)
(do ((x 0 (+ x 1)))
((= x 800) t)
(let ( (xx (/ (- (coerce x 'double-float) 500.0) 700.0) )
(yy (/ (- (coerce y 'double-float) -220.0) 700.0) )
nnn rgb )
(setq nnn (num-iter 0 xx yy xx yy))
(if (= nnn 359) (set-pixel1 bmp 800 600 x y 0 0 0)
(progn
(setq nnn (+ 327 nnn))
(if (> nnn 359) (decf nnn 359))
(setq rgb (hls2rgb (coerce nnn 'double-float) 0.5 1.0))
(set-pixel1 bmp 800 600 x y
(round (* 255.0 (first rgb)))
(round (* 255.0 (second rgb)))
(round (* 255.0 (third rgb))))))
)
)
(print y)
)
(write-bitmap bmp 800 600 "/temp/testlisp.bmp")
)
)
;;; the next line is optional - to suppress some warnings
(setq *warn-on-floating-point-contagion* 'nil)
;;; END
Have fun!