#!/usr/bin/env python3
"""
PROVENANCE: NUMERICAL

Computation 145 -- Is the tension variable uniform under the Bernoulli measure?
================================================================================
Computation 143 found that the threshold-crossing rule of sec:collapse-rule agrees
with the Born rule derived in sec:born-rule precisely when the threshold variable is
UNIFORMLY distributed over substrate configurations.  This computation asks whether
the Bernoulli measure supplies that.

WHAT THE VARIABLE IS
--------------------
sec:tau-derivation fixes tau as a CONSTANT ("not a free parameter of the theory"), so
the run-to-run randomness is not in tau.  It is in the configuration C, through the
tension magnitude |T(C)|, and the threshold condition is |T(C)| > tau.  The scalar
exemplar the book gives is the induced edge count

    T_ex(C) = |{(a,b) : a,b in C, a != b, delta(a,b)}|              (eq:texample)

evaluated on C ~ mu = (x)_{a in D} Bern(1/2).

THE TWO QUESTIONS, AND WHY THE SECOND IS THE REAL ONE
------------------------------------------------------
  (Q1) Is |T(C)| itself uniform under mu?  Almost certainly not, and it is worth
       measuring how far off.
  (Q2) Does that matter?  The book states that every later derivation reads T only
       through the ORDERING of the magnitudes (after eq:texample), and \\compref{140}
       verifies invariance under monotone reparametrisation of those magnitudes.  If
       only the ordering were used, uniformity would be a gauge choice.  WITHDRAWN
       (v32.64): the antecedent is FALSE for the read that matters.  The
       Born-consistency comparison of eq:tau-i-explicit tests tau_0 against Born
       weights carrying their own normalisation sum_i w_i = 1, so it is a CARDINAL
       read; monotone reparametrisation changes every P_i, and even an affine phi
       does (phi(x)=x/2 sends (.5,.3,.2) to (1.0,.6,.4), not a probability vector).
       comp 140 co-transforms BOTH sides, which is the hypothesis being dropped.
       The obstruction is not the shape of the distribution but its ATOMS -- ties in
       |T(C)| cannot be separated by any monotone map, and they bound how close to
       uniform any reparametrisation can get.

CHECKS
------
  (1) The raw distribution of T_ex under mu, against uniform, for several distinction
      graphs.  Kolmogorov-Smirnov distance reported.
  (2) The probability-integral transform of the same data, showing it is uniform up to
      the atoms -- so the shape objection dissolves.
  (3) THE ATOM BOUND, which is the thing that actually limits uniformity.  For any
      random variable X with largest atom of mass m, every monotone reparametrisation
      g satisfies KS(g(X), Uniform) >= m/2.  Measured per graph.
  (4) The regular-graph case, where T_ex degenerates to a function of |C| alone.  Here
      atoms are huge and uniformity is unreachable.  Note this is the same sub-class
      the book already excludes, because T_ex fails complement-asymmetry there.
  (5) Scaling of the largest atom with |D| for generic graphs, which says whether the
      residual non-uniformity vanishes in the large-substrate limit.
"""

import numpy as np
from itertools import combinations

rng = np.random.default_rng(20260926)


def edge_list(D, kind, p=0.5, seed=0):
    r = np.random.default_rng(seed)
    if kind == "complete":
        return list(combinations(range(D), 2))
    if kind == "cycle":
        return [(i, (i + 1) % D) for i in range(D)]
    if kind == "erdos":
        return [e for e in combinations(range(D), 2) if r.random() < p]
    raise ValueError(kind)


def tension_distribution(D, edges):
    """Exact distribution of T_ex over all 2^D configurations, each of mass 2^-D."""
    E = np.array(edges, dtype=int)
    vals = np.empty(1 << D, dtype=np.int64)
    for m in range(1 << D):
        inC = np.array([(m >> a) & 1 for a in range(D)], dtype=bool)
        vals[m] = int(np.sum(inC[E[:, 0]] & inC[E[:, 1]])) if len(E) else 0
    return vals


