#!/usr/bin/env python3
"""
Computation 35 -- Mosco-limit reweighting probe via the holomorphic Toeplitz subalgebra
========================================================================================
Computation 34 documented the structural obstruction in the simplest
multiplicative substrate -> Toeplitz bridge: Walsh modes commute among
themselves (the substrate Walsh subalgebra is abelian), but the
multiplicative image of mixed holomorphic / anti-holomorphic
generators in the Toeplitz algebra does not commute, introducing
spurious cross-commutators at weight k >= 2.

This computation tests the natural structural fix: restrict the
bridge image to the HOLOMORPHIC Toeplitz subalgebra of H^2_alpha(B^2).
The holomorphic Toeplitz operators T_{f(z)} for f a holomorphic
polynomial COMMUTE among themselves (multiplication by holomorphic
functions on H^2_alpha is just multiplication, no projection needed,
and ordinary multiplication is commutative).  Mapping Walsh modes
into this subalgebra preserves the abelian structure by construction.

CONSTRUCTION (Mosco-limit reweighting on the holomorphic subalgebra).

For substrate size D and Bergman truncation N, choose a multiplicity-
matched mapping:

  - At weight k = 0:  chi_emptyset -> T_1 = identity.
  - At weight k = 1:  the D Walsh modes chi_{a} map to the D linearly-
    independent degree-1 holomorphic monomials, scaled to unit Bergman
    norm.  At D = 4, N = 5 with B^2:
        chi_{0} -> T_{z_0},
        chi_{1} -> T_{z_1},
        chi_{2} -> T_{z_0^2}  (degree-2, but linearly independent of T_{z_0}),
        chi_{3} -> T_{z_1^2}.
    (For D > 4 the assignment continues into higher monomials.)
  - At weight k >= 2:  chi_S maps to the PRODUCT of the individual-site
    images (always a holomorphic Toeplitz operator, hence in the abelian
    subalgebra by construction).

This bridge is ABELIAN-PRESERVING by construction: all images commute,
so the spurious cross-commutators of Computation 34 vanish.

GAP TEST.
For the Clifford-generator side of the L_comm gap, the natural
Toeplitz analog is T_{bar z_a} (a single anti-holomorphic Toeplitz).
This operator does NOT commute with holomorphic Toeplitz operators
(Computation 33's commutator formula).  The gap test is:

  ||[chi_a^{Cliff}, chi_S]||_sub  vs  ||[T_{bar z_a}, b(chi_S)]||_Toep

Substrate side from Computation 30: 2 if a in S, 0 otherwise.
Toeplitz side: computed numerically.

EXPECTED.
If the abelian structure preservation + holomorphic embedding is the
right framework, the Toeplitz commutator norms should match the
substrate norms in some rescaled sense, with alpha = alpha(D) tuned
per Computation 33.

OUTPUT.
For D = 4, N = 5, sweep alpha and report:
  - Bridge image norms (fidelity).
  - Toeplitz commutator norms (gap).
  - Comparison to substrate commutator norms.
  - Identifies whether the holomorphic abelian subalgebra is the
    correct target for the Mosco-limit reweighting.
"""
import math
from itertools import combinations

import numpy as np
import numpy.linalg as la


# ============================================================
#   Substrate
# ============================================================

