From d6b615c2c1a5c82e7529a30e07a52e2bbdc6b27b Mon Sep 17 00:00:00 2001
From: wan-may <wnmy@protonmail.com>
Date: Sat, 7 Oct 2023 11:52:55 -0300
Subject: [PATCH] Big commit, basic functionality + plotting.

---
 backend.lua      |   5 ++
 conf.lua         |   3 ++
 enumeration.lua  |   7 +++
 enumerations.lua |  88 +++++++++++++++++++++++++++++++
 main.lua         |  92 +++++++++++++++++++++++++++++++++
 rapport.lua      | 132 +++++++++++++++++++++++++++++++++++++++++++++++
 sequence.lua     |  25 +++++++++
 twintree.lua     | 103 ++++++++++++++++++++++++++++++++++++
 8 files changed, 455 insertions(+)
 create mode 100644 backend.lua
 create mode 100644 conf.lua
 create mode 100644 enumeration.lua
 create mode 100644 enumerations.lua
 create mode 100644 main.lua
 create mode 100644 rapport.lua
 create mode 100644 sequence.lua
 create mode 100644 twintree.lua

diff --git a/backend.lua b/backend.lua
new file mode 100644
index 0000000..eec592e
--- /dev/null
+++ b/backend.lua
@@ -0,0 +1,5 @@
+--Dependencies for displaying things.
+local backend = {}
+
+
+return backend
\ No newline at end of file
diff --git a/conf.lua b/conf.lua
new file mode 100644
index 0000000..948ef2d
--- /dev/null
+++ b/conf.lua
@@ -0,0 +1,3 @@
+function love.conf( t )
+  t.console = true
+end
\ No newline at end of file
diff --git a/enumeration.lua b/enumeration.lua
new file mode 100644
index 0000000..0ffc62e
--- /dev/null
+++ b/enumeration.lua
@@ -0,0 +1,7 @@
+--Lazy enumerator of countable infinite set.
+--Pass a function closure for generating the set element-by-element, in order.
+--This function can return a "simple" value followed by an "exact" value
+--Returns a table including:
+-- an array with all elements generated so far (grows as needed upon access)
+local enum = { __index = function( t, n ) for i = #t + 1, n do t[i] = t.next() end return rawget(t, n) end }
+return function(nextElement) return setmetatable( { ["next"] = nextElement }, enum ) end
\ No newline at end of file
diff --git a/enumerations.lua b/enumerations.lua
new file mode 100644
index 0000000..1c20993
--- /dev/null
+++ b/enumerations.lua
@@ -0,0 +1,88 @@
+local enum = require 'enumeration'
+
+local function CalkinWilf()
+  local n, d = 1, 1
+  return enum( function()
+    repeat n, d = d, ( 2 * d * math.floor( n / d ) - n + d )
+    until n < d
+    return n / d
+  end )
+end
+
+local function Random()
+  return enum( math.random )
+end
+
+--Breadth first traversal of the open left subtree of the Stern-Brocot tree.
+--Formed by inserting the mediants of previous rationals in the tree.
+local function SternBrocot()
+
+  --Rational numbers.
+  local R
+  do
+    local rationals = {}
+    local rmt = {
+      __lt = function(a, b) return a[1]/a[2] < b[1]/b[2] end,
+      __eq = function(a, b) return a[1] == b[1] and a[2] == b[2] end,
+      __pow = function(a, b) return R( a[1] + b[1], a[2] + b[2] ) end, --mediant
+      __call = function(t) return t[1]/t[2] end,
+      __tostring = function(t) return string.format( "%d / %d", t[1], t[2] ) end }
+
+    R = function( n, d )
+      if not rationals[d] then rationals[d] = {} end
+      local den = rationals[d]
+
+      if den[n] then return den[n] end
+
+      local r = setmetatable( { n , d }, rmt )
+      den[n] = r
+      return r
+    end
+  end
+
+  --Flat list of all entries traversed so far.
+  local seq = { R(0, 1), R(1, 2), R(1, 1) }
+  local idx = 0
+  local len = 3
+
+  local function nextRational()
+    --Return every other entry: these are the most recently added mediants.
+    idx = idx + 2
+    if idx < len then return seq[idx] end
+
+    --We have finished iterating over the current level of the tree.
+    --Create a new level by adding all the new mediants.
+    
+    local nSeq = {}
+    for i = 1, len - 1 do
+      nSeq[ 2 * i - 1 ] = seq[i]
+      nSeq[ 2 * i ] = seq[i] ^ seq[i + 1]
+    end
+    nSeq[ 2 * len - 1 ] = seq[ len ]
+    seq = nSeq
+    len = #seq
+    --print( "RATIONALS: NEW LEVEL", len, unpack( nSeq ) )
+    idx = 0
+    return nextRational()
+  end
+
+  return enum( nextRational )
+end
+
+--Enumeration of dyadic rationals in the open unit interval, ascending lexicographic.
+local function Dyad()
+  local n, d = -1, 2
+  return enum( function()
+      n = n + 2
+      if n > d then n, d = 1, d * 2 end
+      return n / d
+    end )
+end
+
+return { 
+  Dyad = Dyad, 
+  SternBrocot = SternBrocot,
+  CalkinWilf = CalkinWilf,
+  Random = Random,
+  QuadIrr = false,
+}
\ No newline at end of file
diff --git a/main.lua b/main.lua
new file mode 100644
index 0000000..8d037d3
--- /dev/null
+++ b/main.lua
@@ -0,0 +1,92 @@
+local enumerations = require 'enumerations'
+local twintree = require 'twintree'
+local step = coroutine.wrap( twintree.buildIncremental )
+local a = enumerations.Dyad()
+local b = enumerations.Random()
+step( a, b )
+
+--Use LOVE for plotting.
+local love = assert( love )
+
+local canvas = love.graphics.newCanvas()
+local numMapped, pointList = 0
+
+local function swap( i, j )
+  b[i], b[j] = b[j], b[i]
+end
+
+local function reset()
+  print( "changing coroutine" )
+  step = coroutine.wrap( twintree.buildIncremental )
+  print( "swapping" )
+  for i = 1, 10000 do
+    swap( i, math.random( 1, 10000 ) )
+  end
+  
+  pointList = nil
+  numMapped = 0
+  print( "step" )
+  step( a, b )
+end
+
+local function new()
+  pointList = step()
+  numMapped = numMapped + 1
+end
+
+local function paint()
+  local type = type
+  --Make sure the whole point list consists of numbers!
+  --We use the __call metamethod in non-number elements to coerce to numbers.
+  --TODO: make this more generic, maybe amenable to FFI types.
+  for i = 1, #pointList do
+    if type( pointList[i] ) ~= 'number' then pointList[i] = pointList[i]() end 
+  end
+  
+  --print( "Points: " )
+  --for i = 1, #pointList - 1, 2 do
+  --  print(pointList[i], pointList[i + 1] )
+  --end
+  
+  love.graphics.setCanvas( canvas )
+  love.graphics.clear()
+  love.graphics.setColor( 1, 1, 1, 0.5 )
+  love.graphics.print( numMapped )
+  love.graphics.push()
+  love.graphics.scale( love.graphics.getWidth(), love.graphics.getHeight() )
+  love.graphics.translate( 0, 1 )
+  love.graphics.scale( 1, -1 )
+  love.graphics.points( pointList )
+  love.graphics.line( pointList )
+  love.graphics.pop()
+  love.graphics.setCanvas()
+end
+
+function love.wheelmoved()
+  for i = 1, 3 do new() end
+  paint()
+ --[[ for i = 1, #sb do
+    print( i, sb[i], dy[i] )
+  end]]
+end
+
+function love.mousepressed()
+  --[[for i = 1, #sb do
+    print( i, sb[i] )
+  end]]
+end
+
+function love.draw()
+  love.graphics.setCanvas()
+  love.graphics.setColor( 1, 1, 1, 1 )
+  love.graphics.draw( canvas )
+end
+
+function love.keypressed()
+  print( "pressed key!" )
+  return reset() and love.wheelmoved()
+end
+
+function love.load()
+  love.graphics.setLineWidth( 0.001 )
+end
\ No newline at end of file
diff --git a/rapport.lua b/rapport.lua
new file mode 100644
index 0000000..dc9261e
--- /dev/null
+++ b/rapport.lua
@@ -0,0 +1,132 @@
+local function Gap(t, i)
+  if not t[i] then return i else return Gap(t, i + 1) end
+end
+
+local function Gaps(t, i)
+  return Gap, t, i
+end
+
+--Data structure supporting fast ordered insertion and search.
+local function Order()
+  local o = {}
+  
+  local list = { val = nil, prev = nil, next = nil }
+  
+  o.Get = function( x )
+    
+  end
+  
+  o.GetAdjacent = function( node )
+    return node.prev.val, node.next.val
+  end
+  
+  o.Insert = function( x )
+    if not list.val then
+      list.val = x
+      return list
+    end
+    
+    local node = list
+    while node and node.val < x do node = node.next end
+    
+    
+  end
+  
+  return o
+end
+
+--Construct an isomorphism between two countable dense total orders
+--these orders are specified by functions which return the next element in some enumeration
+--these elements must be comparable with comparison operators <, >, &c.
+local function Isomorphism( Enumeration, InverseEnumeration, inv )
+
+  --Exported table.
+  local t = {}
+  local enum, order, mapping, inverse
+
+  --Another Isomorphism structure,
+  --which contains the image of this isomorphism and the inverse
+  --isomorphism mapping back to our dense countable linear order.
+  inverse = inv or Isomorphism( InverseEnumeration, nil, t )
+
+  --Array for caching and lazily evaluating the supplied Enumeration,
+  --an iterator that enumerates the elements of a dense countable linear order.
+  enum = { Enumeration() }
+  setmetatable(enum, {__index = function(_, i)
+        for j = #enum + 1, i do
+          enum[j] = Enumeration()
+        end
+        return rawget( enum, i )
+      end})
+  t.enum = enum --Expose as field 
+
+  --Structure supporting fast ordered insertion and search.
+  --Values are i such that mapping[order[i]] is not nil,
+  --sorted according to getmetatable( enum[i] ).__lt
+  order = { 1 }
+
+  --Find the indices i and k such that
+  --enum[i] is the greatest lower bound on enum[j] in the mapping's support
+  --enum[k] is the least upper bound on enum[j] in the mapping's support
+  local function GetBounds( leaf )
+
+  end
+
+  --Array defining isomorphism mapping this dense countable order to the one in inverse.
+  --If mapping[i] = j, then f(enum[i]) = inverse.enum[j]
+  mapping = { 1 }
+  setmetatable( mapping, {__newindex = function(map, k, v)
+    rawset(mapping, k, v)
+    
+    
+      end})
+
+  --Iterates over the natural numbers outside the isomorphism's support,
+  --so returns the least n such that n >= i and mapping[n] == nil.
+  --(Warning, infinite loop! Make sure to break out of this.)
+  local function NextGap(i)
+    if not mapping[i] then return i else return NextGap(i + 1) end --tail-recursive
+  end
+  
+  --Include another point in the mapping's support.
+  --Infinite loop: ))<<>>((
+  function t.ExpandMap()
+    local gapIdx = NextGap(1)
+    local pred, succ = GetBounds( OrderedInsertion( gapIdx ) )
+    mapping[gapIdx] = inverse.AddNextBetween( i, mapping[pred], mapping[succ] )
+    return inverse.ExpandMap()
+  end
+
+  --Find the first j such that mapping[j] == nil,
+  --enum[lower] < enum[j] < enum[upper]
+  --Then map j to i, and add j to the order.
+  function t.AddNextBetween( i, lower, upper )
+    local p, q, r, j = enum[lower], nil, enum[upper], nil
+
+    for k in Gaps( 1 ) do
+      q = enum[k]
+      if  ((not p) or (p < q))
+      and ((not r) or (q < r)) then
+        j = k
+        break
+      end
+    end
+
+    mapping[j] = i
+    OrderedInsertion(j)
+    return j
+  end
+
+  --Get the isomorphism as a dictionary, for plotting.
+  function t.GetMapping( )
+    local map = {}
+    for idx, imageIdx in pairs( mapping ) do
+      map[ enum[ idx ] ] = inverse.enum[ imageIdx ]
+    end
+    return map
+  end
+
+  return t
+end
+
+return setmetatable( { Isomorphism = Isomorphism }, {__call = Isomorphism } )
\ No newline at end of file
diff --git a/sequence.lua b/sequence.lua
new file mode 100644
index 0000000..e4ce1f4
--- /dev/null
+++ b/sequence.lua
@@ -0,0 +1,25 @@
+--Partial sequence s:S->N, S <= N, N = { 1, 2, 3, ... }
+--Consists of:
+  --A sequence stored as a (sparse) array
+  --["cmin"], a key whose value is the smallest element outside the support of the sequence
+local sequence = {}
+sequence.__index = sequence
+function sequence.new( name )
+  return setmetatable( { cmin = 1, name = name }, sequence )
+end
+
+function sequence:insert( x, y )
+  if self.cmin == x then
+    repeat self.cmin = self.cmin + 1 until not self[self.cmin]
+  end
+  self[x] = y
+  --debug:
+  --print( "First unclaimed element:", self.name, self.cmin )
+end
+
+--Maybe leave this for later, we'll inline it for now.
+function sequence:complement()
+  
+end
+
+return sequence
\ No newline at end of file
diff --git a/twintree.lua b/twintree.lua
new file mode 100644
index 0000000..7fd95c3
--- /dev/null
+++ b/twintree.lua
@@ -0,0 +1,103 @@
+local sequence = require 'sequence'
+
+local maxLevels = 20 --Maximum depth of tree.
+local baseIdx = 1 --Base index of implicit tree. Is it a 0-based or 1-based array?
+local tree = {}
+
+local function insert( aTree, aSeq, aVals, bTree, bSeq, bVals )
+
+  --print( "INSERTING" )
+  local treeIdx = 1 --Pointer into implicit binary tree.
+  local aIdx, bIdx = aSeq.cmin, bSeq.cmin --First elements of a and b outside the support of the sequence
+  local aVal, bVal = aVals[aIdx], bVals[bIdx]
+  local aTest, bTest
+
+  --Search for appropriate bIdx.
+  --Sentinel value.
+  local newCandidate = false
+  while aTree[treeIdx] do
+    if newCandidate then
+      newCandidate = false
+      treeIdx = 1
+      repeat bIdx = bIdx + 1 until ( not bSeq[bIdx] )
+      bVal = bVals[bIdx]
+    end
+
+    --invariants guarantee these aren't empty:
+    --loop invariant, that aTree[treeIdx] is nonempty
+    --struct invariant, that bTree and aTree have identical support
+    aTest, bTest = aTree[treeIdx], bTree[treeIdx] 
+    --Go left if bounded above, go right if bounded below,
+    --if there is a discrepancy then start over with the next bIdx
+    --print( "VALUES:", "A[[", aVal, aTest, "B[[", bVal, bTest )
+    if ( aVal < aTest ) then
+      if ( bVal < bTest ) then 
+        --print( "A LEFT B LEFT", aVal, aTest, bVal, bTest ) --debug
+        treeIdx = 2 * treeIdx
+      else
+        --print( "A LEFT B RIGHT", aVal, aTest, bVal, bTest ) --debug
+        newCandidate = true
+      end
+    elseif ( bVal > bTest ) then 
+      --print( "A RIGHT B RIGHT", aVal, aTest, bVal, bTest ) --debug
+      treeIdx = 2 * treeIdx + 1 
+    else 
+      --print( "A RIGHT B LEFT", aVal, aTest, bVal, bTest ) --debug
+      newCandidate = true
+    end
+  end
+
+  --Insert elements into tree.
+  aTree[treeIdx], bTree[treeIdx] = aVal, bVal
+
+  --Insert indices into sequence.
+  aSeq:insert( aIdx, bIdx )
+  bSeq:insert( bIdx, aIdx )
+
+end
+
+--In-order traversal of the tree. 
+--Returns points of twin-tree, sorted and interlaced, for graphing.
+-- { a_min, b_min, a_next, b_next, ...}
+do
+  local a, b, n, result
+  local function flat( idx )
+    if a[ 2 * idx ] then flat( 2 * idx ) end
+    result[ n ], result[ n + 1 ], n = a[idx], b[idx], n + 2
+    if a[ 2 * idx + 1 ] then flat( 2 * idx + 1 ) end
+  end
+
+  function tree.flat( aTree, bTree, t )
+    result = t or {}
+    a, b = aTree, bTree
+    n = 1
+    flat( 1 )
+    for i = n + 1, #result do result[i] = nil end
+    return result
+  end
+end
+
+function tree.build( n, av, bv )
+  local at, bt = {}, {}
+  local as, bs = sequence.new 'a', sequence.new 'b'
+  for i = 1, n do
+    insert( at, as, av, bt, bs, bv ) --Forward map.
+    insert( bt, bs, bv, at, as, av ) --Backward map.
+  end
+  return at, bt
+end
+
+do
+  local yield = assert( coroutine.yield )
+  function tree.buildIncremental( av, bv )
+    local at, bt = {}, {}
+    local as, bs = sequence.new 'a', sequence.new 'b'
+    local flattened = {}
+    repeat
+      insert( at, as, av, bt, bs, bv ) --Forward map.
+      insert( bt, bs, bv, at, as, av ) --Backward map.
+    until yield( tree.flat( at, bt, flattened ) ) and false
+  end
+end
+
+return tree
\ No newline at end of file