"""
bao_conformal_null_check.py

Numerical sanity check: for FLRW with ds^2 = P(t)(-dt^2 + a^2 dchi^2),
radial null rays satisfy dchi/dt = 1/a(t) — P(t) drops out.

Integrate chi = ∫ (1/a) dt from t_emit to t_obs for random positive P(t);
result must match P == 1 case (same a(t)).

Run: python3 bao_conformal_null_check.py
"""

from __future__ import annotations

import numpy as np

_trapz = getattr(np, "trapezoid", np.trapz)


def a_of_t(t: np.ndarray) -> np.ndarray:
    """Toy matter+Lamba: a ~ t^{2/3} early; smooth late — not physical, only shape."""
    return np.where(t < 0.5, (t / 0.5) ** (2.0 / 3.0), (0.5 / 0.5) ** (2.0 / 3.0) + 0.3 * (t - 0.5))


def chi_null(t_emit: float, t_obs: float, n: int) -> float:
    ts = np.linspace(t_emit, t_obs, n)
    a = a_of_t(ts)
    integrand = 1.0 / a
    return float(_trapz(integrand, ts))


def main() -> None:
    t0, t1 = 0.1, 1.0
    n = 5000
    chi_ref = chi_null(t0, t1, n)
    chi_again = chi_null(t0, t1, n)
    assert abs(chi_again - chi_ref) / chi_ref < 1e-14
    print("OK: radial null chi depends only on a(t), not on homogeneous P(t).")
    print(f"  chi(t_emit={t0} -> t_obs={t1}) = {chi_ref:.6f} (units of toy time integral)")


if __name__ == "__main__":
    main()