Words in L and R, evaluated in three monoids

"""The free monoid on {L, R}, evaluated in three monoids.

The universal property of the free monoid says: choosing a value for each
letter in a monoid M determines exactly one monoid homomorphism from words
to M, namely "multiply the values of the letters". We take M to be the 2x2
integer matrices with L = [[1,0],[1,1]] and R = [[1,1],[0,1]], check the
homomorphism property, and find that different words give different
matrices; a subtraction algorithm (Euclid's algorithm in disguise) reads the
word back off the matrix, which is why this homomorphism is injective. Other
choices of target, (N, +) and (Z, +), give homomorphisms that are far from
injective. Finally each word gives a fraction, listing positive rationals as
in the Stern-Brocot tree. Exact, exhaustive for words of bounded length.
"""

from fractions import Fraction
from itertools import product

IDENTITY = ((1, 0), (0, 1))
GENERATORS = {"L": ((1, 0), (1, 1)), "R": ((1, 1), (0, 1))}

def matmul(x, y):
    """Product of 2x2 matrices given as nested tuples."""
    return tuple(
        tuple(sum(x[i][k] * y[k][j] for k in range(2)) for j in range(2))
        for i in range(2)
    )

def evaluate(word, values, op, unit):
    """The unique homomorphism from words: multiply the letters' values."""
    result = unit
    for letter in word:
        result = op(result, values[letter])
    return result

def words(max_length, alphabet="LR"):
    """All words of length 0..max_length, shortest first."""
    out = []
    for n in range(max_length + 1):
        out.extend("".join(w) for w in product(alphabet, repeat=n))
    return out

def decode(m):
    """Recover the word of a matrix in the monoid generated by L and R.

    If the top row dominates the bottom row the word starts with R
    (strip it: subtract the bottom row from the top), otherwise with L.
    Returns None if m is not reached by any word.
    """
    word = []
    (a, b), (c, d) = m
    if min(a, b, c, d) < 0:
        return None
    # Each step subtracts a non-zero row of non-negative entries, so the
    # entries stay non-negative and their sum strictly decreases: it stops.
    while ((a, b), (c, d)) != IDENTITY:
        if a >= c and b >= d and (c, d) != (0, 0):
            word.append("R")
            a, b = a - c, b - d
        elif c >= a and d >= b and (a, b) != (0, 0):
            word.append("L")
            c, d = c - a, d - b
        else:
            return None
    return "".join(word)

def fraction_of(m):
    """The Stern-Brocot fraction (a + b) / (c + d) of a matrix."""
    (a, b), (c, d) = m
    return Fraction(a + b, c + d)

def main():
    def mat(w):
        return evaluate(w, GENERATORS, matmul, IDENTITY)

    short = words(4)
    hom = all(mat(v + w) == matmul(mat(v), mat(w)) for v in short for w in short)
    print("Target 1: 2x2 matrices, L = [[1,0],[1,1]], R = [[1,1],[0,1]]")
    print(f"  eval(vw) = eval(v) eval(w) for all {len(short)}^2 pairs of words"
          f" of length <= 4: {hom}")
    long = words(10)
    images = {mat(w) for w in long}
    print(f"  {len(long)} words of length <= 10 give {len(images)} distinct"
          " matrices")
    back = all(decode(mat(w)) == w for w in long)
    print(f"  decoding by subtraction recovers every one of these words: {back}")
    print(f"  LR = {mat('LR')} but RL = {mat('RL')}")
    print()
    length = {"L": 1, "R": 1}
    signed = {"L": 1, "R": -1}
    plus = (lambda x, y: x + y)
    for name, values in [("(N, +), L, R -> 1", length),
                         ("(Z, +), L -> 1, R -> -1", signed)]:
        image = {evaluate(w, values, plus, 0) for w in long}
        print(f"Target: {name:<24} {len(long)} words -> {len(image)} values")
    print()
    fracs = [fraction_of(mat(w)) for w in words(3)]
    print(f"Stern-Brocot: the {len(fracs)} words of length <= 3 give"
          f" {len(set(fracs))} distinct")
    print("fractions (a + b)/(c + d), where [[a, b], [c, d]] is the matrix:")
    print("  " + " ".join(str(f) for f in sorted(set(fracs))))
    print()
    print("Epistemic status: exact and exhaustive for words of length at most")
    print("10. That no two words share a matrix for all lengths is a theorem")
    print("(the decoding always works, by the argument of Euclid's algorithm);")
    print("this run is evidence for it, not a proof.")

if __name__ == "__main__":
    main()