Add force-first Heawood labelling to the medial tire dual-cut experiment

heawood_labelling(): depth-seeded force-first +/-1 labelling of the
source-dual cut, targeting sum ≡ 0 (mod 3) on each dual face (vertex
link), with bookkeeping for seeds, forced fills, unforceable faces, and
failing faces.

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

papers/medial_tire_cuts/experiments/medial_tire_dual_cut_experiment.py

@@ -259,6 +259,155 @@ def dual_cut_distances(result):
     return nx.single_source_shortest_path_length(dual_cut_graph(result), root), root
 
 
+def _dual_faces_by_vertex(result):
+    """Each *dual face* of the source-dual cut is the link of a primal vertex
+    ``v``: the set of triangular faces (dual vertices) incident to ``v``, of
+    size ``deg(v)``.  Returns ``{v: [face index, ...]}``."""
+    faces = result["faces"]
+    face_of = defaultdict(list)
+    for fi, f in enumerate(faces):
+        for v in f:
+            face_of[v].append(fi)
+    return dict(face_of)
+
+
+def heawood_labelling(result):
+    """Force-first Heawood labelling of the source-dual cut, seeded by walk depth.
+
+    Each dual vertex (triangular face of the source graph) is given a label in
+    ``{+1, -1}``; the target is that every *dual face* -- the link of a primal
+    vertex ``v``, i.e. the ``deg(v)`` triangles around ``v`` -- has label sum
+    ``≡ 0 (mod 3)`` (Heawood's condition, whose existence is equivalent to
+    4-colourability of the source triangulation).
+
+    The procedure interleaves forcing and depth-seeding:
+
+      1. **Force (saturate):** while some dual face has **two or fewer**
+         unlabelled vertices, fill them so that face's sum is ``≡ 0 (mod 3)``.
+
+           * one unlabelled slot ``t``: forced to the unique ``±1`` with
+             ``label(t) ≡ -S (mod 3)`` over the known sum ``S``.  If ``-S ≡ 0``
+             no ``±1`` completes the face -- it is an *unforceable* face,
+             recorded as a violation (the slot is still filled by its own depth
+             parity so the walk can continue).
+           * two unlabelled slots: forced to the same sign when the rule
+             demands it (both ``-1`` or both ``+1``); when the rule only demands
+             *opposite* signs (the known labels already sum to ``≡ 0``), each
+             slot takes its own depth parity -- the even-depth one ``+1``, the
+             odd-depth one ``-1``.
+
+         Repeat until no dual face has ``≤ 2`` unlabelled vertices.
+
+      2. **Seed:** take the unlabelled dual vertex of smallest walk depth (none
+         of whose incident faces is forcing it) and label it by depth parity
+         (``+1`` if even, ``-1`` if odd).  Return to step 1.
+
+    Forcing only ever fills *empty* slots, so no label is overwritten and the
+    only failure mode is a dual face whose final label sum is ``≢ 0 (mod 3)``.
+
+    Returns a dict with ``labels`` (``{dual vertex: ±1}``), ``failing_faces``
+    (primal vertices whose final sum is ``≢ 0 mod 3``), ``success`` (no failing
+    faces), ``seeds`` (depth-seeded vertices, step 2), ``forced`` (vertices
+    filled by step 1, split into ``single``/``pair_same``/``pair_opposite``),
+    ``unforceable`` (one-slot faces with no valid ``±1``), and bookkeeping
+    (``n_dual``, ``max_depth``)."""
+    G = result["G"]
+    dist, root = dual_cut_distances(result)
+    if root is None:
+        raise ValueError("source-dual cut has no entry down tooth to root from")
+    face_of = {v: face_of_v
+               for v, face_of_v in sorted(_dual_faces_by_vertex(result).items())}
+
+    def parity(t):
+        return 1 if dist[t] % 2 == 0 else -1
+
+    labels = {}
+    counts = {"single": 0, "pair_same": 0, "pair_opposite": 0}
+    unforceable = []
+    seeds = 0
+
+    def saturate():
+        changed = True
+        while changed:
+            changed = False
+            for v, tri in face_of.items():
+                unl = [t for t in tri if t not in labels]
+                s = sum(labels[t] for t in tri if t in labels) % 3
+                needed = (-s) % 3      # value (or pair-sum) needed mod 3
+                if len(unl) == 1:
+                    t = unl[0]
+                    if needed == 1:
+                        labels[t] = 1
+                    elif needed == 2:
+                        labels[t] = -1
+                    else:              # -S ≡ 0: no ±1 closes this face
+                        unforceable.append(v)
+                        labels[t] = parity(t)
+                    counts["single"] += 1
+                    changed = True
+                elif len(unl) == 2:
+                    t1, t2 = unl
+                    if needed == 1:
+                        labels[t1] = labels[t2] = -1
+                        counts["pair_same"] += 1
+                    elif needed == 2:
+                        labels[t1] = labels[t2] = 1
+                        counts["pair_same"] += 1
+                    else:              # opposite signs: each by its depth parity
+                        labels[t1], labels[t2] = parity(t1), parity(t2)
+                        counts["pair_opposite"] += 1
+                    changed = True
+
+    order = sorted(dist, key=lambda t: (dist[t], t))
+    saturate()
+    while len(labels) < len(dist):
+        u = next(t for t in order if t not in labels)
+        labels[u] = parity(u)
+        seeds += 1
+        saturate()
+
+    failing_faces = []
+    for v, tri in face_of.items():
+        s = sum(labels[t] for t in tri)
+        if s % 3 != 0:
+            failing_faces.append({
+                "vertex": v, "degree": G.degree(v), "sum": s,
+                "labels": [labels[t] for t in tri]})
+
+    return {
+        "labels": labels, "failing_faces": failing_faces,
+        "seeds": seeds, "forced": counts, "unforceable": unforceable,
+        "root": root, "n_dual": len(result["faces"]),
+        "n_labelled": len(labels),
+        "max_depth": max(dist.values()) if dist else 0,
+        "success": not failing_faces,
+    }
+
+
+def heawood_report(result):
+    """One-line-per-fact summary of ``heawood_labelling`` for printing."""
+    h = heawood_labelling(result)
+    verdict = ("TERMINATES successfully (every dual face sum ≡ 0 mod 3)"
+               if h["success"] else
+               "FAILS (a dual face has label sum ≢ 0 mod 3)")
+    f = h["forced"]
+    lines = [
+        f"Heawood labelling: {verdict}",
+        f"  dual vertices: {h['n_dual']}  labelled: {h['n_labelled']}  "
+        f"max walk depth: {h['max_depth']}",
+        f"  depth-seeded (step 2): {h['seeds']}   forced (step 1): "
+        f"{f['single']} single + {f['pair_same']} same-sign pairs + "
+        f"{f['pair_opposite']} opposite-sign pairs",
+        f"  unforceable faces (1 slot, no valid ±1): {len(h['unforceable'])}",
+        f"  failing dual faces (sum % 3 != 0): {len(h['failing_faces'])}",
+    ]
+    for ff in h["failing_faces"]:
+        lines.append(
+            f"    vertex {ff['vertex']} (deg {ff['degree']}): "
+            f"sum {ff['sum']} from labels {ff['labels']}")
+    return "\n".join(lines)
+
+
 # --------------------------------------------------------------------------- #
 # The four chained entry points.
 # --------------------------------------------------------------------------- #
@@ -881,6 +1030,9 @@ def main():
                         help="render tread 0 (the source cap) to PNG")
     parser.add_argument("--pdf", metavar="PATH",
                         help="render dual, tire tree, and tire cuts in one PDF")
+    parser.add_argument("--heawood", action="store_true",
+                        help="run the walk-depth-seeded Heawood labelling and "
+                             "report whether it terminates successfully")
     args = parser.parse_args()
 
     rng = random.Random(args.seed)
@@ -912,6 +1064,8 @@ def main():
                                  min_degree=args.min_degree)
 
     print(summary(result))
+    if args.heawood:
+        print(heawood_report(result))
     if args.png:
         draw_png(result, args.png)
         print(f"wrote {args.png}")