mismatch.py

1from fractions import Fraction2from itertools import combinations, product34import sympy as sp56DIM = 27BASES = (2, 3)8TOP = 24910def wallis(top):11    even = {2: Fraction(16, 3)}12    odd = {3: Fraction(3, 8)}13    e, o = even[2], odd[3]14    k = 115    while 2 * k + 3 <= top + 1:16        e *= Fraction(4 * k * (k + 1), (2 * k + 1) * (2 * k + 3))17        o *= Fraction((2 * k + 1) * (2 * k + 3), (2 * k + 2) * (2 * k + 4))18        if 2 * k + 2 <= top:19            even[2 * k + 2] = e20        if 2 * k + 3 <= top:21            odd[2 * k + 3] = o22        k += 123    return even, odd2425def squarefree_divisors(q):26    primes = sp.primefactors(q)27    out = []28    for size in range(len(primes) + 1):29        for pick in combinations(primes, size):30            value = 131            for p in pick:32                value *= p33            out.append((value, (-1) ** size))34    return out3536def bracket(design, q):37    total = len(design)38    acc = Fraction(0)39    for e, mu in squarefree_divisors(q):40        count = sum(1 for cell in design if all(c % e == 0 for c in cell))41        acc += mu * Fraction(count, total)42    return acc4344def base_factor(q, dim):45    factor = Fraction(1)46    for p in sp.primefactors(q):47        factor /= 1 - Fraction(1, p**dim)48    return factor4950def designs(q, dim):51    cells = list(product(range(q), repeat=dim))52    for mask in range(1, 1 << len(cells)):53        picked = [cells[i] for i in range(len(cells)) if mask >> i & 1]54        if len(picked) >= 2:55            yield tuple(picked)5657def main():58    even, odd = wallis(TOP)59    print("MISMATCH: design coprime densities against the Version L family")60    print("  bracket B(F) = Sum_{e | rad(q)} mu(e) k_e / k")61    print("  delta = B(F) (1/zeta(D)) Prod_{p | q} (1 - p^-D)^-1")62    print(f"  at D = {DIM} that is delta = numerator/Pi^2 with numerator = 6 B(F) factor")63    print(f"  Version L swept at d = 2..{TOP}")6465    numerators = {}66    for q in BASES:67        factor = base_factor(q, DIM)68        seen = {}69        count = 070        for design in designs(q, DIM):71            count += 172            num = 6 * bracket(design, q) * factor73            seen.setdefault(num, 0)74            seen[num] += 175        numerators[q] = seen76        shown = ", ".join(str(n) for n in sorted(seen))77        print(f"  base {q}: {count} designs, {len(seen)} distinct numerators: {shown}")7879    print("  even d: numerator against numerator")80    hits = []81    for q in BASES:82        for d, r in sorted(even.items()):83            if r in numerators[q]:84                hits.append((q, d, r, numerators[q][r]))85    for q, d, r, mult in hits:86        print(f"  match: base {q}, d = {d}, numerator {r}, carried by {mult} designs")87    print(f"  total even-d matches over d = 2..{TOP}: {len(hits)}")8889    zero_brackets = sum(1 for q in BASES for n in numerators[q] if n == 0)90    print(f"  designs of zero bracket in these bases: {zero_brackets}")91    print("  odd d: Version L is a nonzero rational and every density here is")92    print("  a nonzero rational over Pi^2, so a match would make Pi^2 rational")93    print(f"  odd Version L values swept: {sorted(odd)}")9495    print("  the one match in closed form")96    for q, d, r, _ in hits:97        print(f"  base {q}, d = {d}: ({r})/Pi^2 = {sp.N(sp.Rational(r.numerator, r.denominator) / sp.pi**2, 16)}")9899main()