def ks_vs_uniform(x):
    """KS distance between the empirical law of x (rescaled to [0,1]) and Uniform."""
    x = np.asarray(x, float)
    if np.ptp(x) == 0:
        return 1.0
    u = (x - x.min()) / np.ptp(x)
    u = np.sort(u)
    n = len(u)
    d1 = np.max(np.arange(1, n + 1) / n - u)
    d2 = np.max(u - np.arange(0, n) / n)
    return max(d1, d2)


def pit(x):
    """Probability-integral transform using the empirical CDF (mid-rank for ties)."""
    x = np.asarray(x, float)
    order = np.argsort(x, kind="mergesort")
    ranks = np.empty(len(x), float)
    sx = x[order]
    i = 0
    while i < len(sx):
        j = i
        while j + 1 < len(sx) and sx[j + 1] == sx[i]:
            j += 1
        ranks[order[i:j + 1]] = 0.5 * (i + j + 1) / len(sx)
        i = j + 1
    return ranks


def largest_atom(x):
    _, counts = np.unique(x, return_counts=True)
    return counts.max() / len(x)


print("=" * 78)
print("Computation 145 -- is the tension variable uniform under the Bernoulli measure?")
print("=" * 78)

D = 14
GRAPHS = [("complete K_D", "complete", 0.0), ("cycle C_D", "cycle", 0.0),
          ("Erdos-Renyi p=0.5", "erdos", 0.5), ("Erdos-Renyi p=0.2", "erdos", 0.2)]

print(f"\n|D| = {D}, all 2^{D} = {1<<D} configurations enumerated exactly.\n")
print("(1,2,3) RAW vs REPARAMETRISED, AND THE ATOM BOUND")
print(f"    {'graph':>20} {'#values':>8} {'KS raw':>9} {'KS after PIT':>13}"
      f" {'largest atom':>13} {'atom bound':>11}")
rows = []
for name, kind, p in GRAPHS:
    E = edge_list(D, kind, p, seed=3)
    v = tension_distribution(D, E)
    ks_raw = ks_vs_uniform(v)
    ks_pit = ks_vs_uniform(pit(v))
    m = largest_atom(v)
    rows.append((name, len(np.unique(v)), ks_raw, ks_pit, m))
    print(f"    {name:>20} {len(np.unique(v)):>8} {ks_raw:>9.4f} {ks_pit:>13.4f}"
          f" {m:>13.4f} {m/2:>11.4f}")

print("\n    -> RAW: the tension is far from uniform in its own units, as expected;")
print("       it concentrates near its mean.")
print("    -> AFTER the probability-integral transform -- itself a monotone")
print("       reparametrisation, so inside the invariance \\compref{140} verifies --")
print("       the distribution IS uniform, to within the atom bound.")
print("    -> the atom bound m/2 is the floor no monotone map can beat: ties in")
print("       T_ex cannot be separated by any order-preserving relabelling.")

for name, nv, ks_raw, ks_pit, m in rows:
    assert ks_pit <= m / 2 + 0.02, (name, ks_pit, m)
print("\n    every measured KS-after-PIT respects its own atom bound.")

# ---- CHECK 4 -------------------------------------------------------------
print("\n(4) THE REGULAR-GRAPH DEGENERACY")
for name, kind, p in [("complete K_D", "complete", 0.0), ("cycle C_D", "cycle", 0.0)]:
    E = edge_list(D, kind, p, seed=3)
    v = tension_distribution(D, E)
    print(f"    {name:>20}: {len(np.unique(v)):>4} distinct values over {1<<D} "
          f"configurations, largest atom {largest_atom(v):.4f}")
print("    for a regular graph T_ex depends on C only through |C|, so it collapses to")
print("    ~|D| values with large atoms at the sizes reachable here.")
print("    CORRECTED v32.78 (comp 180).  This block used to add that the book already")
print("    sets these aside 'because T_ex fails complement-asymmetry on regular")
print("    (D, delta)'.  That is wrong twice over.  eq:ca is a condition on v and not")
print("    on delta -- the body says its genuine content is the non-balance")
print("    |sum_a v(a)| > 0 -- so it excludes no graph.  And these atoms DO vanish,")
print("    as |D|^-0.5, slower than generic but vanishing, so regular graphs are not")
print("    the failing class at all.  The class that fails has a BOUNDED ACTIVE PART")
print("    (a star, K_k plus isolated vertices), is admissible on eq:pos and eq:mon,")
print("    and keeps an O(1) atom at every |D| (comp 180).")

