Check in the evidence behind the pd report

The simulator gitignores its runs/ directory, so every table and figure in
reports/blend/pd rested on data that existed only on one machine. This adds the
sweep outputs of record under reports/blend/pd/data -- one directory per study,
1 MB total -- so any number can be checked against its source, or challenged,
without re-running hours of compute.

report_numbers.py comes with them: run it and it prints every value the report
quotes together with its across-topology standard error, straight from these
parquets. It reproduces the report tables exactly.

Kept: default (8000 rounds/cell), redundancy (9600), percolation (6400),
correlated-churn (6400), fullscale (192, the deliberately lighter 1e6 check).
Omitted: the smoke runs, and an earlier 144-rounds/cell redundancy grid whose
sampling error produced a non-monotonic delivery curve -- superseded, and the
reason the kept grid samples 9600.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

reports/blend/pd/README.md

@@ -335,6 +335,8 @@ The simulator, configs, and analytic checks live in [`tools/simulators/blend/pd`
 
 The figures of record for this report are the copies checked in under [`report-figures/`](report-figures); the simulator does not commit its own generated figures. To regenerate: run the sweeps above, then copy `runs/<…>/figures/*.png` into `report-figures/`.
 
+The **evidence** is checked in too: [`data/`](data) holds the sweep outputs behind every table and figure, one directory per study, with [`data/report_numbers.py`](data/report_numbers.py) regenerating every quoted value together with its standard error directly from them. Any number in this report can therefore be checked against its source without re-running the sweeps — see [`data/README.md`](data/README.md) for what each run is and how it was sampled.
+
 ## Figures
 
 All twenty-two rendered figures are versioned in [`report-figures/`](report-figures): `01`–`03` propagation delay (vs degree, vs path length, vs N); `04`–`09` adversary observation and eclipse (vs `f_adv`, vs degree, and heatmaps); `10`–`11` reliability under churn (delivery and coverage); `12`–`15` deanonymization (whole-path and full, vs path length, `f_adv`, and degree); `16`–`18` linkability over time (time to link vs stake, with redundancy, and time to learn stake vs threshold); `19` the redundancy reliability-vs-anonymity trade-off in probability and `21` the same trade in delivery-vs-time-to-link; `20` the churn-percolation threshold; `22` correlated versus uniform outages. Sixteen of the twenty-two are embedded above; the other six (`04`–`06`, `09`, `11`, `13`) are alternative cuts of data already shown — for instance 11 and 20 both plot coverage against churn, and 20 supersedes 11 by walking the churn past every degree's threshold.

reports/blend/pd/data/README.md

@@ -0,0 +1,43 @@
+# Evidence of record
+
+The sweep outputs behind every number in [the report](../README.md). The simulator does not commit
+its own `runs/` directory — these are the copies of record, kept so that any figure or table can be
+re-derived, or challenged, without re-running hours of compute.
+
+Each run directory holds the three tables the simulator writes: `propagation.parquet`,
+`adversary.parquet` and `deanon.parquet`.
+
+| directory | config | sampling | backs |
+|---|---|---|---|
+| `default/` | `configs/default.yaml` | 1 000 rounds × 8 seeds = **8 000/cell** | §3.1–§3.5 — delay, observation, eclipse, deanonymization, delivery, coverage |
+| `redundancy/` | `configs/redundancy.yaml` | 1 200 × 8 = **9 600/cell** | §3.8 — messaging redundancy R = 1…4 |
+| `percolation/` | `configs/percolation.yaml` | 800 × 8 = **6 400/cell** | §3.5 — the churn threshold `u_c = 1 − 1/(degree − 1)` |
+| `correlated-churn/` | `configs/correlated-churn.yaml` | 800 × 8 = **6 400/cell** | §3.9 — correlated AS/region outages vs uniform churn |
+| `fullscale/` | `configs/fullscale.yaml` | 64 × 3 = **192/cell** | §5 — the 10⁶ scaling check (deliberately lighter; not a source of headline numbers) |
+
+The linkability results (§3.6–§3.7) and both deanonymization rates are closed forms over these
+tables rather than separate measurements, so they have no run of their own — `pd.linkability`
+derives them and `make verify` checks them against Monte-Carlo.
+
+## Regenerating the report's numbers
+
+```
+python report_numbers.py
+```
+
+prints every quoted value with its across-topology standard error, straight from the parquets here.
+That is the fastest way to check a table in the report against its evidence. It takes optional
+paths (`report_numbers.py <default> <redundancy> <percolation>`) if you want to point it at fresh
+runs instead.
+
+## Regenerating the data itself
+
+From [`tools/simulators/blend/pd`](../../../../tools/simulators/blend/pd): `make sweep`,
+`make redundancy`, `make percolation`, `make correlated-churn`, `make sweep-fullscale`. Results
+land in that simulator's `runs/<timestamp>_<label>/`. Note that the seed streams depend on the
+configuration, so re-running reproduces the *statistics*, not bit-identical numbers, unless the
+config is unchanged — in which case it does reproduce exactly.
+
+Two runs from the same session are deliberately **not** kept: the smoke runs (throwaway, far too
+noisy to interpret) and an earlier 144-rounds/cell redundancy grid that was superseded because its
+sampling error produced a non-monotonic delivery curve — the reason `redundancy/` samples 9 600.

reports/blend/pd/data/report_numbers.py

@@ -0,0 +1,131 @@
+"""Regenerate every number the report quotes, with its across-topology standard error.
+
+Run from this directory:  python report_numbers.py
+Each printed value is mean +- SEM over the independent topology seeds, computed from the
+parquets checked in beside this script (see README.md for what each run is).
+"""
+import os
+import sys
+
+import numpy as np
+import pandas as pd
+
+_here = os.path.dirname(os.path.abspath(__file__))
+D = sys.argv[1] if len(sys.argv) > 1 else os.path.join(_here, "default")
+R = sys.argv[2] if len(sys.argv) > 2 else os.path.join(_here, "redundancy")
+PC = sys.argv[3] if len(sys.argv) > 3 else os.path.join(_here, "percolation")
+P = pd.read_parquet(D + "/propagation.parquet")
+A = pd.read_parquet(D + "/adversary.parquet")
+Z = pd.read_parquet(D + "/deanon.parquet")
+PR = pd.read_parquet(R + "/propagation.parquet")
+ZR = pd.read_parquet(R + "/deanon.parquet")
+PP = pd.read_parquet(PC + "/propagation.parquet")
+
+
+def sem(s):
+    return s.std(ddof=1) / np.sqrt(s.count()) if s.count() > 1 else np.nan
+
+
+def cell(df, col):
+    return df[col].mean(), sem(df[col])
+
+
+print(f"rounds/cell: default {P.n_rounds.iloc[0]}x{P.graph_seed.nunique()}"
+      f" | redundancy {PR.n_rounds.iloc[0]}x{PR.graph_seed.nunique()}"
+      f" | percolation {PP.n_rounds.iloc[0]}x{PP.graph_seed.nunique()}")
+
+print("\n### 3.1 full delay (s), N=1e5 mbd=3 u=0 : mean+-SEM")
+b = P[(P.n_nodes == 100000) & (P.max_blend_delay == 3) & (P.unresponsive_frac == 0.0)]
+for bh in sorted(b.blend_hops.unique()):
+    out = []
+    for d in sorted(b.degree.unique()):
+        m, e = cell(b[(b.blend_hops == bh) & (b.degree == d)], "full_delay_ms_mean")
+        out.append(f"d{d}:{m/1000:.2f}+-{e/1000:.3f}")
+    print(f" bh={bh}  " + " ".join(out))
+print(" per-hop cost (s) by degree:")
+for d in sorted(b.degree.unique()):
+    g = b[b.degree == d].groupby("blend_hops").full_delay_ms_mean.mean() / 1000
+    print(f"   d={d:<3} 1->2 {g[2]-g[1]:.2f}  2->3 {g[3]-g[2]:.2f}  3->5 {(g[5]-g[3])/2:.2f}")
+
+print("\n### 3.2 composition N=1e5 deg8 bh3 ; and N-scaling")
+g = b[(b.degree == 8) & (b.blend_hops == 3)]
+for c in ("path_delay_ms_mean", "broadcast_delay_ms_mean", "full_delay_ms_mean",
+          "cover50_ms", "cover90_ms", "cover99_ms"):
+    m, e = cell(g, c)
+    print(f"   {c:<24} {m:8.0f} +- {e:.0f} ms")
+for n in sorted(P.n_nodes.unique()):
+    m, e = cell(P[(P.n_nodes == n) & (P.degree == 8) & (P.blend_hops == 3)
+                  & (P.max_blend_delay == 3) & (P.unresponsive_frac == 0)], "full_delay_ms_mean")
+    print(f"   N={n:<8} full {m/1000:.2f}+-{e/1000:.3f} s")
+print("   cover99 by degree (N=1e5,bh3):",
+      {int(d): round(b[(b.degree == d) & (b.blend_hops == 3)].cover99_ms.mean())
+       for d in sorted(b.degree.unique())})
+
+print("\n### 3.3 observed / eclipsed (exact) N=1e5 random")
+ab = A[(A.n_nodes == 100000) & (A.adversary_mode == "random")]
+print(ab.groupby(["f_adv", "degree"]).observed_frac.mean().unstack().round(3).to_string())
+print(ab.groupby(["f_adv", "degree"]).eclipsed_frac.mean().unstack().round(4).to_string())
+print(" random vs worstcase_coverage observed, AT DEGREE 8 (not averaged over degrees):")
+wc = A[(A.n_nodes == 100000) & (A.degree == 8)
+       & A.adversary_mode.isin(["random", "worstcase_coverage"])]
+print(wc.groupby(["f_adv", "adversary_mode"]).observed_frac.mean().unstack().round(3).to_string())
+print(" ...and eclipse random vs worstcase_eclipse at degree 4:")
+we = A[(A.n_nodes == 100000) & (A.degree == 4)
+       & A.adversary_mode.isin(["random", "worstcase_eclipse"])]
+print(we.groupby(["f_adv", "adversary_mode"]).eclipsed_frac.mean().unstack().round(4).to_string())
+
+print("\n### 3.4 deanon (exact) N=1e5 random deg8 : deanon_rate by f_adv x hops")
+zz = Z[(Z.n_nodes == 100000) & (Z.adversary_mode == "random") & (Z.degree == 8)]
+print(zz.groupby(["f_adv", "blend_hops"]).deanon_rate.mean().unstack().round(5).to_string())
+print(" full_deanon vs degree (bh=2):")
+print(Z[(Z.n_nodes == 100000) & (Z.adversary_mode == "random") & (Z.blend_hops == 2)]
+      .groupby(["f_adv", "degree"]).full_deanon_rate.mean().unstack().round(4).to_string())
+
+print("\n### 3.5 delivery (N=1e5 deg8 mbd3): mean+-SEM   [theory (1-u)^bh]")
+dv = P[(P.n_nodes == 100000) & (P.degree == 8) & (P.max_blend_delay == 3)]
+for bh in sorted(dv.blend_hops.unique()):
+    out = []
+    for u in sorted(dv.unresponsive_frac.unique()):
+        if u == 0:
+            continue
+        m, e = cell(dv[(dv.blend_hops == bh) & (dv.unresponsive_frac == u)], "delivery_rate")
+        out.append(f"u{u}:{m:.3f}+-{e:.3f}[{(1-u)**bh:.3f}]")
+    print(f" bh={bh}  " + " ".join(out))
+print("\n### 3.5 coverage (N=1e5 bh1 mbd3): mean+-SEM")
+cc = P[(P.n_nodes == 100000) & (P.blend_hops == 1) & (P.max_blend_delay == 3)]
+for d in sorted(cc.degree.unique()):
+    out = []
+    for u in (0.2, 0.3, 0.5):
+        m, e = cell(cc[(cc.degree == d) & (cc.unresponsive_frac == u)], "frac_reached")
+        out.append(f"u{u}:{m:.4f}+-{e:.4f}")
+    print(f" deg={d:<3} " + " ".join(out))
+
+print("\n### 3.5 PERCOLATION run: coverage vs u (N=1e5, bh=1); u_c = 1-1/(d-1)")
+for d in sorted(PP.degree.unique()):
+    uc = 1 - 1 / (d - 1)
+    g = PP[PP.degree == d].groupby("unresponsive_frac").frac_reached.mean()
+    print(f" deg={d:<3} u_c={uc:.2f}  " + " ".join(f"{u:.1f}:{v:.3f}" for u, v in g.items()))
+
+print("\n### 3.8 REDUNDANCY: delivery vs R (N=20k deg8 bh3): mean+-SEM  [1-(1-p1)^R]")
+pr = PR[(PR.degree == 8) & (PR.blend_hops == 3)]
+for u in sorted(pr.unresponsive_frac.unique()):
+    if u == 0:
+        continue
+    p1 = pr[(pr.unresponsive_frac == u) & (pr.redundancy == 1)].delivery_rate.mean()
+    out, vals = [], []
+    for Rn in (1, 2, 3, 4):
+        m, e = cell(pr[(pr.unresponsive_frac == u) & (pr.redundancy == Rn)], "delivery_rate")
+        vals.append(m)
+        out.append(f"R{Rn}:{m:.3f}+-{e:.3f}[{1-(1-p1)**Rn:.3f}]")
+    mono = all(y >= x - 1e-9 for x, y in zip(vals, vals[1:]))
+    print(f" u={u}  " + " ".join(out) + ("" if mono else "   <<< NON-MONOTONIC"))
+print(" coverage vs R (should be flat -- no union bonus):")
+for d in sorted(PR.degree.unique()):
+    for u in (0.3, 0.5):
+        g = PR[(PR.degree == d) & (PR.blend_hops == 1) & (PR.unresponsive_frac == u)]
+        print(f"   deg={d} u={u}: " +
+              " ".join(f"R{Rn}:{g[g.redundancy==Rn].frac_reached.mean():.4f}" for Rn in (1, 2, 3, 4)))
+print(" deanon vs R (exact, N=20k deg8 bh3 f=0.2 random):")
+zr = ZR[(ZR.degree == 8) & (ZR.blend_hops == 3) & (ZR.f_adv == 0.2)
+        & (ZR.adversary_mode == "random")]
+print(zr.groupby("redundancy")[["deanon_rate", "full_deanon_rate"]].mean().round(4).to_string())