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

Computation 212 -- what a point carries and what a device is, on the theory's own measure:
fibres hold at most seven properties and are singletons almost surely, and no device of
four or more properties is a substrate accident
=========================================================================================
WHAT THIS DOES.  Two open remarks shared a cause.  Backlog B.6: the book leaves unsaid
what the fibres of the realisation map rho are (singletons, or a constant size the
normalisation absorbs).  sec:open-localisation (comp 156): the accident rates that say a
device of six properties is not something an unclustered configuration produces were
computed against a uniform null, because PST fixed no measure over configurations.  Since
comp 141 (the unique G2-invariant law on directions) and comp 211 (maximum entropy is the
product law) PST does fix that measure: directions i.i.d. under the G2-invariant law on
S^6.  This file computes both remarks on it.

  1. A FIBRE IS A CLIQUE OF PAIRWISE-NULL PROPERTIES, AND A CLIQUE IN V_7 IS SMALL.  Two
     properties in one fibre have g = 0, so T({a,b}) = tau_null exactly (eq:metric).  Any
     G2-invariant pair tension is a function of the angle between v(a) and v(b) alone (G2
     is transitive on pairs at a given angle: transitive on S^6, with stabiliser SU(3)
     transitive on the unit sphere of the orthogonal complement), so when that function is
     injective a fibre is a set of unit vectors in R^7 with one pairwise inner product c_0.
     Its Gram matrix is (1 - c_0) I + c_0 J, of rank k when c_0 > -1/(k-1), and a rank-k
     Gram matrix does not embed in R^7 for k > 7.  So |fibre| <= 7 whenever c_0 > -1/6,
     and |fibre| <= 1 + 1/|c_0| (Rankin) when c_0 <= -1/6; two when the null angle is
     antipodal.  Verified by constructing the clique at k_max for nine values of c_0
     (residual 1e-15) and failing at k_max + 1 (residual 0.02 to 1.0).  The seven is the
     dimension of V_7: in R^8 the same construction allows eight.  At the Rankin bound the
     clique's net tension vanishes exactly, |sum v|^2 = k + k(k-1) c_0 = 0: a maximal fibre
     is a balanced sub-configuration.  If the pair tension takes tau_null at s angles rather
     than one, the fibre is an s-distance set and the bound is C(7+s, s) (Delsarte, Goethals
     and Seidel); the point stands that it is finite and set by the dimension.
     Where the null angle sits fixes the bound: at the 5 percent quantile of the pair
     tensions cos theta_0 = -0.62 and a fibre holds at most two properties; at the 20
     percent quantile three; at the median seven.

  2. FIBRES ARE SINGLETONS ALMOST SURELY.  Under i.i.d. directions the pair tension has a
     continuous law, so an exact coincidence T = tau_null has probability zero: over 100
     replicas at N = 50 and 200 and 20 at 800, no pair lands on tau_null, and the nearest
     approach falls like N^-2 (median 7e-4, 5e-5, 2e-6), the spacing of N(N-1)/2 continuous
     values.  A1's |D_n| a_n^4 -> V_M therefore reads |D_n| as the point count, and B.6's
     "constant size absorbed by the normalisation" is the singleton case, with seven the
     most any fibre could ever hold.

  3. NO DEVICE OF THREE OR MORE IS AN ACCIDENT.  comp 156's device criterion (a subset
     whose diameter in the proper distance |T - tau_null|^{1/4} is below its separation
     from the rest, found exhaustively through the nearest-neighbour family) on the
     theory's measure, with tau_null at the 5 percent and the 50 percent quantile, 100
     replicas of 60 properties and 40 of 200: every nearest-neighbour pair qualifies (rate
     1.000 at size two, as on the uniform null); at the low placement no subset of three or
     more ever does, at the median placement a triple appears in 3 to 5 percent of the
     replicas, and nothing of four or more appears at either (the uniform null gave 0.97,
     0.29, 0.04 at sizes three to five).  The uniform null was generous: on the theory's
     own measure a device of any size beyond a pair is a configuration-level object, which
     is where comp 156 placed it, and the margin result there is unchanged.

WHAT THIS CHANGES.  B.6 closes: fibres are singletons mu-almost surely and hold at most
seven properties in every case, the seven being dim V_7.  The caveat on
sec:open-localisation's accident rates closes the other way: computed on PST's own
measure the rates are sharper, not weaker, and the device question is cleanly
configuration-level.  What stays is what comp 156 left: why T has clusters at all is
T(C)-contingent, and why a bound state forms and persists is inherited with the projected
field theory.

