# Copyright Léon Fernando Vlegels -- CC BY 4.0
"""
verify_P337.py
Verifier for Addendum 337: correction of the Monster 3B McKay-Thompson mislabel.

A163/A177 DEFINED the Monster 3B McKay-Thompson series as T_3B(tau) := j(tau)^(1/3)
(real cube root) and reported T_3B(i) = 1728^(1/3) = 12 as a Monster 3B "moonshine
hit". THIS LABELING IS WRONG. j^(1/3) = E4/eta^8 is the E8/G2 cube-root function
(A160 Monodromy-Weyl-A2, A161 E8MonsterBridge), NOT the Monster's 3B series.

The CANONICAL Monster 3B McKay-Thompson series (Conway-Norton / monstrous moonshine,
proved Borcherds 1992) is the genus-zero hauptmodul

      T_3B^std(tau) = (eta(tau) / eta(3 tau))^12 + 12,

whose value at tau = i is NOT 12, but 535.592... . So the "12 = Monster 3B moonshine"
claim is a CONFLATION of two distinct functions both informally called "the cube root
/ 3B". This verifier recomputes EVERY number from scratch (mpmath only -- NO corpus
imports) at dps=40 and proves:

  (a) j(i) = 1728 and j(i)^(1/3) = 12 exactly;
  (b) T_3B^std(i) = (eta(i)/eta(3i))^12 + 12 = 535.592... (NOT 12);
  (c) the CORE correction: j^(1/3)(i) != T_3B^std(i)  (12 != 535.59);
  (d) j^(1/3) = E4/eta^8 has q-expansion q^{-1/3}(1 + 248 q + ...) -- the E8 248 --
      whereas T_3B^std has q-expansion q^{-1} + 0 + 54 q + ... (integer powers):
      two structurally different functions;
  (e) T_1A(i) = j(i) - 744 = 984 -- the surviving correctly-computed value.

What is RETRACTED: the identification T_3B = j^(1/3) AS A MONSTER-3B MOONSHINE
statement, and the "12 = unique 3B McKay-Thompson hit" framing.
What SURVIVES: j(i)^(1/3) = 12 = J_short as E8/G2 ROOT GEOMETRY (A160/A161), the
correctly-computed T_1A/T_2A/T_2B/T_3A values, the perfect-cube lemma, and the
honest A177 flag that T_2B(tau_-163) = -24 is a Ramanujan eta artefact.
"""

import sys
import mpmath as mp

mp.mp.dps = 40

# ---------------------------------------------------------------------------
# Test harness
# ---------------------------------------------------------------------------
_PASS = 0
_FAIL = 0


def check(n, desc, cond):
    global _PASS, _FAIL
    tag = "[PASS]" if cond else "[FAIL]"
    if cond:
        _PASS += 1
    else:
        _FAIL += 1
    print(f"{tag} ({n}) {desc}")


def section(title):
    print()
    print("-- " + title + " " + "-" * max(0, 66 - len(title)))


# ---------------------------------------------------------------------------
# Modular building blocks -- recomputed from scratch
# ---------------------------------------------------------------------------
def eta(tau):
    """Dedekind eta via mpmath's q-series product."""
    return mp.qp(mp.e ** (2j * mp.pi * tau)) * mp.e ** (mp.pi * 1j * tau / 12)


def E4(tau):
    """Eisenstein E4 = 1 + 240 sum sigma_3(n) q^n."""
    q = mp.e ** (2j * mp.pi * tau)
    s = mp.mpf(0)
    for nn in range(1, 400):
        sig3 = sum(d ** 3 for d in range(1, nn + 1) if nn % d == 0)
        term = sig3 * q ** nn
        s += term
        if abs(term) < mp.mpf(10) ** (-50):
            break
    return 1 + 240 * s


def j_invariant(tau):
    """j = E4^3 / eta^24."""
    return E4(tau) ** 3 / eta(tau) ** 24


# CM points
TAU_I = mp.mpc(0, 1)                       # D = -4
TAU_RHO = (1 + 1j * mp.sqrt(3)) / 2        # D = -3, j = 0

print("=" * 72)
print("P337 VERIFICATION: correction of the Monster 3B McKay-Thompson mislabel")
print(f"mpmath precision: {mp.mp.dps} decimal places")
print("=" * 72)

# ---------------------------------------------------------------------------
section("S0  j(i) and its cube root (the E8/G2 cube-root function)")
# ---------------------------------------------------------------------------
ji = j_invariant(TAU_I)
print(f"  j(i)              = {mp.nstr(ji, 25)}")
check(1, "j(i) = 1728", abs(ji - 1728) < mp.mpf(10) ** (-20))

