Koch’s Snowflake
December 23, 2011
It’s Christmas time, and that means winter, and that means snow (sometimes, but probably not this year), so today’s exercise is to draw Koch’s Snowflake. This famous fractal is built from an equilateral triangle, as shown upper left. At each step each side is modified by removing its center third and replacing it with two legs of another equilateral triangle. Thus, at upper right, each of the three legs has been replaced four line segments that form a new shape. Carrying on, the fractal becomes more and more detailed at the edges of the triangle; eventually, the perimeter of the triangle is infinite, though its area is 8/5 the area of the original triangle. Mathematically, the Koch Snowflake can be encoded as a Lindenmayer system with initial string F++F++F, rewriting rule F → F-F++F-F, and angle 60 degrees (in the Lindenmayer notation, F means move forward, + is a right turn and - is a left turn).
Your task is to write a program to draw the Koch Snowflake. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
I noticed there are only six possible directions, and with each axis scaled properly, the arithmetic is very simple. The following computes the coordinates for a given list of moves, F, +, or -, starting at origin and initial angle of 60 degrees. (The x and y coordinates are in separate lists and in the backwards order.)
(define (kochdinates moves)
(let loop ((x 0) (y 0) (xs '(0)) (ys '(0))
(theta 60) (moves moves))
(if (null? moves)
(values xs ys)
(case (car moves)
((F) (call-with-values
(lambda ()
(case theta
(( 0) (values (+ x 2) y))
(( 60) (values (+ x 1) (+ y 1)))
((120) (values (- x 1) (+ y 1)))
((180) (values (- x 2) y))
((240) (values (- x 1) (- y 1)))
((300) (values (+ x 1) (- y 1)))))
(lambda (x1 y1)
(loop x1 y1 (cons x1 xs) (cons y1 ys)
theta (cdr moves)))))
((+) (loop x y xs ys (modulo (- theta 60) 360) (cdr moves)))
((-) (loop x y xs ys (modulo (+ theta 60) 360) (cdr moves)))))))
Here’s the code that produces the snowflake moves for the coordinate computer above. Nest at will.
(define (init) '(F + + F + + F))
(define (rewrite moves)
(apply append
(map (lambda (move)
(case move
((F) '(F - F + + F - F))
((+) '(+))
((-) '(-))))
moves)))
And here’s a compiler that produces an R command to plot a snowflake.
(define (R xs ys)
(display "plot(c(") (display (car xs))
(for-each (lambda (x) (display ",") (display x)) (cdr xs))
(display "), c(") (display (car ys))
(for-each (lambda (y) (display ",") (display y)) (cdr ys))
(display "), type='l', xlab='', ylab='')") (newline))
For example, this:
(call-with-values (lambda () (kochdinates (rewrite (init)))) R)
prints this:
plot(c(0,2,3,4,6,5,6,4,3,2,0,1,0), c(0,0,-1,0,0,1,2,2,3,2,2,1,0), type='l', xlab='', ylab='')
which can be pasted in an R session to open the plot window. For scaling on the screen, grab the corner of the plot window and let the computer do the arithmetic. (I ran Guile and R in one Emacs shell session, so cutting and pasting was easy, but I had to close one interpreter to enter the other. And the R command line editor seemed to choke when there were too many coordinates :)
[…] lovely programming exercise. I took the task one step further and created an interactive explorer of variants of the Koch […]
Use your favourite LOGO interpreter…
(defparameter *koch* '(f + + f + + f)) (defun koch-iterate (koch iterations) (loop while (> iterations 0) do (setq iterations (1- iterations)) (setq koch (apply #'append (mapcar (λ(v) (case v ((f) '(f - f + + f - f)) (t `(,v)))) koch)))) koch) (defun 2logo () (with-output-to-string (s) (mapcar (λ(v) (case v ((f) (format s "fd 5 ")) ((+) (format s "rt 60 ")) ((-) (format s "lt 60 ")))) (koch-iterate *koch* 5))))Here’s a Haskell version, using version 0.4 of the diagrams library. The heart of the program is the
stepfunction. Given an input shape,_, it produces_/\_. This is iterated some number of times (theedgefunction), starting with a single, straight line segment. Lastly, three copies of the final list of segments are glued together in a triangular shape.import Diagrams.Prelude import Diagrams.Backend.Cairo.CmdLine koch n = let e = edge !! n t = mconcat [e, rotateBy (-1/3) e, rotateBy (-2/3) e] in pathFromTrail $ scale (1/3^n) t where edge = iterate step $ fromSegments [straight (1,0)] step p = mconcat [p, rotateBy (1/6) p, rotateBy (-1/6) p, p] main = defaultMain . attrs . stroke $ koch 5 where attrs = lw 0.002 . pad 1.1Here are a few examples of running the program to generate PNG, PDF and SVG output.