TAMPERS.  (a) The dimension: in R^8 the clique at c_0 = 0 reaches eight, so the bound is
V_7's and not an artefact of the construction.  (b) A planted compact device of four
properties (pairwise inner product c_0 to one part in ten thousand) is found by the criterion in 48 of 50
replicas (the two misses have a background point inside the planted separation), so the zero rates in
(3) are the measure's, not the detector's.  (c) The quantile: the
device rates are zero from size four at both tau_null placements, so the result does not
depend on where the null angle sits.

Runs in a few seconds with numpy.
"""
import math
import numpy as np

rng = np.random.default_rng(212)


# ------------------------------------------------------------------ 1. the clique bound
def kmax_theory(c, d=7):
    if c <= -1.0 / (d - 1):
        return 1 + int(math.floor(1.0 / abs(c) + 1e-9))
    return d


def embed_residual(k, c, d=7):
    """best rank-d PSD approximation of the Gram matrix (1-c)I + cJ; zero iff k unit vectors
    with pairwise inner product c exist in R^d"""
    G = (1 - c) * np.eye(k) + c * np.ones((k, k))
    w, V = np.linalg.eigh(G)
    idx = np.argsort(w)[::-1][:d]
    X = V[:, idx] * np.sqrt(np.maximum(w[idx], 0))
    return float(np.abs(X @ X.T - G).max()), X


# ------------------------------------------------------------------ 2 and 3. the measure
def directions(n, d=7):
    v = rng.normal(size=(n, d)); return v / np.linalg.norm(v, axis=1, keepdims=True)


def tension_matrix(v):
    return np.sqrt(np.maximum(2 + 2 * (v @ v.T), 0))        # the exemplar pair tension |v_a + v_b|


def null_cos(q, n=4000):
    v = directions(n); G = (v @ v.T)[np.triu_indices(n, 1)]
    return float(np.quantile(G, q))


def devices(v, tau, kmax=12):
    """comp 156's criterion through its completeness lemma: a compact-and-separated subset is
    exactly the k nearest neighbours of each of its members, so (seed, k) is exhaustive"""
    T = tension_matrix(v); d = np.abs(T - tau) ** 0.25; np.fill_diagonal(d, np.inf)
    N = len(v); found = {}
    order = np.argsort(d, axis=1)
    for s in range(N):
        for k in range(2, kmax + 1):
            S = np.concatenate(([s], order[s, :k - 1]))
            sub = d[np.ix_(S, S)]; diam = sub[np.isfinite(sub)].max()
            mask = np.ones(N, bool); mask[S] = False
            sep = d[np.ix_(S, mask)].min()
            if diam < sep:
                found[k] = found.get(k, 0) + 1
    return found


def device_rates(N, R, tau, kmax=8):
    hit = {k: 0 for k in range(2, kmax + 1)}
    for _ in range(R):
        f = devices(directions(N), tau, kmax)
        for k in hit:
            hit[k] += 1 if f.get(k, 0) > 0 else 0
    return {k: hit[k] / R for k in hit}


def main():
    print("Computation 212 -- what a point carries and what a device is, on the theory's own measure")
    print("=" * 90)
    print("1. a fibre is a clique with one pairwise inner product c0; the largest clique in R^7")
    ok1 = True
    for c in (-0.9, -0.6, -0.4, -0.25, -1 / 6, -0.1, 0.0, 0.3, 0.7):
        km = kmax_theory(c)
        r_ok, X = embed_residual(km, c); r_fail, _ = embed_residual(km + 1, c)
        net = float(np.linalg.norm(X.sum(0)) ** 2)
        print(f"   c0 = {c:+.3f}: k_max = {km}; residual at k_max {r_ok:.1e}, at k_max + 1 {r_fail:.3f}; |sum v|^2 of the clique {net:.3f}")
        ok1 &= r_ok < 1e-9 and r_fail > 0.01
    r8, _ = embed_residual(8, 0.0, d=8)
    print(f"   tamper, the dimension: eight orthonormal vectors in R^8, residual {r8:.1e}; in R^7, {embed_residual(8, 0.0)[0]:.3f}")
    ok1 &= r8 < 1e-9
    print("   where the null angle sits:")
    bounds = {}
    for q in (0.05, 0.2, 0.5):
        c0 = null_cos(q); bounds[q] = kmax_theory(c0)
        print(f"      tau_null at quantile {q:.2f} of the pair tensions: cos theta_0 = {c0:+.3f}, a fibre holds at most {bounds[q]} properties")
    ok1 &= bounds[0.05] == 2 and bounds[0.5] == 7
    print(f"   the bound is seven, V_7's dimension, and smaller for an obtuse null angle: {ok1}")
    print()
    print("2. fibres are singletons almost surely: nearest approach of a pair tension to tau_null (cos theta_0 = -0.6)")
    tau = math.sqrt(2 + 2 * (-0.6)); med = {}
    zeros = 0
    for N, R in ((50, 100), (200, 100), (800, 20)):
        mins = []
        for _ in range(R):
            T = tension_matrix(directions(N))[np.triu_indices(N, 1)]
            mins.append(float(np.abs(T - tau).min()))
        mins = np.array(mins); med[N] = float(np.median(mins)); zeros += int((mins == 0).sum())
        print(f"   N = {N:4d} ({N * (N - 1) // 2:6d} pairs, {R} replicas): median {med[N]:.1e}, smallest {mins.min():.1e}, exact coincidences {(mins == 0).sum()}")
    ok2 = zeros == 0 and med[200] / med[800] > 8 and med[50] / med[200] > 8
    print(f"   no coincidence, and the approach falls like N^-2: {ok2}")
    print()
    print("3. device accident rates on the theory's measure (comp 156's criterion; the uniform null gave 1.00, 0.97, 0.29, 0.04, 0.00 at sizes 2 to 6)")
    ok3 = True
    for q, c0 in ((0.05, -0.62), (0.5, 0.0)):
        tau_q = math.sqrt(2 + 2 * c0)
        for N, R in ((60, 100), (200, 40)):
            rates = device_rates(N, R, tau_q)
            print(f"   tau_null at quantile {q:.2f}, N = {N:3d}, {R} replicas: " + "  ".join(f"size {k}: {rates[k]:.3f}" for k in range(2, 8)))
            ok3 &= rates[2] == 1.0 and rates[3] <= 0.1 and all(rates[k] == 0.0 for k in range(4, 8))
    print(f"   every nearest-neighbour pair qualifies, a triple in at most a few percent of replicas, nothing of four or more, at either placement of the null: {ok3}")
    print()
    print("4. tamper: a planted compact device of four (pairwise inner product c0 to 1e-4) is found")
    hits = 0
    for _ in range(50):
        v = directions(60); _, X = embed_residual(4, -0.3)
        X = np.pad(X, ((0, 0), (0, 7 - X.shape[1]))) + rng.normal(scale=1e-4, size=(4, 7))
        v[:4] = X / np.linalg.norm(X, axis=1, keepdims=True)
        f = devices(v, math.sqrt(2 + 2 * (-0.3)))
        hits += 1 if f.get(4, 0) > 0 else 0
    ok4 = hits >= 45
    print(f"   found in {hits} of 50 replicas: the zero rates above are the measure's, not the detector's: {ok4}")
    print()
    ok = ok1 and ok2 and ok3 and ok4
    print("SUMMARY")
    print(f"  [{'PASS' if ok1 else 'FAIL'}] a fibre of rho holds at most seven properties, fewer for an obtuse null angle")
    print(f"  [{'PASS' if ok2 else 'FAIL'}] fibres are singletons almost surely on the theory's measure")
    print(f"  [{'PASS' if ok3 else 'FAIL'}] no device of four or more properties is a substrate accident, and a triple is rare")
    print(f"  [{'PASS' if ok4 else 'FAIL'}] tamper: a planted device is detected")
    print()
    print("  CONCLUSION: a point of M is a clique of pairwise-null properties, of size at most dim V_7 = 7 and")
    print("  almost surely one; A1 reads |D_n| as the point count (B.6 closes).  On the theory's own measure a")
    print("  device of any size beyond a pair is a configuration-level object, which is where comp 156 placed it;")
    print("  why the configuration has clusters is T(C)-contingent and the bound state's dynamics is inherited.")
    print("  ALL CHECKS PASSED" if ok else "  SOME CHECKS FAILED")


if __name__ == "__main__":
    main()