cube_root_ji = mp.re(ji) ** (mp.mpf(1) / 3)   # real cube root
print(f"  j(i)^(1/3)        = {mp.nstr(cube_root_ji, 25)}")
check(2, "j(i)^(1/3) = 12 exactly", abs(cube_root_ji - 12) < mp.mpf(10) ** (-20))

# ---------------------------------------------------------------------------
section("S1  Canonical Monster 3B series  T_3B^std = (eta/eta3)^12 + 12")
# ---------------------------------------------------------------------------
eta_i = eta(TAU_I)
eta_3i = eta(3 * TAU_I)
ratio12 = (eta_i / eta_3i) ** 12
T3B_std_i = ratio12 + 12
print(f"  (eta(i)/eta(3i))^12 = {mp.nstr(mp.re(ratio12), 25)}")
print(f"  T_3B^std(i)         = {mp.nstr(mp.re(T3B_std_i), 25)}")
# Orchestrator value: 535.592230826...
check(3, "T_3B^std(i) = 535.592... (canonical 3B, NOT 12)",
      abs(mp.re(T3B_std_i) - mp.mpf("535.592230826")) < mp.mpf("1e-6"))
check(4, "T_3B^std(i) is real to high precision",
      abs(mp.im(T3B_std_i)) < mp.mpf(10) ** (-15))
check(5, "T_3B^std(i) != 12 (it is far from 12)",
      abs(mp.re(T3B_std_i) - 12) > mp.mpf(500))

# ---------------------------------------------------------------------------
section("S2  THE CORE CORRECTION:  j^(1/3)(i)  !=  T_3B^std(i)")
# ---------------------------------------------------------------------------
gap = mp.re(T3B_std_i) - cube_root_ji
print(f"  j^(1/3)(i)   = {mp.nstr(cube_root_ji, 18)}   (E8/G2 cube-root fn)")
print(f"  T_3B^std(i)  = {mp.nstr(mp.re(T3B_std_i), 18)}   (canonical Monster 3B)")
print(f"  difference   = {mp.nstr(gap, 18)}")
check(6, "j^(1/3)(i) != T_3B^std(i)  -- the 12 is NOT the Monster 3B value",
      abs(gap) > mp.mpf(500))
check(7, "the two functions disagree by ~523.6 at tau=i",
      abs(gap - mp.mpf("523.592")) < mp.mpf("1e-2"))

# ---------------------------------------------------------------------------
section("S3  q-expansions: two DIFFERENT functions (E8 248 vs integer powers)")
# ---------------------------------------------------------------------------
# j^(1/3) = E4 / eta^8 = q^{-1/3} (1 + 248 q + ...).  Verify the 248 (E8 adjoint dim)
# by computing E4/eta^8 * q^{1/3} as a series in q and reading the q^1 coefficient.
def series_coeff_cuberoot(maxn=6):
    """Coefficients of q^{1/3} * E4(tau)/eta(tau)^8 as a power series in q."""
    # Use a small symbolic-free numerical q-derivative-free approach:
    # f(q) = E4/eta^8 = q^{-1/3} * g(q), g(q) = 1 + a1 q + a2 q^2 + ...
    # Evaluate g at several small real q and least-squares the low coefficients.
    import numpy as np

    qs = [mp.mpf("0.001") * (k + 1) for k in range(maxn + 2)]
    A = []
    b = []
    for q in qs:
        tau = mp.log(q) / (2j * mp.pi)  # q = e^{2 pi i tau}
        f = E4(tau) / eta(tau) ** 8
        g = f * q ** (mp.mpf(1) / 3)     # strip q^{-1/3}
        A.append([float(q) ** k for k in range(maxn)])
        b.append(float(mp.re(g)))
    A = np.array(A)
    b = np.array(b)
    coeffs, *_ = np.linalg.lstsq(A, b, rcond=None)
    return coeffs


cr_coeffs = series_coeff_cuberoot(maxn=4)
print(f"  E4/eta^8 = q^(-1/3) * (a0 + a1 q + ...),  fitted a0={cr_coeffs[0]:.4f}, a1={cr_coeffs[1]:.4f}")
check(8, "j^(1/3) = E4/eta^8 has q^(-1/3) leading term, constant ~ 1",
      abs(cr_coeffs[0] - 1.0) < 1e-2)
check(9, "j^(1/3) = E4/eta^8 has q^1 coefficient = 248 (the E8 adjoint dimension)",
      abs(cr_coeffs[1] - 248.0) < 1.0)

