#!/usr/bin/env python3
"""N22 — arithmetic from D1287 supplied sources.
No network, external packages, fitting operations or shared runtime files.
Printed source cells and unrounded arithmetic are distinct outputs.
"""

import math
from decimal import Decimal, ROUND_HALF_UP, getcontext
from fractions import Fraction

getcontext().prec = 40
PHI = (1.0 + math.sqrt(5.0)) / 2.0


def shown(value, places=9, signed=False):
    """Decimal half-up display; never use display-rounded inputs implicitly."""
    value = Decimal(str(value))
    rounded = value.quantize(Decimal(1).scaleb(-places), rounding=ROUND_HALF_UP)
    return format(rounded, ("+" if signed else "") + "." + str(places) + "f")


def compare(label, theory, reference, lower=None, upper=None):
    """Naive central-value arithmetic, not a likelihood or theory-error model.

    The asymmetric denominator points from the comparator toward the theory:
    lower error below the central value; upper error above it. Any different
    printed convention is audited separately, rather than silently substituted.
    """
    theory, reference = Decimal(str(theory)), Decimal(str(reference))
    if reference == 0:
        raise ValueError("A relative deviation needs a nonzero reference.")
    offset = theory - reference
    percent = 100 * offset / reference
    print(label + ".theory = " + shown(theory))
    print(label + ".reference = " + shown(reference))
    print(label + ".signed_percent = " + shown(percent, signed=True)
          + "%; rounded = " + shown(percent, 2, True) + "%")
    if lower is None or upper is None:
        print(label + ".d = NOT AVAILABLE (no uncertainty supplied for this comparison variable)")
        return
    side = "lower" if offset < 0 else "upper"
    uncertainty = Decimal(str(lower if offset < 0 else upper))
    if uncertainty <= 0:
        raise ValueError("The selected comparator uncertainty must be positive.")
    distance = offset / uncertainty
    print(label + ".uncertainty_used = " + str(uncertainty) + " (" + side + ")")
    print(label + ".d = " + shown(distance, signed=True)
          + "; rounded = " + shown(distance, 2, True))

def main():
    print('=== row 22 ===')
    print('Scorecard (verbatim): Δm² 21 /Δm² 3l φ⁻¹/√5 theory 0.27639 experiment 0.03002 (7.537e−5 / 2.511e−3, NuFIT 6.1) deviation +820% Refuted falsified 9.21× — this is the d = 1 reading; see 22b')
    value = PHI**(-1)/math.sqrt(5.0)
    solar, atmospheric = 7.537e-5, 2.511e-3
    ratio = solar/atmospheric
    print("solar_splitting_eV2 = " + str(solar))
    print("atmospheric_splitting_with_SK_eV2 = " + str(atmospheric))
    print("theory_at_cell_precision = " + shown(value,5))
    print("splitting_ratio = " + shown(ratio,12))
    print("comparator_at_cell_precision = " + shown(ratio,5))
    compare("literal_full_formula", value, ratio)
    compare("printed_cells", "0.27639", "0.03002")
    print("full_ratio_factor = " + shown(value/ratio,9)
          + "; rounded = " + shown(value/ratio,2))
    print("deviation_at_cell_precision = " + shown(100*(value-ratio)/ratio,0,True)+"%")
    print("d_policy = unavailable: no ratio-level uncertainty or joint likelihood supplied")
    print("status = RECORDED NEGATIVE; d=1 is a rung label, not a statistical distance")
    print()
    print('=== row 22b ===')
    print('Scorecard (verbatim): m ν2 /m ν3 = R 2 (d = 2 rung) 1/(φ⁴−1); R 2 ² vs Δm² 21 /Δm² 3l theory 0.02918 experiment 0.03002 (NuFIT 6.1); m 1 =0 NO floor √ratio = 0.17325 vs R 2 = 0.17082 deviation −2.8% · d = −1.44 / −1.79 (below floor) Structural (ladder tier) EDGE TARGET, not a hit — R 2 lies below the m 1 =0 floor and equality is impossible for any m 1 ≥ 0 — i.e. a 1.4–1.8σ tension with the m 1 =0 normal-ordering floor, printed as a tension (Rev32.5); stands only with normal ordering and m 1 ≲ 1 meV; d(ν)=2 is selected, not derived (Paper 5, Rev32.2)')
    rung = 1.0/(PHI**4-1.0)
    solar, u_solar_lower = 7.537e-5, 0.100e-5
    print("selected_rung = " + shown(rung,12))
    print("rung_at_display_precision = " + shown(rung,5))
    print("squared_rung_at_cell_precision = " + shown(rung*rung,5))
    print("solar_splitting_eV2 = " + str(solar))
    print("naive_floor_audit_assumption = independent marginal errors; first-order propagation")
    print("assumption_provenance = NOT IN SUITE; arithmetic diagnostic only, not s1162")
    for name, atmospheric, u_atm_upper in (
            ("without_SK",2.521e-3,0.026e-3),
            ("with_SK",2.511e-3,0.021e-3)):
        ratio = solar/atmospheric
        floor = math.sqrt(ratio)
        sigma = floor/2.0*math.hypot(u_solar_lower/solar,u_atm_upper/atmospheric)
        naive_d = (rung-floor)/sigma
        print(name+".atmospheric_splitting_eV2 = "+str(atmospheric))
        print(name+".solar_lower_error_eV2 = "+str(u_solar_lower))
        print(name+".atmospheric_upper_error_eV2 = "+str(u_atm_upper))
        print(name+".ratio = "+shown(ratio,12))
        print(name+".floor = "+shown(floor,12))
        print(name+".floor_at_display_precision = "+shown(floor,5))
        print(name+".naive_floor_sigma = "+shown(sigma,12))
        print(name+".naive_floor_d = "+shown(naive_d,9,True)
              +"; rounded = "+shown(naive_d,2,True))
        print(name+".squared_rung_signed_percent = "
              +shown(100*(rung*rung-ratio)/ratio,9,True)+"%"
              +"; rounded_1dp = "+shown(100*(rung*rung-ratio)/ratio,1,True)+"%")
        assert rung < floor
    print("printed_squared_cells_percent = "
          +shown(100*(Decimal('0.02918')-Decimal('0.03002'))/Decimal('0.03002'),9,True)+"%")
    print("source_floor_distances = -1.44 / -1.79 (without / with SK); reported s1162")
    print("source_conditions = normal ordering; lightest mass approximately at most 1 meV")
    print("status = EDGE TARGET, not a hit; rung assignment selected, not derived")
    print()


if __name__ == "__main__":
    main()