[Lars H], 2008-08-27: I stumbled upon a [Wikipedia] article [http://en.wikipedia.org/wiki/Hashlife] about this curious little algorithm, and once I had read enough to understand how it works, I couldn't help implementing it (especially since Wikipedia claims ''"While an amateur programmer might be able to write a simple Life player over an afternoon, few Hashlife implementations exist"''). So here it is, as a weekend fun project. Basic principles of the algorithm: 1. Use a quadtree to store patterns. This simplifies supporting unlimited grid sizes (just use a tree with enough levels). 2. Only ever create one node for each pattern (so the tree is really a [DAG]). This requires introducing the [hash] to map pattern represented by a node to the actual node index. The hash key is formed from the indices of the four subnodes. 3. In [Life], cells are affected by all their neighbours, so a given ''n''×''n'' grid only contains enough information to compute the next state for the middle (''n''–2)×(''n''–2) cells. If ''n'' = 2**''k'' (for a level ''k'' node in the quadtree), then the largest size of a quadtree node that fits within the inner square is 2**(''k''–1)× 2**(''k''–1) — one level further down. Hence this is the basic mode of operation for stepping in [HashLife]: given a tree node, compute the one step smaller node for the next state of the middle part of the given grid. 4. Use [memoizing] to speed things up when the same node occurs in several places — cache (in the node itself) the ''next'' value for a node once you've computed it. 5. Make the step size for a node proportional to the side of the node! In the code, this amounts to making '''next::bignode''' advance its components two steps instead of one, but the result is an exponential speedup in the simulated time. **Short demo** ====== % encode {{1 1 1} {1 0 1} {1 0 1} {} {} {} {1 1} {1 1}} 10 % decode 10 {1 1 1 0 0 0 0 0} {1 0 1 0 0 0 0 0} {1 0 1 0 0 0 0 0} {0 0 0 0 0 0 0 0} {0 0 0 0 0 0 0 0} {0 0 0 0 0 0 0 0} {1 1 0 0 0 0 0 0} {1 1 0 0 0 0 0 0} % chardecode 10 OOO O O O O OO OO % grow 10 14 % chardecode 14 OOO O O O O OO OO % step 14 137 «» chardecode 137 OOO O O O O O O OOO OOO OO OO % set generations 8 ====== Since node 14 is a 16×16 grid, '''step'''ping it forward advances time 8 generations. '''grow'''ing larger nodes also speeds up the processing. ====== % step [grow 137] 743 % step [grow 743] 2769 % set generations 56 ====== Larger nodes are difficult to display in character graphics. '''node2xbm" can be used to make a bitmap instead, like so: % package require Tk % image create bitmap -background yellow image1 % pack [label .l -image image1] % image1 configure -data [node2xbm 2769 3] ; # 3 for making a 3x3 block of pixels for each cell. High level commands defined above: '''encode''': Compute node (number) for a matrix of 0s and 1s. '''decode''': Given node (number), return corresponding matrix. '''chardecode''': Given node (number), return corresponding block of chars. '''grow''': Double the side of a node, by padding with zeroes. Return new node number. '''side''': Given node (number), return its side length. '''step''': Given node (number), advance it ''side''/2 generations, return new node number. (Combination of grow and next, so new node has the same size as the old one.) '''step1''': Given node (number), advance it 1 generation, return new node number. (Combination of grow and next1, so new node has the same size as the old one.) '''node2xbm''': Given node and magnification factor (number of pixels in side of one cell), return [XBM] data for a corresponding bitmap. **The code** ====== # HashLife -- Hashed quadtree-based Life algorithm with memoisation # The basic data structure is a heap (implemented as a list) and a hash # table (array) that stores the positions in the heap of particular nodes. # Heap item contents: # # cell <0 or 1> # 2node # 4node # bignode # cached # cached1 # # Evaluating one of them in the next namespace computes the node # and mutates the heap item to a [cached] node. Evaluating one of them # in the next1 namespace computes the node and mutates the heap # item to a [cached1] node. # # The , , , , and are positions of other nodes. # The is the position of the node itself, which is needed for # nodes that mutate themselves into [cached] nodes. # In a [2node], they're [cell]s. In a [4node], they're [2node]s. # In a [bignode], they're [4node]s, [bignode]s, or [cached]s. # In a [cached], they're [2node]s or [cached]s. # # set heap {{cell 0} {cell 1}} # Must be items 0 and 1 -- relied upon by several procs below. array unset hash ## # The following procedure finds the heap position (pointer) for a # particular node, given its four corners. It will *create* the heap # item for the node if it doesn't already exist. ## proc find {nw ne sw se} { variable hash set key [list $nw $ne $sw $se] if {![info exists hash($key)]} then { variable heap set hash($key) [llength $heap] switch -- [lindex $heap $nw 0] "cell" { lappend heap [list 2node $nw $ne $sw $se] } "2node" { lappend heap [list 4node $nw $ne $sw $se $hash($key)] } default { lappend heap [list bignode $nw $ne $sw $se $hash($key)] } } return $hash($key) } ## # side $node -- Compute and return the side of the $node. ## proc side {node} {expr { [lindex $::heap $node 0] eq "cell" ? 1 : 2*[side [lindex $::heap $node 1]] }} ## # encode $matrix # Encode a matrix of booleans as a node in the heap. # The matrix is roughly centered within the node's grid. # # encode $matrix $side $col $row # Encode the $side-by-$side submatrix of $matrix whose # NW corner is column $col and row $row. This submatrix # need not be contained within the $matrix. The $side # must be a power of 2. ## proc encode {matrix {side 0} {x 0} {y 0}} { if {!$side} then { set min [llength $matrix] foreach row $matrix { if {[llength $row]>$min} then {set min [llength $row]} } set side 1 while {$side<$min} {incr side $side} set x [expr {($min-$side)/2}] set y [expr {([llength $matrix]-$side)/2}] } if {$side==1} then { return [lindex {0 1} [string is true -strict [lindex $matrix $y $x]]] } set half [expr {$side/2}] set nw [encode $matrix $half $x $y] set ne [encode $matrix $half [expr {$x+$half}] $y] set sw [encode $matrix $half $x [expr {$y+$half}]] set se [encode $matrix $half [expr {$x+$half}] [expr {$y+$half}]] return [find $nw $ne $sw $se] } ## # Decode a heap node as a matrix of 0s and 1s. ## proc decode {node} { if {[lindex $::heap $node 0] == "cell"} then { return [list [list [lindex $::heap $node 1]]] } set res {} foreach low [decode [lindex $::heap $node 1]]\ high [decode [lindex $::heap $node 2]] { lappend res [concat $low $high] } foreach low [decode [lindex $::heap $node 3]]\ high [decode [lindex $::heap $node 4]] { lappend res [concat $low $high] } return $res } ## # Decode a heap node as a character block. Dead cells # are spaces, live cells some char (defaults to O). ## proc chardecode {node {cell O}} { string map [list { } {} 0 { } 1 $cell] [ join [decode $node] \n ] } ## # Expand a heap node $level times, into a 2**$level by 2**$level # matrix of nodes. ## proc expand {node level} { if {$level<=0} then {return [list [list $node]]} incr level -1 set res {} foreach low [expand [lindex $::heap $node 1] $level]\ high [expand [lindex $::heap $node 2] $level] { lappend res [concat $low $high] } foreach low [expand [lindex $::heap $node 3] $level]\ high [expand [lindex $::heap $node 4] $level] { lappend res [concat $low $high] } return $res } ## # Contract the $side-by-$side submatrix, whose NW corner is at # ($row,$col), of the $matrix of nodes into a single node. ## proc contract {matrix side row col} { if {$side<=1} then {return [lindex $matrix $row $col]} set half [expr {$side/2}] set nw [contract $matrix $half $row $col] set ne [contract $matrix $half $row [expr {$col+$half}]] set sw [contract $matrix $half [expr {$row+$half}] $col] set se [contract $matrix $half [expr {$row+$half}] [expr {$col+$half}]] return [find $nw $ne $sw $se] } namespace eval next { proc cached {nw ne sw se next} {return $next} proc cell {c} {error "Can't advance a single cell"} proc 2node {nw ne sw se} {error "Can't advance a 2x2 node"} proc cached1 {nw ne sw se next} { set self [find $nw $ne $sw $se] if {[lindex $::heap $nw 0] == "2node"} then { 4node $nw $ne $sw $se $self } else { bignode $nw $ne $sw $se $self } } } ## # The following procedure computes the next generation of the middle # [2node] in a [4node]. It does not mutate the node itself. ## proc 4node {nw ne sw se} { set grid {} foreach low [decode $nw] high [decode $ne] { lappend grid [concat $low $high] } foreach low [decode $sw] high [decode $se] { lappend grid [concat $low $high] } set call [list find] foreach y1 {1 2} { foreach x1 {1 2} { set sum 0 foreach y2 {-1 0 1} { foreach x2 {-1 0 1} { incr sum [lindex $grid [expr {$y1+$y2}] [expr {$x1+$x2}]] } } set sum [expr {2*$sum - [lindex $grid $y1 $x1]}] lappend call [expr {$sum>=5 && $sum<=7}] } } eval $call } ## # The following procedure computes the next generation of the middle # [2node] in a [4node], and mutates the node into a [cached] node. ## proc next::4node {nw ne sw se self} { set next [::4node $nw $ne $sw $se] lset ::heap $self 0 cached lset ::heap $self 5 $next return $next } proc next::bignode {nw ne sw se self} { set grid {} foreach low [expand $nw 1] high [expand $ne 1] { lappend grid [concat $low $high] } foreach low [expand $sw 1] high [expand $se 1] { lappend grid [concat $low $high] } set n [contract $grid 2 0 1] set w [contract $grid 2 1 0] set mid [contract $grid 2 1 1] set e [contract $grid 2 1 2] set s [contract $grid 2 2 1] set grid2 {} foreach row [list [list $nw $n $ne] [list $w $mid $e] [list $sw $s $se]] { set row2 {} foreach node $row { lappend row2 [eval [lindex $::heap $node]] } lappend grid2 $row2 } set call [list find] foreach row {0 1} { foreach col {0 1} { set node [contract $grid2 2 $row $col] lappend call [eval [lindex $::heap $node]] } } set next [eval $call] lset ::heap $self 0 cached lset ::heap $self 5 $next return $next } namespace eval next1 { proc cached1 {nw ne sw se next} {return $next} proc cell {c} {error "Can't advance a single cell"} proc 2node {nw ne sw se} {error "Can't advance a 2x2 node"} proc cached {nw ne sw se next} { set self [find $nw $ne $sw $se] if {[lindex $::heap $nw 0] == "2node"} then { lset ::heap $self 0 cached1 set next } else { bignode $nw $ne $sw $se $self } } } ## # The following procedure computes the next generation of the middle # [2node] in a [4node], and mutates the node into a [cached1] node. ## proc next1::4node {nw ne sw se self} { set next [::4node $nw $ne $sw $se] lset ::heap $self 0 cached1 lset ::heap $self 5 $next return $next } proc next1::bignode {nw ne sw se self} { set M [expand $self 3] set call [list find] foreach {row col} {1 1 1 3 3 1 3 3} { set node [contract $M 4 $row $col] lappend call [eval [lindex $::heap $node]] } set next [eval $call] lset ::heap $self 0 cached1 lset ::heap $self 5 $next return $next } ## # Commands in the zeroed namespace return a node with the same size as # the node which they are, but with all cells 0. ## namespace eval zeroed { proc cell {c} {return 0} proc 2node {a b c d} { find 0 0 0 0 } proc cached {a b c d next} { set sub [eval [lindex $::heap $a]] find $sub $sub $sub $sub } } interp alias {} zeroed::bignode {} zeroed::cached interp alias {} zeroed::4node {} zeroed::cached interp alias {} zeroed::cached1 {} zeroed::cached ## # The [grow] procedure return a node whose side is twice # that of the argument, by padding it with a frame of all # zeroes. ## proc grow {node} { if {[lindex $::heap $node 0] == "cell"} then { error "You can't grow a single cell" } foreach {type nw ne sw se} [lindex $::heap $node] break set clean [namespace eval ::zeroed [lindex $::heap $nw]] set call [list find] lappend call [find $clean $clean $clean $nw] lappend call [find $clean $clean $ne $clean] lappend call [find $clean $sw $clean $clean] lappend call [find $se $clean $clean $clean] eval $call } set generations 0 proc step {node} { variable generations set side [side $node] set node2 [grow $node] set res [namespace eval ::next [lindex $::heap $node2]] incr generations [expr {$side/2}] return $res } proc step1 {node} { variable generations set node2 [grow $node] set res [namespace eval ::next1 [lindex $::heap $node2]] incr generations return $res } # A [bitmap] actually seems more appropriate than a [photo] for # representing Life patterns, so we need some way of converting # nodes to XBM data. Intermediate representation: # # $side $matrix # # If the $side is divisible by 8, then the $matrix is a list of # $side rows, where each row is a string matching the regexp: # (0x[[:xdigit:]]{2}, )+ # # If the $side is not divisible by 8, then the $matrix is a list # of $side rows, where each row is a string of $side binary digits. proc node2bm {node mag} { if {$node <= 1} then { set row [string repeat $node $mag] if {!($mag % 8)} then { binary scan [binary format b* $row] H* row regsub -all {..} $row {0x&, } row } set mat {} for {set n 1} {$n<=$mag} {incr n} {lappend mat $row} return [list $mag $mat] } set mat {} foreach {side leftM} [node2bm [lindex $::heap $node 1] $mag] break foreach l $leftM r [ lindex [node2bm [lindex $::heap $node 2] $mag] 1 ] {lappend mat $l$r} foreach l [ lindex [node2bm [lindex $::heap $node 3] $mag] 1 ] r [ lindex [node2bm [lindex $::heap $node 4] $mag] 1 ] {lappend mat $l$r} if {($side%8) && !($side%4)} then { set mat2 {} foreach row $mat { binary scan [binary format b* $row] H* row regsub -all {..} $row {0x&, } row lappend mat2 $row } set mat $mat2 } return [list [expr {2*$side}] $mat] } proc node2xbm {node mag} { foreach {side mat} [node2bm $node $mag] break if {$side % 8} then { set width [expr {8*(1 + $side/8)}] set mat2 {} foreach row $mat { binary scan [binary format b* $row] H* row regsub -all {..} $row {0x&, } row lappend mat2 $row } set mat $mat2 } else { set width $side } format { #define data_width %d #define data_height %d static unsigned char data_bits[] = { %s }; } $width $side [string trimright [join $mat "\n "] ", "] } ====== ---- !!!!!! %| [Category Algorithm] | [Category Toys] | [Category Compression] |% !!!!!!