# T_3B^std = (eta/eta3)^12 + 12 = q^{-1} + 0 + 54 q - 76 q^2 + ...  (integer powers)
# Verify leading q^{-1} and the +12 constant and the q^1 coefficient 54.
def series_coeff_t3bstd(maxn=4):
    import numpy as np

    qs = [mp.mpf("0.0008") * (k + 1) for k in range(maxn + 3)]
    A = []
    b = []
    for q in qs:
        tau = mp.log(q) / (2j * mp.pi)
        f = (eta(tau) / eta(3 * tau)) ** 12 + 12  # = q^{-1} + c0 + c1 q + ...
        h = f * q  # strip q^{-1}: h = 1 + c0 q + c1 q^2 + ...
        A.append([float(q) ** k for k in range(maxn)])
        b.append(float(mp.re(h)))
    A = np.array(A)
    b = np.array(b)
    coeffs, *_ = np.linalg.lstsq(A, b, rcond=None)
    return coeffs


t3b_coeffs = series_coeff_t3bstd(maxn=4)
print(f"  T_3B^std = q^(-1)*(b0 + b1 q + b2 q^2 + ...), fitted b0={t3b_coeffs[0]:.4f}, b1={t3b_coeffs[1]:.4f}, b2={t3b_coeffs[2]:.4f}")
check(10, "T_3B^std has q^(-1) leading term (coeff 1)",
      abs(t3b_coeffs[0] - 1.0) < 1e-2)
check(11, "T_3B^std constant term c0 = 0",
      abs(t3b_coeffs[1] - 0.0) < 1.0)
check(12, "T_3B^std q^1 coefficient = 54 (Conway-Norton 3B coefficient)",
      abs(t3b_coeffs[2] - 54.0) < 1.5)
check(13, "the two q-expansions differ structurally (248 != 54; cube-root vs integer powers)",
      abs(cr_coeffs[1] - t3b_coeffs[2]) > 100)

# ---------------------------------------------------------------------------
section("S4  T_1A(i) = j(i) - 744 = 984  (SURVIVES -- correct McKay-Thompson value)")
# ---------------------------------------------------------------------------
T1A_i = mp.re(ji) - 744
print(f"  T_1A(i) = j(i) - 744 = {mp.nstr(T1A_i, 20)}")
check(14, "T_1A(i) = j(i) - 744 = 984 (1A is the standard hauptmodul; correct)",
      abs(T1A_i - 984) < mp.mpf(10) ** (-15))

# ---------------------------------------------------------------------------
section("S5  Surviving facts: perfect-cube lemma + the 12 under E8/G2 provenance")
# ---------------------------------------------------------------------------
# Perfect-cube lemma sample: j(rho) = 0 = 0^3, j(i) = 1728 = 12^3.
j_rho = j_invariant(TAU_RHO)
print(f"  j(rho) = {mp.nstr(mp.re(j_rho), 8)}  (expect 0)")
check(15, "perfect-cube lemma sample: j(rho) = 0 = 0^3",
      abs(j_rho) < mp.mpf(10) ** (-10))
check(16, "perfect-cube lemma sample: j(i) = 1728 = 12^3 (stands, classical CM)",
      abs(mp.re(ji) - 12 ** 3) < mp.mpf(10) ** (-15))
# The 12 stands under its correct provenance: 12 = j^(1/3)(i) = E4(i)/eta(i)^8,
# the E8/G2 cube-root reading (A160/A161), NOT the Monster 3B series value.
e4_over_eta8 = mp.re(E4(TAU_I) / eta(TAU_I) ** 8)
print(f"  E4(i)/eta(i)^8 = {mp.nstr(e4_over_eta8, 20)}  (= j^(1/3)(i), the E8/G2 fn)")
check(17, "the 12 stands as E8/G2 root geometry: E4(i)/eta(i)^8 = 12 = j^(1/3)(i)",
      abs(e4_over_eta8 - 12) < mp.mpf(10) ** (-12))
check(18, "12 (E8/G2) is distinct from T_3B^std(i)=535.59 (Monster 3B): the conflation",
      abs(e4_over_eta8 - mp.re(T3B_std_i)) > mp.mpf(500))

# ---------------------------------------------------------------------------
# Summary
# ---------------------------------------------------------------------------
print()
print("=" * 72)
print(f"RESULT: {_PASS} PASS / {_FAIL} FAIL")
if _FAIL == 0:
    print("VERDICT: CORRECTION CONFIRMED. j^(1/3)(i) = 12 is the E8/G2 cube-root "
          "function\n         (E4/eta^8, A160/A161), NOT the Monster 3B McKay-Thompson "
          "series\n         (canonical 3B = (eta/eta3)^12 + 12 = 535.592 at i). The "
          "'12 = Monster\n         3B moonshine hit' claim of A163/A177 is a conflation "
          "and is RETRACTED.\n         The 12 = J_short fact and T_1A(i)=984 SURVIVE under "
          "correct provenance.")
else:
    print("VERDICT: FAIL -- a recomputed number did not reproduce; do NOT ship A337.")
print("=" * 72)

sys.exit(0 if _FAIL == 0 else 1)
