chriswarbo-net: dd4187b169c22d8ce8cb24307d6390b9105f407e

     1: ---
     2: title: Computing Numbers
     3: packages: [ "ghcWithOmega" ]
     4: ---
     5: 
     6: This is a sequel to [an earlier blog post about computable numbers](
     7: /blog/2017-02-27-listing_things.html). In that post I gave a simple suggestion
     8: for computing numbers to unlimited precision using a Turing machine. Personally,
     9: I much prefer linguistic/algebraic approaches to computability, so I also wrote
    10: up an alternative language-based approach that I find cleaner:
    11: 
    12:  - Similar to before, we will represent a number as two never-ending "streams"
    13:    of bits, one on each side of the radix point.
    14:  - Rather than fiddling around with binary tapes, interleaving, etc. we will
    15:    keep everything abstract by using [free variables](
    16:    https://en.wikipedia.org/wiki/Free_variables_and_bound_variables). We will
    17:    use 4 variables `l0`, `l1`, `r0` and `r1`.
    18:  - Now we choose a universal programming language; I'm quite fond of
    19:    [SK combinatory logic](
    20:    https://en.wikipedia.org/wiki/SKI_combinator_calculus). SK programs are
    21:    binary trees with leaves labelled either `S` or `K`. Computation is performed
    22:    by rewriting parts of the tree which match one of the following rules:
    23: 
    24: ```
    25:   /\               x
    26:  /\ y   becomes
    27: K  x
    28: 
    29:    /\              /\
    30:   /\ z  becomes   /  \
    31:  /\ y            /\  /\
    32: S  x            x z  y z
    33: ```
    34: 
    35:  - We can list all SK trees by starting with `S` and `K`, then using a similar
    36:    procedure as Colin used to list the rationals: we list all pairs of natural
    37:    numbers $(x, y)$ such that $x + y = 0$, followed by all pairs where
    38:    $x + y = 1$, then $2$, and so on. For each $(x, y)$ pair, rather than turning
    39:    them into rationals $\frac{x}{y}$ as Colin did, we instead construct a tree
    40:    where the left and right children are the $x$th and $y$th trees from our
    41:    list. We can see that this recursive definition is well-founded: we begin
    42:    with the trees `S` and `K`, which don't require recursion, then each time we
    43:    look one-step-further in the list, we end up adding many trees to the end; so
    44:    our lookups can never 'overtake' our definitions.
    45:  - For each tree `t` in our list, we replace it with a tree of the form:
    46: 
    47: ```
    48:     /\
    49:    /\ r1
    50:   /\ r0
    51:  /\ l1
    52: t  l0
    53: ```
    54: 
    55:  - We evaluate each of these trees to (weak head) normal form by following the
    56:    above rewrite rules (this may not halt).
    57:  - We interpret the results as follows:
    58:   - We begin with a radix point "."; on its own this may be thought of as the
    59:     number 0 (equivalent to …000.000…).
    60:   - If the tree is a just a leaf, our interpretation is complete.
    61:   - Otherwise, the tree is a node with two children `x` and `y`.
    62:   - If `x` is the variable `l0`, this tree represents the same number as `y`,
    63:     but with a 0 inserted to the left of the radix point. For example, if `y`
    64:     represents a number like …00111.110011…, then this tree represents
    65:     …001110.110011….
    66:   - Likewise if `x` is the variable `l1`, a 1 is insert to the left of the radix
    67:     point.
    68:   - Similarly, `r0` and `r1` insert a 0 or 1 to the *right* of the radix point,
    69:     respectively.
    70:   - If `x` has any other form, it is ignored; this tree represents the same
    71:     number as `y`.
    72: 
    73: I've made a quick implementation of this in Haskell; here's the beginning of the
    74: list:
    75: 
    76: ```{.haskell pipe="runhaskell"}
    77: {-# LANGUAGE BangPatterns, MonadComprehensions #-}
    78: module X where
    79: 
    80: import Control.Monad
    81: import Control.Monad.Omega
    82: import Data.List
    83: import Data.Maybe
    84: 
    85: data Comb = S | K | Node Comb Comb | V Var deriving (Eq, Show)
    86: 
    87: data Var = Pair | Zero | One deriving (Eq, Show)
    88: 
    89: step (Node (Node K x) y)          = (True, x)
    90: step (Node (Node (Node S x) y) z) = (True, Node (Node x z) (Node y z))
    91: step (Node x y)                   = case step x of
    92:                                          (True, x') -> (True, Node x' y)
    93:                                          (False, _) -> Node x <$> step y
    94: step x                            = (False, x)
    95: 
    96: reduce 0 x = Nothing
    97: reduce n x = case step x of
    98:                  (True, x') -> reduce (n-1) x'
    99:                  (False, _) -> Just x
   100: 
   101: trees = S : K : runOmega [Node x y | x <- each trees, y <- each trees]
   102: 
   103: readTree n t = do t' <- reduce n (Node t (V Pair))
   104:                   case t' of
   105:                        Node (Node (V Pair) l) r -> readStreams n l r
   106:                        _                        -> Nothing
   107: 
   108: readStreams n l r = do l' <- readStream n l
   109:                        r' <- readStream n r
   110:                        pure (l', r')
   111: 
   112: readStream n = go [] 8
   113:   where go !acc 0 _ = Just acc
   114:         go !acc m x = do t' <- reduce n (Node (Node x (V Zero)) (V One))
   115:                          case t' of
   116:                               Node (V Zero) rest -> go (acc ++ [0]) (m-1) rest
   117:                               Node (V One)  rest -> go (acc ++ [1]) (m-1) rest
   118:                               _                  -> Just acc
   119: 
   120: evaledTrees = map (readTree 10) trees
   121: 
   122: bitsToNum :: ([Int], [Int]) -> Double
   123: bitsToNum (l, r) = parseWhole 0 1 (reverse l) + parseFrac 0 0.5 r
   124:   where parseWhole !n !power []     = n
   125:         parseWhole !n !power (0:xs) = parseWhole  n          (power * 2) xs
   126:         parseWhole !n !power (1:xs) = parseWhole (n + power) (power * 2) xs
   127: 
   128:         parseFrac  !n !power []     = n
   129:         parseFrac  !n !power (0:xs) = parseFrac   n          (power / 2) xs
   130:         parseFrac  !n !power (1:xs) = parseFrac  (n + power) (power / 2) xs
   131: 
   132: nubble = go []
   133:   where go seen []     = []
   134:         go seen (x:xs) = if x `elem` seen
   135:                             then     go    seen  xs
   136:                             else x : go (x:seen) xs
   137: 
   138: -- FIXME: We're only taking 1 since even taking 5 ends up blowing the memory
   139: -- after a while. The problem is that we're throwing away so many programs; I'd
   140: -- like to think of a more direct encoding, but haven't come up with one that's
   141: -- functional (e.g. we could read the S and K occurrences in post-order, but
   142: -- that violates things like beta-equality).
   143: render = intercalate ", " . take 1 . nubble .
   144:          -- Turn into a list of strings of numbers
   145:          map (show . bitsToNum) .
   146:          -- Discard anything that failed to normalise
   147:          catMaybes
   148: 
   149: main = putStr . render $ evaledTrees
   150: ```

Generated by git2html.