# ---- CHECK 5 -------------------------------------------------------------
print("\n(5) DOES THE RESIDUAL VANISH AS THE SUBSTRATE GROWS?  (Erdos-Renyi p=0.5)")
scal = []
print(f"    {'|D|':>5} {'#values':>9} {'largest atom':>13} {'KS after PIT':>13}")
prev = None
for Dn in [8, 10, 12, 14, 16]:
    E = edge_list(Dn, "erdos", 0.5, seed=3)
    v = tension_distribution(Dn, E)
    m = largest_atom(v)
    k = ks_vs_uniform(pit(v))
    print(f"    {Dn:>5} {len(np.unique(v)):>9} {m:>13.5f} {k:>13.5f}")
    scal.append((Dn, m, k))

Dv = np.array([r[0] for r in scal], float)
for label, idx in (("largest atom", 1), ("KS after PIT", 2)):
    y = np.array([r[idx] for r in scal], float)
    q = np.polyfit(np.log(Dv), np.log(y), 1)
    ss = 1 - np.sum((np.log(y) - np.polyval(q, np.log(Dv))) ** 2) / \
             np.sum((np.log(y) - np.log(y).mean()) ** 2)
    print(f"    fit: {label} ~ |D|^({q[0]:.2f}), R^2 = {ss:.4f};"
          f"  extrapolated to |D|=256: {np.exp(np.polyval(q, np.log(256))):.1e}")
    assert q[0] < -1.0 and ss > 0.99
print("    -> both decay faster than 1/|D|, so the residual non-uniformity is")
print("       negligible for any substrate large enough to carry the emergent")
print("       physics (|D| = 6 is the matched-scaling value used elsewhere, but the")
print("       measurement setting involves an apparatus, i.e. many more properties).")

print("\n" + "=" * 78)
print("""SUMMARY

  The question "does the Bernoulli measure give a uniform tau_0" has a two-part
  answer, and the first part is misleading on its own.

  NO, not in raw units.  The tension magnitude concentrates near its mean under
  mu; KS distance from uniform is large for every distinction graph tried.

  BUT THAT IS THE WRONG TEST.  The book states, right after eq:texample, that
  every later derivation reads T only through the ORDERING of the magnitudes,
  and comp 140 verifies invariance under monotone reparametrisation -- but only with
  BOTH sides co-transformed.  The crossing rule compares tau_0 against Born weights
  that are NOT reparametrised, so that invariance does not reach it.  Uniformity is
  therefore a SUBSTANTIVE CONSTRAINT on the law of T over substrate configurations,
  not a gauge choice.  What this computation establishes is that it is SATISFIABLE.

  THE REAL OBSTRUCTION IS ATOMS, and it is small.  Ties in T_ex survive every
  monotone map and put a floor m/2 on the achievable KS distance.  For generic
  (Erdos-Renyi) distinction graphs that floor is already small at |D| = 14 and
  shrinks as |D| grows.  For regular graphs it is larger at these sizes but still
  shrinks.  CORRECTED v32.78 (comp 180): this used to say regular graphs are
  excluded by complement-asymmetry.  eq:ca constrains v rather than delta and
  excludes no graph, and regular graphs are not the failing class anyway.  The
  graphs that keep an O(1) atom at every |D| have a bounded active part, and they
  are admissible.

  CONSEQUENCE.  The Born-compatibility condition found in Computation 143 is
  SATISFIABLE, and satisfied in the large-substrate limit for generic
  distinction structures.  It is NOT satisfied for every admissible one: comp 180
  exhibits admissible graphs whose atom floor never vanishes, so discharging the
  item needs a genericity constraint on delta the theory does not supply.  It should be stated as a CONDITION, and NOT discharged by
  the monotone-reparametrisation argument, which does not reach a cardinal read.
""")
