#!/usr/bin/env python3
"""
walsh_higher_commutators.py
============================
Higher-order Dirac commutators of Walsh modes.

Tests whether ||ad_D^k(chi_S)||_op grows like (2 sqrt(|S|))^k, so that the
Frechet smooth Lip-norm L_smooth(T) := sup_k ||ad_D^k(T)|| / k! goes like
e^{2 sqrt(|S|)} -- still subexponential in |S|, hence NOT sufficient to
close L_round.

D = sum_a chi_a is the discrete Dirac on Boolean-CAR substrate.
ad_D(T) := [D, T].
"""
import math
import numpy as np
import numpy.linalg as la

sz = np.array([[1, 0], [0, -1]], dtype=complex)
sx = np.array([[0, 1], [1, 0]], dtype=complex)
I2 = np.eye(2, dtype=complex)


def kron_chain(ops):
    out = ops[0]
    for op in ops[1:]:
        out = np.kron(out, op)
    return out


def chi_clifford(D, a):
    return kron_chain([sz] * a + [sx] + [I2] * (D - 1 - a))


def chi_walsh(D, S):
    return kron_chain([sz if a in S else I2 for a in range(D)])


def dirac(D):
    N = 1 << D
    out = np.zeros((N, N), dtype=complex)
    for a in range(D):
        out = out + chi_clifford(D, a)
    return out


def comm(A, B):
    return A @ B - B @ A


def op_norm(M):
    return float(la.norm(M, ord=2))


# Sweep ad_D^k(chi_S) for various D, S, k.
KMAX = 10

print("=" * 92)
print("  walsh_higher_commutators.py  --  higher-order Dirac commutators of Walsh modes")
print("=" * 92)
print()
print("  ad_D^k(chi_S) := [D, [D, ..., [D, chi_S]...]]   k-fold,  D = sum_a chi_a")
print()

for D in (4, 6, 8):
    D_op = dirac(D)
    for k in (1, 2, 3, D - 1, D):
        S = tuple(range(k))
        chi_S_op = chi_walsh(D, S)
        T = chi_S_op
        print(f"  D = {D:>2},  |S| = {k}:  ad_D^j norms for j = 1..{KMAX} (and (2sqrt|S|)^j for reference)")
        norms = []
        for j in range(1, KMAX + 1):
            T = comm(D_op, T)
            n = op_norm(T)
            norms.append(n)
        # Reference  (2 sqrt|S|)^j
        refs = [(2 * math.sqrt(k)) ** j for j in range(1, KMAX + 1)]
        print(f"           {'j':>3} {'||ad_D^j chi_S||':>20} {'(2 sqrt|S|)^j':>16} {'ratio':>10} {'log_ratio /j':>14}")
        for j in range(1, KMAX + 1):
            n = norms[j - 1]
            r = refs[j - 1]
            if r > 0:
                ratio = n / r
                # log_growth = log(n) / j -- effective exponential rate
                log_growth = math.log(n) / j if n > 0 else float('-inf')
                print(f"           {j:>3} {n:>20.4f} {r:>16.4f} {ratio:>10.4f} {log_growth:>14.4f}")
        print()
    print()

print("=" * 92)
print("  Interpretation")
print("=" * 92)
print()
print("  If ||ad_D^j chi_S|| ~ (2 sqrt|S|)^j, then the Frechet smooth Lip-norm")
print("    L_smooth(chi_S) := sup_j ||ad_D^j chi_S|| / j!  ~  e^{2 sqrt|S|}")
print("  by Stirling.")
print()
print("  For L_round to close under L_smooth, we'd need 2^D / L_smooth(T_const_tail) -> 0.")
print("  L_smooth(T_const_tail) <= D-many terms with sup_S L_smooth ~ e^{2 sqrt D}.")
print("  Ratio: 2^D / e^{2 sqrt D}  =  e^{D log 2  -  2 sqrt D}.")
print("  For large D,  D log 2 >> 2 sqrt D,  so ratio BLOWS UP.")
print()
print("  Hence even the Frechet smooth Lip-norm is NOT sufficient to close L_round")
print("  for arbitrary T in A_D.  The exponential Lip-norm e^{c|S|} (c > log 2)")
print("  remains the structural minimum.")