#!/usr/bin/env python3 """ Test 1, change the medium, not the state (October 14, 2026): the exploratory pilot planned in plan.md (revision 2). In the oscillator model the pull within each group, K, is raised or lowered by the same fraction for 200 time units and set back exactly. The arms are compared with the untouched run, and with pushes matched to the same coherence at the restoration time R, each push with its own sham. nudge_test.py's model and step code and the cause test's helpers are imported unchanged; the earlier folders are only read. A sandbox script kept with the project. Nothing here is served by the site. Usage (each command runs in the foreground and ends by itself; nothing is left running) python3 -B medium_test.py checks # the code checks on old start 1: equality and mechanics only; checks.txt python3 -B medium_test.py dryreport # the pilot's report on made-up numbers; writes only to /tmp python3 -B medium_test.py calibrate # one piece of the calibration on old starts 16 and 17 (coherence only) python3 -B medium_test.py freeze # the calibration's report and the freeze record; refused once any pilot # output exists python3 -B medium_test.py pilot # one piece of the pilot on starts 103 to 114 (needs the freeze, and # Colson's go) python3 -B medium_test.py report # the pilot's readings, from the pilot's own starts only calibrate and pilot take an optional last argument: the seconds into a piece after which no task begins (default 120; operational only, it changes no number). Each is run with a hard limit: timeout -k 10 285 python3 -B medium_test.py """ import json import shutil import sys import time import warnings from pathlib import Path HERE = Path(__file__).resolve().parent ROOT = HERE.parents[2] ANALYSIS = ROOT / "analysis" sys.dont_write_bytecode = True # leave every folder as it is sys.path.insert(0, str(ANALYSIS / "cause-prelock-2026-10-13")) from prelock_common import (ct, nt, np, stamp, shown, sha256, write_json, plain, save_npz, run_piece, # noqa: E402 guard, lock_changes, machine) N, K, EPS, DT, GAP, EVERY, RUN = ct.N, ct.K, ct.EPS, ct.DT, ct.GAP, ct.EVERY, ct.RUN SEED = 20261014 # draws: (SEED, 1, seed) a start's moment offset; (SEED, 2, seed, moment, size, direction) the # calibration's sham signs; (SEED, 3, seed, moment, direction) the pilot's; (SEED, 9, ...) checks CAL_STARTS = (15, 16) # old starts 16 and 17 (seeds) PILOT_STARTS = tuple(range(102, 114)) # fresh starts 103 to 114 (seeds 102 to 113) CHECK_START = 0 # old start 1: the code checks only MOMENTS, FIRST, SPACING, OFFSETS = 10, 2000, 330, 66 # T = 2,000 + 5 x (0 to 65) + 330 m, m = 0 to 9 WINDOW, AFTER, KEPT_FOR = 200, 500, 700 # the changed pull; followed after R; the kept-changed reference SIZES = (0.001, 0.0025, 0.005, 0.01, 0.02, 0.05) # fractions of K, the same size raised and lowered DIRECTIONS = ("raised", "lowered") MARKS = dict(raise_low=1.5, raise_high=3.0, lower_high=0.67, built=18) LATER, EARLY, BLOCK = (200, 500), (0, 200), 100 # growth windows after R; blocks for the reported curves BENCHMARK, HALF, CLEAR, SIDE = 0.1602, 0.5, 3.0, 2.0 RESAMPLES, RESAMPLE_SEED, MOST_MISSING = 400, 20261016, 12 SHAM_BAND = (0.95, 1.05) SCAN, LARGEST, BISECT, SOLVED, MATCH, KEPT = ct.SCAN, ct.LARGEST, ct.BISECT, ct.SOLVED, ct.MATCH, ct.KEPT BUDGET, START_BY = 190.0, 120.0 # seconds into a piece: keep a base run's state; begin no task STATES = ("untouched", "medium raised", "medium lowered", "push raised", "push lowered", "sham raised", "sham lowered") REC_KEYS = ("R", "target", "sign", "amp", "status", "push_ratio", "sham_ratio", "sumsq", "gaprel", "match", "same", "built") DESC_KEYS = ("R", "inner", "outer", "z2", "z3", "z4", "C", "ab") WORK = ROOT / "generated-artifacts" / "medium-test" # saved states and parts (kept out of the project's history) CAL_PARTS, PILOT_PARTS = WORK / "calibration", WORK / "pilot" PLAN, DRAFT = HERE / "plan.md", HERE / "plan-draft-1.md" CHECKS_TXT = HERE / "checks.txt" CAL_LOG, CAL_TXT, CAL_JSON = (HERE / "calibration_runlog.json", HERE / "calibration_results.txt", HERE / "calibration_results.json") FREEZE = HERE / "medium_freeze.json" PILOT_LOG, PILOT_TXT, PILOT_JSON = HERE / "pilot_runlog.json", HERE / "pilot_results.txt", HERE / "pilot_results.json" DRY = Path("/tmp/medium-test-dry") PROTECTED_FOLDERS = tuple(ANALYSIS / name for name in ( "cause-test-2026-10-12", "cause-explore-2026-10-12", "cause-prelock-2026-10-13", "nudge-test-2026-10-10", "coherence-growth-2026-10-11")) SETTINGS = dict( n=N, k=K, eps=EPS, delta=nt.DELTA, dt=DT, settle=nt.SETTLE, run=RUN, gap=GAP, put_back_every=EVERY, gap_seed=nt.LYAP_SEED, seed=SEED, calibration_starts=[s + 1 for s in CAL_STARTS], pilot_starts=[s + 1 for s in PILOT_STARTS], moments=MOMENTS, first=FIRST, spacing=SPACING, offsets=OFFSETS, window=WINDOW, after=AFTER, kept_for=KEPT_FOR, sizes=list(SIZES), marks=MARKS, later=list(LATER), early=list(EARLY), block=BLOCK, benchmark=BENCHMARK, half=HALF, clear=CLEAR, side=SIDE, resamples=RESAMPLES, resample_seed=RESAMPLE_SEED, most_missing=MOST_MISSING, sham_band=list(SHAM_BAND), scan=SCAN, largest=LARGEST, bisect=BISECT, solved=SOLVED, match=MATCH, gap_kept=KEPT) PLAN_NUMBERS = dict( # the plan's numbers, which the code must match n=4000, k=1.9, eps=0.1, dt=0.025, settle=300.0, run=6000.0, gap=1e-7, put_back_every=5.0, calibration_starts=[16, 17], pilot_starts=list(range(103, 115)), moments=10, first=2000, spacing=330, offsets=66, window=200, after=500, kept_for=700, sizes=[0.001, 0.0025, 0.005, 0.01, 0.02, 0.05], marks=dict(raise_low=1.5, raise_high=3.0, lower_high=0.67, built=18), later=[200, 500], benchmark=0.1602, half=0.5, clear=3.0, side=2.0, resamples=400, resample_seed=20261016, most_missing=12, sham_band=[0.95, 1.05], largest=np.pi / 2, solved=1e-6, match=1e-9, gap_kept=1e-6, seed=20261014) MEANING = { "differs": ("The medium differs from the push: the two intervention histories left states that behave differently," " inside the model. This is a difference between histories and the states they leave, not evidence of" " an independent memory variable."), "like push": ("The medium is like the push, a useful null: a brief change to the medium did no more than a push to" " the same coherence. This sharpens Test 2's question without committing to it."), "non-monotonic": ("Changing the pull at all moved later growth, either way. For harm reduction: judge any" " intervention against an untouched comparison."), "opposite direction": ("A result in its own right, not a failure: the natural correlation's direction doesn't carry" " over to this lever."), "same direction": "Later growth moved the natural correlation's way, against untouched runs.", "lasting switch": ("The difference stopped fading. First confirm the pull really was restored (see the checks), then" " describe what holds the difference."), "ruled out": "No direction-dependent effect as large as half the benchmark, in either direction.", "unresolved": "The error bars can't tell half the benchmark from zero: inconclusive about size, not negative.", "push inconclusive": ("More than 1 moment in 10 lacked a matched push or sham: the medium-vs-push readings are" " inconclusive."), } # --------------------------------------------------------------------------------------------------------------------- # small helpers def moment_times(seed): """A start's ten moments T: 2,000 + one offset (a multiple of 5 from 0 to 325, drawn from (SEED, 1, seed)) + 330 m. The last ends its arms at T + 700 <= 5,995, inside the 6,000-unit base run.""" offset = 5 * int(np.random.default_rng((SEED, 1, int(seed))).integers(0, OFFSETS)) return [FIRST + offset + SPACING * m for m in range(MOMENTS)] def keep_times(seed): return sorted({t for T in moment_times(seed) for t in (T, T + WINDOW)}) _FORMS = {} def form(copies): """Arrays made exactly as the base run makes its two copies (copies 2) or its settling state (copies 1). Sums over the free turners follow the memory layout to the last digit, so a state read back from a file, or a branch's copy, is laid out like these.""" if not _FORMS: slow, _ = nt.model(N, K, EPS) u, f = nt.starting_phases(N, slow, 0, 1) _FORMS[1] = dict(u=u, f=f) _FORMS[2] = dict(u=np.concatenate([u, u * np.exp(1j * np.zeros((1,) + u.shape[1:]))]), f=np.concatenate([f, f])) return _FORMS[copies] def shaped(copies, name, values): out = np.empty_like(form(copies)[name]) out[...] = values return out def coherence(us, uf): """Each row's and group's coherence, exactly as the base run measures it.""" return (us.sum(-1) + uf.sum(-1)) / N def stacked(items, empty_shape): return np.stack(items) if items else np.zeros(empty_shape, complex) def growth(logs, lo, hi): """Gap growth per time unit over (lo, hi], from the log growth per 5-unit interval counted from the arm's start.""" return float(np.sum(logs[lo // 5:hi // 5]) / (hi - lo)) def protected_fingerprints(): """Every file in the earlier folders and in cancer-research (outside this test's folder), with the first draft.""" paths = [] for folder in PROTECTED_FOLDERS: paths += [p for p in sorted(folder.rglob("*")) if p.is_file() and "__pycache__" not in p.parts] paths += [p for p in sorted(HERE.parent.rglob("*")) if p.is_file() and HERE not in p.parents] paths.append(DRAFT) return {str(p.relative_to(ROOT)): sha256(p) for p in paths} def current_fingerprints(): return dict(plan=sha256(PLAN), code=sha256(Path(__file__).resolve()), calibration_json=sha256(CAL_JSON), calibration_txt=sha256(CAL_TXT)) def shown_path(p): return str(p.relative_to(ROOT)) if p.is_relative_to(ROOT) else str(p) def pilot_outputs(): found = [p for p in sorted(PILOT_PARTS.rglob("*")) if p.is_file()] if PILOT_PARTS.exists() else [] return [shown_path(p) for p in found + [q for q in (PILOT_LOG, PILOT_TXT, PILOT_JSON) if q.exists()]] # --------------------------------------------------------------------------------------------------------------------- # the model runs def base_run(s, deadline, path, keep, run=RUN): """Start s (seed), step for step as the cause test's base run: settling, the second copy with its tiny gap, then `run` time units with the gap measured and put back every 5. The full state (both copies, just after the put-back) is kept at each time in `keep`. At the first put-back after `deadline` (time.time()) the exact state is kept at `path` and None is returned; the next call goes on from it, laid out as the base run's own arrays, so the result is the same to the last digit as one unbroken run (checked).""" clock = time.time() slow, stepper = nt.model(N, K, EPS) run_steps, reset_every, sample_every = (int(round(x / DT)) for x in (run, EVERY, 1.0)) keep_steps = {int(round(t / DT)): t for t in keep} if path.exists(): with np.load(path) as kept: i, spent, pieces = int(kept["i"]), float(kept["spent"]), int(kept["pieces"]) + 1 us, uf = shaped(2, "u", kept["us"]), shaped(2, "f", kept["uf"]) trace = kept["trace"].copy() samples, logs = list(kept["samples"].copy()), list(kept["logs"].copy()) kept_t, kept_u, kept_f = list(kept["kept_t"]), list(kept["kept_u"].copy()), list(kept["kept_f"].copy()) else: us, uf = nt.starting_phases(N, slow, s, 1) settle_steps, trace_every = int(round(nt.SETTLE / DT)), int(round(5.0 / DT)) trace = [np.abs(coherence(us, uf))] step, none = stepper(DT), np.zeros(1) for j in range(1, settle_steps + 1): step(us, uf, none) if j % trace_every == 0: trace.append(np.abs(coherence(us, uf))) trace = np.stack(trace, -1).astype(np.float32) offset = np.stack([np.random.default_rng((nt.LYAP_SEED, s)).standard_normal(us.shape[1:])]) offset -= offset.mean(axis=(-2, -1), keepdims=True) offset *= GAP / np.sqrt((offset * offset).mean(axis=(-2, -1), keepdims=True)) us, uf = np.concatenate([us, us * np.exp(1j * offset)]), np.concatenate([uf, uf]) samples, logs, kept_t, kept_u, kept_f = [coherence(us[:1], uf[:1])], [], [], [], [] i, spent, pieces = 0, 0.0, 1 if us.strides != form(2)["u"].strides or uf.strides != form(2)["f"].strides: raise RuntimeError("the state is not laid out as the base run's") count, free = us.shape[-1], uf.shape[-1] step, none = stepper(DT), np.zeros(2) while i < run_steps: i += 1 step(us, uf, none) if i % sample_every == 0: samples.append(coherence(us[:1], uf[:1])) if i % reset_every == 0: shift = np.angle(us[1:] * np.conj(us[:1])) shift -= shift.mean(axis=(-2, -1), keepdims=True) size = np.sqrt((shift * shift).mean(axis=(-2, -1))) logs.append(np.log(size / GAP)) us[1:] = us[:1] * np.exp(1j * shift * (GAP / size)[:, None, None]) if i in keep_steps: kept_t.append(keep_steps[i]) kept_u.append(us.copy()) kept_f.append(uf.copy()) if time.time() > deadline and i < run_steps: save_npz(path, i=i, spent=spent + time.time() - clock, pieces=pieces, us=us, uf=uf, trace=trace, samples=np.stack(samples), logs=np.stack(logs), kept_t=np.array(kept_t, float), kept_u=stacked(kept_u, (0, 2, 2, count)), kept_f=stacked(kept_f, (0, 2, 2, free))) return None return dict(z=np.stack(samples, -1)[0], logs=np.stack(logs, -1)[0], trace=trace[0], kept_t=np.array(kept_t, float), kept_u=stacked(kept_u, (0, 2, 2, count)), kept_f=stacked(kept_f, (0, 2, 2, free)), spent=spent + time.time() - clock, pieces=pieces) def branch(us, uf, schedule, copies, change=None): """One arm from a kept state (us: (copy, group, pulled), uf: (copy, group, free); the base run's two copies), laid out as the base run's arrays. change (group, pulled), if given, is added to every copy's pulled phases first (a push or a sham: the same change on both copies, so the gap is as it was). schedule: stretches (k, units) run one after another, each with nudge_test's step for the pull k within each group (eps unchanged). copies 2: the gap is measured and put back every 5 units exactly as in the base run; copies 1: the first copy only, with no gap. Returns the first copy's coherence every time unit from the start, (group, samples); the gap's log growth in each 5-unit interval (copies 2); the pull in use in each time unit; and the full state at the end of each stretch.""" U = shaped(copies, "u", us[:copies] if change is None else us[:copies] * np.exp(1j * change)) F = shaped(copies, "f", uf[:copies]) reset_every, sample_every = (int(round(x / DT)) for x in (EVERY, 1.0)) samples, logs, pulls, ends = [coherence(U[:1], F[:1])[0]], [], [], [] none = np.zeros(copies) for k, units in schedule: step = nt.model(N, k, EPS)[1](DT) for i in range(1, int(round(units / DT)) + 1): step(U, F, none) if i % sample_every == 0: samples.append(coherence(U[:1], F[:1])[0]) pulls.append(k) if copies == 2 and i % reset_every == 0: shift = np.angle(U[1:] * np.conj(U[:1])) shift -= shift.mean(axis=(-2, -1), keepdims=True) size = np.sqrt((shift * shift).mean(axis=(-2, -1))) logs.append(np.log(size / GAP)[0]) U[1:] = U[:1] * np.exp(1j * shift * (GAP / size)[:, None, None]) ends.append((U.copy(), F.copy())) return dict(z=np.stack(samples, -1), logs=np.array(logs), pulls=np.array(pulls), ends=ends) # --------------------------------------------------------------------------------------------------------------------- # the matched push and its sham def amplitude(u, along, z, target, sign): """The smallest a in (0, pi/2] with R(a) = target, where R(a) = |z + (1/N) sum_j u_j (e^{i sign a along_j} - 1)| is the group's coherence after moving every pulled oscillator (u: their e^{i theta}) by sign * a * along_j, along_j = sin(psi - theta_j): toward the group's rhythm for sign +1, away for -1. Scanned at pi/512, 2 pi/512, ..., pi/2; the first value reaching the target and the one before are halved between until less than 1e-12 apart (the cause test's solver, one-sided). Returns (a, 0); (nan, 1) if no move up to a quarter turn reaches it (infeasible); (nan, 2) if a number is not finite or the move misses by more than a relative 1e-6 (a solver failure).""" real, imag = u.real.copy(), u.imag.copy() def R(a): angle = sign * np.multiply.outer(np.atleast_1d(a), along) c, s = np.cos(angle) - 1, np.sin(angle) cr, ci, sr, si = c @ real, c @ imag, s @ real, s @ imag return np.abs(z + ((cr - si) + 1j * (ci + sr)) / N) if not (np.isfinite(target) and target > 0 and np.all(np.isfinite(along))): return np.nan, 2 def reached(r): return r >= target if sign > 0 else r <= target grid = LARGEST * np.arange(1, SCAN + 1) / SCAN values = R(grid) if not np.all(np.isfinite(values)): return np.nan, 2 hits = np.flatnonzero(reached(values)) if not len(hits): return np.nan, 1 i = int(hits[0]) low, high = (float(grid[i - 1]) if i else 0.0), float(grid[i]) while high - low > BISECT: middle = (low + high) / 2 if reached(R(middle)[0]): high = middle else: low = middle value = R(high)[0] if not (np.isfinite(value) and abs(value / target - 1) <= SOLVED): return np.nan, 2 return high, 0 def matched(us, uf, targets, key): """The push matched to a medium arm, and its sham, at R from the untouched state (us, uf: both copies), worked out from the first copy. Push: each group taken to targets[g] (the medium arm's coherence at R) by the smallest dense move, toward the group's rhythm if the target is above its coherence, away if below. Sham: the same oscillators moved by the same sizes, with the cause test's paired signs (draw (SEED, *key), group A before B). Built only if both groups are solved, the push lands on target (relative 1e-6), the sham keeps each group within 0.95 to 1.05 of where it was, the sizes are identical (sums of squares within 1e-9) and the gap is as it was (relative 1e-6). Nothing is clipped or replaced; the reason is kept. Returns the record and the two changes (push, sham).""" u, f = us[0], uf[0] z = ct.harmonic(u[None], f[None])[0] R, psi = np.abs(z), np.angle(z) along = np.sin(psi[:, None] - np.angle(u)) targets = np.asarray(targets, float) nan = np.nan rec = dict(R=R, target=targets, sign=np.zeros(2), amp=np.full(2, nan), status=np.zeros(2, int), push_ratio=np.full(2, nan), sham_ratio=np.full(2, nan), sumsq=np.full(2, nan), gaprel=np.full(2, nan), match=nan, same=False, built=False, reason="") dtheta = np.zeros((2,) + u.shape) reasons = [] for g in range(2): sign = 1.0 if targets[g] >= R[g] else -1.0 rec["sign"][g] = sign if np.isfinite(targets[g]) and abs(targets[g] / R[g] - 1) <= SOLVED: a, status = 0.0, 0 else: a, status = amplitude(u[g], along[g], z[g], targets[g], sign) rec["amp"][g], rec["status"][g] = a, status if status: reasons.append(f"group {'AB'[g]}: {'no move up to a quarter turn reaches' if status == 1 else 'solver failure at'}" f" {targets[g] / R[g]:.3f} x R") if (rec["status"] == 0).all(): draw = np.random.default_rng((SEED, *key)) signs = [ct.paired_signs(along[g], draw) for g in range(2)] for g in range(2): move = rec["sign"][g] * rec["amp"][g] * along[g] dtheta[0, g], dtheta[1, g] = move, signs[g] * move d0 = ct.gap_of(us[0], us[1]) for b in range(2): turn = np.exp(1j * dtheta[b]) after = ct.state_numbers(u * turn, f)["R"] if b == 0: rec["push_ratio"] = after / targets else: rec["sham_ratio"] = after / R rec["sumsq"][b] = float(np.square(dtheta[b]).sum()) rec["gaprel"][b] = ct.rms(ct.gap_of(us[0] * turn, us[1] * turn) - d0) / ct.rms(d0) rec["same"] = bool(np.array_equal(np.abs(dtheta[0]), np.abs(dtheta[1]))) rec["match"] = abs(rec["sumsq"][1] - rec["sumsq"][0]) / rec["sumsq"][0] if rec["sumsq"][0] > 0 else 0.0 hit = bool(np.all(np.abs(rec["push_ratio"] - 1) <= SOLVED)) inband = bool(np.all((SHAM_BAND[0] <= rec["sham_ratio"]) & (rec["sham_ratio"] <= SHAM_BAND[1]))) kept = bool(np.all(rec["gaprel"] < KEPT)) if not hit: reasons.append("the push missed its target") if not inband: reasons.append("sham outside its band (" + ", ".join(f"{x:.4f}" for x in rec["sham_ratio"]) + " x R)") if not (rec["same"] and rec["match"] <= MATCH): reasons.append("the sham's sizes differ from the push's") if not kept: reasons.append("the change moved the gap") rec["built"] = hit and inband and rec["same"] and rec["match"] <= MATCH and kept rec["reason"] = "; ".join(reasons) return rec, dtheta # --------------------------------------------------------------------------------------------------------------------- # descriptions of a state (reported only) _INNER = [] def inner_half(): """The slower-turning half of the pulled oscillators (natural rates closest to the group's): 1,936 of 3,872.""" if not _INNER: omega = nt.natural_rates(N) w = np.abs(omega[np.abs(omega) <= nt.FAST]) _INNER.append(w <= np.median(w)) return _INNER[0] def describe(u, f): """One copy's state (u: (group, pulled), f: (group, free)): R_g, |z2| to |z4|, C, the A-B timing, and the coherence within each group's slower-turning and faster-turning halves of the pulled oscillators.""" out = ct.state_numbers(u, f) inner = inner_half() out["inner"], out["outer"] = np.abs(u[:, inner].mean(-1)), np.abs(u[:, ~inner].mean(-1)) return out def distance(u1, u2): """How far apart two states' pulled phases are, oscillator by oscillator, after lining up each group's overall timing: the rms of the wrapped differences (radians).""" d = np.angle(u1 * np.conj(u2)) d = ct.wrap(d - np.angle(np.exp(1j * d).sum(-1, keepdims=True))) return float(np.sqrt(np.mean(d * d))) # --------------------------------------------------------------------------------------------------------------------- # the tasks def base_task(task, deadline, parts): state, part = parts / f"{task['name']}.state.npz", parts / f"{task['name']}.npz" out = base_run(task["seed"], deadline, state, task["keep"]) if out is None: return dict(status="saved") arrays = dict(start=task["seed"] + 1, z=out["z"], kept_t=out["kept_t"], kept_u=out["kept_u"], kept_f=out["kept_f"], trace=out["trace"], spent=out["spent"], pieces=out["pieces"]) if task["logs"]: # the calibration keeps no gap growth arrays["logs"] = out["logs"] save_npz(part, **arrays) state.unlink(missing_ok=True) return dict(status="done", sha256=sha256(part)) def kept_state(base, t): times = [float(x) for x in base["kept_t"]] i = times.index(float(t)) return base["kept_u"][i], base["kept_f"][i] def calibration_moment(task, parts, window=WINDOW): """Every size, raised and lowered, for `window` units from T on the first copy only (no gap): each group's coherence up to R, and whether the matched push and its sham can be built at R. Nothing after R is run.""" s, m = task["seed"], task["moment"] T = task.get("T", moment_times(s)[m]) with np.load(parts / f"base_{s + 1:03d}.npz") as base: us_T, uf_T = kept_state(base, T) us_R, uf_R = kept_state(base, T + window) none = np.abs(base["z"][:, T:T + window + 1]) paths = np.zeros((len(SIZES), 2, 2, window + 1)) records = [] for a, f in enumerate(SIZES): for d, sign in enumerate((1, -1)): paths[a, d] = np.abs(branch(us_T, uf_T, [(K * (1 + sign * f), window)], 1)["z"]) records.append(matched(us_R, uf_R, paths[a, d, :, -1], (2, s, m, a, d))[0]) shape = (len(SIZES), 2) arrays = {f"rec_{key}": np.array([r[key] for r in records]).reshape(shape + np.shape(records[0][key])) for key in REC_KEYS} part = parts / f"{task['name']}.npz" save_npz(part, start=s + 1, moment=m + 1, T=T, none=none, R_none=none[:, -1], paths=paths, R_end=paths[..., -1], reason=np.array([r["reason"] for r in records]).reshape(shape), **arrays) return part, "" def reference_moment(task, parts): """The chosen size kept changed for 700 units from T, raised and lowered, first copy only: coherence only.""" s, m, f = task["seed"], task["moment"], task["size"] T = moment_times(s)[m] with np.load(parts / f"base_{s + 1:03d}.npz") as base: us_T, uf_T = kept_state(base, T) none = np.abs(base["z"][:, T:T + KEPT_FOR + 1]) paths = np.stack([np.abs(branch(us_T, uf_T, [(K * (1 + sign * f), KEPT_FOR)], 1)["z"]) for sign in (1, -1)]) part = parts / f"{task['name']}.npz" save_npz(part, start=s + 1, moment=m + 1, T=T, size=f, none=none, paths=paths) return part, "" def pilot_moment(task, parts, window=WINDOW, after=AFTER): """One moment of the pilot: the two medium arms (the pull raised or lowered by the frozen size for `window` units from T, restored exactly at R, followed `after` more), then from the untouched state at R the push matched to each and its sham, each followed `after` units; both copies throughout, the gap put back every 5 units.""" s, m, f = task["seed"], task["moment"], task["size"] T = task.get("T", moment_times(s)[m]) R = T + window with np.load(parts / f"base_{s + 1:03d}.npz") as base: us_T, uf_T = kept_state(base, T) us_R, uf_R = kept_state(base, R) z_none = base["z"][:, T:R + after + 1] logs_none = base["logs"][T // 5:(R + after) // 5] medium = [branch(us_T, uf_T, [(K * (1 + sign * f), window), (K, after)], 2) for sign in (1, -1)] records, changes, pushes, shams = [], [], [], [] for d in range(2): rec, dtheta = matched(us_R, uf_R, np.abs(medium[d]["z"][:, window]), (3, s, m, d)) records.append(rec) changes.append(dtheta) pushes.append(branch(us_R, uf_R, [(K, after)], 2, change=dtheta[0]) if rec["built"] else None) shams.append(branch(us_R, uf_R, [(K, after)], 2, change=dtheta[1]) if rec["built"] else None) states = [(us_R[0], uf_R[0])] + [(arm["ends"][0][0][0], arm["ends"][0][1][0]) for arm in medium] for b in range(2): # the pushes', then the shams', states at R states += [(us_R[0] * np.exp(1j * changes[d][b]), uf_R[0]) if records[d]["built"] else None for d in range(2)] desc = {key: np.full((len(STATES), 2), np.nan) for key in ("R", "inner", "outer", "z2", "z3", "z4")} desc.update(C=np.full(len(STATES), np.nan), ab=np.full(len(STATES), np.nan)) for i, state in enumerate(states): if state is not None: numbers = describe(*state) for key in DESC_KEYS: desc[key][i] = numbers[key] apart = [distance(states[1 + d][0], states[3 + d][0]) if records[d]["built"] else np.nan for d in range(2)] delivered = np.array([f * K * np.abs(arm["z"][:, 1:window + 1]).sum(-1) for arm in medium]) restored = np.array([bool(len(arm["pulls"]) == window + after and np.all(arm["pulls"][window:] == K) and np.all(arm["pulls"][:window] == K * (1 + sign * f))) for arm, sign in zip(medium, (1, -1))]) blank_z, blank_logs = np.full((2, after + 1), np.nan + 0j), np.full(after // 5, np.nan) arrays = {f"rec_{key}": np.array([r[key] for r in records]) for key in REC_KEYS} part = parts / f"{task['name']}.npz" save_npz(part, start=s + 1, moment=m + 1, T=T, R=R, size=f, z_none=z_none, logs_none=logs_none, z_medium=np.stack([arm["z"] for arm in medium]), logs_medium=np.stack([arm["logs"] for arm in medium]), z_push=np.stack([p["z"] if p else blank_z for p in pushes]), logs_push=np.stack([p["logs"] if p else blank_logs for p in pushes]), z_sham=np.stack([p["z"] if p else blank_z for p in shams]), logs_sham=np.stack([p["logs"] if p else blank_logs for p in shams]), built=np.array([r["built"] for r in records]), reason=np.array([r["reason"] for r in records]), distance=np.array(apart), delivered=delivered, restored=restored, **arrays, **{f"desc_{key}": value for key, value in desc.items()}) medium_finite = all(np.isfinite(arm["z"]).all() and np.isfinite(arm["logs"]).all() for arm in medium) problem = "" if restored.all() and medium_finite else f"restored {restored.tolist()}, finite {medium_finite}" return part, problem def work(task, deadline): parts = Path(task["parts"]) parts.mkdir(parents=True, exist_ok=True) if task["kind"] == "base": return base_task(task, deadline, parts) run = dict(cal=calibration_moment, ref=reference_moment, pilot=pilot_moment)[task["kind"]] part, problem = run(task, parts) return dict(status="done", sha256=sha256(part), error=problem) def job(args): return guard(work, args) def unfinished(tasks): return [t for t in tasks if not (Path(t["parts"]) / f"{t['name']}.npz").exists()] def paced_piece(log_path, tasks, start_by): """One foreground piece (prelock_common.run_piece), then the plan's and this file's fingerprints noted with it.""" run_piece(log_path, tasks, job, BUDGET, start_by) log = json.loads(log_path.read_text()) log["pieces"][-1].update(plan_sha256=sha256(PLAN), code_sha256=sha256(Path(__file__).resolve()), start_by=start_by) write_json(log_path, log) # --------------------------------------------------------------------------------------------------------------------- # the calibration def calibration_tasks(): base = [dict(name=f"base_{s + 1:03d}", kind="base", seed=s, keep=keep_times(s), logs=False, parts=str(CAL_PARTS)) for s in CAL_STARTS] moments = [dict(name=f"cal_{s + 1:03d}_{m + 1:02d}", kind="cal", seed=s, moment=m, parts=str(CAL_PARTS)) for s in CAL_STARTS for m in range(MOMENTS)] return base, moments def reference_tasks(size): return [dict(name=f"ref_{s + 1:03d}_{m + 1:02d}", kind="ref", seed=s, moment=m, size=size, parts=str(CAL_PARTS)) for s in CAL_STARTS for m in range(MOMENTS)] def load_parts(tasks): rows = [] for t in tasks: with np.load(Path(t["parts"]) / f"{t['name']}.npz") as p: rows.append({key: p[key] for key in p.files}) return rows def calibration_table(rows): """The rule of plan.md: for each size, the ratio of the average coherence at R (C, the two groups' average) over the 20 moments to the untouched run's, raised and lowered, and at how many moments the push and sham are built. The smallest size meeting every mark is chosen.""" C_none = np.array([r["R_none"].mean() for r in rows]) # (moment,) C = np.array([r["R_end"].mean(-1) for r in rows]) # (moment, size, direction) built = np.array([r["rec_built"] for r in rows]) # (moment, size, direction) ratio = C.mean(0) / C_none.mean() counts = built.sum(0) meets = dict(raised=(ratio[:, 0] >= MARKS["raise_low"]) & (ratio[:, 0] <= MARKS["raise_high"]), lowered=ratio[:, 1] <= MARKS["lower_high"], built=(counts >= MARKS["built"]).all(-1)) qualifies = meets["raised"] & meets["lowered"] & meets["built"] chosen = int(np.flatnonzero(qualifies)[0]) if qualifies.any() else None return dict(C_none=C_none, C=C, built=built, ratio=ratio, counts=counts, meets=meets, qualifies=qualifies, chosen=chosen) def table_lines(table): yes = {True: "yes", False: "no"} lines = [" size of K raised: C(R) / untouched lowered: C(R) / untouched pushes built (raised, lowered)" " marks met (raised, lowered, pushes)"] for a, f in enumerate(SIZES): lines.append(f" {100 * f:>8.2f}% {table['ratio'][a, 0]:>24.3f} {table['ratio'][a, 1]:>25.3f} " f"{table['counts'][a, 0]:>12d} and {table['counts'][a, 1]:<2d} of 20 " f"{yes[bool(table['meets']['raised'][a])]}, {yes[bool(table['meets']['lowered'][a])]}, " f"{yes[bool(table['meets']['built'][a])]}") chosen = table["chosen"] lines.append(f" chosen: {'none: no size meets every mark' if chosen is None else f'{100 * SIZES[chosen]:g}% of K'}") return lines def calibrate(start_by): if FREEZE.exists(): sys.exit("The calibration is frozen (medium_freeze.json): nothing more runs in it.") if pilot_outputs(): sys.exit("Pilot output exists: the calibration can't run once the pilot has begun.") base, moments = calibration_tasks() for stage in (base, moments): todo = unfinished(stage) if todo: paced_piece(CAL_LOG, todo, start_by) return table = calibration_table(load_parts(moments)) print("\n".join(table_lines(table))) if table["chosen"] is None: print("No size meets every mark. The calibration stops here (plan.md): report what each size did and propose a" " revised rule. Nothing more runs until Colson approves it.") return todo = unfinished(reference_tasks(SIZES[table["chosen"]])) if todo: paced_piece(CAL_LOG, todo, start_by) return print("The calibration is complete: run freeze.") def calibration_report(table, refs, log): a = table["chosen"] f = SIZES[a] rows_C, C_none = table["C"], table["C_none"] per_moment = rows_C[:, a, :] / C_none[:, None] # (moment, direction), shown only recs = load_parts(calibration_tasks()[1]) amp = np.array([r["rec_amp"][a] for r in recs]) # (moment, direction, group) sham = np.array([r["rec_sham_ratio"][a] for r in recs]) built = np.array([r["rec_built"][a] for r in recs]) reasons = [(f"start {int(r['start'])}, moment {int(r['moment'])}, {DIRECTIONS[d]}: {r['reason'][a, d]}") for r in recs for d in range(2) if not r["rec_built"][a, d]] window_C = np.array([r["paths"][a].mean(-2) for r in recs]) # (moment, direction, time) none_C = np.array([r["none"].mean(-2) for r in recs]) # (moment, time) ticks = list(range(0, WINDOW + 1, 50)) ref_C = np.array([r["paths"].mean(-2) for r in refs]) # (moment, direction, time) ref_none = np.array([r["none"].mean(-2) for r in refs]) blocks = [(b, b + BLOCK) for b in range(0, KEPT_FOR, BLOCK)] def block_mean(x, lo, hi): return float(x[..., lo + 1:hi + 1].mean()) seconds = sum(p.get("seconds", 0) for p in log["pieces"]) versions = dict(plan=sorted({str(p.get("plan_sha256")) for p in log["pieces"]}), code=sorted({str(p.get("code_sha256")) for p in log["pieces"]})) lines = [ f"Test 1, change the medium, not the state: calibration ({shown(stamp())}).", f"Plan: plan.md, revision 2 (SHA-256 {sha256(PLAN)}).", "Old starts 16 and 17 (seeds 15 and 16), ten moments each: 20 moments. These starts stay out of the pilot's" " results.", "Coherence only, up to the restoration point R = T + 200, and whether the matched pushes and shams can be built." " No gap growth is saved or looked at, and no restored arm is run after R.", "C is the average of the two groups' coherence. The marks use the ratio of the averages over the 20 moments.", f"Marks: raising {MARKS['raise_low']} to {MARKS['raise_high']} times the untouched run's C at R; lowering at" f" most {MARKS['lower_high']} times; the push and its sham built at {MARKS['built']} or more of the 20 moments," " raised and lowered separately. The smallest size meeting every mark is used unchanged in the pilot.", "", *table_lines(table), "", f"At the chosen size ({100 * f:g}% of K, so K = {K * (1 + f):.6g} raised and {K * (1 - f):.6g} lowered for 200" " units):", f" C(R) / untouched at single moments (shown only; the mark uses the ratio of averages): raised" f" {per_moment[:, 0].min():.2f} to {per_moment[:, 0].max():.2f} (median {np.median(per_moment[:, 0]):.2f});" f" lowered {per_moment[:, 1].min():.2f} to {per_moment[:, 1].max():.2f}" f" (median {np.median(per_moment[:, 1]):.2f})", f" matched pushes: amplitude (quarter turn = {np.pi / 2:.3f}) raised median {np.nanmedian(amp[:, 0]):.4f}," f" largest {np.nanmax(amp[:, 0]):.4f}; lowered median {np.nanmedian(amp[:, 1]):.4f}, largest" f" {np.nanmax(amp[:, 1]):.4f}", f" shams: each group's coherence {np.nanmin(sham):.4f} to {np.nanmax(sham):.4f} times where it was", f" built: raised {int(built[:, 0].sum())} of 20, lowered {int(built[:, 1].sum())} of 20", *(f" not built: {x}" for x in reasons), " average C during the window (20 moments), at T + " + ", ".join(str(t) for t in ticks) + ":", " raised " + " ".join(f"{window_C[:, 0, t].mean():.4f}" for t in ticks), " lowered " + " ".join(f"{window_C[:, 1, t].mean():.4f}" for t in ticks), " untouched " + " ".join(f"{none_C[:, t].mean():.4f}" for t in ticks), "", "Kept-changed reference (the chosen size never restored, 700 units from T; coherence only). Reported here only;" " it isn't used in any pilot reading.", " average C in 100-unit blocks after T: " + ", ".join(f"({lo}, {hi}]" for lo, hi in blocks), " raised " + " ".join(f"{block_mean(ref_C[:, 0], lo, hi):.4f}" for lo, hi in blocks), " lowered " + " ".join(f"{block_mean(ref_C[:, 1], lo, hi):.4f}" for lo, hi in blocks), " untouched " + " ".join(f"{block_mean(ref_none, lo, hi):.4f}" for lo, hi in blocks), "", "Moments: " + "; ".join(f"start {s + 1}: T = " + ", ".join(str(t) for t in moment_times(s)) for s in CAL_STARTS), f"Compute: {len(log['pieces'])} piece{'s' * (len(log['pieces']) != 1)}, {seconds / 60:.1f} minutes of" " foreground time (7 workers at a time). Machine:" f" {log.get('machine', {}).get('cpu', 'unrecorded')}.", f"The cause test's locked files unchanged: {not lock_changes()}.", "Every calibration piece ran under this plan and this code: " + ("yes" if versions == dict(plan=[sha256(PLAN)], code=[sha256(Path(__file__).resolve())]) else f"NO (plans {versions['plan']}, code {versions['code']})"), ] result = dict(made=stamp(), sizes=list(SIZES), ratio=table["ratio"], counts=table["counts"], meets=table["meets"], qualifies=table["qualifies"], chosen=dict(index=a, size=f), C_none=C_none, C=rows_C, per_moment_at_chosen=per_moment, amp_at_chosen=amp, sham_at_chosen=sham, built_at_chosen=built, not_built=reasons, window_ticks=ticks, window_C=dict(raised=window_C[:, 0, ticks].mean(0), lowered=window_C[:, 1, ticks].mean(0), untouched=none_C[:, ticks].mean(0)), reference=dict(blocks=blocks, raised=[block_mean(ref_C[:, 0], lo, hi) for lo, hi in blocks], lowered=[block_mean(ref_C[:, 1], lo, hi) for lo, hi in blocks], untouched=[block_mean(ref_none, lo, hi) for lo, hi in blocks]), moments={str(s + 1): moment_times(s) for s in CAL_STARTS}, pieces=len(log["pieces"]), seconds=seconds, machine=log.get("machine"), plan_sha256=sha256(PLAN), pieces_ran_under=versions) return result, lines def freeze(): if FREEZE.exists(): made = json.loads(FREEZE.read_text())["made"] sys.exit(f"Already frozen ({shown(made)}). A new freeze needs Colson's approval and is labelled as a change.") found = pilot_outputs() if found: sys.exit(f"Pilot output already exists ({found[0]} ...): the freeze is refused (plan.md).") base, moments = calibration_tasks() if unfinished(base) or unfinished(moments): sys.exit("The calibration is incomplete: the freeze is refused.") table = calibration_table(load_parts(moments)) if table["chosen"] is None: sys.exit("No size qualified: the freeze is refused. The calibration's numbers and a proposal go to Colson.") size = SIZES[table["chosen"]] refs = reference_tasks(size) if unfinished(refs): sys.exit("The kept-changed reference is incomplete: the freeze is refused.") if lock_changes(): sys.exit("A file fixed by the cause test's lock has changed. Stopping.") log = json.loads(CAL_LOG.read_text()) result, lines = calibration_report(table, load_parts(refs), log) write_json(CAL_JSON, result) CAL_TXT.write_text("\n".join(lines) + "\n") parts = sorted(CAL_PARTS.glob("*.npz")) record = dict( made=stamp(), what=("The calibration's choice for Test 1, frozen before any pilot output exists (plan.md, revision 2). The" " pilot takes the size only from here, and checks every fingerprint below before each piece and before" " its report. Nothing here changes afterwards; any change would be a new, dated freeze labelled as a" " change, with this record kept."), chosen=dict(size=size, index=table["chosen"], raised_k=K * (1 + size), lowered_k=K * (1 - size)), table=dict(sizes=list(SIZES), ratio=table["ratio"], counts=table["counts"], meets=table["meets"], qualifies=table["qualifies"]), sha256=current_fingerprints(), settings=SETTINGS, pilot=dict(starts=[s + 1 for s in PILOT_STARTS], moments={str(s + 1): moment_times(s) for s in PILOT_STARTS}), calibration_parts={p.name: sha256(p) for p in parts}, machine=machine(), protected=protected_fingerprints(), cause_lock_unchanged=True, pilot_outputs_at_freeze=found) write_json(FREEZE, record) print("\n".join(lines)) print(f"\nFrozen {shown(record['made'])}: size {100 * size:g}% of K; fingerprints of the plan, the code and the" " calibration's results kept in medium_freeze.json. The pilot waits for Colson's go.") # --------------------------------------------------------------------------------------------------------------------- # the pilot def frozen(): """The freeze record, after checking that the plan, the code, the calibration's results and the settings are exactly as frozen. The pilot reads nothing else from the calibration.""" if not FREEZE.exists(): sys.exit("No freeze record: the pilot can't run before the calibration is frozen.") record = json.loads(FREEZE.read_text()) now = current_fingerprints() changed = [name for name in now if record["sha256"].get(name) != now[name]] if record["settings"] != json.loads(json.dumps(plain(SETTINGS))): changed.append("settings") if changed: sys.exit(f"Changed since the freeze ({shown(record['made'])}): {', '.join(changed)}. Stopping.") return record def pilot_tasks(size): base = [dict(name=f"base_{s + 1:03d}", kind="base", seed=s, keep=keep_times(s), logs=True, parts=str(PILOT_PARTS)) for s in PILOT_STARTS] moments = [dict(name=f"moment_{s + 1:03d}_{m + 1:02d}", kind="pilot", seed=s, moment=m, size=size, parts=str(PILOT_PARTS)) for s in PILOT_STARTS for m in range(MOMENTS)] return base, moments def pilot(start_by): record = frozen() if protected_fingerprints() != record["protected"]: sys.exit("A protected file has changed since the freeze. Stopping.") base, moments = pilot_tasks(record["chosen"]["size"]) todo = unfinished(base) ready = [t for t in unfinished(moments) if (PILOT_PARTS / f"base_{t['seed'] + 1:03d}.npz").exists()] if not todo and not ready: print("All 132 tasks are done: run report.") return paced_piece(PILOT_LOG, todo + ready, start_by) def by_start(x): """Each start's average over its moments (moments without a value left out).""" with warnings.catch_warnings(): warnings.simplefilter("ignore", category=RuntimeWarning) return np.nanmean(np.asarray(x, float), axis=1) def summary(per_start): """The average over the starts, with its error bar: the larger of the resampled one (starts drawn with replacement, 400 times) and SD / sqrt(starts).""" v = np.asarray(per_start, float) v = v[np.isfinite(v)] n = len(v) if n < 2: return dict(mean=float(v.mean()) if n else np.nan, bar=np.nan, n=n, above=int((v > 0).sum()), clear=False) rng = np.random.default_rng(RESAMPLE_SEED) boot = v[rng.integers(0, n, (RESAMPLES, n))].mean(-1) mean = float(v.mean()) bar = max(float(boot.std(ddof=1)), float(v.std(ddof=1) / np.sqrt(n))) return dict(mean=mean, bar=bar, n=n, above=int((v > 0).sum()), clear=bool(abs(mean) > CLEAR * bar)) def ratio_summary(numerator, denominator): """A ratio of averages over the starts, with the resampled error bar.""" a, b = np.asarray(numerator, float), np.asarray(denominator, float) keep = np.isfinite(a) & np.isfinite(b) a, b = a[keep], b[keep] n = len(a) if n < 2: return dict(mean=np.nan, bar=np.nan, n=n, clear=False) idx = np.random.default_rng(RESAMPLE_SEED).integers(0, n, (RESAMPLES, n)) boot = a[idx].mean(-1) / b[idx].mean(-1) mean, bar = float(a.mean() / b.mean()), float(boot.std(ddof=1)) return dict(mean=mean, bar=bar, n=n, clear=bool(abs(mean) > CLEAR * bar)) def shape_of(S, P): if S["clear"] and abs(S["mean"]) > abs(P["mean"]): return "same direction" if S["mean"] > 0 else "opposite direction" if P["clear"] and abs(P["mean"]) > abs(S["mean"]): return "non-monotonic" if not S["clear"] and not P["clear"]: return "no clear effect" return "mixed" def ruled_out(S, half): return bool(np.isfinite(S["bar"]) and S["mean"] - SIDE * S["bar"] > -half and S["mean"] + SIDE * S["bar"] < half) def load_pilot(parts): """Every pilot moment, from the pilot's own starts only.""" rows = [] for s in PILOT_STARTS: for m in range(MOMENTS): path = parts / f"moment_{s + 1:03d}_{m + 1:02d}.npz" if not path.exists(): sys.exit(f"{path.name} is missing: the report needs all 120 moments.") with np.load(path) as p: rows.append({key: p[key] for key in p.files}) return rows def analyse(rows, size): """The readings fixed in plan.md (revision 2), from the pilot's starts 103 to 114 only.""" S_, M_ = len(PILOT_STARTS), MOMENTS if len(rows) != S_ * M_ or sorted({int(r["start"]) for r in rows}) != [s + 1 for s in PILOT_STARTS]: sys.exit("The report reads only the pilot's own starts, 103 to 114, with all their moments. Stopping.") w, (lo, hi) = WINDOW, LATER def grid(fn): return np.array([[fn(rows[i * M_ + m]) for m in range(M_)] for i in range(S_)], float) def C(z): return np.abs(z).mean(-2) built = np.array([[rows[i * M_ + m]["built"] for m in range(M_)] for i in range(S_)], bool) # (start, moment, d) both = built.all(-1) mask = np.where(both, 1.0, np.nan) missing = int((~both).sum()) push_ok = missing <= MOST_MISSING def growths(lo_, hi_, origin_shift): """Growth over (R + lo_, R + hi_] for every arm; medium and untouched logs count from T, the others from R.""" a, b = origin_shift + lo_, origin_shift + hi_ return dict( none=grid(lambda r: growth(r["logs_none"], a, b)), raised=grid(lambda r: growth(r["logs_medium"][0], a, b)), lowered=grid(lambda r: growth(r["logs_medium"][1], a, b)), push_raised=grid(lambda r: growth(r["logs_push"][0], lo_, hi_)), push_lowered=grid(lambda r: growth(r["logs_push"][1], lo_, hi_)), sham_raised=grid(lambda r: growth(r["logs_sham"][0], lo_, hi_)), sham_lowered=grid(lambda r: growth(r["logs_sham"][1], lo_, hi_))) dC = grid(lambda r: (C(r["z_medium"][0])[w] - C(r["z_medium"][1])[w]) / 2) dC_mean = float(np.mean(by_start(dC))) natural = 1.0 if dC_mean >= 0 else -1.0 # S > 0: the arm higher at R grew more slowly later benchmark = BENCHMARK * abs(dC_mean) half = HALF * benchmark def parts_of(G): """S and P per moment for the medium arms (against untouched), the pushes against their shams, the pushes and shams against untouched, and medium minus push like for like (moments with both pushes built).""" Dr, Dl = G["raised"] - G["none"], G["lowered"] - G["none"] out = dict(S_medium=natural * (Dl - Dr) / 2, P_medium=(Dr + Dl) / 2) Pr, Pl = (G["push_raised"] - G["sham_raised"]) * mask, (G["push_lowered"] - G["sham_lowered"]) * mask out.update(S_push_sham=natural * (Pl - Pr) / 2, P_push_sham=(Pr + Pl) / 2) Ur, Ul = (G["push_raised"] - G["none"]) * mask, (G["push_lowered"] - G["none"]) * mask out.update(S_push=natural * (Ul - Ur) / 2, P_push=(Ur + Ul) / 2) Hr, Hl = (G["sham_raised"] - G["none"]) * mask, (G["sham_lowered"] - G["none"]) * mask out.update(sham_raised=Hr, sham_lowered=Hl, sham_both=(Hr + Hl) / 2, D_raised=Dr, D_lowered=Dl) out.update(S_difference=out["S_medium"] * mask - out["S_push"], P_difference=out["P_medium"] * mask - out["P_push"]) return {key: summary(by_start(value)) for key, value in out.items()} later = parts_of(growths(lo, hi, w)) windows = {f"({lo_}, {hi_}]": parts_of(growths(lo_, hi_, w)) for lo_, hi_ in [EARLY] + [(b, b + BLOCK) for b in range(0, AFTER, BLOCK)]} during = dict(S_medium=summary(by_start(natural * (grid(lambda r: growth(r["logs_medium"][1], 0, w)) - grid(lambda r: growth(r["logs_medium"][0], 0, w))) / 2)), P_medium=summary(by_start((grid(lambda r: growth(r["logs_medium"][0], 0, w)) + grid(lambda r: growth(r["logs_medium"][1], 0, w))) / 2 - grid(lambda r: growth(r["logs_none"], 0, w))))) diffR = grid(lambda r: C(r["z_medium"][0])[w] - C(r["z_medium"][1])[w]) wm = grid(lambda r: (C(r["z_medium"][0]) - C(r["z_medium"][1]))[w + lo + 1:w + hi + 1].mean()) wp = grid(lambda r: (C(r["z_push"][0]) - C(r["z_push"][1]))[lo + 1:hi + 1].mean()) persistence = dict( medium=ratio_summary(by_start(wm * mask), by_start(diffR * mask)), push=ratio_summary(by_start(wp * mask), by_start(diffR * mask)), difference=ratio_summary(by_start((wm - wp) * mask), by_start(diffR * mask)), medium_all_moments=ratio_summary(by_start(wm), by_start(diffR))) sgn = 1.0 if np.nanmean(by_start(diffR * mask)) >= 0 else -1.0 def block(z, a, b): return (C(z[0]) - C(z[1]))[a + 1:b + 1].mean() last_m = grid(lambda r: block(r["z_medium"], w + 400, w + 500)) * sgn * mask decline_m = grid(lambda r: block(r["z_medium"], w + 200, w + 300) - block(r["z_medium"], w + 400, w + 500)) \ * sgn * mask decline_p = grid(lambda r: block(r["z_push"], 200, 300) - block(r["z_push"], 400, 500)) * sgn * mask lasting = dict(medium_last=summary(by_start(last_m)), medium_decline=summary(by_start(decline_m)), push_decline=summary(by_start(decline_p))) lasting["switch"] = bool(lasting["medium_last"]["clear"] and lasting["medium_last"]["mean"] > 0 and not (lasting["medium_decline"]["clear"] and lasting["medium_decline"]["mean"] > 0) and lasting["push_decline"]["clear"] and lasting["push_decline"]["mean"] > 0) curves = {} for d, name in enumerate(DIRECTIONS): curves[f"medium {name}"] = [float(np.mean(by_start(grid( lambda r, b=b, d=d: (C(r["z_medium"][d]) - C(r["z_none"]))[b + 1:b + BLOCK + 1].mean())))) for b in range(0, w + AFTER, BLOCK)] for arm in ("push", "sham"): curves[f"{arm} {name}"] = [float(np.mean(by_start(grid( lambda r, b=b, d=d, arm=arm: (C(r[f"z_{arm}"][d]) - C(r["z_none"])[w:])[b + 1:b + BLOCK + 1].mean()) * mask))) for b in range(0, AFTER, BLOCK)] shapes = dict(medium=shape_of(later["S_medium"], later["P_medium"]), push_sham=shape_of(later["S_push_sham"], later["P_push_sham"]), push_untouched=shape_of(later["S_push"], later["P_push"])) out_ruled = ruled_out(later["S_medium"], half) differs = bool(persistence["difference"]["clear"] or later["S_difference"]["clear"] or later["P_difference"]["clear"]) readings = dict( medium_shape=shapes["medium"], medium_ruled_out=out_ruled, medium_unresolved=bool(not later["S_medium"]["clear"] and not out_ruled), push_shape_vs_sham=shapes["push_sham"], push_shape_vs_untouched=shapes["push_untouched"], push_inconclusive=not push_ok, differs=bool(push_ok and differs), like_push=bool(push_ok and not differs and out_ruled), lasting_switch=bool(push_ok and lasting["switch"])) patterns = [key for key, holds in ( ("differs", readings["differs"]), ("like push", readings["like_push"]), ("non-monotonic", shapes["medium"] == "non-monotonic"), ("opposite direction", shapes["medium"] == "opposite direction"), ("same direction", shapes["medium"] == "same direction"), ("lasting switch", readings["lasting_switch"]), ("ruled out", out_ruled), ("unresolved", readings["medium_unresolved"]), ("push inconclusive", readings["push_inconclusive"])) if holds] descriptions = dict( distance={name: summary(by_start(grid(lambda r, d=d: r["distance"][d]) * mask)) for d, name in enumerate(DIRECTIONS)}, delivered={name: float(np.mean(by_start(grid(lambda r, d=d: r["delivered"][d].mean())))) for d, name in enumerate(DIRECTIONS)}, states={state: {key: float(np.mean(by_start(grid( lambda r, i=i, key=key: np.mean(np.abs(r[f"desc_{key}"][i]))) * (mask if i >= 3 else 1.0)))) for key in DESC_KEYS} for i, state in enumerate(STATES)}) checks = dict( restored=bool(all(r["restored"].all() for r in rows)), medium_finite=bool(all(np.isfinite(r["z_medium"]).all() and np.isfinite(r["logs_medium"]).all() for r in rows)), built={name: int(built[..., d].sum()) for d, name in enumerate(DIRECTIONS)}, both=int(both.sum()), missing=missing, push_ok=push_ok, built_checks=bool(all(np.all(np.abs(r["rec_push_ratio"][d] - 1) <= SOLVED) and np.all((SHAM_BAND[0] <= r["rec_sham_ratio"][d]) & (r["rec_sham_ratio"][d] <= SHAM_BAND[1])) and r["rec_same"][d] and r["rec_match"][d] <= MATCH and np.all(r["rec_gaprel"][d] < KEPT) for r in rows for d in range(2) if r["built"][d])), not_built=[f"start {int(r['start'])}, moment {int(r['moment'])}, {DIRECTIONS[d]}: {r['reason'][d]}" for r in rows for d in range(2) if not r["built"][d]]) return dict(size=size, natural=natural, dC=dC_mean, benchmark=benchmark, half=half, later=later, windows=windows, during=during, persistence=persistence, lasting=lasting, curves=curves, shapes=shapes, readings=readings, patterns=patterns, descriptions=descriptions, checks=checks) def report_lines(res, record, title): def est(x, digits=5): if not np.isfinite(x["mean"]): return "n/a" if not np.isfinite(x["bar"]): return f"{x['mean']:+.{digits}f}" out = f"{x['mean']:+.{digits}f} +- {x['bar']:.{digits}f} ({x['mean'] / x['bar']:+.1f} bars" return out + (f"; {x['above']} of {x['n']} starts above zero)" if "above" in x else ")") yes = {True: "yes", False: "no"} later, ch = res["later"], res["checks"] lines = [ title, f"Frozen {shown(record['made'])} (medium_freeze.json): the pull within each group changed by" f" {100 * res['size']:g}% of K for 200 time units, raised and lowered, then set back exactly to K = {K}.", "Calibration data used: only the chosen size, read from the freeze record. Every number below comes from the" " pilot's own starts, 103 to 114 (12 starts, 10 moments each); each start's moments are averaged first, then" " the 12 starts. Error bars: the larger of the resampled one and SD/sqrt(12). Clear: more than 3 bars.", "A toy model: structure inside the model, not evidence about cells. Every medium-vs-push difference below is a" " difference between two intervention histories and the states they leave, not evidence of an independent" " memory variable.", "", "Checks", f" medium arms at all 120 moments, every number present: {yes[ch['medium_finite']]}; the pull back at exactly K" f" after R in every medium arm: {yes[ch['restored']]}", f" matched push and sham built: raised {ch['built']['raised']} of 120, lowered {ch['built']['lowered']} of" f" 120, both {ch['both']} ({ch['missing']} moments without; at most {MOST_MISSING} allowed): medium-vs-push" f" readings {'conclusive' if ch['push_ok'] else 'INCONCLUSIVE'}", f" every built push on target, its sham in band, the same sizes, the gap kept: {yes[ch['built_checks']]}", *(f" not built: {x}" for x in ch["not_built"]), "", f"Benchmark: 0.1602 x half the raised-lowered coherence difference at R ({res['dC']:+.5f}) = {res['benchmark']:.5f};" f" half the benchmark = {res['half']:.5f}. A yardstick for size, in either direction.", "", "Later growth over (R + 200, R + 500]; S > 0 is the natural correlation's direction (the arm with higher" " coherence at R grew more slowly)", f" medium against untouched: S = {est(later['S_medium'])}; P = {est(later['P_medium'])}", f" D_raised = {est(later['D_raised'])}; D_lowered = {est(later['D_lowered'])}", f" shape: {res['shapes']['medium']}; S ruled out at half the benchmark: {yes[res['readings']['medium_ruled_out']]}", f" pushes against their shams: S = {est(later['S_push_sham'])}; P = {est(later['P_push_sham'])};" f" shape: {res['shapes']['push_sham']}", f" pushes against untouched: S = {est(later['S_push'])}; P = {est(later['P_push'])};" f" shape: {res['shapes']['push_untouched']}", f" shams against untouched (pushing itself): raised {est(later['sham_raised'])}; lowered" f" {est(later['sham_lowered'])}; average {est(later['sham_both'])}", "", "Medium against matched push, like for like (both against untouched, moments with both pushes built)", f" later growth: S difference {est(later['S_difference'])}; P difference {est(later['P_difference'])}", f" persistence, the share of the raised-lowered coherence difference at R still there over (R + 200, R + 500]:" f" medium {est(res['persistence']['medium'], 3)}; push {est(res['persistence']['push'], 3)}; difference" f" {est(res['persistence']['difference'], 3)}; medium at all moments {est(res['persistence']['medium_all_moments'], 3)}", f" lasting switch: the medium's difference over (R + 400, R + 500] {est(res['lasting']['medium_last'])};" f" its decline from (R + 200, R + 300] {est(res['lasting']['medium_decline'])}; the push's decline" f" {est(res['lasting']['push_decline'])}: {yes[res['readings']['lasting_switch']]}", "", "What holds (the readings fixed in plan.md):", *(f" {key}: {MEANING[key]}" for key in res["patterns"]), *([" none of the named patterns: reported as it stands."] if not res["patterns"] else []), "", "Also reported (descriptions; not used for any reading)", f" growth during the window (T, R], medium: S = {est(res['during']['S_medium'])};" f" P = {est(res['during']['P_medium'])}", " growth after R by window: S and P for the medium, S for the pushes against their shams and against untouched", *(f" {name}: S {est(v['S_medium'])}; P {est(v['P_medium'])}; push-sham S {est(v['S_push_sham'])};" f" push S {est(v['S_push'])}" for name, v in res["windows"].items()), " coherence difference from untouched, 100-unit blocks (medium from T, the others from R):", *(f" {name}: " + " ".join(f"{x:+.4f}" for x in values) for name, values in res["curves"].items()), " delivered extra pull through the window (|change of K| x the group's coherence, summed over time units):" f" raised {res['descriptions']['delivered']['raised']:.4f}; lowered {res['descriptions']['delivered']['lowered']:.4f}", " states at R: medium against its matched push, rms phase distance after lining up the timing: raised" f" {est(res['descriptions']['distance']['raised'], 3)}; lowered {est(res['descriptions']['distance']['lowered'], 3)}", " states at R, averages over the groups (R, slower half, faster half, |z2|, |z3|, |z4|, |A-B timing|):", *(f" {state}: " + " ".join(f"{v[key]:.4f}" for key in ("R", "inner", "outer", "z2", "z3", "z4", "ab")) for state, v in res["descriptions"]["states"].items()), ] return lines def report(): record = frozen() if protected_fingerprints() != record["protected"]: sys.exit("A protected file has changed since the freeze. Stopping.") if lock_changes(): sys.exit("A file fixed by the cause test's lock has changed. Stopping.") rows = load_pilot(PILOT_PARTS) res = analyse(rows, record["chosen"]["size"]) log = json.loads(PILOT_LOG.read_text()) lines = report_lines(res, record, f"Test 1, change the medium, not the state: pilot ({shown(stamp())}).") seconds = sum(p.get("seconds", 0) for p in log["pieces"]) lines += ["", f"Compute: {len(log['pieces'])} pieces, {seconds / 60:.1f} minutes of foreground time.", f"Plan, code and calibration results as frozen; protected folders and files unchanged since the freeze;" f" the cause test's locked files unchanged: yes."] res.update(made=stamp(), freeze=dict(made=record["made"], sha256=record["sha256"]), moment_sha256={p.name: sha256(p) for p in sorted(PILOT_PARTS.glob("moment_*.npz"))}) write_json(PILOT_JSON, res) PILOT_TXT.write_text("\n".join(lines) + "\n") print("\n".join(lines)) # --------------------------------------------------------------------------------------------------------------------- # tries before the calibration: the report on made-up numbers, and the code checks def made_up_rows(size=0.0025): """Numbers with the pilot's shapes, drawn at random (no model run): for trying the report only.""" rng = np.random.default_rng(7) rows = [] for s in PILOT_STARTS: for m, T in enumerate(moment_times(s)): def z(shape): return 0.026 * (1 + 0.3 * rng.standard_normal(shape)) * np.exp(1j * rng.uniform(-np.pi, np.pi, shape)) medium = z((2, 2, WINDOW + AFTER + 1)) push, sham = z((2, 2, AFTER + 1)), z((2, 2, AFTER + 1)) built = rng.random(2) > 0.04 push[:, :, 0] = np.abs(medium[:, :, WINDOW]) * np.exp(1j * np.angle(push[:, :, 0])) for d in range(2): if not built[d]: push[d], sham[d] = np.nan, np.nan logs = lambda shape: 0.07 + 0.02 * rng.standard_normal(shape) # noqa: E731 row = dict(start=s + 1, moment=m + 1, T=T, R=T + WINDOW, size=size, z_none=z((2, WINDOW + AFTER + 1)), logs_none=logs((WINDOW + AFTER) // 5), z_medium=medium, logs_medium=logs((2, (WINDOW + AFTER) // 5)), z_push=push, logs_push=np.where(built[:, None], logs((2, AFTER // 5)), np.nan), z_sham=sham, logs_sham=np.where(built[:, None], logs((2, AFTER // 5)), np.nan), built=built, reason=np.array(["" if b else "sham outside its band (made up)" for b in built]), distance=np.where(built, rng.uniform(0.1, 1.0, 2), np.nan), delivered=rng.uniform(0, 1, (2, 2)), restored=np.ones(2, bool), rec_push_ratio=np.ones((2, 2)), rec_sham_ratio=np.ones((2, 2)), rec_same=np.ones(2, bool), rec_match=np.zeros(2), rec_gaprel=np.zeros((2, 2))) for key in DESC_KEYS: row[f"desc_{key}"] = rng.uniform(0, 0.1, (len(STATES), 2) if key not in ("C", "ab") else len(STATES)) rows.append(row) return rows def dryreport(): if DRY.exists(): shutil.rmtree(DRY) DRY.mkdir(parents=True) record = dict(made=stamp()) res = analyse(made_up_rows(), 0.0025) lines = report_lines(res, record, "DRY RUN ON MADE-UP NUMBERS (no model run; nothing here is a result).") (DRY / "dryreport.txt").write_text("\n".join(lines) + "\n") write_json(DRY / "dryreport.json", res) print(f"The report ran on made-up numbers: {len(lines)} lines, written to {DRY} only. Patterns found in the noise:" f" {res['patterns'] or 'none'} (meaningless; random numbers).") def checks(): if DRY.exists(): shutil.rmtree(DRY) DRY.mkdir(parents=True) s, run, T0, short = CHECK_START, 400, 100, 20 keep = (T0, T0 + short) lines = [f"Code checks, {shown(stamp())}. Old start 1 (seed 0), short runs. Equality and mechanics only: no outcome" " of any change is read, shown or kept (plan.md, revision 2).", ""] passed = [] now = json.loads(json.dumps(plain(SETTINGS))) wrong = [key for key, value in json.loads(json.dumps(plain(PLAN_NUMBERS))).items() if now.get(key) != value] passed.append(not wrong) lines.append(f"0. The code's settings match the plan's numbers: {'yes' if not wrong else 'NO: ' + ', '.join(wrong)}") clock = time.time() reference = ct.base_run(s, keep=keep, run=run) reference_seconds = time.time() - clock mine = base_run(s, np.inf, DRY / "never.state.npz", keep, run) same = dict(coherence=np.array_equal(reference["z"][0], mine["z"]), gap=np.array_equal(reference["logs"][0], mine["logs"]), settling=np.array_equal(reference["trace"][0], mine["trace"]), states=len(reference["states"]) == len(mine["kept_u"]) and all( np.array_equal(a, b) and np.array_equal(c, d) for (a, c), b, d in zip(reference["states"], mine["kept_u"], mine["kept_f"]))) passed.append(all(same.values())) lines.append(f"1. This file's base run against the cause test's unchanged base run, both unbroken, {run} units:" f" identical to the last digit: {same}") stops, out = 0, None while out is None: out = base_run(s, time.time() + 2.0, DRY / "pieces.state.npz", keep, run) stops += out is None resumed = {key: bool(np.array_equal(mine[key], out[key])) for key in ("z", "logs", "trace", "kept_u", "kept_f")} passed.append(all(resumed.values()) and stops > 0) lines.append(f"2. The same run stopped every 2 seconds ({stops} times) and gone on from its kept state: identical" f" to one go: {resumed}") us, uf = kept_state(dict(kept_t=mine["kept_t"], kept_u=mine["kept_u"], kept_f=mine["kept_f"]), T0) b = branch(us, uf, [(K, 50), (K, 50)], 2) zero = dict(coherence=bool(np.array_equal(b["z"], mine["z"][:, T0:T0 + 101])), gap=bool(np.array_equal(b["logs"], mine["logs"][T0 // 5:T0 // 5 + 20]))) passed.append(all(zero.values())) lines.append(f"3. A branch from t = {T0} with the pull changed by zero, in two 50-unit stretches as a medium arm's" f" change and restoration: reproduces the untouched run to the last digit: {zero}") c = branch(us, uf, [(K * 1.05, 50), (K, 50)], 2) d = branch(c["ends"][0][0], c["ends"][0][1], [(K, 50)], 2) restore = dict(restored_is_unchanged_model=bool(np.array_equal(c["z"][:, 50:], d["z"]) and np.array_equal(c["logs"][10:], d["logs"])), change_took_effect=not np.array_equal(c["z"][:, 1:51], mine["z"][:, T0 + 1:T0 + 51]), pull_as_scheduled=bool(np.all(c["pulls"][:50] == K * 1.05) and np.all(c["pulls"][50:] == K))) passed.append(all(restore.values())) lines.append(f"4. The pull raised 5% for 50 units, then restored: the restored stretch is the unchanged model run" f" from the same state, to the last digit; the change took effect; the pull in use is as scheduled:" f" {restore}") usR, ufR = kept_state(dict(kept_t=mine["kept_t"], kept_u=mine["kept_u"], kept_f=mine["kept_f"]), T0 + short) R0 = np.abs(coherence(usR[:1], ufR[:1])[0]) for name, factor, code in (("raising", MARKS["raise_low"], 0), ("lowering", MARKS["lower_high"], 1)): rec, dtheta = matched(usR, ufR, factor * R0, (9, s, 0, code)) mech = dict(built=rec["built"], on_target=bool(np.all(np.abs(rec["push_ratio"] - 1) <= SOLVED)), sham_in_band=bool(np.all((SHAM_BAND[0] <= rec["sham_ratio"]) & (rec["sham_ratio"] <= SHAM_BAND[1]))), same_sizes=bool(rec["same"] and rec["match"] <= MATCH), gap_kept=bool(np.all(rec["gaprel"] < KEPT))) passed.append(all(mech.values())) lines.append(f"5{'ab'[code]}. A push {name} each group's coherence to {factor} times, and its sham, on the state" f" at t = {T0 + short}: {mech}{' (' + rec['reason'] + ')' if rec['reason'] else ''}") save_npz(DRY / "base_001.npz", start=1, z=mine["z"], logs=mine["logs"], kept_t=mine["kept_t"], kept_u=mine["kept_u"], kept_f=mine["kept_f"]) clock = time.time() part, problem = pilot_moment(dict(name="check_pilot", seed=s, moment=0, size=0.05, T=T0), DRY, short, 50) pilot_seconds = time.time() - clock with np.load(part) as p: shapes_ok = (p["z_medium"].shape == (2, 2, short + 51) and p["logs_medium"].shape == (2, (short + 50) // 5) and p["z_push"].shape == (2, 2, 51) and p["logs_push"].shape == (2, 10) and p["z_none"].shape == (2, short + 51) and p["logs_none"].shape == ((short + 50) // 5,)) built = p["built"].tolist() finite = bool(all(np.isfinite(p[key][d]).all() for key in ("z_push", "logs_push", "z_sham", "logs_sham") for d in range(2) if p["built"][d])) mech = dict(shapes=shapes_ok, restored=bool(p["restored"].all()), problem=problem or "none", pushes_built=built, built_arms_complete=finite, descriptions_present=bool(np.isfinite(p["desc_R"][:3]).all())) passed.append(shapes_ok and mech["restored"] and not problem and finite) lines.append(f"6. The pilot's moment task on a stand-in (T = {T0}, the pull changed 5% for {short} units, 50 after):" f" runs end to end: {mech}") clock = time.time() part, _ = calibration_moment(dict(name="check_calibration", seed=s, moment=0, T=T0), DRY, short) calibration_seconds = time.time() - clock with np.load(part) as p: mech = dict(shapes=p["paths"].shape == (len(SIZES), 2, 2, short + 1) and p["rec_built"].shape == (len(SIZES), 2), every_number_present=bool(np.isfinite(p["paths"]).all() and np.isfinite(p["R_none"]).all()), records=int(p["rec_built"].size)) passed.append(mech["shapes"] and mech["every_number_present"]) lines.append(f"7. The calibration's moment task on the stand-in ({short} units, every size both ways, first copy" f" only): runs end to end: {mech}") passed.append(not lock_changes()) lines.append(f"8. The cause test's locked files unchanged: {not lock_changes()}") full_pilot = pilot_seconds * (2 * (WINDOW + AFTER) + 4 * AFTER) / (2 * (short + 50) + 4 * 50) full_calibration = calibration_seconds * WINDOW / short lines += ["", f"Time: the cause test's base run took {reference_seconds:.0f} s for {run} units after settling;" f" scaled up, a pilot moment takes about {full_pilot:.0f} s and a calibration moment about" f" {full_calibration:.0f} s on one worker (more with seven at once).", f"All checks passed: {'yes' if all(passed) else 'NO'}", f"Run under plan.md SHA-256 {sha256(PLAN)} and medium_test.py SHA-256 {sha256(Path(__file__).resolve())}."] shutil.rmtree(DRY) CHECKS_TXT.write_text("\n".join(lines) + "\n") print("\n".join(lines)) if __name__ == "__main__": commands = dict(checks=checks, dryreport=dryreport, freeze=freeze, report=report) paced = dict(calibrate=calibrate, pilot=pilot) if len(sys.argv) >= 2 and sys.argv[1] in commands and len(sys.argv) == 2: commands[sys.argv[1]]() elif len(sys.argv) in (2, 3) and sys.argv[1] in paced: paced[sys.argv[1]](float(sys.argv[2]) if len(sys.argv) == 3 else START_BY) else: sys.exit(__doc__)