Listing the integers, pairs and rationals

"""Listing the integers, the grid N x N and the positive rationals.

A set is countable when its elements can be listed a_0, a_1, a_2, ... with
every element appearing. This experiment runs three explicit listings:

* the zigzag 0, -1, 1, -2, 2, ... of the integers;
* Cantor's pairing function pi(m, n) = (m + n)(m + n + 1)/2 + n, which walks
  the grid N x N diagonal by diagonal, together with its inverse;
* the Calkin-Wilf sequence 1, 1/2, 2, 1/3, 3/2, 2/3, 3, ..., generated by
  Newman's rule x -> 1 / (2 floor(x) - x + 1), which lists every positive
  rational number exactly once.

Runs unchanged in the browser (Pyodide) using only the standard library.
"""

from fractions import Fraction
from math import gcd, isqrt

def zigzag(n):
    """The n-th integer in the list 0, -1, 1, -2, 2, ..."""
    if n < 0:
        raise ValueError("n must be a natural number")
    return n // 2 if n % 2 == 0 else -(n // 2) - 1

def cantor_pair(m, n):
    """Cantor's pairing function: the position of (m, n) in the diagonal walk."""
    d = m + n
    return d * (d + 1) // 2 + n

def cantor_unpair(k):
    """Inverse of ``cantor_pair``: find the diagonal d with
    d(d + 1)/2 <= k < (d + 1)(d + 2)/2, then the position n along it."""
    if k < 0:
        raise ValueError("k must be a natural number")
    d = (isqrt(8 * k + 1) - 1) // 2
    n = k - d * (d + 1) // 2
    return d - n, n

def calkin_wilf(count):
    """The first ``count`` terms of the Calkin-Wilf sequence (Newman's rule)."""
    terms = []
    x = Fraction(1)
    for _ in range(count):
        terms.append(x)
        x = 1 / (2 * (x.numerator // x.denominator) - x + 1)
    return terms

def first_time_all_seen(terms, bound):
    """Number of terms needed before every reduced p/q with 1 <= p, q <= bound
    has appeared, or None if some never appears in ``terms``."""
    wanted = {Fraction(p, q) for p in range(1, bound + 1)
              for q in range(1, bound + 1) if gcd(p, q) == 1}
    for i, x in enumerate(terms):
        wanted.discard(x)
        if not wanted:
            return i + 1
    return None

def main():
    print("Zigzag listing of Z:", ", ".join(str(zigzag(n)) for n in range(11)),
          "...")
    print()
    print("Cantor's pairing pi(m, n), rows m = 0..4, columns n = 0..5:")
    for m in range(5):
        print("  m=%d: " % m + " ".join(f"{cantor_pair(m, n):3d}"
                                       for n in range(6)))
    limit = 100_000
    ok = all(cantor_pair(*cantor_unpair(k)) == k for k in range(limit))
    ok2 = all(cantor_unpair(cantor_pair(m, n)) == (m, n)
              for m in range(200) for n in range(200))
    print(f"pi(unpair(k)) = k for all k < {limit}: {ok}")
    print(f"unpair(pi(m, n)) = (m, n) for all m, n < 200: {ok2}")
    print()
    terms = calkin_wilf(2047)
    print("Calkin-Wilf:", ", ".join(str(x) for x in terms[:13]), "...")
    distinct = len(set(terms)) == len(terms)
    print(f"First {len(terms)} terms all different: {distinct}")
    n_all = first_time_all_seen(terms, 10)
    print(f"All 63 reduced p/q with p, q <= 10 appear within {n_all} terms.")
    print(f"Where is 3/5? Term number {terms.index(Fraction(3, 5)) + 1}.")
    print()
    print("Epistemic status: exact computation on finitely many values. It")
    print("illustrates the listings and checks them as far as it goes; that")
    print("they list every integer, pair and positive rational exactly once is")
    print("proved in the lesson (and for Calkin-Wilf in Calkin and Wilf, 2000).")

if __name__ == "__main__":
    main()