#!/usr/bin/env python3 """ verify_P002.py -- Paper 02: Geometric Spectral Theory. This verifier checks the main numerical chain in toe/02_Paper_GeometricSpectralTheory.tex: cubic-density moments, self-lensing energy, dynamic beta comparison, QED beta ratio, the critical point, and the mass-hierarchy conjecture. Most basic moment and energy arithmetic is reproducible. The audit flags: * the QED beta ratio 70,205.483 is obtained only if the experimental alpha is used, while the displayed formula uses the geometric m0 value; * the comparator-density E_self values in the optimal-density table are not the values produced by the stated normalization; * the wave equation cannot have rho_cubic as a Neumann equilibrium as stated; * the critical-point proof has the wrong s -> 0 limit; * the zeta-prime numeric value is right, but the displayed log formula is not; * M_Pl / exp(C) with C ~= 1362 is not 2.7e16 GeV. """ import math import sys from pathlib import Path sys.path.insert(0, str(Path(__file__).resolve().parent)) PASS = FAIL = 0 _N = 0 def check(n, desc, cond): global PASS, FAIL ok = bool(cond) PASS += ok FAIL += not ok print(f" [{'PASS' if ok else 'FAIL'}] {n:>2}. {desc}") return ok class Verifier: """Output adapter: identical check semantics, modern [PASS]/[FAIL] format.""" def __init__(self, name): print(name) def check(self, label, computed, claimed, *, rel=1e-3, abs_tol=None, detail=""): if abs_tol is not None: ok = abs(computed - claimed) <= abs_tol err_detail = f"abs err={abs(computed - claimed):.6g}, tol={abs_tol:.6g}" else: if claimed == 0: ok = abs(computed) <= (rel or 1e-12) err_detail = f"abs value={abs(computed):.6g}, tol={rel:.6g}" else: err = (computed - claimed) / abs(claimed) ok = abs(err) <= (rel or 0) err_detail = f"rel err={100 * err:+.6g}%, tol={100 * (rel or 0):.6g}%" return self.record(label, ok, computed, claimed, detail, err_detail) def record(self, label, ok, computed="", claimed="", detail="", err_detail=""): global _N _N += 1 desc = label + (f" -- {detail}" if detail else "") check(_N, desc, ok) if computed != "" or claimed != "": print(f" computed: {computed}") print(f" claimed : {claimed}") if err_detail: print(f" {err_detail}") return ok def summary(self): print(f"\n{'='*60}\nRESULT: {PASS} PASS / {FAIL} FAIL") return 1 if FAIL else 0 v = Verifier("P002 -- Geometric Spectral Theory") PI = math.pi GAMMA = 0.5772156649015329 PHI_GOLDEN = (1 + math.sqrt(5)) / 2 A3 = 16 * PI**3 A2 = 3 * PI**2 A1 = 2 * PI M0 = 4 * PI**3 + PI**2 + PI M1 = A3 / 5 + A2 / 4 + A1 / 3 ALPHA_GEOM = 1 / M0 ALPHA_EXP_INV = 137.035999084 ALPHA_EXP = 1 / ALPHA_EXP_INV def rho_prime(x: float) -> float: return 48 * PI**3 * x**2 + 6 * PI**2 * x + 2 * PI def rho_second(x: float) -> float: return 96 * PI**3 * x + 6 * PI**2 def phi_laplace(s: float) -> float: e = math.exp(-s) return ( A3 * (6 - e * (6 + 6 * s + 3 * s**2 + s**3)) / s**4 + A2 * (2 - e * (2 + 2 * s + s**2)) / s**3 + A1 * (1 - e * (1 + s)) / s**2 ) def solve_critical() -> float: target = M0**2 / M1 lo, hi = 0.52, 0.53 for _ in range(100): mid = (lo + hi) / 2 if phi_laplace(mid) / mid > target: lo = mid else: hi = mid return (lo + hi) / 2 energy = 0.5 * ( (48 * PI**3) ** 2 / 5 + (6 * PI**2) ** 2 / 3 + (2 * PI) ** 2 + 2 * (48 * PI**3) * (6 * PI**2) / 4 + 2 * (48 * PI**3) * (2 * PI) / 3 + 2 * (6 * PI**2) * (2 * PI) / 2 ) e_self = energy / M0**2 delta_e = e_self - 4 * PI beta_geom = M1 / M0 beta_dyn = 2 * PI / (3 * PHI_GOLDEN**2) - GAMMA / 100 delta_beta_over_beta = abs(beta_geom - beta_dyn) / beta_dyn kappa_beta = 2 * delta_beta_over_beta kappa_alpha = ALPHA_GEOM ** (5 / 4) v.check("m0 = 4*pi^3 + pi^2 + pi", M0, 137.036303776, rel=1e-12) v.check( "geometric alpha inverse vs CODATA relative percent", 100 * (M0 - ALPHA_EXP_INV) / ALPHA_EXP_INV, 0.000222, rel=2e-3, ) v.check("m1 moment", M1, 108.716683780, rel=1e-11) v.check("beta_geom = m1/m0", beta_geom, 0.793342207849, rel=2e-10) v.check("Dirichlet energy E[rho]", energy, 247444.81, rel=5e-8) v.check("self-lensing energy E_self", e_self, 13.17671, rel=5e-7) v.check("floor 4*pi", 4 * PI, 12.56637, rel=5e-7) v.check("gap Delta E", delta_e, 0.61034, rel=1e-4) v.check("fractional gap Delta E / E_self", delta_e / e_self, 0.04633, rel=3e-4) v.check("beta_dyn = 2*pi/(3 phi^2) - gamma/100", beta_dyn, 0.794215778, rel=3e-7) v.check("static-dynamic relative difference", delta_beta_over_beta, 0.0011, rel=5e-4) v.check("kappa from beta difference", kappa_beta, 0.0022, rel=5e-4) v.check("kappa = alpha^(5/4) comparison ratio", kappa_beta / kappa_alpha, 1.031, rel=1e-3) v.check( "Delta E / E_self agreement with 22*kappa is within 1.3 percent", 100 * abs(delta_e / e_self - 22 * kappa_beta) / (22 * kappa_beta), 1.3, rel=1e-2, detail="Expected fail: the stated quantities differ by about 4.27%, as the proof later acknowledges.", ) v.check( "Delta E / E_self agreement with 22*kappa is within 5 percent", abs(delta_e / e_self - 22 * kappa_beta) / (22 * kappa_beta), 0.0, abs_tol=0.05, ) # The optimal-density table uses normalized densities with integral m0. linear_e_self = 2.0 quadratic_e_self = 6.0 v.check("optimal-density table: uniform E_self", 0.0, 0.0, abs_tol=1e-12) v.check( "optimal-density table: linear E_self", linear_e_self, 0.171, rel=1e-3, detail="Expected fail: for rho=2*m0*x, E/m0^2 = 2, not 0.171.", ) v.check( "optimal-density table: quadratic E_self", quadratic_e_self, 0.383, rel=1e-3, detail="Expected fail: for rho=3*m0*x^2, E/m0^2 = 6, not 0.383.", ) v.record( "optimal-density proof exhausts polynomial densities", False, computed="four monomial test cases", claimed="minimum among polynomial densities", detail="Expected fail: checking uniform/linear/quadratic/cubic monomials is not an optimization proof over polynomial coefficients.", ) beta_qed_geom = 2 * ALPHA_GEOM**2 / (3 * PI) ratio_geom_alpha = beta_geom / beta_qed_geom beta_qed_exp = 2 * ALPHA_EXP**2 / (3 * PI) ratio_exp_alpha = beta_geom / beta_qed_exp log_planck_e = math.log(1.220890e19 / (5.10999e-4)) c_emp_exp = ratio_exp_alpha / log_planck_e c_formula = 10 * M0**3 / (M0**2 + M1) v.check("QED beta using displayed geometric m0 formula", beta_qed_geom, 1.130024e-5, rel=5e-7) v.check( "beta_geom/beta_QED using displayed geometric m0 formula", ratio_geom_alpha, 70205.483, rel=1e-6, detail="Expected fail: 70,205.483 uses experimental alpha, not beta_QED=2/(3*pi*m0^2).", ) v.check("beta_geom/beta_QED using experimental alpha", ratio_exp_alpha, 70205.483, rel=1e-6) v.check("log(M_Pl/m_e) from stated masses", log_planck_e, 51.527840, rel=1e-8) v.check("empirical C from experimental-alpha ratio", c_emp_exp, 1362.477, rel=5e-7) v.check( "exact C formula value", c_formula, 1362.482800, abs_tol=0.001, detail="Expected fail: exact m0,m1 give about 1362.47528; the quoted value is a small arithmetic mismatch.", ) s_star = solve_critical() v.check("critical point s*", s_star, 0.525159133415, rel=1e-11) v.check("critical point proximity to pi/6 percent", 100 * (s_star - PI / 6) / s_star, 0.297, rel=2e-3) v.record( "critical-point proof limit as s -> 0+", False, computed="Phi(s)/s -> +infinity", claimed="Phi(s)/s -> Phi'(0) = -m1", detail="Expected fail: the root value is correct, but this proof step is not.", ) rho_prime_0 = rho_prime(0.0) rho_prime_1 = rho_prime(1.0) v.check( "rho_cubic satisfies Neumann boundary at x=0", rho_prime_0, 0.0, abs_tol=1e-12, detail="Expected fail: rho'(0)=2*pi.", ) v.check( "rho_cubic satisfies Neumann boundary at x=1", rho_prime_1, 0.0, abs_tol=1e-12, detail="Expected fail: rho'(1) is large and nonzero.", ) v.check( "rho_cubic is static equilibrium of stated wave equation at x=1/2", rho_second(0.5), 0.0, abs_tol=1e-12, detail="Expected fail: at rho=rho_cubic, the RHS leaves rho_cubic'' instead of zero.", ) v.record( "linearization reduces to eta_tt = (1-kappa) eta_xx", False, computed="eta_xx - kappa*eta*rho_cubic'' plus a nonzero zeroth-order residual", claimed="(1-kappa) eta_xx", detail="Expected fail: differentiating -kappa(rho-rho_cubic)rho'' gives a multiplicative rho_cubic'' term, not -kappa eta''.", ) rhs_log = (3 * PI / 20) * M1 / (1 - M1 * ALPHA_GEOM**2) v.check("geometric-series mass-log RHS", rhs_log, 51.529851, rel=5e-9) v.check("geometric-series relative error percent", 100 * abs(rhs_log - log_planck_e) / log_planck_e, 0.0039, rel=2e-3) v.check("series convergence parameter m1*alpha^2", M1 * ALPHA_GEOM**2, 0.005789285, rel=2e-7) v.check("leading mass-log term", (3 * PI / 20) * M1, 51.231530, rel=5e-8) v.check("alpha^2 correction term", (3 * PI / 20) * M1**2 * ALPHA_GEOM**2, 0.2966, rel=1e-4) zeta_prime_integral = A3 / 16 + A2 / 9 + A1 / 4 zeta_prime_printed_log_formula = A3 * math.log(1 / 4) + A2 * math.log(1 / 3) + A1 * math.log(1 / 2) v.check("zeta'(0) from integral definition", zeta_prime_integral, 35.867, rel=2e-6) v.check( "displayed zeta'(0) log formula", zeta_prime_printed_log_formula, 35.867, rel=1e-3, detail="Expected fail: the printed sum A_k log(1/k) is negative and not the derivative of the integral zeta function.", ) v.check("zeta'(0)/zeta(0)", zeta_prime_integral / M0, 0.262, rel=2e-3) v.check("one-loop action -1/2 zeta'(0)", -0.5 * zeta_prime_integral, -17.93, rel=2e-4) log10_e_geom = math.log10(1.220890e19) - 1362.477 / math.log(10) v.check( "log10(M_Pl / exp(C)) with C=1362.477", log10_e_geom, math.log10(2.7e16), abs_tol=0.1, detail="Expected fail: exp(1362) is enormous; the scale is around 10^-573 GeV, not 2.7e16 GeV.", ) sys.exit(v.summary())