The failure is the optimiser, not the information: a 257x faster descent and a benchmark

The gate settles which failure mode each field is in, and refutes the
symmetry hypothesis I proposed. On the caption-omitted field the truth is
a STRICT local minimum -- descent started at the truth does not move at
all -- the anchor bound says the information is 99.7% intact, and our
solver stops 0.44 above it at 4.9% accuracy. That is a pure optimiser
failure. Natural data is the opposite: descent from the truth falls a
further 0.167, so the truth is not even locally optimal, which is the
information-deficit signature the 0.291 bound predicted.

Steepest descent was brute-forcing all 32,640 candidate permutations
through the full energy every step, including a batched cube trace with
the triangle term active -- 203 seconds per descent, which is why the
gates were hopeless. The pairwise term needs one matrix product for the
whole table: swapping p,q changes the alignment sum by
2(C_pq + C_qp - C_pp - C_qq + 2 A_pq B_pq) with C = A @ B. Verified
against brute force to 1e-9 before use, and the fast descent reaches the
same optimum. 203s -> 0.79s.

Adds a matching benchmark with known-reachable answers and the solver
families never tried on these fields: Gromov-Wasserstein, entropic GW
with an annealed regulariser, BAPG.

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

12 files changed, 1060 insertions, 0 deletions

diff --git a/artifacts/synth_v1/degeneracy_gate.json b/artifacts/synth_v1/degeneracy_gate.json
new file mode 100644
index 0000000..6a910d0
--- /dev/null
+++ b/artifacts/synth_v1/degeneracy_gate.json
@@ -0,0 +1,45 @@
+{
+  "protocol": "One spectral-start descent per trial, plus one descent started at the truth. A positive gap means the searcher stopped above the truth; zero means it found something the energy cannot distinguish from the truth; negative means the truth is beaten.",
+  "rows": [
+    {
+      "label": "omit-size",
+      "size": 256,
+      "energy_gap_found_minus_truth": 0.43996140951295937,
+      "truth_descends_further_by": 0.0,
+      "blind_accuracy": 0.049479166666666664,
+      "verdict": "truth is DEEPER than what is found -- optimiser failure"
+    },
+    {
+      "label": "noise035",
+      "size": 256,
+      "energy_gap_found_minus_truth": 0.04761915721057547,
+      "truth_descends_further_by": 0.0020415746981433114,
+      "blind_accuracy": 0.88671875,
+      "verdict": "truth is DEEPER than what is found -- optimiser failure"
+    },
+    {
+      "label": "rank8",
+      "size": 256,
+      "energy_gap_found_minus_truth": 0.364986204015837,
+      "truth_descends_further_by": 0.006232799331109179,
+      "blind_accuracy": 0.06510416666666667,
+      "verdict": "truth is DEEPER than what is found -- optimiser failure"
+    },
+    {
+      "label": "synth-full",
+      "size": 256,
+      "energy_gap_found_minus_truth": 0.07695912511043244,
+      "truth_descends_further_by": 0.0005474964148536876,
+      "blind_accuracy": 0.8828125,
+      "verdict": "truth is DEEPER than what is found -- optimiser failure"
+    },
+    {
+      "label": "natural-best",
+      "size": 256,
+      "energy_gap_found_minus_truth": 0.08529010447794294,
+      "truth_descends_further_by": 0.16736025860425965,
+      "blind_accuracy": 0.02734375,
+      "verdict": "truth is DEEPER than what is found -- optimiser failure"
+    }
+  ]
+}
\ No newline at end of file
diff --git a/artifacts/vg_5k/anchor_bound_reprice.json b/artifacts/vg_5k/anchor_bound_reprice.json
new file mode 100644
index 0000000..4571545
--- /dev/null
+++ b/artifacts/vg_5k/anchor_bound_reprice.json
@@ -0,0 +1,70 @@
+{
+  "protocol": "Half the scenes are anchors with their correspondence given. The rest are described by their field rows against the anchors and matched by Hungarian assignment. Upper-bounds blind recovery, which must additionally discover the anchor correspondence.",
+  "rows": [
+    {
+      "label": "baseline-keep8",
+      "correlation": 0.6558829559111277,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.1296875,
+      "anchor_bound": 0.1953125
+    },
+    {
+      "label": "base-6seg",
+      "correlation": 0.7309253746041517,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.1671875,
+      "anchor_bound": 0.290625
+    },
+    {
+      "label": "large-6seg",
+      "correlation": 0.7089623334444304,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.1609375,
+      "anchor_bound": 0.284375
+    },
+    {
+      "label": "base-16seg",
+      "correlation": 0.6794444345566255,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.15625,
+      "anchor_bound": 0.24375
+    },
+    {
+      "label": "oracle-boxes",
+      "correlation": 0.7027488907202334,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.165625,
+      "anchor_bound": 0.2734375
+    },
+    {
+      "label": "structured-phrase",
+      "correlation": 0.7252267991153302,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.19375,
+      "anchor_bound": 0.321875
+    },
+    {
+      "label": "third-order",
+      "correlation": 0.7171326909156865,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.196875,
+      "anchor_bound": 0.325
+    },
+    {
+      "label": "bow-text",
+      "correlation": 0.4699873859667931,
+      "probe_scenes": 128,
+      "chance": 0.0078125,
+      "nearest_neighbour": 0.14375,
+      "anchor_bound": 0.3296875
+    }
+  ],
+  "reading": "A low bound means the representation does not identify scenes and no solver will help. A high bound with low blind recovery means the information is present and the search or the landscape is at fault."
+}
\ No newline at end of file
diff --git a/artifacts/vg_5k/field_geometry.json b/artifacts/vg_5k/field_geometry.json
new file mode 100644
index 0000000..7d36682
--- /dev/null
+++ b/artifacts/vg_5k/field_geometry.json
@@ -0,0 +1,54 @@
+{
+  "protocol": "Each geometry is applied to each modality independently, so no transform crosses the gauge or uses the pairing. The bound is the anchor upper bound; the last column re-computes it with a rank comparison operator instead of a linear one.",
+  "label": "natural-best",
+  "rows": [
+    {
+      "geometry": "euclidean (current)",
+      "correlation": 0.7157944505313568,
+      "anchor_bound": 0.3359375,
+      "anchor_bound_rank_operator": 0.1453125
+    },
+    {
+      "geometry": "global rank (copula)",
+      "correlation": 0.6453667782045476,
+      "anchor_bound": 0.171875,
+      "anchor_bound_rank_operator": 0.1453125
+    },
+    {
+      "geometry": "row rank",
+      "correlation": 0.6374571651517471,
+      "anchor_bound": 0.165625,
+      "anchor_bound_rank_operator": 0.1484375
+    },
+    {
+      "geometry": "local scaling",
+      "correlation": 0.7162509779848556,
+      "anchor_bound": 0.3625,
+      "anchor_bound_rank_operator": 0.1703125
+    },
+    {
+      "geometry": "diffusion t=3",
+      "correlation": 0.5659469014331384,
+      "anchor_bound": 0.0546875,
+      "anchor_bound_rank_operator": 0.0453125
+    },
+    {
+      "geometry": "heat kernel",
+      "correlation": 0.7381033036202107,
+      "anchor_bound": 0.31875,
+      "anchor_bound_rank_operator": 0.1546875
+    },
+    {
+      "geometry": "hyperbolic d=16",
+      "correlation": 0.7135662664726405,
+      "anchor_bound": 0.246875,
+      "anchor_bound_rank_operator": 0.1703125
+    }
+  ],
+  "best": {
+    "geometry": "local scaling",
+    "correlation": 0.7162509779848556,
+    "anchor_bound": 0.3625,
+    "anchor_bound_rank_operator": 0.1703125
+  }
+}
\ No newline at end of file
diff --git a/artifacts/vg_5k/split_half.json b/artifacts/vg_5k/split_half.json
new file mode 100644
index 0000000..5140582
--- /dev/null
+++ b/artifacts/vg_5k/split_half.json
@@ -0,0 +1,18 @@
+{
+  "protocol": "Each scene's parts split into disjoint halves within one modality; a relation field is built from each half and the two are matched by the anchor bound. Both halves describe the same scenes, so this upper-bounds any cross-modal result.",
+  "results": {
+    "text": {
+      "scenes": 256,
+      "mean_parts": 15.9765625,
+      "split_half_correlation": 0.7933697484313135,
+      "split_half_anchor_bound": 0.490625
+    },
+    "vision": {
+      "scenes": 256,
+      "mean_parts": 6.0,
+      "split_half_correlation": 0.9836462733134468,
+      "split_half_anchor_bound": 0.99375
+    }
+  },
+  "cross_modal_bound_for_reference": 0.291
+}
\ No newline at end of file
diff --git a/logs/battery.log b/logs/battery.log
new file mode 100644
index 0000000..6030612
--- /dev/null
+++ b/logs/battery.log
@@ -0,0 +1,7 @@
+
+=== omit-size-BOUND0.997  (N=256, chance=0.0039)
+  grampa (current)       raw=0.0169  refined=0.0443  (0.0s)
+  umeyama                raw=0.0052  refined=0.0143  (0.0s)
+  gromov-wasserstein     raw=0.0182  refined=0.0221  (0.1s)
+/home/yurenh2/.local/lib/python3.13/site-packages/ot/bregman/_sinkhorn.py:666: UserWarning: Sinkhorn did not converge. You might want to increase the number of iterations `numItermax` or the regularization parameter `reg`.
+  warnings.warn(
diff --git a/logs/degeneracy.log b/logs/degeneracy.log
new file mode 100644
index 0000000..ff31564
--- /dev/null
+++ b/logs/degeneracy.log
@@ -0,0 +1,11 @@
+omit-size            gap(found-truth)=+0.4400  truth-descends=+0.0000  acc=0.049
+                     -> truth is DEEPER than what is found -- optimiser failure
+noise035             gap(found-truth)=+0.0476  truth-descends=+0.0020  acc=0.887
+                     -> truth is DEEPER than what is found -- optimiser failure
+rank8                gap(found-truth)=+0.3650  truth-descends=+0.0062  acc=0.065
+                     -> truth is DEEPER than what is found -- optimiser failure
+synth-full           gap(found-truth)=+0.0770  truth-descends=+0.0005  acc=0.883
+                     -> truth is DEEPER than what is found -- optimiser failure
+natural-best         gap(found-truth)=+0.0853  truth-descends=+0.1674  acc=0.027
+                     -> truth is DEEPER than what is found -- optimiser failure
+{"done": true}
diff --git a/tests/test_core.py b/tests/test_core.py
index 09e58f1..429eed0 100644
--- a/tests/test_core.py
+++ b/tests/test_core.py
@@ -323,3 +323,60 @@ def test_structured_phrase_encoding_separates_head_from_modifier():
         np.mean([vectors[t] for t in ["red", "bus"]], axis=0),
         np.mean([vectors[t] for t in ["bus", "red"]], axis=0),
     )
