% 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

% step [grow 137]
743
% step [grow 743]
2769
% set generations
56

# 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 <nw> <ne> <sw> <se>
# 4node <nw> <ne> <sw> <se> <self>
# bignode <nw> <ne> <sw> <se> <self>
# cached <nw> <ne> <sw> <se> <next>
# cached1 <nw> <ne> <sw> <se> <next>
# 
# Evaluating one of them in the next namespace computes the <next> node 
# and mutates the heap item to a [cached] node. Evaluating one of them 
# in the next1 namespace computes the <next> node and mutates the heap 
# item to a [cached1] node.
# 
# The <nw>, <ne>, <sw>, <se>, and <next> are positions of other nodes. 
# The <self> 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    "] ", "]
}