#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
SPACE CUBION — annexe reproductible (RS-1, 2026-08-11)

Recalcule les 7 blocs de calcul du papier resources/papers/space-cubion.md
(version v1 corrigée : table carré-cube ×3, cellule cubion ~1,2×10⁻⁷ s, n=3
divisions pour r/2) et imprime un tableau attendu-vs-calculé avec ✓/✗.

Régime : tout ce qui est vérifié ici est RÉEL (physique établie, arithmétique).
Aucune dépendance hors stdlib. Usage : python3 space-cubion-verify.py
Sortie 0 si tout est ✓, 1 sinon.
"""

import math
import sys

# ── Constantes ────────────────────────────────────────────────────────────────
G      = 6.674e-11        # m³·kg⁻¹·s⁻² (constante gravitationnelle)
c_ms   = 2.998e8          # m/s
g0     = 9.81             # m/s² (1 g)
AL     = 9.4607e15        # m par année-lumière
AN     = 3.15576e7        # s par an (année julienne)
L_SUN  = 3.8e26           # W (luminosité solaire totale, valeur du papier)

CHECKS = []


def check(bloc, label, attendu, calcule, tol_rel=0.02, unite=""):
    """Compare attendu vs calculé à tol_rel près, mémorise ✓/✗."""
    ok = attendu == calcule if attendu == 0 else abs(calcule - attendu) <= tol_rel * abs(attendu)
    CHECKS.append((bloc, label, attendu, calcule, unite, ok))
    return ok


def fmt(x):
    if x == 0:
        return "0"
    ax = abs(x)
    if ax >= 1e5 or ax < 1e-3:
        return f"{x:.3g}"
    return f"{x:.4g}"


# ── Bloc 1 · Temps propres à 1 g avec retournement (§7.1, RÉEL) ──────────────
# τ = (2c/a)·acosh(1 + a·d/(2c²)) — unités naturelles c=1 al/an, a en al/an².
def bloc1():
    a = g0 * AN**2 / AL           # 1 g ≈ 1,0326 al/an²
    tau = lambda d_al: (2.0 / a) * math.acosh(1.0 + a * d_al / 2.0)
    # (cible, distance al — table canonique Baez pour Vega/centre/Andromède, attendu ans)
    cibles = [
        ("Proxima Centauri (4,24 al)",      4.24,      3.54),
        ("Alpha du Centaure (4,37 al)",     4.37,      3.58),
        ("Tau Ceti (11,9 al)",              11.9,      5.14),
        ("Vega (27 al)",                    27.0,      6.58),
        ("TRAPPIST-1 (40,7 al)",            40.7,      7.33),
        ("Centre galactique (30 000 al)",   30000.0,   20.0),
        ("Andromède (2×10⁶ al)",            2.0e6,     28.2),
    ]
    for nom, d, att in cibles:
        check("1 · τ à 1 g", f"τ {nom}", att, tau(d), tol_rel=0.01, unite="ans")


# ── Bloc 2 · Table carré-cube CORRIGÉE (§5.1, RÉEL) ──────────────────────────
# a = 3·F_s/(ρ·r), F_s = 10⁵ N/m², ρ = 2 000 kg/m³  →  1,5 / 0,15 / 0,015 / 0,0015.
def bloc2():
    Fs, rho = 1e5, 2000.0
    a = lambda r: 3.0 * Fs / (rho * r)
    for r, att in [(100.0, 1.5), (1e3, 0.15), (1e4, 0.015), (1e5, 0.0015)]:
        check("2 · carré-cube", f"a(r={fmt(r)} m)", att, a(r), tol_rel=1e-9, unite="m/s²")


# ── Bloc 3 · Tsiolkovsky m₀/m_f = exp(Δv/v_e) (§7.3 + §15, RÉEL) ─────────────
def bloc3():
    dv1 = 3.0e7        # m/s — 0,1 c, valeur explicite du papier (§7.3, c arrondi à 3×10⁸)
    check("3 · Tsiolkovsky", "exposant chimique (v_e=4,5 km/s)", 6667.0, dv1 / 4.5e3, 1e-3)
    check("3 · Tsiolkovsky", "exposant ionique (v_e=50 km/s)",   600.0,  dv1 / 5.0e4, 1e-9)
    check("3 · Tsiolkovsky", "m₀/m_f fusion (v_e=10⁴ km/s), 0,1 c", 20.09, math.exp(dv1 / 1e7), 1e-3)
    check("3 · Tsiolkovsky", "m₀/m_f fusion, 0,2 c", 403.4, math.exp(2.0 * dv1 / 1e7), 1e-3)
    check("3 · Tsiolkovsky", "m₀/m_f rendez-vous 0,01 c (§15)", 1.35, math.exp(0.1 * dv1 / 1e7), 5e-3)


# ── Bloc 4 · Budgets ½mv² vs L_☉ (§7.2, RÉEL — cellule cubion CORRIGÉE) ──────
# E = ½·m·(0,1c)² ; temps équivalent = E / L_☉.  Cubion 100 t → ~1,2×10⁻⁷ s
# (~0,1 µs — l'ancien « ~0,1 s » du papier v0 était faux d'un facteur 10⁶).
def bloc4():
    v = 3e7                     # m/s — 0,1 c (arrondi du papier, c≈3×10⁸)
    t_eq = lambda m: 0.5 * m * v**2 / L_SUN
    lignes = [
        ("cubion 100 t (10⁵ kg)",        1e5,     1.18e-7, "s"),
        ("vaisseau 1 km (10¹² kg)",      1e12,    1.18,    "s"),
        ("vaisseau 10 km (10¹⁵ kg)",     1e15,    19.7,    "min"),
        ("vaisseau 100 km (10¹⁸ kg)",    1e18,    13.7,    "jours"),
        ("masse lunaire (7,35×10²² kg)", 7.35e22, 2758.0,  "ans"),
    ]
    div = {"s": 1.0, "min": 60.0, "jours": 86400.0, "ans": AN}
    for nom, m, att, u in lignes:
        check("4 · ½mv²/L_☉", f"E/L_☉ {nom}", att, t_eq(m) / div[u], tol_rel=0.01, unite=u)


# ── Bloc 5 · Vitesse de libération (§2, RÉEL) ────────────────────────────────
# g_s = (4/3)·π·G·ρ·r ; v_esc = √(2·g_s·r).
def bloc5():
    g_s   = lambda rho, r: (4.0 / 3.0) * math.pi * G * rho * r
    v_esc = lambda rho, r: math.sqrt(2.0 * g_s(rho, r) * r)
    check("5 · v_esc", "g_s Lune (ρ=3340, r=1737 km)", 1.62, g_s(3340.0, 1.737e6), 0.01, "m/s²")
    check("5 · v_esc", "v_esc Lune", 2.37, v_esc(3340.0, 1.737e6) / 1e3, 0.01, "km/s")
    check("5 · v_esc", "g_s agrégat 100 km (ρ=2500)", 0.070, g_s(2500.0, 1e5), 0.01, "m/s²")
    check("5 · v_esc", "v_esc agrégat 100 km", 0.118, v_esc(2500.0, 1e5) / 1e3, 0.01, "km/s")


# ── Bloc 6 · Ceinture d'astéroïdes et cubions (§8, RÉEL) ─────────────────────
# ⚠ ERRATUM v1 : le papier v0 écrit « ~0,03 % de la ceinture » pour 8×10⁹
# cubions ; le vrai ratio est 8×10¹⁴/3×10²¹ ≈ 2,7×10⁻⁷ soit ~0,00003 %
# (facteur ~10³) — l'argument s'en trouve RENFORCÉ, mais le chiffre est à
# corriger dans le papier au même titre que les 3 erreurs déjà relevées.
def bloc6():
    M_ceinture, m_cubion, humains = 3e21, 1e5, 8e9
    check("6 · ceinture", "cubions fabricables (3×10²¹/10⁵)", 3e16, M_ceinture / m_cubion, 1e-9)
    check("6 · ceinture", "fraction pour 8×10⁹ humains", 2.67e-7, humains * m_cubion / M_ceinture, 0.01)
    check("6 · ceinture", "doublements 1 → 10¹⁶", 53.2, math.log2(1e16), 0.01)


# ── Bloc 7 · Mitose 2^(1/3) et échelle 2ⁿ (§5.2 + §8.1, RÉEL — CORRIGÉ n=3) ──
# Couper en 2 → rayon × 2^(-1/3) ≈ 0,794. Il faut n = 3 divisions par 2
# successives pour diviser le rayon par 2 (0,794³ = 0,5) — pas « huit » ;
# 8 = le nombre de morceaux d'UN octa-découpage (2³). Et l'arête décidée :
# 4 m = 2⁸ bions de 16 mm (16 mm × 2⁸ = 4,096 m).
def bloc7():
    f = 2.0 ** (1.0 / 3.0)
    check("7 · mitose", "facteur de rayon 2^(1/3)", 1.26, f, 0.01)
    check("7 · mitose", "rayon d'une moitié (r×2^(-1/3))", 0.794, 1.0 / f, 0.01)
    check("7 · mitose", "n divisions pour r/2 (ln2/ln2^(1/3))", 3.0, math.log(2.0) / math.log(f), 1e-9)
    check("7 · mitose", "morceaux d'un octa-découpage (2³)", 8.0, 2.0**3, 1e-9)
    check("7 · mitose", "arête cubion 16 mm × 2⁸", 4.096, 0.016 * 2**8, 1e-9, "m")


# ── Exécution + tableau ──────────────────────────────────────────────────────
def main():
    for bloc in (bloc1, bloc2, bloc3, bloc4, bloc5, bloc6, bloc7):
        bloc()

    print("SPACE CUBION — annexe reproductible : attendu vs calculé (papier v1 corrigé)")
    print("=" * 96)
    print(f"{'bloc':<16} {'grandeur':<44} {'attendu':>12} {'calculé':>12}  {'ok':<2}")
    print("-" * 96)
    dernier_bloc = None
    for bloc, label, att, cal, unite, ok in CHECKS:
        if bloc != dernier_bloc and dernier_bloc is not None:
            print("-" * 96)
        dernier_bloc = bloc
        u = f" {unite}" if unite else ""
        print(f"{bloc:<16} {label:<44} {fmt(att):>12} {fmt(cal):>12}  {'✓' if ok else '✗'}{u}")
    print("=" * 96)

    n_ok = sum(1 for chk in CHECKS if chk[-1])
    n = len(CHECKS)
    tout_bon = n_ok == n
    print(f"{n_ok}/{n} vérifications ✓ — {'papier v1 corrigé COHÉRENT' if tout_bon else 'DIVERGENCE : corriger le papier ou le script'}")
    return 0 if tout_bon else 1


if __name__ == "__main__":
    sys.exit(main())