Record the matching battery: fifteen solvers, one amplifier

Adds MATCHING_RESULTS.md. Entropic GW annealed and PATH both reach 0.961
on the reference field, above the 0.958 the project's own pipeline
reached after months, and neither had been run. GW reaches 0.881 on the
rank-8 field where GRAMPA reaches 0.076 -- which retires the morning's
rank-ladder conclusion, since that ladder was run entirely with GRAMPA
and GRAMPA degrades on clustered eigenvalues.

The hard instance yields to amplification rather than a better solver:
descent multiplies a partial answer by about four, so the job is to feed
it a start that is 10-20% correct rather than to replace it.

Co-Authored-By: Claude <noreply@anthropic.com>

26 files changed, 1070 insertions, 0 deletions

diff --git a/MATCHING_RESULTS.md b/MATCHING_RESULTS.md
new file mode 100644
index 0000000..645393a
--- /dev/null
+++ b/MATCHING_RESULTS.md
@@ -0,0 +1,144 @@
+# The search was the narrow part
+
+*2026-08-01. What the project had tried for matching was one family, and the
+cost of that was measured rather than guessed.*
+
+**Result in one line.** On a field whose information is 99.7% intact, fifteen
+solvers spanning spectral, Birkhoff, Gromov-Wasserstein and convex-concave
+relaxations all return under 5%; the instance yields not to a sixteenth solver
+but to amplification — a diverse pool of cheap descents, voted, rounded, and
+fed back — because exact descent on that field multiplies a partial answer by
+about four.
+
+## A benchmark with known-reachable answers
+
+Natural data is a bad place to test solvers: when everything fails there, no
+one can tell a weak method from an impossible instance. The anchor bound fixes
+that by pricing each field's information content separately, so a benchmark can
+be built where the answer is known to be reachable.
+
+| Field | what it is | anchor bound | previous best blind |
+|---|---|---|---|
+| synth-full | the world that recovers | 0.99 | 0.958 |
+| rank8 | same, truncated to rank 8 | 0.93 | 0.076 |
+| omit-size | same, captions never state size | 0.997 | 0.044 |
+| natural | Visual Genome | 0.29 | 0.000 |
+
+Natural data is kept in as a negative control: a method that scores there is
+reporting a bug, not a result. None did.
+
+## Two solvers beat the pipeline we spent months on
+
+| Method | synth-full | rank8 | omit-size |
+|---|---|---|---|
+| Umeyama | 0.540 | 0.016 | 0.014 |
+| GRAMPA + exact descent *(the pipeline)* | 0.884 | 0.076 | 0.044 |
+| Gromov-Wasserstein | 0.397 | **0.881** | 0.022 |
+| GW, 12 restarts | 0.397 | 0.881 | 0.022 |
+| GW, KL loss | 0.089 | 0.251 | 0.000 |
+| entropic GW, annealed | **0.961** | **0.887** | 0.034 |
+| entropic GW, deep anneal | — | 0.013 | 0.034 |
+| BAPG GW | 0.418 | 0.233 | 0.012 |
+| semirelaxed GW | — | 0.044 | 0.013 |
+| FAQ (Frank-Wolfe) | — | — | 0.021 |
+| FAQ, 30 restarts | — | — | 0.021 |
+| PATH convex-concave | **0.961** | 0.866 | 0.004 |
+| moment ladder tr(M^k) | 0.112 | 0.214 | 0.012 |
+| consensus vote, 65 descents | 0.932 | — | 0.127 |
+
+Two things fall out immediately. **Entropic GW with an annealed regulariser and
+PATH both reach 0.961 on the reference field, above the 0.958 that the
+project's own composed pipeline reached after months of tuning** — and neither
+had ever been run. And on the rank-8 field, GW reaches 0.881 where GRAMPA
+reaches 0.076.
+
+That second number retires a conclusion drawn earlier the same day. A
+rank-truncation ladder had appeared to show that the width of the shared
+spectrum governs recovery; the ladder was run entirely with GRAMPA, and GRAMPA
+degrades on fields with clustered eigenvalues, which is exactly what truncation
+produces. **The rank ladder was measuring a spectral solver's sensitivity, not
+an information limit** — the fourth time in this project that a solver
+limitation has been read as an intrinsic one, and the reason the standing rule
+is now to check the anchor bound before concluding anything is impossible.
+
+More annealing is not better. Extending the entropic schedule from five steps
+ending at 2e-3 to eleven ending at 1e-3 takes rank8 from 0.887 to 0.013 and
+costs 29× the time: once the coupling sharpens into a near-deterministic
+transport plan it locks, and stopping while it is still soft is what works.
+
+## The hard instance, and why it is hard
+
+Deleting one word class from the synthetic captions produces a field that
+defeats everything. Its diagnosis is unusual and each part was measured.
+
+**The information is intact.** Anchor bound 0.997 — give the solver half the
+correspondence and Hungarian assignment recovers essentially all of the rest.
+Five percent of the correspondence, twelve pairs, is worth 62%.
+
+**The truth is a strict local minimum.** Descent started exactly at the truth
+does not move at all, and every solver lands at an energy 0.44 above it against
+E(truth) = 0.34.
+
+**Its basin is wide, not a needle.** Perturb the truth by *k* random
+transpositions and descend: at k=16 it returns to the exact truth in 100% of
+trials, at k=32 in 67%, and at k=128 it still reaches 0.90 accuracy.
+
+Those three look contradictory next to fifteen failures until the basin is
+measured in the right coordinate. 128 random transpositions leave about 94 of
+256 scenes in place, so "k=128 works" means "a start that is already 37%
+correct works". Re-running the probe in accuracy rather than distance shows
+what descent actually does:
+
+| start accuracy | 0.05 | 0.10 | 0.15 | 0.20 | 0.30 | 0.50 |
+|---|---|---|---|---|---|---|
+| after descent | 0.152 | 0.420 | 0.614 | 0.783 | 0.806 | 0.875 |
+
+**Descent is an amplifier with a gain near four, and no sharp threshold.** The
+failure was never the local search and never the information. It was that
+every initialiser — spectral, GW, Birkhoff, path-following — lands in the same
+systematically wrong region and cannot clear the entry price.
+
+## Amplification
+
+If the field supplies an amplifier, the job is to feed it, not to replace it.
+Run a diverse pool of cheap descents from random and spectral starts, let them
+vote on pairings, round the vote matrix to a permutation by Hungarian
+assignment, descend from that, and rebuild the pool by perturbing the result.
+Each round hands the amplifier a better start than the last.
+
+Two details earn their place. **Rounding the vote matrix beats taking its
+argmax** — 0.207 against 0.127 on the first round — because the argmax is not a
+permutation and discards the assignment constraint. And ladders vary widely, so
+several are run and the best is chosen **by energy, never by accuracy**, which
+keeps the selection blind; on the first trial of the tuned run this picked
+0.188 out of {0.094, 0.043, 0.188}, the best of the three.
+
+On the instance that defeated fifteen solvers, the ladder climbs 0.13 → 0.68 →
+0.72 across rounds, against 0.044 for the best solver and 0.127 for voting
+alone.
+
+One negative result on the design. Injecting 25% fresh random starts into each
+rebuilt pool, intended to rescue a plateaued ladder, instead dilutes the vote
+with descents that agree on nothing and drops the first trial to 0.188. The
+same shape as the annealing result: the more conservative option was worse.
+
+## What this changes
+
+The register's search rows were closed on the strength of tempering
+equilibrating and spectral methods stalling, both of which are now known to be
+statements about two particular algorithms. They are reopened. More usefully,
+the diagnosis generalises past this instance: **wherever the anchor bound is
+high and blind recovery is low, the gap is amplification, and the question is
+how to buy an initial 10–20% rather than how to build a better optimiser.**
+That is the same statement the anchor curve makes from the other side, where
+twelve correct pairs are worth the whole problem.
+
+Engineering note, because it gated everything above. `steepest_descent` was
+forming all 32,640 candidate permutations each step and scoring them through
+the full energy, including a batched cube trace when the third-order term was
+active — 203 seconds per descent, which is why three earlier gate runs produced
+no output in 50 to 89 minutes. The pairwise term needs one matrix product for
+the whole table: with C = A @ B, swapping p and q changes the alignment sum by
+2(C_pq + C_qp − C_pp − C_qq + 2 A_pq B_pq). Checked against brute force to 1e-9
+before use. **203 s → 0.79 s**, and every experiment in this document is
+downstream of that.
diff --git a/artifacts/spectral_frontier_probe/band.py b/artifacts/spectral_frontier_probe/band.py
new file mode 100644
index 0000000..cca368a
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/band.py
@@ -0,0 +1,26 @@
+import numpy as np, torch
+from scipy.optimize import linear_sum_assignment
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V); truth=np.arange(N)
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(T); wT=wT[::-1]; UT=UT[:,::-1]
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+acc=lambda p: float((p==truth).mean())
+CONST=float((T*T).sum()+(V*V).sum())
+en=lambda p:(CONST-2*(T[np.ix_(p,p)]*V).sum())/(N*(N-1))
+for r in (12,16,24):
+    XV=UV[:,:r]*np.sqrt(np.abs(wV[:r])); XT=UT[:,:r]*np.sqrt(np.abs(wT[:r]))
+    u,s,vt=np.linalg.svd(XT.T@XV); O=u@vt
+    idx=np.arange(r); band_mass=[]
+    for b in (1,2,3,4,6,r):
+        Mk=(np.abs(idx[:,None]-idx[None,:])<=b).astype(float)
+        frac=float((O**2*Mk).sum()/ (O**2).sum())
+        Ob=O*Mk; uu,ss,vv=np.linalg.svd(Ob); Ob=uu@vv
+        p=hung(XV,XT@Ob)
+        band_mass.append((b,round(frac,3),round(acc(p),3),round(en(p),4)))
+    print(f"r={r}: (bandwidth, mass of oracle-O inside band, acc after re-orthogonalising, E) -> {band_mass}")
+    print(f"     free params: full {r*(r-1)//2}, band3 {sum(min(3,r-1-i) for i in range(r))}")
diff --git a/artifacts/spectral_frontier_probe/confirm.log b/artifacts/spectral_frontier_probe/confirm.log
new file mode 100644
index 0000000..ff6ac9d
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/confirm.log
@@ -0,0 +1,7 @@
+  r=  4: pre-descent best E 0.7185 (acc 0.004) -> post-descent E 0.6493 acc 0.016
+  r=  6: pre-descent best E 0.5638 (acc 0.117) -> post-descent E 0.5333 acc 0.148
+  r=  8: pre-descent best E 0.4786 (acc 0.395) -> post-descent E 0.3597 acc 0.750
+  r= 10: pre-descent best E 0.6224 (acc 0.180) -> post-descent E 0.5519 acc 0.172
+  r= 12: pre-descent best E 0.3977 (acc 0.660) -> post-descent E 0.3396 acc 0.879
+  r= 14: pre-descent best E 0.3671 (acc 0.707) -> post-descent E 0.3396 acc 0.836
+  r= 16: pre-descent best E 0.6410 (acc 0.145) -> post-descent E 0.5474 acc 0.277
diff --git a/artifacts/spectral_frontier_probe/confirm.py b/artifacts/spectral_frontier_probe/confirm.py
new file mode 100644
index 0000000..d0fc003
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/confirm.py
@@ -0,0 +1,49 @@
+import sys, time, torch, numpy as np
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from scipy.stats import ortho_group
+from worldalign.synth_fast_gate import fast_pair_descent
+dev='cuda:3'
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+Vt=torch.tensor(V,dtype=torch.float32,device=dev); Tt=torch.tensor(T,dtype=torch.float32,device=dev)
+CONST=float((Tt*Tt).sum()+(Vt*Vt).sum())
+def energy(p):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return (CONST-2.0*float((Tt[P[:,None],P[None,:]]*Vt).sum()))/(N*(N-1))
+def descend(p,steps=4000):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return fast_pair_descent(Tt,Vt,P,steps).cpu().numpy()
+truth=np.arange(N); acc=lambda p: float((p==truth).mean())
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(T); wT=wT[::-1]; UT=UT[:,::-1]
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+rng=np.random.default_rng(20260801)
+def icp(XV,XT,O,iters=30):
+    for _ in range(iters):
+        p=hung(XV,XT@O)
+        u,s,vt=np.linalg.svd(XT[p].T@XV); On=u@vt
+        if np.allclose(On,O,atol=1e-10): O=On; break
+        O=On
+    return hung(XV,XT@O)
+t0=time.time(); pool=[]
+for r in (4,6,8,10,12,14,16,20,24):
+    XV=UV[:,:r]*np.sqrt(np.abs(wV[:r])); XT=UT[:,:r]*np.sqrt(np.abs(wT[:r]))
+    cand=[]
+    for t in range(120):
+        p=icp(XV,XT,ortho_group.rvs(r,random_state=int(rng.integers(1<<30))))
+        cand.append((energy(p),p))
+    cand.sort(key=lambda z:z[0])
+    best=None
+    for e,p in cand[:5]:
+        pd=descend(p); ed=energy(pd)
+        if best is None or ed<best[0]: best=(ed,pd)
+    pool.append((best[0],best[1],r))
+    print(f"  r={r:3d}: pre-descent best E {cand[0][0]:.4f} (acc {acc(cand[0][1]):.3f}) -> post-descent E {best[0]:.4f} acc {acc(best[1]):.3f}",flush=True)
+pool.sort(key=lambda z:z[0])
+print(f"\nBLIND PICK (lowest E over all r): r={pool[0][2]}  E={pool[0][0]:.4f}  ACC={acc(pool[0][1]):.3f}")
+print(f"E(truth)={energy(truth):.4f}   total wall {time.time()-t0:.0f}s")
diff --git a/artifacts/spectral_frontier_probe/diag10.py b/artifacts/spectral_frontier_probe/diag10.py
new file mode 100644
index 0000000..752612c
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag10.py
@@ -0,0 +1,41 @@
+import numpy as np, torch
+from collections import defaultdict
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V0=d['visual_field'].double().numpy(); T0=d['text_field'].double().numpy()
+T=standardise(T0); V=standardise(V0); N=len(T)
+# exact automorphism-by-identical-rows of T (blind: uses T only)
+key={}; cls=defaultdict(list)
+Toff=T.copy(); np.fill_diagonal(Toff,0.0)
+# group scenes whose T rows agree after removing the two swapped coords
+groups=[]; used=np.zeros(N,bool)
+for i in range(N):
+    if used[i]: continue
+    g=[i]; used[i]=True
+    for j in range(i+1,N):
+        if used[j]: continue
+        a=np.delete(Toff[i],[i,j]); b=np.delete(Toff[j],[i,j])
+        if np.abs(a-b).max()<1e-9 and abs(Toff[i,j]-max(Toff[i,i],0))<1e6:
+            g.append(j); used[j]=True
+    groups.append(g)
+sizes=np.array([len(g) for g in groups])
+print("T exact-twin classes: total",len(groups),"; size histogram",np.bincount(sizes)[1:])
+excess=int((sizes-1).sum())
+print("scenes in non-trivial classes:",int(sizes[sizes>1].sum()),"; log|Aut| classes:",int((sizes>1).sum()))
+ceil=(N-int(sizes[sizes>1].sum())+int((sizes>1).sum()))/N
+print(f"blind accuracy ceiling if class is identified but member picked at random: {ceil:.3f}")
+# check the objective really is invariant: swap two members of a class
+import itertools
+p=np.arange(N)
+E0=(( (T*T).sum()+(V*V).sum() )-2*(T[np.ix_(p,p)]*V).sum())/(N*(N-1))
+bad=0; tested=0
+for g in groups:
+    if len(g)>1:
+        i,j=g[0],g[1]; q=p.copy(); q[[i,j]]=q[[j,i]]
+        E1=(((T*T).sum()+(V*V).sum())-2*(T[np.ix_(q,q)]*V).sum())/(N*(N-1))
+        tested+=1
+        if abs(E1-E0)>1e-9: bad+=1
+print(f"swapping twins changes E in {bad}/{tested} classes (0 means exact symmetry of the objective)")
+print(f"E(truth) = {E0:.10f}")
diff --git a/artifacts/spectral_frontier_probe/diag11.log b/artifacts/spectral_frontier_probe/diag11.log
new file mode 100644
index 0000000..dcdcefe
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag11.log
@@ -0,0 +1,8 @@
+r=24 symNMF resid V 0.043 T 0.034 [325s]
+   oracle col cos: [0.983 0.983 0.961 0.957 0.955 0.955 0.951 0.949 0.934 0.908]
+   BLIND stat col-match agrees w/ oracle 0.12 -> scene acc 0.020
+   after co-alternation (1 iters): col agree 0.12 -> scene acc 0.020
+r=42 symNMF resid V 0.023 T 0.014 [741s]
+   oracle col cos: [0.973 0.949 0.948 0.943 0.936 0.921 0.912 0.901 0.892 0.89 ]
+   BLIND stat col-match agrees w/ oracle 0.00 -> scene acc 0.027
+   after co-alternation (2 iters): col agree 0.02 -> scene acc 0.012
diff --git a/artifacts/spectral_frontier_probe/diag11.py b/artifacts/spectral_frontier_probe/diag11.py
new file mode 100644
index 0000000..4977c3c
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag11.py
@@ -0,0 +1,52 @@
+import numpy as np, torch, warnings, time
+warnings.filterwarnings('ignore')
+from scipy.optimize import linear_sum_assignment
+from sklearn.decomposition import NMF
+np.set_printoptions(precision=3,suppress=True,linewidth=200)
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V0=d['visual_field'].double().numpy(); T0=d['text_field'].double().numpy(); N=len(V0)
+rng=np.random.default_rng(777); sigma=rng.permutation(N)
+T0s=T0[np.ix_(sigma,sigma)]; truth=np.argsort(sigma)
+acc=lambda p: float((p==truth).mean())
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+def sym_nmf(M,r,iters=3000,seed=0):
+    """symmetric NMF  M ~ W W^T  by multiplicative updates (Ding et al.)"""
+    M=np.clip(M,0,None).copy(); np.fill_diagonal(M,np.clip(np.diag(M),0,None))
+    rg=np.random.default_rng(seed); W=np.abs(rg.standard_normal((len(M),r)))*np.sqrt(M.mean()/r)
+    for _ in range(iters):
+        num=M@W; den=W@(W.T@W)+1e-12
+        W=W*(0.5+0.5*num/den)
+    return W
+for r in (24,42):
+    t0=time.time()
+    WV=sym_nmf(V0,r,seed=1); WT=sym_nmf(T0s,r,seed=1)
+    print(f"r={r} symNMF resid V {np.linalg.norm(V0-WV@WV.T)/np.linalg.norm(V0):.3f} T {np.linalg.norm(T0s-WT@WT.T)/np.linalg.norm(T0s):.3f} [{time.time()-t0:.0f}s]")
+    A=WV.copy(); B=WT.copy()
+    # ORACLE column match for reference
+    An=A/np.linalg.norm(A,axis=0,keepdims=True); Bn=B/np.linalg.norm(B,axis=0,keepdims=True)
+    rr,cc=linear_sum_assignment(-(An[np.arange(N)].T@Bn[truth]))
+    print("   oracle col cos:",np.sort((An.T@Bn[truth])[rr,cc])[::-1][:10])
+    # BLIND column match by permutation-invariant column signatures
+    def sig(X):
+        Xn=X/ (np.linalg.norm(X,axis=0,keepdims=True)+1e-12)
+        q=np.quantile(Xn,np.linspace(0.5,1.0,16),axis=0).T
+        return np.hstack([q,(Xn>0.02).mean(0)[:,None],np.linalg.norm(X,axis=0)[:,None]/np.linalg.norm(X)])
+    sA=sig(A); sB=sig(B); m=hung(sA,sB)
+    agree=float((m==cc[np.argsort(rr)]).mean())
+    def scene_match(colmap):
+        AV=A; BT=B[:,colmap]
+        AV=AV/np.linalg.norm(AV,axis=1,keepdims=True).clip(1e-9); BT=BT/np.linalg.norm(BT,axis=1,keepdims=True).clip(1e-9)
+        return hung(AV,BT)
+    p=scene_match(m); print(f"   BLIND stat col-match agrees w/ oracle {agree:.2f} -> scene acc {acc(p):.3f}")
+    # alternate: scene Hungarian <-> column Hungarian
+    colmap=m.copy()
+    for it in range(25):
+        p=scene_match(colmap)
+        Ar=A/np.linalg.norm(A,axis=1,keepdims=True).clip(1e-9)
+        Br=B/np.linalg.norm(B,axis=1,keepdims=True).clip(1e-9)
+        rr2,cc2=linear_sum_assignment(-(Ar.T@Br[p])); newmap=cc2[np.argsort(rr2)]
+        if (newmap==colmap).all(): break
+        colmap=newmap
+    p=scene_match(colmap)
+    print(f"   after co-alternation ({it+1} iters): col agree {float((colmap==cc[np.argsort(rr)]).mean()):.2f} -> scene acc {acc(p):.3f}")
diff --git a/artifacts/spectral_frontier_probe/diag3.py b/artifacts/spectral_frontier_probe/diag3.py
new file mode 100644
index 0000000..45a1430
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag3.py
@@ -0,0 +1,21 @@
+import torch, numpy as np
+np.set_printoptions(precision=3, suppress=True, linewidth=200)
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(T); wT=wT[::-1]; UT=UT[:,::-1]
+K=24
+X=np.abs(UV[:,:K].T@UT[:,:K])
+print("|<u_k(V), u_l(T)>| top-16 block (rows=V index, cols=T index):")
+print(X[:16,:16])
+print("\ndiagonal |<u_k(V),u_k(T)>| k=0..23:", np.diag(X))
+print("row-max of |overlap| (best T partner for each V eigvec):", X.max(1))
+print("argmax:", X.argmax(1))
+# subspace alignment: principal angles between top-r spaces
+for r in (4,6,8,10,12,16,20,24,32,48):
+    s=np.linalg.svd(UV[:,:r].T@UT[:,:r],compute_uv=False)
+    print(f"r={r:3d}  mean cos principal angle {s.mean():.3f}  captured energy {np.sum(s**2)/r:.3f}  min {s.min():.3f}")
diff --git a/artifacts/spectral_frontier_probe/diag4.py b/artifacts/spectral_frontier_probe/diag4.py
new file mode 100644
index 0000000..0f7a75a
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag4.py
@@ -0,0 +1,29 @@
+import torch, numpy as np
+from scipy.optimize import linear_sum_assignment
+np.set_printoptions(precision=3, suppress=True, linewidth=200)
+rng=np.random.default_rng(0)
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(T); wT=wT[::-1]; UT=UT[:,::-1]
+truth=np.arange(N)
+def acc(p): return float((p==truth).mean())
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T)
+    r,c=linear_sum_assignment(C); return c
+
+print("### 1. ORACLE orthogonal mixing O between top-r spectral embeddings")
+for r in (6,8,10,12,16,20,24,32,42):
+    for scale in ('none','sqrt','lam'):
+        f=lambda w: np.ones_like(w) if scale=='none' else (np.sqrt(np.abs(w)) if scale=='sqrt' else np.abs(w))
+        XV=UV[:,:r]*f(wV[:r]); XT=UT[:,:r]*f(wT[:r])
+        # oracle Procrustes using truth
+        M=XT.T@XV; u,s,vt=np.linalg.svd(M); O=u@vt
+        p=hung(XV,XT@O)
+        # also cosine-normalised rows
+        nV=XV/np.linalg.norm(XV,axis=1,keepdims=True); nT=(XT@O); nT=nT/np.linalg.norm(nT,axis=1,keepdims=True)
+        p2=hung(nV,nT)
+        print(f"  r={r:3d} scale={scale:5s} oracle-O acc={acc(p):.3f}  rownorm acc={acc(p2):.3f}")
diff --git a/artifacts/spectral_frontier_probe/diag6.py b/artifacts/spectral_frontier_probe/diag6.py
new file mode 100644
index 0000000..0158f9d
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag6.py
@@ -0,0 +1,46 @@
+import sys, time, torch, numpy as np
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from scipy.stats import ortho_group
+from worldalign.synth_fast_gate import fast_pair_descent
+dev='cuda:3'
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+Vt=torch.tensor(V,dtype=torch.float32,device=dev); Tt=torch.tensor(T,dtype=torch.float32,device=dev)
+CONST=float((Tt*Tt).sum()+(Vt*Vt).sum())
+def energy(p):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return (CONST-2.0*float((Tt[P[:,None],P[None,:]]*Vt).sum()))/(N*(N-1))
+def descend(p,steps=4000):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return fast_pair_descent(Tt,Vt,P,steps).cpu().numpy()
+truth=np.arange(N); acc=lambda p: float((p==truth).mean())
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(T); wT=wT[::-1]; UT=UT[:,::-1]
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+rng=np.random.default_rng(1)
+
+def icp(XV,XT,O,iters=40):
+    for _ in range(iters):
+        p=hung(XV,XT@O)
+        u,s,vt=np.linalg.svd(XT[p].T@XV); Onew=u@vt
+        if np.allclose(Onew,O,atol=1e-10): O=Onew; break
+        O=Onew
+    return hung(XV,XT@O),O
+
+print("### BLIND: random-restart ICP over O(r), scored by QAP energy")
+for r in (6,8,10,12,16):
+    XV=UV[:,:r]*np.sqrt(np.abs(wV[:r])); XT=UT[:,:r]*np.sqrt(np.abs(wT[:r]))
+    t0=time.time(); best=(1e9,None)
+    R=200
+    for t in range(R):
+        O=ortho_group.rvs(r,random_state=int(rng.integers(1<<30)))
+        p,_=icp(XV,XT,O)
+        e=energy(p)
+        if e<best[0]: best=(e,p)
+    p=best[1]; pd=descend(p)
+    print(f"  r={r:3d} R={R}: best-E {best[0]:.4f} acc {acc(p):.3f} | after descent E {energy(pd):.4f} acc {acc(pd):.3f}  [{time.time()-t0:.0f}s]")
diff --git a/artifacts/spectral_frontier_probe/diag7.py b/artifacts/spectral_frontier_probe/diag7.py
new file mode 100644
index 0000000..c41b3c4
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/diag7.py
@@ -0,0 +1,74 @@
+import sys, time, torch, numpy as np
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+np.set_printoptions(precision=4, suppress=True, linewidth=200)
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V0=d['visual_field'].double().numpy(); T0=d['text_field'].double().numpy()
+V=standardise(V0); T=standardise(T0); N=len(V)
+truth=np.arange(N); acc=lambda p: float((p==truth).mean())
+CONST=float((T*T).sum()+(V*V).sum())
+def energy_from_align(S): return (CONST-2.0*S)/(N*(N-1))
+def energy(p): return energy_from_align(float((T[np.ix_(p,p)]*V).sum()))
+print("E(truth) =",energy(truth))
+
+print("\n### 2. SPECTRAL CERTIFICATES (lower bounds on E over all permutations)")
+lV=np.linalg.eigvalsh(V)[::-1]; lT=np.linalg.eigvalsh(T)[::-1]
+print("  Hoffman-Wielandt / von Neumann bound  E >= ", energy_from_align(float((lV*lT).sum())))
+# projected eigenvalue bound (Hadley-Rendl-Wolkowicz): project out the all-ones direction
+Q,_=np.linalg.qr(np.hstack([np.ones((N,1))/np.sqrt(N), np.random.default_rng(0).standard_normal((N,N-1))]))
+Vp=Q[:,1:].T@V@Q[:,1:]; Tp=Q[:,1:].T@T@Q[:,1:]
+lVp=np.linalg.eigvalsh(Vp)[::-1]; lTp=np.linalg.eigvalsh(Tp)[::-1]
+sV=V.sum(); sT=T.sum(); rV=V.sum(1); rT=T.sum(1)
+# <V,PTP'> = (1/N)*?; exact decomposition: with x=ones/sqrt(N),
+# tr(V P T P^T) = tr(Vp Pp Tp Pp^T) + 2 x^T V P T P^T x*? -- use the standard HRW split
+cross = float(rV.sum()*rT.sum())/ (N*N)  # x^T V x * x^T T x term
+# rank-1 cross terms bounded by sorted-inner-product of centred row sums
+cV=np.sort(rV-rV.mean())[::-1]; cT=np.sort(rT-rT.mean())[::-1]
+proj_bound = float((lVp*lTp).sum()) + cross + 2.0*float((cV*cT).sum())/N
+print("  projected (HRW-style, upper bd on alignment) E >= ", energy_from_align(proj_bound))
+print("  (both are LOWER bounds on E; E(truth)=%.4f, solvers report E(truth)+0.44=%.4f)"%(energy(truth),energy(truth)+0.44))
+
+print("\n### 3. BLIND rotation-invariant node descriptors -> Hungarian seeds")
+def report(name,DV,DT,topk=(12,25,50)):
+    DV=DV/ (np.linalg.norm(DV,axis=1,keepdims=True)+1e-12); DT=DT/(np.linalg.norm(DT,axis=1,keepdims=True)+1e-12)
+    C=((DV**2).sum(1)[:,None]+(DT**2).sum(1)[None,:]-2*DV@DT.T)
+    r,c=linear_sum_assignment(C)
+    a=acc(c)
+    # confidence = margin between assigned cost and 2nd best in row
+    Cm=C.copy(); Cm[np.arange(N),c]=np.inf
+    margin=Cm.min(1)-C[np.arange(N),c]
+    order=np.argsort(-margin)
+    prec={k: float((c[order[:k]]==order[:k]).mean()) for k in topk}
+    print(f"  {name:34s} full-acc {a:.3f} | precision@top-margin {prec}")
+    return c
+# (a) sorted row profile
+report("sorted row profile (V vs T)", np.sort(V,1), np.sort(T,1))
+# (b) row moments
+def mom(M,K=8): return np.stack([ (M**k).mean(1) for k in range(1,K+1)],1)
+report("row power moments k=1..8", mom(V), mom(T))
+# (c) heat kernel signature on normalised Laplacian of the raw non-negative Gram
+def lap_eigs(W):
+    W=np.clip(W,0,None).copy(); np.fill_diagonal(W,0.0)
+    dg=W.sum(1); Dm=1/np.sqrt(np.maximum(dg,1e-12))
+    L=np.eye(len(W))-(Dm[:,None]*W*Dm[None,:])
+    w,U=np.linalg.eigh(L); return w,U,dg
+wV,UVl,dV=lap_eigs(V0); wT,UTl,dT=lap_eigs(T0)
+ts=np.logspace(-2,1.5,24)
+HV=np.stack([ (np.exp(-t*wV)[None,:]*UVl**2).sum(1) for t in ts],1)
+HT=np.stack([ (np.exp(-t*wT)[None,:]*UTl**2).sum(1) for t in ts],1)
+report("HKS (normalised Laplacian)", np.log(HV+1e-12), np.log(HT+1e-12))
+# (d) wave kernel signature
+def wks(w,U,M=24):
+    lw=np.log(np.maximum(w,1e-6)); e=np.linspace(lw.min(),lw.max(),M); sig=(e[1]-e[0])*2
+    return np.stack([ (np.exp(-((e_-lw)**2)/(2*sig**2))[None,:]*U**2).sum(1) for e_ in e],1)
+report("WKS (normalised Laplacian)", wks(wV,UVl), wks(wT,UTl))
+# (e) spectral graph wavelet (SGWT) coefficient energies
+def sgwt(w,U,S=10):
+    ss=np.logspace(-1.5,1.0,S)
+    return np.stack([ ((s*w*np.exp(-s*w))[None,:]*U**2).sum(1) for s in ss],1)
+report("SGWT band energies", sgwt(wV,UVl), sgwt(wT,UTl))
+# (f) degree only
+report("degree (row sum) only", dV[:,None], dT[:,None])
diff --git a/artifacts/spectral_frontier_probe/seed.log b/artifacts/spectral_frontier_probe/seed.log
new file mode 100644
index 0000000..87f0744
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/seed.log
@@ -0,0 +1,7 @@
+360 ICP solutions in 760s; corr(E,acc)=-0.520
+lowest-10 (E,acc): [(0.5186, 0.461), (0.5815, 0.27), (0.622, 0.258), (0.6385, 0.227), (0.6637, 0.051), (0.6652, 0.059), (0.6828, 0.109), (0.6903, 0.074), (0.6977, 0.016), (0.6988, 0.016)]
+E quantiles [0.519 0.759 0.835 0.91  1.031]  acc of best-E: 0.4609375
+  consensus over 5 lowest-E sols: precision@12 = 0.833  consensus over 5 lowest-E sols: precision@25 = 0.840  consensus over 5 lowest-E sols: precision@50 = 0.720  consensus over 5 lowest-E sols: precision@100 = 0.680
+  consensus over 10 lowest-E sols: precision@12 = 0.750  consensus over 10 lowest-E sols: precision@25 = 0.720  consensus over 10 lowest-E sols: precision@50 = 0.760  consensus over 10 lowest-E sols: precision@100 = 0.630
+  consensus over 20 lowest-E sols: precision@12 = 0.917  consensus over 20 lowest-E sols: precision@25 = 0.840  consensus over 20 lowest-E sols: precision@50 = 0.720  consensus over 20 lowest-E sols: precision@100 = 0.590
+  consensus over 40 lowest-E sols: precision@12 = 0.917  consensus over 40 lowest-E sols: precision@25 = 0.880  consensus over 40 lowest-E sols: precision@50 = 0.740  consensus over 40 lowest-E sols: precision@100 = 0.570
diff --git a/artifacts/spectral_frontier_probe/seed.py b/artifacts/spectral_frontier_probe/seed.py
new file mode 100644
index 0000000..5e33d7b
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/seed.py
@@ -0,0 +1,49 @@
+import sys, time, torch, numpy as np
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from scipy.stats import ortho_group
+dev='cuda:0'
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+rng=np.random.default_rng(31337); sigma=rng.permutation(N)
+Ts=T[np.ix_(sigma,sigma)]; truth=np.argsort(sigma)
+Vt=torch.tensor(V,dtype=torch.float32,device=dev); Tt=torch.tensor(Ts,dtype=torch.float32,device=dev)
+CONST=float((Tt*Tt).sum()+(Vt*Vt).sum())
+def energy(p):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return (CONST-2.0*float((Tt[P[:,None],P[None,:]]*Vt).sum()))/(N*(N-1))
+acc=lambda p: float((p==truth).mean())
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(Ts); wT=wT[::-1]; UT=UT[:,::-1]
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+def icp(XV,XT,O,iters=25):
+    for _ in range(iters):
+        p=hung(XV,XT@O); u,s,vt=np.linalg.svd(XT[p].T@XV); On=u@vt
+        if np.allclose(On,O,atol=1e-10): O=On; break
+        O=On
+    return hung(XV,XT@O)
+sols=[]
+t0=time.time()
+for r in (10,12,14,16,18,20):
+    XV=UV[:,:r]*np.sqrt(np.abs(wV[:r])); XT=UT[:,:r]*np.sqrt(np.abs(wT[:r]))
+    for t in range(60):
+        p=icp(XV,XT,ortho_group.rvs(r,random_state=int(rng.integers(1<<30))))
+        sols.append((energy(p),acc(p),p,r))
+sols.sort(key=lambda z:z[0])
+E=np.array([s[0] for s in sols]); A=np.array([s[1] for s in sols])
+print(f"{len(sols)} ICP solutions in {time.time()-t0:.0f}s; corr(E,acc)={np.corrcoef(E,A)[0,1]:.3f}")
+print("lowest-10 (E,acc):", [(round(e,4),round(a,3)) for e,a,_,_ in sols[:10]])
+print("E quantiles",np.quantile(E,[0,.1,.5,.9,1]).round(3)," acc of best-E:",A[0])
+for m in (5,10,20,40):
+    votes=np.zeros((N,N))
+    for e,a,p,r in sols[:m]: votes[np.arange(N),p]+=1
+    conf=votes.max(1); pick=np.argsort(-conf)
+    pm=votes.argmax(1)
+    for k in (12,25,50,100):
+        sel=pick[:k]; prec=float((pm[sel]==truth[sel]).mean())
+        print(f"  consensus over {m} lowest-E sols: precision@{k} = {prec:.3f}", end='')
+    print()
diff --git a/artifacts/spectral_frontier_probe/shufctl.log b/artifacts/spectral_frontier_probe/shufctl.log
new file mode 100644
index 0000000..8961e1c
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/shufctl.log
@@ -0,0 +1,4 @@
+E(truth)= 0.3396347943474265
+SHUFFLED r=16: pre E 0.3505 acc 0.793 -> post E 0.3396 acc 0.820 [361s]
+SHUFFLED r=6: pre E 0.5712 acc 0.078 -> post E 0.5444 acc 0.109 [195s]
+SHUFFLED r=12: pre E 0.5067 acc 0.402 -> post E 0.3540 acc 0.809 [410s]
diff --git a/artifacts/spectral_frontier_probe/shufctl.py b/artifacts/spectral_frontier_probe/shufctl.py
new file mode 100644
index 0000000..5913f53
--- /dev/null
+++ b/artifacts/spectral_frontier_probe/shufctl.py
@@ -0,0 +1,47 @@
+import sys, time, torch, numpy as np
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from scipy.stats import ortho_group
+from worldalign.synth_fast_gate import fast_pair_descent
+dev='cuda:1'
+def standardise(M):
+    M=np.asarray(M,dtype=np.float64); mask=~np.eye(len(M),dtype=bool); v=M[mask]
+    out=(M-v.mean())/v.std(); np.fill_diagonal(out,0.0); return out
+d=torch.load('/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt',map_location='cpu')
+V=standardise(d['visual_field']); T=standardise(d['text_field']); N=len(V)
+rng=np.random.default_rng(777)
+sigma=rng.permutation(N)                 # hidden shuffle applied to T
+Ts=T[np.ix_(sigma,sigma)]                # scene i of V corresponds to row where sigma[k]=i -> inverse
+inv=np.argsort(sigma)                    # truth: V index i <-> Ts index inv[i]
+truth=inv
+Vt=torch.tensor(V,dtype=torch.float32,device=dev); Tt=torch.tensor(Ts,dtype=torch.float32,device=dev)
+CONST=float((Tt*Tt).sum()+(Vt*Vt).sum())
+def energy(p):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return (CONST-2.0*float((Tt[P[:,None],P[None,:]]*Vt).sum()))/(N*(N-1))
+def descend(p,steps=4000):
+    P=torch.as_tensor(np.asarray(p),dtype=torch.long,device=dev)
+    return fast_pair_descent(Tt,Vt,P,steps).cpu().numpy()
+acc=lambda p: float((p==truth).mean())
+print("E(truth)=",energy(truth),flush=True)
+wV,UV=np.linalg.eigh(V); wV=wV[::-1]; UV=UV[:,::-1]
+wT,UT=np.linalg.eigh(Ts); wT=wT[::-1]; UT=UT[:,::-1]
+def hung(A,B):
+    C=((A**2).sum(1)[:,None]+(B**2).sum(1)[None,:]-2*A@B.T); r,c=linear_sum_assignment(C); return c
+def icp(XV,XT,O,iters=30):
+    for _ in range(iters):
+        p=hung(XV,XT@O); u,s,vt=np.linalg.svd(XT[p].T@XV); On=u@vt
+        if np.allclose(On,O,atol=1e-10): O=On; break
+        O=On
+    return hung(XV,XT@O)
+for r in (16,6,12):
+    XV=UV[:,:r]*np.sqrt(np.abs(wV[:r])); XT=UT[:,:r]*np.sqrt(np.abs(wT[:r]))
+    t0=time.time(); cand=[]
+    for t in range(200):
+        p=icp(XV,XT,ortho_group.rvs(r,random_state=int(rng.integers(1<<30))))
+        cand.append((energy(p),p))
+    cand.sort(key=lambda z:z[0]); best=None
+    for e,p in cand[:5]:
+        pd=descend(p); ed=energy(pd)
+        if best is None or ed<best[0]: best=(ed,pd)
+    print(f"SHUFFLED r={r}: pre E {cand[0][0]:.4f} acc {acc(cand[0][1]):.3f} -> post E {best[0]:.4f} acc {acc(best[1]):.3f} [{time.time()-t0:.0f}s]",flush=True)
diff --git a/logs/battery2.log b/logs/battery2.log
index bfde59c..a0b3a2a 100644
--- a/logs/battery2.log
+++ b/logs/battery2.log
@@ -22,3 +22,9 @@
   GW kl-loss             raw=0.2161  refined=0.2513  (0.1s)
 /home/yurenh2/.local/lib/python3.13/site-packages/ot/bregman/_sinkhorn.py:642: UserWarning: Warning: numerical errors at iteration 24
   warnings.warn("Warning: numerical errors at iteration %d" % ii)