+
+
+def test_closed_form_pair_swap_deltas_match_brute_force():
+    """Every transposition's energy change, against the definition."""
+    import numpy as np
+    from worldalign.synth_fast_gate import (
+        ClosedFormEnergy, all_pair_swap_deltas, all_swaps, apply_swaps,
+    )
+
+    size = 24
+    generator = torch.Generator().manual_seed(0)
+    def field():
+        raw = torch.randn(size, size, generator=generator, dtype=torch.float64)
+        out = (raw + raw.T) / 2
+        out.fill_diagonal_(0.0)
+        return out
+
+    text, visual = field(), field()
+    energy = ClosedFormEnergy(text, visual, 1.0, 0.0, 64)
+    permutation = torch.randperm(size, generator=generator)
+
+    base = float(energy.energy(permutation[None])[0])
+    swaps = all_swaps(size, torch.device("cpu"))
+    brute = energy.energy(apply_swaps(permutation, swaps)) - base
+
+    permuted = text[permutation[:, None], permutation[None, :]]
+    table = all_pair_swap_deltas(permuted, visual)
+    closed = table[swaps[:, 0], swaps[:, 1]]
+
+    assert torch.allclose(closed, brute, atol=1e-9), (
+        closed[:4].tolist(), brute[:4].tolist()
+    )
+
+
+def test_fast_pair_descent_matches_brute_force_descent():
+    """The fast descent must reach the same energy as the slow one."""
+    from worldalign.synth_fast_gate import (
+        ClosedFormEnergy, all_swaps, fast_pair_descent, steepest_descent,
+    )
+
+    size = 24
+    generator = torch.Generator().manual_seed(1)
+    def field():
+        raw = torch.randn(size, size, generator=generator, dtype=torch.float64)
+        out = (raw + raw.T) / 2
+        out.fill_diagonal_(0.0)
+        return out
+
+    text, visual = field(), field()
+    energy = ClosedFormEnergy(text, visual, 1.0, 0.0, 64)
+    start = torch.randperm(size, generator=generator)
+    slow, slow_value = steepest_descent(
+        energy, start.clone(), all_swaps(size, torch.device("cpu")), 500
+    )
+    fast = fast_pair_descent(text, visual, start.clone(), 500)
+    fast_value = float(energy.energy(fast[None])[0])
+    assert abs(fast_value - slow_value) < 1e-9, (fast_value, slow_value)
diff --git a/worldalign/degeneracy_gate.py b/worldalign/degeneracy_gate.py
new file mode 100644
index 0000000..4ba3dc1
--- /dev/null
+++ b/worldalign/degeneracy_gate.py
@@ -0,0 +1,130 @@
+"""Is the truth the optimum, tied with it, or beaten? Three descents, not seven hundred.
+
+The question that decides the next move is narrow, and the first attempt at it
+was not. Suppressing a caption factor leaves the information intact -- the
+anchor bound is 0.997 -- while blind matching returns 5.6%, and there are three
+possible reasons with three different responses:
+
+    E(truth) <  E(best found)   the searcher simply missed it   -> better optimiser
+    E(truth) == E(best found)   the truth is one of many minima  -> quotient the symmetry
+    E(truth) >  E(best found)   the truth is not the optimum     -> information is gone
+
+Answering it does not need a strong searcher. It needs the energy of the truth
+and the energy of whatever the standard pipeline actually converges to, which
+is one descent per trial rather than the sixty-one the first version ran -- at
+N=256 a single descent over all 32,640 transpositions costs minutes, so the
+difference is hours against seconds.
+"""
+
+from __future__ import annotations
+
+import argparse
+import json
+
+import numpy as np
+import torch
+
+from .common import write_json
+from .spectral_match import grampa
+from .synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from .synth_triangle_gate import standardized
+
+
+def parse_args() -> argparse.Namespace:
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--fields", nargs="+", required=True)
+    parser.add_argument("--labels", nargs="+", default=None)
+    parser.add_argument("--trials", type=int, default=3)
+    parser.add_argument("--iterations", type=int, default=600)
+    parser.add_argument("--device", default="cuda:3")
+    parser.add_argument("--tolerance", type=float, default=1e-3)
+    parser.add_argument("--output", default="artifacts/synth_v1/degeneracy_gate.json")
+    return parser.parse_args()
+
+
+def standardise(matrix: np.ndarray) -> np.ndarray:
+    mask = ~np.eye(len(matrix), dtype=bool)
+    values = matrix[mask]
+    out = (matrix - values.mean()) / values.std()
+    np.fill_diagonal(out, 0.0)
+    return out
+
+
+def analyse(path: str, label: str, args: argparse.Namespace) -> dict:
+    state = torch.load(path, map_location="cpu", weights_only=False)
+    visual = standardise(state["visual_field"].double().numpy())
+    text = standardise(state["text_field"].double().numpy())
+    size = len(visual)
+    device = torch.device(args.device)
+
+    gaps, accuracies, truths = [], [], []
+    for trial in range(args.trials):
+        generator = np.random.default_rng(trial)
+        hidden = generator.permutation(size)
+        shuffled = text[np.ix_(hidden, hidden)]
+        shuffled_gpu = standardized(torch.from_numpy(shuffled).to(device)).double()
+        visual_gpu = standardized(torch.from_numpy(visual).to(device)).double()
+        energy = ClosedFormEnergy(shuffled_gpu, visual_gpu, 1.0, 0.0, 256)
+        truth = torch.from_numpy(np.argsort(hidden).copy()).to(device)
+        truth_energy = float(energy.energy(truth[None])[0])
+
+        start = torch.from_numpy(grampa(visual, shuffled, 1.0).copy()).to(device)
+        final = fast_pair_descent(shuffled_gpu, visual_gpu, start, args.iterations)
+        found_energy = float(energy.energy(final[None])[0])
+        # descent from the truth itself: if it moves, the truth is not a local
+        # minimum either, which separates "not optimal" from "merely not found"
+        from_truth = fast_pair_descent(
+            shuffled_gpu, visual_gpu, truth.clone(), args.iterations
+        )
+        truth_local = float(energy.energy(from_truth[None])[0])
+
+        gaps.append(found_energy - truth_energy)
+        truths.append(truth_energy - truth_local)
+        accuracies.append(float((hidden[final.cpu().numpy()] == np.arange(size)).mean()))
+
+    gap = float(np.mean(gaps))
+    scale = abs(float(np.mean([g for g in gaps]))) + 1.0
+    if gap < -args.tolerance:
+        verdict = "truth is NOT the optimum -- information deficit"
+    elif abs(gap) <= args.tolerance:
+        verdict = "truth is TIED with what is found -- degenerate, quotient the symmetry"
+    else:
+        verdict = "truth is DEEPER than what is found -- optimiser failure"
+    return {
+        "label": label,
+        "size": size,
+        "energy_gap_found_minus_truth": gap,
+        "truth_descends_further_by": float(np.mean(truths)),
+        "blind_accuracy": float(np.mean(accuracies)),
+        "verdict": verdict,
+    }
+
+
+def main() -> None:
+    args = parse_args()
+    labels = args.labels or [p.split("/")[-1] for p in args.fields]
+    rows = []
+    for path, label in zip(args.fields, labels):
+        row = analyse(path, label, args)
+        rows.append(row)
+        print(
+            f"{row['label']:<20} gap(found-truth)={row['energy_gap_found_minus_truth']:+.4f}  "
+            f"truth-descends={row['truth_descends_further_by']:+.4f}  "
+            f"acc={row['blind_accuracy']:.3f}\n"
+            f"{'':<20} -> {row['verdict']}",
+            flush=True,
+        )
+    write_json(args.output, {
+        "protocol": (
+            "One spectral-start descent per trial, plus one descent started at "
+            "the truth. A positive gap means the searcher stopped above the "
+            "truth; zero means it found something the energy cannot "
+            "distinguish from the truth; negative means the truth is beaten."
+        ),
+        "rows": rows,
+    })
+    print(json.dumps({"done": True}))
+
+
+if __name__ == "__main__":
+    main()
diff --git a/worldalign/field_geometry.py b/worldalign/field_geometry.py
new file mode 100644
index 0000000..9360427
--- /dev/null
+++ b/worldalign/field_geometry.py
@@ -0,0 +1,222 @@
+"""The space and the operators, which we never varied.
+
+Everything measured so far fixed one geometry without saying so. Scene states
+are compared by Euclidean inner product, the relation field holds those raw
+values, and two fields are compared by Pearson correlation -- three
+commitments to flatness made by default rather than by argument. The
+measured deficit is that the two modalities disagree about which scenes are
+similar, and "similar" is exactly what those three choices define.
+
+The mismatch that costs the most is likely the cheapest to remove. Two
+independently trained encoders have no reason to put their similarities on
+the same scale: if vision similarity is any monotone but nonlinear function
+of text similarity, the two fields describe identical relational structure
+and Pearson correlation still reports disagreement. A rank transform is
+invariant to every monotone reparametrisation and costs one sort.
+
+Density is the second. A scene sitting in a crowded region of one modality
+and a sparse region of the other has systematically shifted similarities to
+everything, which local scaling and the diffusion kernel both remove -- the
+latter by replacing similarity with how probability spreads, a quantity
+defined by the graph rather than by whatever units the encoder happened to
+emit.
+
+Each transform is applied to each modality independently, so nothing here
+uses the pairing or crosses the gauge.
+"""
+
+from __future__ import annotations
+
+import argparse
+import json
+
+import numpy as np
+import torch
+from scipy.optimize import linear_sum_assignment
+from scipy.stats import rankdata
+
+from .common import write_json
+
+
+def parse_args() -> argparse.Namespace:
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--fields", default="artifacts/vg_5k/np_n256.pt")
+    parser.add_argument("--label", default="natural-best")
+    parser.add_argument("--repeats", type=int, default=5)
+    parser.add_argument("--output", default="artifacts/vg_5k/field_geometry.json")
+    return parser.parse_args()
+
+
+def offdiagonal_mask(size: int) -> np.ndarray:
+    return ~np.eye(size, dtype=bool)
+
+
+def standardise(matrix: np.ndarray) -> np.ndarray:
+    mask = offdiagonal_mask(len(matrix))
+    values = matrix[mask]
+    out = (matrix - values.mean()) / values.std()
+    np.fill_diagonal(out, 0.0)
+    return out
+
+
+# --- geometries, each applied within one modality -------------------------
+
+def identity(field: np.ndarray) -> np.ndarray:
+    return standardise(field)
+
+
+def global_rank(field: np.ndarray) -> np.ndarray:
+    """Copula transform: keep the ordering of similarities, discard the scale."""
+    size = len(field)
+    mask = offdiagonal_mask(size)
+    out = np.zeros_like(field)
+    out[mask] = rankdata(field[mask]) / mask.sum()
+    return standardise(out)
+
+
+def row_rank(field: np.ndarray) -> np.ndarray:
+    """Rank within each scene's own neighbourhood, removing per-scene scale."""
+    size = len(field)
+    out = np.zeros_like(field)
+    for row in range(size):
+        others = np.delete(np.arange(size), row)
+        out[row, others] = rankdata(field[row, others]) / len(others)
+    out = (out + out.T) / 2.0
+    return standardise(out)
+
+
+def local_scaling(field: np.ndarray, neighbours: int = 10) -> np.ndarray:
+    """Self-tuning: divide by each scene's own neighbourhood scale (CSLS-like)."""
+    size = len(field)
+    work = field.copy()
+    np.fill_diagonal(work, -np.inf)
+    order = np.sort(work, axis=1)[:, ::-1]
+    scale = order[:, :neighbours].mean(1)
+    out = field - 0.5 * (scale[:, None] + scale[None, :])
+    np.fill_diagonal(out, 0.0)
+    return standardise(out)
+
+
+def diffusion(field: np.ndarray, steps: int = 3) -> np.ndarray:
+    """Replace similarity by how probability spreads on the similarity graph."""
+    size = len(field)
+    affinity = field - field.min()
+    np.fill_diagonal(affinity, 0.0)
+    degree = affinity.sum(1).clip(1e-9)
+    walk = affinity / degree[:, None]
+    powered = np.linalg.matrix_power(walk, steps)
+    out = powered @ powered.T
+    np.fill_diagonal(out, 0.0)
+    return standardise(out)
+
+
+def heat_kernel(field: np.ndarray, time: float = 1.0) -> np.ndarray:
+    """Heat flow on the normalised Laplacian: a metric the graph defines itself."""
+    size = len(field)
+    affinity = field - field.min()
+    np.fill_diagonal(affinity, 0.0)
+    degree = affinity.sum(1).clip(1e-9)
+    normalised = affinity / np.sqrt(np.outer(degree, degree))
+    values, vectors = np.linalg.eigh((normalised + normalised.T) / 2.0)
+    out = (vectors * np.exp(time * (values - values.max()))) @ vectors.T
+    np.fill_diagonal(out, 0.0)
+    return standardise(out)
+
+
+def hyperbolic(field: np.ndarray, dimension: int = 16) -> np.ndarray:
+    """Embed each field in hyperbolic space and re-derive distances there.
+
+    If scene similarity is hierarchical -- broad categories containing finer
+    ones -- a flat inner product cannot hold it without distortion, while
+    negative curvature can. The field is embedded by eigenmap and read back
+    as Lorentzian distance.
+    """
+    size = len(field)
+    values, vectors = np.linalg.eigh((field + field.T) / 2.0)
+    order = np.argsort(values)[::-1][:dimension]
+    coordinates = vectors[:, order] * np.sqrt(np.abs(values[order]))
+    scale = np.abs(coordinates).max().clip(1e-9)
+    spatial = coordinates / scale * 0.9
+    time = np.sqrt(1.0 + (spatial ** 2).sum(1))
+    product = np.outer(time, time) - spatial @ spatial.T
+    out = -np.arccosh(np.clip(product, 1.0, None))
+    np.fill_diagonal(out, 0.0)
+    return standardise(out)
+
+
+GEOMETRIES = {
+    "euclidean (current)": identity,
+    "global rank (copula)": global_rank,
+    "row rank": row_rank,
+    "local scaling": local_scaling,
+    "diffusion t=3": diffusion,
+    "heat kernel": heat_kernel,
+    "hyperbolic d=16": hyperbolic,
+}
+
+
+# --- comparison operators -------------------------------------------------
+
+def anchor_bound(first: np.ndarray, second: np.ndarray, repeats: int,
+                 spearman: bool = False) -> float:
+    size = len(first)
+    scores = []
+    for repeat in range(repeats):
+        generator = np.random.default_rng(repeat)
+        order = generator.permutation(size)
+        half = size // 2
+        anchors, probe = order[:half], order[half:]
+
+        def profile(field):
+            block = field[np.ix_(probe, anchors)]
+            if spearman:
+                block = np.apply_along_axis(rankdata, 1, block)
+            block = block - block.mean(1, keepdims=True)
+            return block / block.std(1, keepdims=True).clip(1e-9)
+
+        similarity = profile(first) @ profile(second).T / half
+        _, columns = linear_sum_assignment(-similarity)
+        scores.append(float((columns == np.arange(len(probe))).mean()))
+    return float(np.mean(scores))
+
+
+def main() -> None:
+    args = parse_args()
+    state = torch.load(args.fields, map_location="cpu", weights_only=False)
+    visual = state["visual_field"].double().numpy()
+    text = state["text_field"].double().numpy()
+    mask = offdiagonal_mask(len(visual))
+
+    rows = []
+    print("%-24s %8s %10s %12s" % ("geometry", "rho", "bound", "bound+rank-op"))
+    for name, transform in GEOMETRIES.items():
+        try:
+            first, second = transform(visual), transform(text)
+        except Exception as error:  # a geometry that cannot be formed is a result
+            print("%-24s  failed: %s" % (name, error))
+            continue
+        rho = float(np.corrcoef(first[mask], second[mask])[0, 1])
+        plain = anchor_bound(first, second, args.repeats)
+        ranked = anchor_bound(first, second, args.repeats, spearman=True)
+        rows.append({"geometry": name, "correlation": rho,
+                     "anchor_bound": plain, "anchor_bound_rank_operator": ranked})
+        print("%-24s %8.3f %10.4f %12.4f" % (name, rho, plain, ranked))
+
+    best = max(rows, key=lambda row: max(row["anchor_bound"],
+                                         row["anchor_bound_rank_operator"]))
+    write_json(args.output, {
+        "protocol": (
+            "Each geometry is applied to each modality independently, so no "
+            "transform crosses the gauge or uses the pairing. The bound is the "
+            "anchor upper bound; the last column re-computes it with a rank "
+            "comparison operator instead of a linear one."
+        ),
+        "label": args.label,
+        "rows": rows,
+        "best": best,
+    })
+    print(json.dumps(best))
+
+
+if __name__ == "__main__":
+    main()
diff --git a/worldalign/match_battery.py b/worldalign/match_battery.py
new file mode 100644
index 0000000..45fa26e
--- /dev/null
+++ b/worldalign/match_battery.py
@@ -0,0 +1,205 @@
+"""A matching benchmark, and the solvers we never tried on it.
+
+Everything the project has used to search is one narrow family: spectral
+initialisation on the pairwise field, exact steepest descent over
+transpositions, tempering, and Sinkhorn. All of them act on the permutation
+group with a pairwise objective in flat space. The measured cost of that
+narrowness is large and specific -- on a field whose anchor bound is 0.997,
+meaning the information is intact almost perfectly, blind matching returns
+5.6%.
+
+That makes a clean benchmark, because unlike natural data the answer is known
+to be reachable. Three synthetic fields span the range: one that already
+recovers, one whose captions omit a factor, one truncated to rank eight.
+Natural data is included as the case where no method should be expected to
+work, so a method that scores there is reporting a bug rather than a result.
+
+Gromov-Wasserstein is the canonical formulation for exactly our problem --
+align two metric spaces with no shared embedding, by matching their internal
+distance structures -- and it has never been run on these fields. Its entropic
+form with an annealed regulariser is the landscape-deformation method that the
+rank ladder argued for and that we never reached. Each solver is scored blind
+and then re-scored after exact refinement, since composition is what made the
+spectral solver work.
+"""
+
+from __future__ import annotations
+
+import argparse
+import json
+import time
+
+import numpy as np
+import torch
+
+from .common import write_json
+from .spectral_match import grampa, umeyama
+from .synth_fast_gate import ClosedFormEnergy, fast_pair_descent
+from .synth_triangle_gate import standardized
+
+
+def parse_args() -> argparse.Namespace:
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--fields", nargs="+", required=True)
+    parser.add_argument("--labels", nargs="+", default=None)
+    parser.add_argument("--trials", type=int, default=3)
+    parser.add_argument("--iterations", type=int, default=600)
+    parser.add_argument("--device", default="cuda:3")
+    parser.add_argument("--methods", nargs="+", default=None)
+    parser.add_argument("--output", default="artifacts/synth_v1/match_battery.json")
+    return parser.parse_args()
+
+
+def standardise(matrix: np.ndarray) -> np.ndarray:
+    mask = ~np.eye(len(matrix), dtype=bool)
+    values = matrix[mask]
+    out = (matrix - values.mean()) / values.std()
+    np.fill_diagonal(out, 0.0)
+    return out
+
+
+def coupling_to_permutation(coupling: np.ndarray) -> np.ndarray:
+    from scipy.optimize import linear_sum_assignment
+
+    _, columns = linear_sum_assignment(-coupling)
+    return columns
+
+
+# --- solvers: each returns column assignment for rows of `visual` ---------
+
+def solve_grampa(visual, shuffled, **_):
+    return grampa(visual, shuffled, 1.0)
+
+
+def solve_umeyama(visual, shuffled, **_):
+    return umeyama(visual, shuffled)
+
+
+def solve_gw(visual, shuffled, **_):
+    import ot
+
+    size = len(visual)
+    weights = ot.unif(size)
+    coupling = ot.gromov.gromov_wasserstein(
+        visual, shuffled, weights, weights, "square_loss",
+        max_iter=200, tol_rel=1e-9,
+    )
+    return coupling_to_permutation(coupling)
+
+
+def solve_gw_annealed(visual, shuffled, **_):
+    """Entropic GW with the regulariser annealed from smooth to sharp."""
+    import ot
+
+    size = len(visual)
+    weights = ot.unif(size)
+    coupling = None
+    for epsilon in (5e-2, 2e-2, 1e-2, 5e-3, 2e-3):
+        try:
+            coupling = ot.gromov.entropic_gromov_wasserstein(
+                visual, shuffled, weights, weights, "square_loss",
+                epsilon=epsilon, G0=coupling, max_iter=200, tol=1e-9,
+                solver="PGD",
+            )
+        except Exception:
+            break
+        if not np.isfinite(coupling).all():
+            break
+    if coupling is None or not np.isfinite(coupling).all():
+        raise RuntimeError("entropic GW diverged at every regulariser")
+    return coupling_to_permutation(coupling)
+
+
+def solve_bapg(visual, shuffled, **_):
+    """Bregman alternating projected gradient: a newer GW solver."""
+    import ot
+
+    size = len(visual)
+    weights = ot.unif(size)
+    coupling = ot.gromov.BAPG_gromov_wasserstein(
+        visual, shuffled, weights, weights, "square_loss",
+        epsilon=1e-1, max_iter=1000, tol=1e-9,
+    )
+    return coupling_to_permutation(coupling)
+
+
+SOLVERS = {
+    "grampa (current)": solve_grampa,
+    "umeyama": solve_umeyama,
+    "gromov-wasserstein": solve_gw,
+    "entropic GW annealed": solve_gw_annealed,
+    "BAPG GW": solve_bapg,
+}
+
+
+def main() -> None:
+    args = parse_args()
+    labels = args.labels or [p.split("/")[-1] for p in args.fields]
+    chosen = args.methods or list(SOLVERS)
+    device = torch.device(args.device)
+    rows = []
+
+    for path, label in zip(args.fields, labels):
+        state = torch.load(path, map_location="cpu", weights_only=False)
+        visual = standardise(state["visual_field"].double().numpy())
+        text = standardise(state["text_field"].double().numpy())
+        size = len(visual)
+        text_gpu = standardized(torch.from_numpy(text).to(device)).double()
+        visual_gpu = standardized(torch.from_numpy(visual).to(device)).double()
+        print(f"\n=== {label}  (N={size}, chance={1/size:.4f})", flush=True)
+
+        for name in chosen:
+            solver = SOLVERS[name]
+            raw, refined, seconds, failures = [], [], [], []
+            for trial in range(args.trials):
+                generator = np.random.default_rng(trial)
+                hidden = generator.permutation(size)
+                shuffled = text[np.ix_(hidden, hidden)]
+                start = time.time()
+                try:
+                    columns = solver(visual=visual, shuffled=shuffled)
+                except Exception as error:
+                    failures.append(str(error)[:80])
+                    continue
+                seconds.append(time.time() - start)
+                raw.append(float((hidden[columns] == np.arange(size)).mean()))
+
+                shuffled_gpu = standardized(
+                    torch.from_numpy(shuffled).to(device)
+                ).double()
+                final = fast_pair_descent(
+                    shuffled_gpu, visual_gpu,
+                    torch.from_numpy(np.ascontiguousarray(columns)).to(device),
+                    args.iterations,
+                )
+                refined.append(
+                    float((hidden[final.cpu().numpy()] == np.arange(size)).mean())
+                )
+            row = {
+                "field": label, "method": name,
+                "raw": float(np.mean(raw)) if raw else None,
+                "refined": float(np.mean(refined)) if refined else None,
+                "seconds": float(np.mean(seconds)) if seconds else None,
+                "failures": failures,
+            }
+            rows.append(row)
+            if raw:
+                print("  %-22s raw=%.4f  refined=%.4f  (%.1fs)"
+                      % (name, row["raw"], row["refined"], row["seconds"]), flush=True)
+            else:
+                print("  %-22s FAILED: %s" % (name, failures[0] if failures else "?"),
+                      flush=True)
+
+    write_json(args.output, {
+        "protocol": (
+            "Blind: the hidden permutation is generated per trial and read only "
+            "for scoring. Each solver is scored as returned and again after "
+            "exact steepest-descent refinement."
+        ),
+        "rows": rows,
+    })
+    print(json.dumps({"done": True}))
+
+
+if __name__ == "__main__":
+    main()
diff --git a/worldalign/split_half.py b/worldalign/split_half.py
new file mode 100644
index 0000000..ee3f273
--- /dev/null
+++ b/worldalign/split_half.py
@@ -0,0 +1,190 @@
+"""The ceiling nobody measured: can a modality identify a scene against itself?
+
+Every deficit so far has been read as a cross-modal one, and that assumed each
+side describes its scenes well enough to identify them. Nothing tested the
+assumption. Split each scene's parts into two disjoint halves -- odd and even
+regions, odd and even segments -- build a relation field from each half, and
+match the halves against each other with the anchor bound. Both halves come
+from the same modality and the same scene, so any failure is that modality
+failing to pin down its own scenes.
+
+This upper-bounds anything cross-modal. If one half of an image's description
+cannot pick that image out from among the others given the other half, no
+vision encoder can be blamed for failing to, and no better cross-modal recipe
+exists to find. It costs one pass over data already on disk, and it should have
+been the first measurement on natural data rather than one of the last.
+"""
+
+from __future__ import annotations
+
+import argparse
+import json
+import re
+from pathlib import Path
+
+import numpy as np
+import torch
+import torch.nn.functional as F
+from scipy.optimize import linear_sum_assignment
+
+from .common import seed_everything, write_json
+from .natural_families import load_phrases
+from .natural_pipeline import content_directions, word_vectors
+from .synth_cc_battery import moment_field
+
+
+def parse_args() -> argparse.Namespace:
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--vg-dir", default="artifacts/vg_5k")
+    parser.add_argument("--objects", default="artifacts/vg_5k/nobj_base_gpu.pt")
+    parser.add_argument("--samples", type=int, default=256)
+    parser.add_argument("--fit-scenes", type=int, default=1200)
+    parser.add_argument("--word-vectors", type=int, default=128)
+    parser.add_argument("--min-count", type=int, default=60)
+    parser.add_argument("--max-vocabulary", type=int, default=600)
+    parser.add_argument("--context-window", type=int, default=2)
+    parser.add_argument("--keep", type=int, default=64)
+    parser.add_argument("--shrinkage", type=float, default=0.05)
+    parser.add_argument("--weight-power", type=float, default=0.5)
+    parser.add_argument("--seed", type=int, default=0)
+    parser.add_argument("--output", default="artifacts/vg_5k/split_half.json")
+    return parser.parse_args()
+
+
+def standardise(matrix: np.ndarray) -> np.ndarray:
+    values = matrix[~np.eye(len(matrix), dtype=bool)]
+    out = (matrix - values.mean()) / values.std()
+    np.fill_diagonal(out, 0.0)
+    return out
+
+
+def anchor_bound(first: np.ndarray, second: np.ndarray, repeats: int = 5) -> float:
+    size = len(first)
+    scores = []
+    for repeat in range(repeats):
+        generator = np.random.default_rng(repeat)
+        order = generator.permutation(size)
+        half = size // 2
+        anchors, probe = order[:half], order[half:]
+
+        def profile(field):
+            block = field[np.ix_(probe, anchors)]
+            block = block - block.mean(1, keepdims=True)
+            return block / block.std(1, keepdims=True).clip(1e-9)
+
+        similarity = profile(first) @ profile(second).T / half
+        _, columns = linear_sum_assignment(-similarity)
+        scores.append(float((columns == np.arange(len(probe))).mean()))
+    return float(np.mean(scores))
+
+
+def build_field(sets, fitted, keep, power) -> torch.Tensor:
+    centre, basis, values = fitted
+    width = min(keep, basis.shape[1])
+    scale = values[:width] ** power if power else 1.0
+    projected = [
+        F.normalize(
+            torch.tensor((item - centre) @ basis[:, :width] * scale, dtype=torch.float32),
+            dim=-1,
+        )
+        for item in sets
+    ]
+    return moment_field(projected)
+
+
+def main() -> None:
+    args = parse_args()
+    seed_everything(args.seed)
+    vg_dir = Path(args.vg_dir)
+    truth = [
+        json.loads(line)
+        for line in (vg_dir / "ground_truth.private.jsonl").read_text(encoding="utf-8").splitlines()
+        if line.strip()
+    ]
+    records = {
+        json.loads(line)["node_id"]: json.loads(line)
+        for line in (vg_dir / "text_nodes.jsonl").read_text(encoding="utf-8").splitlines()
+        if line.strip()
+    }
+    state = torch.load(args.objects, map_location="cpu", weights_only=False)
+    index = {node: position for position, node in enumerate(state["node_ids"])}
+    vectors = word_vectors(load_phrases(vg_dir / "text_nodes.jsonl", "region_closed"), args)
+
+    def text_parts(node):
+        rows = []
+        for phrase in records[node]["region_closed"]:
+            hits = [vectors[t] for t in re.findall(r"[a-z]+", phrase.lower()) if t in vectors]
+            if hits:
+                rows.append(np.mean(hits, axis=0))
+        return rows
+
+    def vision_parts(node):
+        return [s["feature"].astype(np.float64) for s in state["segments"][index[node]]]
+
+    pairs = [p for p in truth if p["vision_node_id"] in index and p["text_node_id"] in records]
+
+    def collect(getter, key, rows, minimum):
+        out = []
+        for pair in rows:
+            parts = getter(pair[key])
+            out.append(np.stack(parts) if len(parts) >= minimum else None)
+        return out
+
+    results = {}
+    for name, getter, key in (
+        ("text", text_parts, "text_node_id"),
+        ("vision", vision_parts, "vision_node_id"),
+    ):
+        fit_rows = pairs[args.samples : args.samples + args.fit_scenes]
+        fit = [x for x in collect(getter, key, fit_rows, 2) if x is not None]
+        fitted = content_directions(fit, args.shrinkage)
+
+        evaluated, first_half, second_half = [], [], []
+        for pair in pairs[: args.samples]:
+            parts = getter(pair[key])
+            if len(parts) < 4:
+                continue
+            stacked = np.stack(parts)
+            evaluated.append(pair)
+            first_half.append(stacked[0::2])
+            second_half.append(stacked[1::2])
+
+        field_a = standardise(
+            build_field(first_half, fitted, args.keep, args.weight_power).double().numpy()
+        )
+        field_b = standardise(
+            build_field(second_half, fitted, args.keep, args.weight_power).double().numpy()
+        )
+        mask = ~np.eye(len(evaluated), dtype=bool)
+        results[name] = {
+            "scenes": len(evaluated),
+            "mean_parts": float(np.mean([len(getter(p[key])) for p in evaluated])),
+            "split_half_correlation": float(
+                np.corrcoef(field_a[mask], field_b[mask])[0, 1]
+            ),
+            "split_half_anchor_bound": anchor_bound(field_a, field_b),
+        }
+        print(
+            f"{name:<7} scenes={results[name]['scenes']:4d} "
+            f"parts={results[name]['mean_parts']:5.1f}  "
+            f"rho={results[name]['split_half_correlation']:.3f}  "
+            f"self-identification bound={results[name]['split_half_anchor_bound']:.3f}",
+            flush=True,
+        )
+
+    summary = {
+        "protocol": (
+            "Each scene's parts split into disjoint halves within one modality; "
+            "a relation field is built from each half and the two are matched by "
+            "the anchor bound. Both halves describe the same scenes, so this "
+            "upper-bounds any cross-modal result."
+        ),
+        "results": results,
+        "cross_modal_bound_for_reference": 0.291,
+    }
+    print(json.dumps(summary["results"]))
+    write_json(args.output, summary)
+
+
+if __name__ == "__main__":
+    main()
diff --git a/worldalign/synth_fast_gate.py b/worldalign/synth_fast_gate.py
index 54b6891..01fb18f 100644
--- a/worldalign/synth_fast_gate.py
+++ b/worldalign/synth_fast_gate.py
@@ -343,3 +343,54 @@ def main() -> None:
 
 if __name__ == "__main__":
     main()