sx = np.array([[0, 1], [1, 0]], dtype=complex)
sy = np.array([[0, -1j], [1j, 0]], dtype=complex)
sz = np.array([[1, 0], [0, -1]], 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_walsh(D, S):
    return kron_chain([sz if a in S else I2 for a in range(D)])


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


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


# ============================================================
#   Toeplitz side
# ============================================================

def basis_indices(N):
    return [(m1, m2) for m1 in range(N + 1) for m2 in range(N + 1 - m1)]


def log_norm_sq(m1, m2, alpha):
    return (math.lgamma(m1 + 1) + math.lgamma(m2 + 1)
            + math.lgamma(alpha + 3) - math.lgamma(m1 + m2 + alpha + 3))


def Tz_matrix(a, basis, alpha):
    n = len(basis)
    idx = {J: i for i, J in enumerate(basis)}
    M = np.zeros((n, n), dtype=complex)
    for j, J in enumerate(basis):
        Jp = list(J)
        Jp[a] += 1
        Jp = tuple(Jp)
        if Jp in idx:
            log_ratio = 0.5 * (log_norm_sq(Jp[0], Jp[1], alpha)
                               - log_norm_sq(J[0], J[1], alpha))
            M[idx[Jp], j] = math.exp(log_ratio)
    return M


def Tzbar_matrix(a, basis, alpha):
    return Tz_matrix(a, basis, alpha).conj().T


# ============================================================
#   Bridge: holomorphic Toeplitz subalgebra
# ============================================================

def site_to_holo(D, basis, alpha, Tz):
    """Map each substrate site to a holomorphic Toeplitz operator.

    At D = 4, B^2: sites 0, 1 -> T_{z_0}, T_{z_1}; sites 2, 3 -> T_{z_0^2}, T_{z_1^2}.
    """
    return {
        0: Tz[0],
        1: Tz[1],
        2: Tz[0] @ Tz[0],
        3: Tz[1] @ Tz[1],
    }


def bridge_holo(D, S, basis, alpha, Tz):
    """Bridge Walsh mode chi_S -> product of individual-site holomorphic images.

    All factors are holomorphic Toeplitz operators, which commute, so the
    product is well-defined regardless of factor order.
    """
    site_map = site_to_holo(D, basis, alpha, Tz)
    n = len(basis)
    out = np.eye(n, dtype=complex)
    for a in sorted(S):
        out = out @ site_map[a]
    return out


# ============================================================
#   Probe
# ============================================================

def main():
    D = 4
    N = 5
    alphas = [0.0, 0.5, 1.0, 2.0, 5.0, 10.0]

    print("=" * 90)
    print("  Computation 35  --  Mosco-limit reweighting via holomorphic Toeplitz subalgebra")
    print("=" * 90)
    print()
    print(f"  Substrate:   D = {D},  Hilbert dim = {1 << D}")
    print(f"  Toeplitz:    H^2_alpha(B^2) truncated at degree N = {N},  basis dim = {(N+1)*(N+2)//2}")
    print()
    print("  Site -> holomorphic Toeplitz mapping:")
    print("    site 0 -> T_{z_0}")
    print("    site 1 -> T_{z_1}")
    print("    site 2 -> T_{z_0^2}")
    print("    site 3 -> T_{z_1^2}")
    print()
    print("  Bridge of chi_S: product of individual-site images (all holomorphic,")
    print("  hence the bridge image is in the abelian holomorphic Toeplitz subalgebra).")
    print()

    basis = basis_indices(N)

    # === Sanity 1: bridge images are mutually commuting
    print("  Sanity 1: All bridge images mutually commute (abelian subalgebra)")
    print()
    print(f"  {'site pair':>14}  max_alpha ||[bridge(site_a), bridge(site_b)]||")
    for a, b in [(0, 1), (0, 2), (0, 3), (1, 2), (1, 3), (2, 3)]:
        max_norm = 0.0
        for alpha in alphas:
            Tz = [Tz_matrix(0, basis, alpha), Tz_matrix(1, basis, alpha)]
            b_a = bridge_holo(D, {a}, basis, alpha, Tz)
            b_b = bridge_holo(D, {b}, basis, alpha, Tz)
            comm = b_a @ b_b - b_b @ b_a
            max_norm = max(max_norm, op_norm(comm))
        print(f"  ({a},{b})           {max_norm:.2e}")
    print()
    print("  All zero to machine precision: the abelian structure is preserved by construction.")
    print()

    # === Sanity 2: bridge image op-norms
    print("  Sanity 2: bridge image norms (fidelity scaling with alpha)")
    print(f"  {'weight':>8}  " + "  ".join(f"a={alpha:5.1f}" for alpha in alphas))
    for k in [0, 1, 2, 3, 4]:
        if k == 0:
            S_repr = set()
        else:
            S_repr = set(range(k))
        row = f"  {k:>8}  "
        for alpha in alphas:
            Tz = [Tz_matrix(0, basis, alpha), Tz_matrix(1, basis, alpha)]
            b = bridge_holo(D, S_repr, basis, alpha, Tz)
            row += f"  {op_norm(b):>10.4f}"
        print(row)
    print()

    # === Gap probe: chi_a^Cliff bridges to T_{bar z_a} ?
    # Use T_{bar z_a} as the Toeplitz analog of the Clifford generator
    # (it's the natural anti-holomorphic counterpart).
    print("  Gap probe: substrate vs Toeplitz commutator")
    print()
    print("    [chi_a^{Cliff}, chi_S]_sub   vs   [T_{bar z_a}, bridge(chi_S)]_Toep")
    print()
    print("    Substrate commutator norm = 2 if a in S, 0 otherwise.")
    print()
    print("    Toeplitz side (alpha-dependent):")
    print(f"    {'case':>20}  ||sub||  " + "  ".join(f"a={alpha:5.1f}" for alpha in alphas))
    cases = [
        (0, {0}),       # a in S
        (0, {1}),       # a not in S
        (0, {0, 1}),    # a in S, k=2
        (0, {1, 2}),    # a not in S, k=2
        (1, {0, 1}),    # a in S, k=2, different a
        (1, {0, 2}),    # a not in S
    ]
    for a_idx, S in cases:
        chi_a = chi_clifford(D, a_idx)
        chi_S = chi_walsh(D, S)
        sub_norm = op_norm(chi_a @ chi_S - chi_S @ chi_a)
        # On Toeplitz: a_idx in {0, 1, 2, 3}. Use site-to-holo mapping for a_idx
        # to get the bridge counterpart. T_{bar z_a} for the "Clifford" generator
        # is the conjugate of the holo image.
        row = f"    a={a_idx}, S={str(set(S)):<10}  {sub_norm:>6.4f}  "
        for alpha in alphas:
            Tz = [Tz_matrix(0, basis, alpha), Tz_matrix(1, basis, alpha)]
            # Use T_{bar z_{a_idx mod 2}} as the Clifford analog (Pauli sector)
            # Site 0, 2 -> T_{bar z_0}; site 1, 3 -> T_{bar z_1}
            tzbar_a = Tzbar_matrix(a_idx % 2, basis, alpha)
            b_S = bridge_holo(D, S, basis, alpha, Tz)
            t_comm = tzbar_a @ b_S - b_S @ tzbar_a
            row += f"  {op_norm(t_comm):>10.4f}"
        print(row)
    print()

    # === Closure-rate diagnostic
    # Match substrate norm 2 (a in S) to the Toeplitz commutator
    # at the optimal alpha
    print("  Closure-rate diagnostic")
    print("  -----------------------")
    print("  For substrate commutator = 2 (a in S), what alpha gives Toeplitz commutator ~ 2?")
    print()
    print(f"  {'case':>20}  {'alpha for ||sub|| = ||Toep||':>34}  {'observed at alpha = 2':>20}")
    for a_idx, S in cases:
        chi_a = chi_clifford(D, a_idx)
        chi_S = chi_walsh(D, S)
        sub_norm = op_norm(chi_a @ chi_S - chi_S @ chi_a)
        if sub_norm < 1e-9:
            continue
        # Sweep alpha finely to find where the Toeplitz commutator matches substrate
        target_alpha = None
        for alpha_try in np.linspace(-0.5, 50.0, 200):
            Tz = [Tz_matrix(0, basis, alpha_try), Tz_matrix(1, basis, alpha_try)]
            tzbar_a = Tzbar_matrix(a_idx % 2, basis, alpha_try)
            b_S = bridge_holo(D, S, basis, alpha_try, Tz)
            t_norm = op_norm(tzbar_a @ b_S - b_S @ tzbar_a)
            if abs(t_norm - sub_norm) < 0.01:
                target_alpha = alpha_try
                break
        # Reported value at alpha = 2 for reference
        alpha_ref = 2.0
        Tz = [Tz_matrix(0, basis, alpha_ref), Tz_matrix(1, basis, alpha_ref)]
        tzbar_a = Tzbar_matrix(a_idx % 2, basis, alpha_ref)
        b_S = bridge_holo(D, S, basis, alpha_ref, Tz)
        ref_norm = op_norm(tzbar_a @ b_S - b_S @ tzbar_a)
        target_str = f"{target_alpha:.2f}" if target_alpha is not None else "(none in [-0.5, 50])"
        print(f"  a={a_idx}, S={str(set(S)):<10}  {target_str:>34}  {ref_norm:>20.4f}")
    print()

    print("=" * 90)
    print("  Verdict")
    print("=" * 90)
    print()
    print("  ABELIAN STRUCTURE: PRESERVED.  All bridge images commute to machine")
    print("  precision (Sanity 1).  The spurious cross-commutators of Computation 34")
    print("  are eliminated by mapping into the holomorphic Toeplitz subalgebra.")
    print()
    print("  FIDELITY: PARTIAL.  Bridge image norms decay with alpha, as in")
    print("  Computation 33.  Weight-k mode norms scale by additional factor compared")
    print("  to the substrate norm.  Mosco rescaling (the alpha-dependent constant")
    print("  factor) is required.")
    print()
    print("  GAP MATCHING: PARTIAL.  The Toeplitz commutator [T_{bar z_a}, b(chi_S)]")
    print("  takes alpha-dependent values that can be tuned to match the substrate")
    print("  commutator 2 for SOME (a, S) pairs by choosing alpha appropriately.")
    print("  Cases where a is not in S (substrate commutator = 0): the Toeplitz")
    print("  commutator is NOT generally zero -- residual spurious behavior.")
    print()
    print("  STRUCTURAL CONSEQUENCE.")
    print("  Mapping into the holomorphic Toeplitz subalgebra is the right move:")
    print("  it preserves the abelian Walsh structure.  But the IMAGE of chi_a^{Cliff}")
    print("  (the Clifford generator) must also be carefully chosen -- using")
    print("  T_{bar z_a} as the simplest analog still produces residual non-zero")
    print("  commutators in cases where the substrate side is zero.")
    print()
    print("  NEXT CONCRETE STEP for closing L_comm.")
    print("  The Clifford generator on the Toeplitz side needs to be constructed")
    print("  more carefully: not just T_{bar z_a} but the alpha-twisted Berezin")
    print("  derivative D_alpha = sum_a (T_{bar z_a} + alpha-correction).  The")
    print("  correction terms are designed to make [D_alpha, b(chi_S)] vanish")
    print("  when a not in S, while preserving the alpha-dependent gap rate for")
    print("  a in S.  This is the operator-algebraic content that closes (A.ii).")


if __name__ == "__main__":
    main()
