% 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
}


## 
 # Check if a node is a node of all zeroes.
 ##

proc iszero {node} {
    while {$node>=2} {
        foreach {type nw ne sw se} [lindex $::heap $node] break
        if {$nw!=$ne || $nw!=$sw || $nw!=$se} then {return 0}
        set node $nw
    }
    return [expr {$node==0}]
}

## 
 # 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
}

## 
 # The [safe] procedure returns $node if the border of the level 
 # 2 expansion of the $node is all zeroes, and [grow $node] 
 # otherwise. This ensures [step [safe $node]] hasn't lost any 
 # cells at the edge.
 ##

proc safe {node} {
    set grid [expand $node 2]
    set clean [namespace eval ::zeroed [lindex $::heap\
      [lindex $grid 0 0]]]
    foreach row $grid rc {1 0 0 1} {
        foreach item $row cc {1 0 0 1} {
            if {($rc||$cc) && $item!=$clean} then {
                return [grow $node]
            }
        }
    }
    return $node
}

##
 # Garbage collection
 ##

namespace eval gc {
    
    ## 
     # Node commands in the gc namespace create a new heap node 
     # for this old heap element and return the number of that 
     # new node.
     ##
    
    variable newheap
    variable tonew
    
    ## 
     # The gc::asnew command generally returns the number of a 
     # new node, given the number of the corresponding old node.
     # It may as a side-effect modify the newheap.
     ##
    
    proc asnew {node} {
        variable tonew
        if {![info exists tonew($node)]} then {
            set tonew($node) [eval [lindex $::heap $node]]
        }
        return $tonew($node)
    }
    
    ## 
     # The gc::find command generally returns the number of a 
     # new node, given the number of the corresponding old node.
     ##
    
    proc cell {bit} {return $bit}
    
    proc 2node {nw ne sw se} {
        variable newheap
        lappend newheap\
          [list 2node [asnew $nw] [asnew $ne] [asnew $sw] [asnew $se]]
        return [expr {[llength $newheap]-1}]
    }
    
    proc 4node {nw ne sw se self} {
        variable newheap
        lappend newheap\
          [list 4node [asnew $nw] [asnew $ne] [asnew $sw] [asnew $se] ""]
        set num [expr {[llength $newheap]-1}]
        lset newheap end end $num
        return $num
    }
    
    proc bignode {nw ne sw se self} {
        variable newheap
        lappend newheap\
          [list bignode [asnew $nw] [asnew $ne] [asnew $sw] [asnew $se] ""]
        set num [expr {[llength $newheap]-1}]
        lset newheap end end $num
        return $num
    }
    
    proc cached {nw ne sw se next} {
        variable newheap
        lappend newheap [list cached [asnew $nw] [asnew $ne]\
          [asnew $sw] [asnew $se] [asnew $next]]
        return [expr {[llength $newheap]-1}]
    }
    
    proc cached1 {nw ne sw se next} {
        variable newheap
        lappend newheap [list cached1 [asnew $nw] [asnew $ne]\
          [asnew $sw] [asnew $se] [asnew $next]]
        return [expr {[llength $newheap]-1}]
    }
}

proc garbagecollect {args} {
    set gc::newheap {{cell 0} {cell 1}}
    array unset gc::tonew
    
    set res {}
    foreach node $args {
        lappend res [gc::asnew $node]
    }
    
    array unset gc::tonew {[01]}
    variable hash
    array unset hash
    foreach {old node} [array get gc::tonew] {
        set hash([lrange [lindex $gc::newheap $node] 1 4]) $node
    }
    array unset gc::tonew
    
    set ::heap $gc::newheap
    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    "] ", "]
}


package require Tk 8.5

namespace eval ::tophoto {
    image create photo ::tophoto::dead -width 1 -height 1
    image create photo ::tophoto::live -width 1 -height 1
    
    variable where
    variable ms
    
    proc idleupdate {period} {
        variable ms
        if {$ms + $period < [clock milliseconds]} then {
            update idletasks
            set ms [clock milliseconds]
        }
    }
    
    proc init {fg bg period} {
        dead put $bg
        live put $fg
        variable where
        unset -nocomplain where
        if {$period < infinity} then {
            interp alias {} ::tophoto::maybe_idleupdate\
              {} ::tophoto::idleupdate $period
            variable ms [clock milliseconds]
        } else {
            proc ::tophoto::maybe_idleupdate {} {}
        }
    }
    
    proc pixelnode {node img side x y w h} {
        variable where
        if {![info exists where($node)]} then {
           if {[iszero $node]} then {
               set where($node) [list [namespace which dead]]
           } else {
               set where($node) [list [namespace which live]]
           }
        }
        $img copy {*}$where($node) -to\
          [expr {$x<0 ? 0 : $x}] [expr {$y<0 ? 0 : $y}]\
          [expr {$x+$side>=$w ? $w : $x+$side}]\
          [expr {$y+$side>=$h ? $h : $y+$side}]
    }
    
    proc anynode {node img side unit x y w h} {
        maybe_idleupdate
        variable where
        if {$side <= $unit} then {
            pixelnode $node $img $side $x $y $w $h
        } elseif {$x>=0 && $y>=0 && $x+$side<=$w && $y+$side<=$h} then {
            if {[info exists where($node)]} then {
                $img copy {*}$where($node) -to $x $y
            } else {
                set half [expr {$side/2}]
                anynode [lindex $::heap $node 1] $img $half $unit\
                  $x $y $w $h
                anynode [lindex $::heap $node 2] $img $half $unit\
                  [expr {$x+$half}] $y $w $h
                anynode [lindex $::heap $node 3] $img $half $unit\
                  $x [expr {$y+$half}] $w $h
                anynode [lindex $::heap $node 4] $img $half $unit\
                  [expr {$x+$half}] [expr {$y+$half}] $w $h
                set where($node) [list $img -from $x $y\
                  [expr {$x+$side}] [expr {$y+$side}]]
            }
        } elseif {$x>-$side && $y>-$side && $x<$w && $y<$h} then {
            set half [expr {$side/2}]
            anynode [lindex $::heap $node 1] $img $half $unit $x $y $w $h
            anynode [lindex $::heap $node 2] $img $half $unit\
              [expr {$x+$half}] $y $w $h
            anynode [lindex $::heap $node 3] $img $half $unit\
              $x [expr {$y+$half}] $w $h
            anynode [lindex $::heap $node 4] $img $half $unit\
              [expr {$x+$half}] [expr {$y+$half}] $w $h
        }
    }
}

proc node2photo {node args} {
    array set Opt {
        -fg black
        -bg white
        -unzoom 0
        -center {0.5 0.5}
        -shift {0 0}
        -opts {}
        -block 1
        -period infinity
    }
    array set Opt $args
    if {[info exists Opt(-image)]} then {
        set img $Opt(-image)
        $img configure {*}$Opt(-opts)
    } else {
        set img [image create photo {*}$Opt(-opts)]
    }
    
    $img blank
    set w [image width $img]
    set h [image width $img]
    tophoto::init $Opt(-fg) $Opt(-bg) $Opt(-period)
    
    set side [expr {[side $node]*$Opt(-block)}]
    for {} {$Opt(-unzoom)>0} {incr Opt(-unzoom) -1} {
        set side [expr {$side/2}]
    }
    foreach {rx ry} $Opt(-center) break
    foreach {dx dy} $Opt(-shift) break
    
    tophoto::anynode $node $img $side $Opt(-block) [expr {
        round(0.5*$w - $side*$rx) + $dx
    }] [expr {
        round(0.5*$h - $side*$ry) + $dy
    }] $w $h
    return $img
}