+
+
+def all_pair_swap_deltas(
+    permuted_text: torch.Tensor, visual: torch.Tensor
+) -> torch.Tensor:
+    """Energy change of every transposition at once, for the pairwise term.
+
+    The brute-force descent forms all N(N-1)/2 candidate permutations and
+    scores each through the full energy, which at 256 scenes means 32,640
+    matrices per step and, with the triangle term active, a batched cube
+    trace over every one of them -- 203 seconds for a single descent. The
+    pairwise term does not need any of that. Writing S for the alignment sum
+    and C for `permuted_text @ visual`, swapping positions p and q changes S
+    by 2(C_pq + C_qp - C_pp - C_qq + 2 A_pq B_pq), so one matrix product
+    yields the entire table.
+
+    Both fields must be symmetric with zero diagonal, which standardisation
+    already guarantees. Returns a dense [N, N] table of energy deltas; the
+    sign convention matches `ClosedFormEnergy.energy`, so negative is an
+    improvement.
+    """
+    size = permuted_text.shape[-1]
+    product = permuted_text @ visual
+    diagonal = torch.diagonal(product)
+    gain = 2.0 * (
+        product + product.T
+        - diagonal[:, None] - diagonal[None, :]
+        + 2.0 * permuted_text * visual
+    )
+    return -2.0 * gain / (size * (size - 1))
+
+
+def fast_pair_descent(
+    text: torch.Tensor, visual: torch.Tensor, start: torch.Tensor, max_steps: int
+) -> torch.Tensor:
+    """Steepest descent on the pairwise energy using the closed-form table."""
+    size = len(visual)
+    current = start.clone()
+    triangle = torch.triu(
+        torch.ones(size, size, dtype=torch.bool, device=visual.device), diagonal=1
+    )
+    for _ in range(max_steps):
+        permuted = text[current[:, None], current[None, :]]
+        deltas = all_pair_swap_deltas(permuted, visual)
+        deltas = torch.where(triangle, deltas, torch.inf)
+        best = int(deltas.argmin())
+        if float(deltas.view(-1)[best]) >= -1e-12:
+            break
+        p, q = best // size, best % size
+        current[[p, q]] = current[[q, p]]
+    return current