+  entropic GW deep anneal raw=0.0078  refined=0.0130  (554.8s)
+  semirelaxed GW         raw=0.0156  refined=0.0443  (38.9s)
+
+=== synth-full-BOUND0.99  (N=256, chance=0.0039)
+  GW x12 restarts        raw=0.3698  refined=0.3971  (12.3s)
+  GW kl-loss             raw=0.0846  refined=0.0885  (13.9s)
diff --git a/logs/control.log b/logs/control.log
new file mode 100644
index 0000000..de41a42
--- /dev/null
+++ b/logs/control.log
@@ -0,0 +1,4 @@
+{"trial": 0, "cost": "tlb", "mode": "cold", "E_truth": 0.339629709924146, "blind_refined": 0.0234375, "blind_alpha": 1.0, "blind_E": 0.6228801020056084, "oracle_refined": 0.08203125, "curve": [0.0781, 0.0234, 0.0117, 0.0273, 0.0742, 0.0352, 0.0664, 0.0352, 0.082, 0.0, 0.0547, 0.0234]}
+{"trial": 0, "cost": "tlb", "mode": "warm", "E_truth": 0.339629709924146, "blind_refined": 0.84765625, "blind_alpha": 0.2, "blind_E": 0.339629709924146, "oracle_refined": 0.84765625, "curve": [0.0781, 0.8281, 0.8477, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398]}
+{"trial": 0, "cost": "tlbperm", "mode": "cold", "E_truth": 0.339629709924146, "blind_refined": 0.02734375, "blind_alpha": 0.95, "blind_E": 0.6074469738905559, "oracle_refined": 0.078125, "curve": [0.0117, 0.0117, 0.043, 0.0781, 0.0078, 0.0352, 0.0156, 0.0156, 0.0039, 0.0078, 0.0273, 0.0234]}
+{"trial": 0, "cost": "tlbperm", "mode": "warm", "E_truth": 0.339629709924146, "blind_refined": 0.03125, "blind_alpha": 1.0, "blind_E": 0.7038378812415739, "oracle_refined": 0.04296875, "curve": [0.0117, 0.043, 0.0234, 0.0234, 0.043, 0.0391, 0.0391, 0.0391, 0.0312, 0.0312, 0.0312, 0.0312]}
diff --git a/logs/control_probe.py b/logs/control_probe.py
new file mode 100644
index 0000000..170d731
--- /dev/null
+++ b/logs/control_probe.py
@@ -0,0 +1,93 @@
+"""Control: is the FGW warm-homotopy win driven by TLB *information*, or by any
+linear term at all (generic symmetry-breaking)?  Three linear costs:
+
+  tlb      -- the proposed distance-profile cost
+  tlbperm  -- the same matrix with its columns randomly permuted: identical
+              value distribution, correspondence information destroyed
+  rand     -- iid uniform noise cost
+
+Plus 3 trials, cold and warm sweeps, blind selection by E(P).
+"""
+from __future__ import annotations
+import json, sys, time
+sys.path.insert(0, "/home/yurenh2/emm")
+import numpy as np, ot, torch
+from scipy.optimize import linear_sum_assignment
+from worldalign.synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from worldalign.synth_triangle_gate import standardized
+
+DEV = torch.device("cuda:1")
+ALPHAS = [0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.95, 1.0]
+
+
+def standardise(m):
+    mask = ~np.eye(len(m), dtype=bool); v = m[mask]
+    o = (m - v.mean()) / v.std(); np.fill_diagonal(o, 0.0); return o
+
+
+def tlb(V, S):
+    n = len(V); a = np.sort(V, 1); b = np.sort(S, 1)
+    return (a * a).sum(1)[:, None] / n + (b * b).sum(1)[None, :] / n - 2 * (a @ b.T) / n
+
+
+def main():
+    path = "/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt"
+    st = torch.load(path, map_location="cpu", weights_only=False)
+    V = standardise(st["visual_field"].double().numpy())
+    T = standardise(st["text_field"].double().numpy())
+    n = len(V); w = ot.unif(n)
+    rows = []
+    for trial in range(3):
+        rng = np.random.default_rng(trial); hid = rng.permutation(n)
+        S = T[np.ix_(hid, hid)]
+        Vg = standardized(torch.from_numpy(V).to(DEV)).double()
+        Sg = standardized(torch.from_numpy(S).to(DEV)).double()
+        en = ClosedFormEnergy(Sg, Vg, 1.0, 0.0, 64)
+        et = float(en.energy(torch.from_numpy(
+            np.ascontiguousarray(np.argsort(hid))).to(DEV)[None])[0])
+
+        def acc(c): return float((hid[c] == np.arange(n)).mean())
+
+        def refine(c):
+            f = fast_pair_descent(Sg, Vg,
+                                  torch.from_numpy(np.ascontiguousarray(c)).to(DEV),
+                                  600).cpu().numpy()
+            e = float(en.energy(torch.from_numpy(np.ascontiguousarray(f)).to(DEV)[None])[0])
+            return e, acc(f)
+
+        base = tlb(V, S)
+        shuf = np.random.default_rng(100 + trial).permutation(n)
+        variants = {
+            "tlb": base / base.max(),
+            "tlbperm": base[:, shuf] / base.max(),
+            "rand": np.random.default_rng(200 + trial).random((n, n)),
+        }
+        for vname, M in variants.items():
+            for mode in ("cold", "warm"):
+                G = np.outer(w, w); recs = []
+                for a in ALPHAS:
+                    g0 = None if mode == "cold" else G
+                    G2 = ot.gromov.fused_gromov_wasserstein(
+                        M, V, S, w, w, "square_loss", alpha=a, G0=g0,
+                        max_iter=200, tol_rel=1e-9)
+                    if mode == "warm":
+                        G = G2
+                    _, cols = linear_sum_assignment(-G2)
+                    e, ar = refine(cols)
+                    recs.append({"alpha": a, "raw": acc(cols), "refined": ar, "E": e})
+                best = min(recs, key=lambda x: x["E"])
+                row = {"trial": trial, "cost": vname, "mode": mode,
+                       "E_truth": et, "blind_refined": best["refined"],
+                       "blind_alpha": best["alpha"], "blind_E": best["E"],
+                       "oracle_refined": max(x["refined"] for x in recs),
+                       "curve": [round(x["refined"], 4) for x in recs]}
+                rows.append(row)
+                print(json.dumps(row), flush=True)
+    scratch = ("/tmp/claude-1273071/-home-yurenh2-emm/"
+               "bb97201d-d9b5-46e6-82fc-981d37e7ab98/scratchpad/")
+    json.dump(rows, open(scratch + "control.json", "w"))
+    print(json.dumps({"done": True}))
+
+
+if __name__ == "__main__":
+    main()
diff --git a/logs/fgw_alpha_control.log b/logs/fgw_alpha_control.log
new file mode 100644
index 0000000..de41a42
--- /dev/null
+++ b/logs/fgw_alpha_control.log
@@ -0,0 +1,4 @@
+{"trial": 0, "cost": "tlb", "mode": "cold", "E_truth": 0.339629709924146, "blind_refined": 0.0234375, "blind_alpha": 1.0, "blind_E": 0.6228801020056084, "oracle_refined": 0.08203125, "curve": [0.0781, 0.0234, 0.0117, 0.0273, 0.0742, 0.0352, 0.0664, 0.0352, 0.082, 0.0, 0.0547, 0.0234]}
+{"trial": 0, "cost": "tlb", "mode": "warm", "E_truth": 0.339629709924146, "blind_refined": 0.84765625, "blind_alpha": 0.2, "blind_E": 0.339629709924146, "oracle_refined": 0.84765625, "curve": [0.0781, 0.8281, 0.8477, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398, 0.8398]}
+{"trial": 0, "cost": "tlbperm", "mode": "cold", "E_truth": 0.339629709924146, "blind_refined": 0.02734375, "blind_alpha": 0.95, "blind_E": 0.6074469738905559, "oracle_refined": 0.078125, "curve": [0.0117, 0.0117, 0.043, 0.0781, 0.0078, 0.0352, 0.0156, 0.0156, 0.0039, 0.0078, 0.0273, 0.0234]}
+{"trial": 0, "cost": "tlbperm", "mode": "warm", "E_truth": 0.339629709924146, "blind_refined": 0.03125, "blind_alpha": 1.0, "blind_E": 0.7038378812415739, "oracle_refined": 0.04296875, "curve": [0.0117, 0.043, 0.0234, 0.0234, 0.043, 0.0391, 0.0391, 0.0391, 0.0312, 0.0312, 0.0312, 0.0312]}
diff --git a/logs/fgw_alpha_lean.log b/logs/fgw_alpha_lean.log
new file mode 100644
index 0000000..0787289
--- /dev/null
+++ b/logs/fgw_alpha_lean.log
@@ -0,0 +1,71 @@
+# omit-size trial0  E_truth=0.3396
+  cold a=0.0   raw=0.0625 ref=0.0781 E=0.6734  (8.4s)
+  cold a=0.2   raw=0.0117 ref=0.0117 E=0.6397  (9.5s)
+  cold a=0.4   raw=0.0703 ref=0.0742 E=0.6363  (2.8s)
+  cold a=0.6   raw=0.0586 ref=0.0664 E=0.6366  (1.2s)
+  cold a=0.75  raw=0.0508 ref=0.0430 E=0.6678  (5.7s)
+  cold a=0.85  raw=0.0312 ref=0.0234 E=0.6678  (1.6s)
+  cold a=0.88  raw=0.0586 ref=0.0508 E=0.6401  (2.4s)
+  cold a=0.9   raw=0.0000 ref=0.0000 E=0.6654  (2.2s)
+  cold a=0.92  raw=0.0352 ref=0.0312 E=0.7672  (0.5s)
+  cold a=0.94  raw=0.1680 ref=0.1875 E=0.6087  (0.5s)
+  cold a=0.95  raw=0.0625 ref=0.0547 E=0.6642  (2.1s)
+  cold a=0.96  raw=0.0664 ref=0.0664 E=0.6551  (0.3s)
+  cold a=0.97  raw=0.8203 ref=0.8203 E=0.3396  (0.4s)
+  cold a=0.98  raw=0.0039 ref=0.0039 E=0.6459  (1.0s)
+  cold a=0.99  raw=0.0156 ref=0.0156 E=0.6409  (1.0s)
+  cold a=1.0   raw=0.0195 ref=0.0234 E=0.6229  (0.7s)
+  == cold: blind-by-E refined=0.8203 (alpha=0.97), oracle-best=0.8203
+  warm a=0.0   raw=0.0625 ref=0.0781 E=0.6734  (0.6s)
+  warm a=0.2   raw=0.8242 ref=0.8320 E=0.3396  (1.2s)
+  warm a=0.4   raw=0.8281 ref=0.8320 E=0.3396  (0.4s)
+  warm a=0.6   raw=0.8320 ref=0.8320 E=0.3396  (0.3s)
+  warm a=0.75  raw=0.8320 ref=0.8320 E=0.3396  (0.2s)
+  warm a=0.85  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.88  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.9   raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.92  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.94  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.95  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.3s)
+  warm a=0.96  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.97  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.98  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.3s)
+  warm a=0.99  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=1.0   raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.0s)
+  == warm: blind-by-E refined=0.8320 (alpha=0.2), oracle-best=0.8320
+# synth-full trial0  E_truth=0.1430
+  cold a=0.0   raw=0.2578 ref=0.8633 E=0.2199  (2.1s)
+  cold a=0.2   raw=0.9531 ref=0.9609 E=0.1425  (2.6s)
+  cold a=0.4   raw=0.9531 ref=0.9688 E=0.1425  (1.7s)
+  cold a=0.6   raw=0.9570 ref=0.9766 E=0.1425  (2.2s)
+  cold a=0.75  raw=0.9375 ref=0.9609 E=0.1425  (3.0s)
+  cold a=0.85  raw=0.9531 ref=0.9609 E=0.1425  (4.3s)
+  cold a=0.88  raw=0.9453 ref=0.9609 E=0.1425  (3.8s)
+  cold a=0.9   raw=0.9766 ref=0.9844 E=0.1425  (3.5s)
+  cold a=0.92  raw=0.9453 ref=0.9531 E=0.1425  (0.9s)
+  cold a=0.94  raw=0.9531 ref=0.9766 E=0.1425  (2.0s)
+  cold a=0.95  raw=0.9453 ref=0.9531 E=0.1425  (1.1s)
+  cold a=0.96  raw=0.9453 ref=0.9531 E=0.1425  (0.8s)
+  cold a=0.97  raw=0.9375 ref=0.9531 E=0.1425  (6.2s)
+  cold a=0.98  raw=0.3789 ref=0.3789 E=0.5018  (7.1s)
+  cold a=0.99  raw=0.4219 ref=0.4141 E=0.4734  (10.8s)
+  cold a=1.0   raw=0.3633 ref=0.3906 E=0.5082  (12.6s)
+  == cold: blind-by-E refined=0.9609 (alpha=0.2), oracle-best=0.9844
+  warm a=0.0   raw=0.2578 ref=0.8633 E=0.2199  (4.5s)
+  warm a=0.2   raw=0.9375 ref=0.9453 E=0.1425  (1.9s)
+  warm a=0.4   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (2.5s)
+  warm a=0.6   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.3s)
+  warm a=0.75  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (1.4s)
+  warm a=0.85  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.9s)
+  warm a=0.88  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.6s)
+  warm a=0.9   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.8s)
+  warm a=0.92  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (1.2s)
+  warm a=0.94  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.8s)
+  warm a=0.95  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.4s)
+  warm a=0.96  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.6s)
+  warm a=0.97  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.7s)
+  warm a=0.98  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.4s)
+  warm a=0.99  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.1s)
+  warm a=1.0   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.2s)
+  == warm: blind-by-E refined=0.9453 (alpha=0.2), oracle-best=0.9453
+{"done": true}
diff --git a/logs/fgw_alpha_mech.log b/logs/fgw_alpha_mech.log
new file mode 100644
index 0000000..0c3f85f
--- /dev/null
+++ b/logs/fgw_alpha_mech.log
@@ -0,0 +1,12 @@
+--- trial 0  E_truth=0.339630
+  tlb-vertex     start=0.0586 | 1 power step ->0.0312 (E=0.7114) | 5 steps ->0.1719 (E=0.6409) | warm FGW ->0.8555 (E=0.3396)
+  random-vertex  start=0.0000 | 1 power step ->0.0312 (E=0.7279) | 5 steps ->0.0781 (E=0.6819) | warm FGW ->0.8125 (E=0.3396)
+  grampa-vertex  start=0.0156 | 1 power step ->0.0938 (E=0.7821) | 5 steps ->0.0469 (E=0.7049) | warm FGW ->0.8203 (E=0.3396)
+--- trial 1  E_truth=0.339630
+  tlb-vertex     start=0.0508 | 1 power step ->0.0352 (E=0.7114) | 5 steps ->0.1836 (E=0.6409) | warm FGW ->0.8164 (E=0.3396)
+  random-vertex  start=0.0039 | 1 power step ->0.0156 (E=0.7051) | 5 steps ->0.0195 (E=0.6692) | warm FGW ->0.0195 (E=0.6554)
+  grampa-vertex  start=0.0195 | 1 power step ->0.0742 (E=0.7821) | 5 steps ->0.0547 (E=0.7049) | warm FGW ->0.8164 (E=0.3396)
+--- trial 2  E_truth=0.339630
+  tlb-vertex     start=0.0430 | 1 power step ->0.0312 (E=0.7114) | 5 steps ->0.1797 (E=0.6409) | warm FGW ->0.8359 (E=0.3396)
+  random-vertex  start=0.0000 | 1 power step ->0.0195 (E=0.7149) | 5 steps ->0.0117 (E=0.6700) | warm FGW ->0.0391 (E=0.6387)
+  grampa-vertex  start=0.0156 | 1 power step ->0.0820 (E=0.7821) | 5 steps ->0.0508 (E=0.7049) | warm FGW ->0.8203 (E=0.3396)
diff --git a/logs/lean.log b/logs/lean.log
new file mode 100644
index 0000000..0787289
--- /dev/null
+++ b/logs/lean.log
@@ -0,0 +1,71 @@
+# omit-size trial0  E_truth=0.3396
+  cold a=0.0   raw=0.0625 ref=0.0781 E=0.6734  (8.4s)
+  cold a=0.2   raw=0.0117 ref=0.0117 E=0.6397  (9.5s)
+  cold a=0.4   raw=0.0703 ref=0.0742 E=0.6363  (2.8s)
+  cold a=0.6   raw=0.0586 ref=0.0664 E=0.6366  (1.2s)
+  cold a=0.75  raw=0.0508 ref=0.0430 E=0.6678  (5.7s)
+  cold a=0.85  raw=0.0312 ref=0.0234 E=0.6678  (1.6s)
+  cold a=0.88  raw=0.0586 ref=0.0508 E=0.6401  (2.4s)
+  cold a=0.9   raw=0.0000 ref=0.0000 E=0.6654  (2.2s)
+  cold a=0.92  raw=0.0352 ref=0.0312 E=0.7672  (0.5s)
+  cold a=0.94  raw=0.1680 ref=0.1875 E=0.6087  (0.5s)
+  cold a=0.95  raw=0.0625 ref=0.0547 E=0.6642  (2.1s)
+  cold a=0.96  raw=0.0664 ref=0.0664 E=0.6551  (0.3s)
+  cold a=0.97  raw=0.8203 ref=0.8203 E=0.3396  (0.4s)
+  cold a=0.98  raw=0.0039 ref=0.0039 E=0.6459  (1.0s)
+  cold a=0.99  raw=0.0156 ref=0.0156 E=0.6409  (1.0s)
+  cold a=1.0   raw=0.0195 ref=0.0234 E=0.6229  (0.7s)
+  == cold: blind-by-E refined=0.8203 (alpha=0.97), oracle-best=0.8203
+  warm a=0.0   raw=0.0625 ref=0.0781 E=0.6734  (0.6s)
+  warm a=0.2   raw=0.8242 ref=0.8320 E=0.3396  (1.2s)
+  warm a=0.4   raw=0.8281 ref=0.8320 E=0.3396  (0.4s)
+  warm a=0.6   raw=0.8320 ref=0.8320 E=0.3396  (0.3s)
+  warm a=0.75  raw=0.8320 ref=0.8320 E=0.3396  (0.2s)
+  warm a=0.85  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.88  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.9   raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.92  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.94  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.95  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.3s)
+  warm a=0.96  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.2s)
+  warm a=0.97  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=0.98  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.3s)
+  warm a=0.99  raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.1s)
+  warm a=1.0   raw=0.8320 ref=0.8320 E=0.3396 FROZEN (0.0s)
+  == warm: blind-by-E refined=0.8320 (alpha=0.2), oracle-best=0.8320
+# synth-full trial0  E_truth=0.1430
+  cold a=0.0   raw=0.2578 ref=0.8633 E=0.2199  (2.1s)
+  cold a=0.2   raw=0.9531 ref=0.9609 E=0.1425  (2.6s)
+  cold a=0.4   raw=0.9531 ref=0.9688 E=0.1425  (1.7s)
+  cold a=0.6   raw=0.9570 ref=0.9766 E=0.1425  (2.2s)
+  cold a=0.75  raw=0.9375 ref=0.9609 E=0.1425  (3.0s)
+  cold a=0.85  raw=0.9531 ref=0.9609 E=0.1425  (4.3s)
+  cold a=0.88  raw=0.9453 ref=0.9609 E=0.1425  (3.8s)
+  cold a=0.9   raw=0.9766 ref=0.9844 E=0.1425  (3.5s)
+  cold a=0.92  raw=0.9453 ref=0.9531 E=0.1425  (0.9s)
+  cold a=0.94  raw=0.9531 ref=0.9766 E=0.1425  (2.0s)
+  cold a=0.95  raw=0.9453 ref=0.9531 E=0.1425  (1.1s)
+  cold a=0.96  raw=0.9453 ref=0.9531 E=0.1425  (0.8s)
+  cold a=0.97  raw=0.9375 ref=0.9531 E=0.1425  (6.2s)
+  cold a=0.98  raw=0.3789 ref=0.3789 E=0.5018  (7.1s)
+  cold a=0.99  raw=0.4219 ref=0.4141 E=0.4734  (10.8s)
+  cold a=1.0   raw=0.3633 ref=0.3906 E=0.5082  (12.6s)
+  == cold: blind-by-E refined=0.9609 (alpha=0.2), oracle-best=0.9844
+  warm a=0.0   raw=0.2578 ref=0.8633 E=0.2199  (4.5s)
+  warm a=0.2   raw=0.9375 ref=0.9453 E=0.1425  (1.9s)
+  warm a=0.4   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (2.5s)
+  warm a=0.6   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.3s)
+  warm a=0.75  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (1.4s)
+  warm a=0.85  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.9s)
+  warm a=0.88  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.6s)
+  warm a=0.9   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.8s)
+  warm a=0.92  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (1.2s)
+  warm a=0.94  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.8s)
+  warm a=0.95  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.4s)
+  warm a=0.96  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.6s)
+  warm a=0.97  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.7s)
+  warm a=0.98  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.4s)
+  warm a=0.99  raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.1s)
+  warm a=1.0   raw=0.9375 ref=0.9453 E=0.1425 FROZEN (0.2s)
+  == warm: blind-by-E refined=0.9453 (alpha=0.2), oracle-best=0.9453
+{"done": true}
diff --git a/logs/lean_probe.py b/logs/lean_probe.py
new file mode 100644
index 0000000..f3b394c
--- /dev/null
+++ b/logs/lean_probe.py
@@ -0,0 +1,84 @@
+"""Lean version: 15-alpha grid, TLB cost, cold vs warm, omit-size then synth-full."""
+from __future__ import annotations
+import json, sys, time
+sys.path.insert(0, "/home/yurenh2/emm")
+import numpy as np, ot, torch
+from scipy.optimize import linear_sum_assignment
+from worldalign.synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from worldalign.synth_triangle_gate import standardized
+
+DEV = torch.device("cuda:1")
+ALPHAS = [0.0, 0.2, 0.4, 0.6, 0.75, 0.85, 0.88, 0.90, 0.92, 0.94, 0.95,
+          0.96, 0.97, 0.98, 0.99, 1.0]
+
+
+def standardise(m):
+    mask = ~np.eye(len(m), dtype=bool); v = m[mask]
+    o = (m - v.mean()) / v.std(); np.fill_diagonal(o, 0.0); return o
+
+
+def tlb(V, S):
+    n = len(V); a = np.sort(V, 1); b = np.sort(S, 1)
+    return (a * a).sum(1)[:, None] / n + (b * b).sum(1)[None, :] / n - 2 * (a @ b.T) / n
+
+
+def run(path, label, trial=0):
+    st = torch.load(path, map_location="cpu", weights_only=False)
+    V = standardise(st["visual_field"].double().numpy())
+    T = standardise(st["text_field"].double().numpy())
+    n = len(V); w = ot.unif(n)
+    rng = np.random.default_rng(trial); hid = rng.permutation(n)
+    S = T[np.ix_(hid, hid)]
+    Vg = standardized(torch.from_numpy(V).to(DEV)).double()
+    Sg = standardized(torch.from_numpy(S).to(DEV)).double()
+    en = ClosedFormEnergy(Sg, Vg, 1.0, 0.0, 64)
+    M = tlb(V, S); M = M / M.max()
+
+    def acc(c): return float((hid[c] == np.arange(n)).mean())
+
+    def refine(c):
+        f = fast_pair_descent(Sg, Vg, torch.from_numpy(np.ascontiguousarray(c)).to(DEV),
+                              600).cpu().numpy()
+        e = float(en.energy(torch.from_numpy(np.ascontiguousarray(f)).to(DEV)[None])[0])
+        return e, acc(f)
+
+    et = float(en.energy(torch.from_numpy(np.ascontiguousarray(np.argsort(hid))).to(DEV)[None])[0])
+    print(f"# {label} trial{trial}  E_truth={et:.4f}", flush=True)
+
+    out = {"label": label, "E_truth": et, "cold": [], "warm": []}
+    for mode in ("cold", "warm"):
+        G = np.outer(w, w)
+        prev = None
+        for a in ALPHAS:
+            t0 = time.time()
+            g0 = None if mode == "cold" else G
+            G2 = ot.gromov.fused_gromov_wasserstein(M, V, S, w, w, "square_loss",
+                                                    alpha=a, G0=g0, max_iter=200,
+                                                    tol_rel=1e-9)
+            if mode == "warm":
+                G = G2
+            _, cols = linear_sum_assignment(-G2)
+            e, ar = refine(cols)
+            same = (prev is not None and np.array_equal(cols, prev))
+            prev = cols
+            rec = {"alpha": a, "raw": acc(cols), "refined": ar, "E": e,
+                   "sec": round(time.time() - t0, 1), "same_as_prev": bool(same)}
+            out[mode].append(rec)
+            print(f"  {mode} a={a:<5} raw={rec['raw']:.4f} ref={ar:.4f} E={e:.4f}"
+                  f" {'FROZEN' if same else ''} ({rec['sec']}s)", flush=True)
+        good = out[mode]
+        best = min(good, key=lambda x: x["E"])
+        print(f"  == {mode}: blind-by-E refined={best['refined']:.4f} (alpha={best['alpha']}),"
+              f" oracle-best={max(x['refined'] for x in good):.4f}", flush=True)
+    return out
+
+
+if __name__ == "__main__":
+    base = "/home/yurenh2/emm/artifacts/synth_v1/"
+    res = []
+    for p, l in [("omit_size.pt", "omit-size"), ("omit_none.pt", "synth-full")]:
+        res.append(run(base + p, l))
+    scratch = ("/tmp/claude-1273071/-home-yurenh2-emm/"
+               "bb97201d-d9b5-46e6-82fc-981d37e7ab98/scratchpad/")
+    json.dump(res, open(scratch + "lean.json", "w"))
+    print(json.dumps({"done": True}))
diff --git a/logs/mech.log b/logs/mech.log
new file mode 100644
index 0000000..817501a
--- /dev/null
+++ b/logs/mech.log
@@ -0,0 +1,11 @@
+--- trial 0  E_truth=0.339630
+  tlb-vertex     start=0.0586 | 1 power step ->0.0312 (E=0.7114) | 5 steps ->0.1719 (E=0.6409) | warm FGW ->0.8555 (E=0.3396)
+  random-vertex  start=0.0000 | 1 power step ->0.0312 (E=0.7279) | 5 steps ->0.0781 (E=0.6819) | warm FGW ->0.8125 (E=0.3396)
+  grampa-vertex  start=0.0156 | 1 power step ->0.0938 (E=0.7821) | 5 steps ->0.0469 (E=0.7049) | warm FGW ->0.8203 (E=0.3396)
+--- trial 1  E_truth=0.339630
+  tlb-vertex     start=0.0508 | 1 power step ->0.0352 (E=0.7114) | 5 steps ->0.1836 (E=0.6409) | warm FGW ->0.8164 (E=0.3396)
+  random-vertex  start=0.0039 | 1 power step ->0.0156 (E=0.7051) | 5 steps ->0.0195 (E=0.6692) | warm FGW ->0.0195 (E=0.6554)
+  grampa-vertex  start=0.0195 | 1 power step ->0.0742 (E=0.7821) | 5 steps ->0.0547 (E=0.7049) | warm FGW ->0.8164 (E=0.3396)
+--- trial 2  E_truth=0.339630
+  tlb-vertex     start=0.0430 | 1 power step ->0.0312 (E=0.7114) | 5 steps ->0.1797 (E=0.6409) | warm FGW ->0.8359 (E=0.3396)
+  random-vertex  start=0.0000 | 1 power step ->0.0195 (E=0.7149) | 5 steps ->0.0117 (E=0.6700) | warm FGW ->0.0391 (E=0.6387)
diff --git a/logs/mech.py b/logs/mech.py
new file mode 100644
index 0000000..6762777
--- /dev/null
+++ b/logs/mech.py
@@ -0,0 +1,51 @@
+import sys, numpy as np, ot, torch
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from worldalign.synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from worldalign.synth_triangle_gate import standardized
+DEV=torch.device("cuda:1")
+def standardise(m):
+    mask=~np.eye(len(m),dtype=bool); v=m[mask]
+    o=(m-v.mean())/v.std(); np.fill_diagonal(o,0.0); return o
+def tlb(V,S):
+    n=len(V); a=np.sort(V,1); b=np.sort(S,1)
+    return (a*a).sum(1)[:,None]/n+(b*b).sum(1)[None,:]/n-2*(a@b.T)/n
+st=torch.load("/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt",map_location='cpu',weights_only=False)
+V=standardise(st["visual_field"].double().numpy()); T=standardise(st["text_field"].double().numpy())
+n=len(V); w=ot.unif(n)
+for trial in range(3):
+    rng=np.random.default_rng(trial); hid=rng.permutation(n); S=T[np.ix_(hid,hid)]
+    Vg=standardized(torch.from_numpy(V).to(DEV)).double(); Sg=standardized(torch.from_numpy(S).to(DEV)).double()
+    en=ClosedFormEnergy(Sg,Vg,1.0,0.0,64)
+    Et=float(en.energy(torch.from_numpy(np.ascontiguousarray(np.argsort(hid))).to(DEV)[None])[0])
+    acc=lambda c: float((hid[c]==np.arange(n)).mean())
+    def refine(c,it=2000):
+        f=fast_pair_descent(Sg,Vg,torch.from_numpy(np.ascontiguousarray(c)).to(DEV),it).cpu().numpy()
+        return float(en.energy(torch.from_numpy(np.ascontiguousarray(f)).to(DEV)[None])[0]), acc(f)
+    M=tlb(V,S); M=M/M.max()
+    _,tlbcols=linear_sum_assignment(M)
+    def vert(cols):
+        P=np.zeros((n,n)); P[np.arange(n),cols]=1.0/n; return P
+    starts={"tlb-vertex":tlbcols,
+            "random-vertex":np.random.default_rng(50+trial).permutation(n),
+            "grampa-vertex":None}
+    from worldalign.spectral_match import grampa
+    starts["grampa-vertex"]=grampa(V,S,1.0)
+    print(f"--- trial {trial}  E_truth={Et:.6f}")
+    for name,cols in starts.items():
+        # (i) pure power step: LSAP on gradient V P S  (anchor-propagation, one round)
+        P=vert(cols); grad=V@P@S
+        _,c1=linear_sum_assignment(-grad)
+        e1,a1=refine(c1)
+        # (ii) 5 rounds of power step
+        c=cols.copy()
+        for _ in range(5):
+            P=vert(c); _,c=linear_sum_assignment(-(V@P@S))
+        e2,a2=refine(c)
+        # (iii) warm FGW at alpha=0.1 then 0.2 from this vertex
+        G=vert(cols)
+        for a in (0.1,0.2,0.5):
+            G=ot.gromov.fused_gromov_wasserstein(M,V,S,w,w,"square_loss",alpha=a,G0=G,max_iter=200,tol_rel=1e-9)
+        _,c3=linear_sum_assignment(-G); e3,a3=refine(c3)
+        print("  %-14s start=%.4f | 1 power step ->%.4f (E=%.4f) | 5 steps ->%.4f (E=%.4f) | warm FGW ->%.4f (E=%.4f)"
+              %(name,acc(cols),a1,e1,a2,e2,a3,e3), flush=True)
diff --git a/logs/verify.py b/logs/verify.py
new file mode 100644
index 0000000..17fd643
--- /dev/null
+++ b/logs/verify.py
@@ -0,0 +1,59 @@
+import sys, numpy as np, ot, torch
+sys.path.insert(0,'/home/yurenh2/emm')
+from scipy.optimize import linear_sum_assignment
+from worldalign.synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from worldalign.synth_triangle_gate import standardized
+DEV=torch.device("cuda:1")
+def standardise(m):
+    mask=~np.eye(len(m),dtype=bool); v=m[mask]
+    o=(m-v.mean())/v.std(); np.fill_diagonal(o,0.0); return o
+def tlb(V,S):
+    n=len(V); a=np.sort(V,1); b=np.sort(S,1)
+    return (a*a).sum(1)[:,None]/n+(b*b).sum(1)[None,:]/n-2*(a@b.T)/n
+st=torch.load("/home/yurenh2/emm/artifacts/synth_v1/omit_size.pt",map_location='cpu',weights_only=False)
+V=standardise(st["visual_field"].double().numpy()); T=standardise(st["text_field"].double().numpy())
+n=len(V); w=ot.unif(n)
+rng=np.random.default_rng(0); hid=rng.permutation(n); S=T[np.ix_(hid,hid)]
+Vg=standardized(torch.from_numpy(V).to(DEV)).double(); Sg=standardized(torch.from_numpy(S).to(DEV)).double()
+en=ClosedFormEnergy(Sg,Vg,1.0,0.0,64)
+M=tlb(V,S); M=M/M.max()
+G=np.outer(w,w)
+for a in (0.0,0.2):
+    G=ot.gromov.fused_gromov_wasserstein(M,V,S,w,w,"square_loss",alpha=a,G0=G,max_iter=200,tol_rel=1e-9)
+_,cols=linear_sum_assignment(-G)
+fin=fast_pair_descent(Sg,Vg,torch.from_numpy(np.ascontiguousarray(cols)).to(DEV),2000).cpu().numpy()
+truth=np.argsort(hid)
+Ef=float(en.energy(torch.from_numpy(np.ascontiguousarray(fin)).to(DEV)[None])[0])
+Et=float(en.energy(torch.from_numpy(np.ascontiguousarray(truth)).to(DEV)[None])[0])
+acc=float((hid[fin]==np.arange(n)).mean())
+print("accuracy=%.4f  E_found=%.10f  E_truth=%.10f  diff=%.3e"%(acc,Ef,Et,Ef-Et))
+wrong=np.where(hid[fin]!=np.arange(n))[0]
+print("n_wrong=%d"%len(wrong))
+# for each wrong scene i: correlation between S-row of assigned node and S-row of true node
+Sal=S  # in shuffled index space
+corr=[]
+vcorr=[]
+for i in wrong:
+    a_idx=fin[i]; t_idx=truth[i]
+    x=np.delete(Sal[a_idx],[a_idx,t_idx]); y=np.delete(Sal[t_idx],[a_idx,t_idx])
+    corr.append(np.corrcoef(x,y)[0,1])
+    u=np.delete(V[i],[i]); 
+    # visual side: is the visually-assigned scene similar to i? compare V rows of i and of the scene truly at a_idx
+    j=hid[a_idx]
+    v1=np.delete(V[i],[i,j]); v2=np.delete(V[j],[i,j])
+    vcorr.append(np.corrcoef(v1,v2)[0,1])
+rng2=np.random.default_rng(7); pairs=rng2.integers(0,n,(len(wrong),2))
+base=[np.corrcoef(np.delete(Sal[p],[p,q]),np.delete(Sal[q],[p,q]))[0,1] for p,q in pairs if p!=q]
+print("text-row corr of confused partners: median=%.3f mean=%.3f  | random pair baseline median=%.3f"%(np.median(corr),np.mean(corr),np.median(base)))
+print("visual-row corr of confused partners: median=%.3f"%np.median(vcorr))
+# cycle structure of the error
+perm=hid[fin]
+seen=set(); cyc=[]
+for i in range(n):
+    if i in seen or perm[i]==i: continue
+    c=0; j=i
+    while j not in seen:
+        seen.add(j); j=perm[j]; c+=1
+    cyc.append(c)
+import collections
+print("error cycle lengths:", dict(collections.Counter(cyc)))