feat(numeric): add gale_shapley_traced with proposal stats and selector hook

Expose proposal counts and a pluggable proposer-selection rule on the
numpy-backed M-W loop. Motivated by downstream RL work that needs to
benchmark "fewer proposals" strategies against deferred-acceptance baselines.

Added:
- GSStats: frozen dataclass with match, total proposals, per-proposer counts.
- Selector type alias + lifo_selector / fifo_selector / random_selector helpers.
- gale_shapley_traced: same M-W loop, returns GSStats, accepts a selector.
- men_optimal_traced / women_optimal_traced: convenience wrappers that mirror
  the existing men_optimal_gs / women_optimal_gs orientation conventions.
- gale_shapley refactored to a thin wrapper: returns gale_shapley_traced(...).match.

Tests cover Knuth invariance of the match (and proposal count) under
selector choice, the canonical n*(n+1)/2 worst case for shared
preferences, women_optimal_traced orientation, and validation parity.

Coverage on numeric/gs.py: 100%.

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

Author: oedokumaci

src/gale_shapley_algorithm/numeric/__init__.py

@@ -30,7 +30,19 @@
     >>> lattice = enumerate_stable_matchings(proposer_rank, responder_rank)
 """
 
-from gale_shapley_algorithm.numeric.gs import gale_shapley, men_optimal_gs, women_optimal_gs
+from gale_shapley_algorithm.numeric.gs import (
+    GSStats,
+    Selector,
+    fifo_selector,
+    gale_shapley,
+    gale_shapley_traced,
+    lifo_selector,
+    men_optimal_gs,
+    men_optimal_traced,
+    random_selector,
+    women_optimal_gs,
+    women_optimal_traced,
+)
 from gale_shapley_algorithm.numeric.lattice import (
     apply_rotation,
     enumerate_stable_matchings,
@@ -39,12 +51,20 @@
 from gale_shapley_algorithm.numeric.stability import find_blocking_pairs, is_stable
 
 __all__ = [
+    "GSStats",
+    "Selector",
     "apply_rotation",
     "enumerate_stable_matchings",
     "exposed_rotations",
+    "fifo_selector",
     "find_blocking_pairs",
     "gale_shapley",
+    "gale_shapley_traced",
     "is_stable",
+    "lifo_selector",
     "men_optimal_gs",
+    "men_optimal_traced",
+    "random_selector",
     "women_optimal_gs",
+    "women_optimal_traced",
 ]

src/gale_shapley_algorithm/numeric/gs.py

@@ -2,6 +2,8 @@
 
 from __future__ import annotations
 
+from collections.abc import Callable, Sequence
+from dataclasses import dataclass
 from typing import TYPE_CHECKING
 
 try:
@@ -16,6 +18,76 @@
     from numpy.typing import NDArray
 
 
+type Selector = Callable[[Sequence[int]], int]
+"""Picks which free proposer goes next.
+
+Receives the current free list of proposer ids and returns an index into
+that list (``0..len(free)-1``; negative indices supported per ``list.pop``
+semantics). The default is :func:`lifo_selector`, which matches the
+historical behavior of :func:`gale_shapley` (``free.pop()``).
+"""
+
+
+def lifo_selector(free: Sequence[int]) -> int:
+    """Return the last index — equivalent to ``list.pop()`` (LIFO)."""
+    del free
+    return -1
+
+
+def fifo_selector(free: Sequence[int]) -> int:
+    """Return the first index — pops in arrival order (FIFO)."""
+    del free
+    return 0
+
+
+def random_selector(rng: np.random.Generator) -> Selector:
+    """Build a uniform-random selector backed by a numpy ``Generator``.
+
+    Roth-Vande Vate-style: pick a free proposer uniformly at random each step.
+    Final matching is invariant to selector choice (Knuth's order independence)
+    but the proposal count is not — random order is a useful baseline against
+    LIFO/FIFO for amortized analysis.
+
+    Args:
+        rng: a numpy ``Generator`` (use ``np.random.default_rng(seed)`` for
+            determinism).
+
+    Returns:
+        A :data:`Selector` closure that draws ``rng.integers(len(free))`` per call.
+    """
+
+    def _selector(free: Sequence[int]) -> int:
+        return int(rng.integers(len(free)))
+
+    return _selector
+
+
+@dataclass(slots=True, frozen=True, eq=False)
+class GSStats:
+    """Result of a traced Gale-Shapley run.
+
+    ``eq=False`` because :attr:`match` and :attr:`proposals_per_proposer` are
+    numpy arrays whose ``__eq__`` returns an array, not a bool — so the
+    dataclass-generated ``__eq__`` would raise. Compare fields directly.
+
+    Attributes:
+        match: shape ``(n,)``, dtype ``int16``. Proposer-indexed:
+            ``match[p]`` is the responder paired with proposer ``p``. The
+            wrapper :func:`women_optimal_traced` re-indexes this to be
+            men-indexed; see its docstring.
+        proposals: total number of proposal events. Equals
+            ``proposals_per_proposer.sum()``.
+        proposals_per_proposer: shape ``(n,)``, dtype ``int16``.
+            ``proposals_per_proposer[p]`` is the number of proposals made by
+            proposer ``p``, which also equals the rank (1-indexed) of ``p``'s
+            final match in ``p``'s preference list.
+    """
+
+    match: NDArray[np.int16]
+    proposals: int
+    proposals_per_proposer: NDArray[np.int16]
+
+
 def _validate_rank_matrices(
     proposer_rank: NDArray[np.integer],
     responder_rank: NDArray[np.integer],
@@ -47,31 +119,44 @@ def _validate_rank_matrices(
             )
 
 
-def gale_shapley(proposer_rank: NDArray[np.integer], responder_rank: NDArray[np.integer]) -> NDArray[np.int16]:
-    """Run proposer-optimal deferred acceptance.
+def gale_shapley_traced(
+    proposer_rank: NDArray[np.integer],
+    responder_rank: NDArray[np.integer],
+    *,
+    selector: Selector = lifo_selector,
+) -> GSStats:
+    """Run proposer-optimal deferred acceptance and return per-proposer statistics.
+
+    Mechanically identical to :func:`gale_shapley` (sequential McVitie-Wilson
+    on the free-proposer pool), but exposes the proposal count and a
+    pluggable proposer-selection rule. The final matching is invariant under
+    selector choice (Knuth); the proposal count is not.
 
     Args:
         proposer_rank: ``(n, n)`` array where ``proposer_rank[i, j]`` is
             the 1-indexed position of responder ``j`` in proposer ``i``'s
             preference list. Each row must be a permutation of ``1..n``.
         responder_rank: ``(n, n)`` array with the symmetric convention.
+        selector: picks which free proposer acts next; defaults to
+            :func:`lifo_selector`. See :data:`Selector`.
 
     Returns:
-        ``match`` of shape ``(n,)``, dtype ``int16``, where ``match[i]`` is
-        the responder paired with proposer ``i``.
+        :class:`GSStats` with the matching, total proposal count, and
+        per-proposer counts.
 
     Raises:
         ValueError: if the two arrays don't have the same square shape, or
             any row is not a permutation of ``1..n``.
+        IndexError: if ``selector`` returns an index outside the free list.
     """
     _validate_rank_matrices(proposer_rank, responder_rank)
     n = proposer_rank.shape[0]
     next_proposal = np.zeros(n, dtype=np.int16)
     responder_match = np.full(n, -1, dtype=np.int16)
     proposer_pref = np.argsort(proposer_rank, axis=1).astype(np.int16)
-    free = list(range(n))
+    free: list[int] = list(range(n))
     while free:
-        p = free.pop()
+        p = free.pop(selector(free))
         r = int(proposer_pref[p, next_proposal[p]])
         next_proposal[p] += 1
         current = int(responder_match[r])
@@ -82,7 +167,31 @@ def gale_shapley(proposer_rank: NDArray[np.integer], responder_rank: NDArray[np.
             free.append(current)
         else:
             free.append(p)
-    return np.argsort(responder_match).astype(np.int16)
+    return GSStats(
+        match=np.argsort(responder_match).astype(np.int16),
+        proposals=int(next_proposal.sum()),
+        proposals_per_proposer=next_proposal.copy(),
+    )
+
+
+def gale_shapley(proposer_rank: NDArray[np.integer], responder_rank: NDArray[np.integer]) -> NDArray[np.int16]:
+    """Run proposer-optimal deferred acceptance.
+
+    Args:
+        proposer_rank: ``(n, n)`` array where ``proposer_rank[i, j]`` is
+            the 1-indexed position of responder ``j`` in proposer ``i``'s
+            preference list. Each row must be a permutation of ``1..n``.
+        responder_rank: ``(n, n)`` array with the symmetric convention.
+
+    Returns:
+        ``match`` of shape ``(n,)``, dtype ``int16``, where ``match[i]`` is
+        the responder paired with proposer ``i``.
+
+    Raises:
+        ValueError: if the two arrays don't have the same square shape, or
+            any row is not a permutation of ``1..n``.
+    """
+    return gale_shapley_traced(proposer_rank, responder_rank).match
 
 
 def men_optimal_gs(men_rank: NDArray[np.integer], women_rank: NDArray[np.integer]) -> NDArray[np.int16]:
@@ -93,3 +202,57 @@ def men_optimal_gs(men_rank: NDArray[np.integer], women_rank: NDArray[np.integer
 def women_optimal_gs(men_rank: NDArray[np.integer], women_rank: NDArray[np.integer]) -> NDArray[np.int16]:
     """Return the women-optimal stable matching, still in men-indexed form ``match[m] = w``."""
     return np.argsort(gale_shapley(women_rank, men_rank)).astype(np.int16)
+
+
+def men_optimal_traced(
+    men_rank: NDArray[np.integer],
+    women_rank: NDArray[np.integer],
+    *,
+    selector: Selector = lifo_selector,
+) -> GSStats:
+    """Men-optimal stable matching with per-man proposal stats.
+
+    Equivalent to ``gale_shapley_traced(men_rank, women_rank, selector=...)``;
+    provided for symmetry with :func:`men_optimal_gs`.
+
+    Args:
+        men_rank: ``(n, n)`` rank matrix for men (proposers).
+        women_rank: ``(n, n)`` rank matrix for women (responders).
+        selector: see :data:`Selector`.
+
+    Returns:
+        :class:`GSStats` where ``match`` and ``proposals_per_proposer`` are
+        both men-indexed.
+    """
+    return gale_shapley_traced(men_rank, women_rank, selector=selector)
+
+
+def women_optimal_traced(
+    men_rank: NDArray[np.integer],
+    women_rank: NDArray[np.integer],
+    *,
+    selector: Selector = lifo_selector,
+) -> GSStats:
+    """Women-optimal stable matching with per-woman proposal stats.
+
+    Runs :func:`gale_shapley_traced` with women as the proposing side. The
+    returned ``match`` is re-indexed to men-indexed form for consistency with
+    :func:`women_optimal_gs`. ``proposals_per_proposer`` remains
+    **women-indexed**, since women are the proposers in this run.
+
+    Args:
+        men_rank: ``(n, n)`` rank matrix for men.
+        women_rank: ``(n, n)`` rank matrix for women.
+        selector: see :data:`Selector`.
+
+    Returns:
+        :class:`GSStats` with ``match[m] = w`` (men-indexed) and
+        ``proposals_per_proposer[w] = number of proposals woman w made``
+        (women-indexed). ``proposals`` is the unambiguous total.
+    """
+    stats = gale_shapley_traced(women_rank, men_rank, selector=selector)
+    return GSStats(
+        match=np.argsort(stats.match).astype(np.int16),
+        proposals=stats.proposals,
+        proposals_per_proposer=stats.proposals_per_proposer,
+    )