beam: add lis (O(n log n) patience sorting) + tests

wizard-work/beam/lis.py

@@ -0,0 +1,128 @@
+"""Longest increasing subsequence (LIS) utilities.
+
+Provides:
+  - lis_length(seq): length of a longest strictly increasing subsequence.
+  - lis(seq): one such subsequence (list of items).
+  - lis_indices(seq): indices into seq of one such subsequence.
+
+Algorithm: patience sorting with binary search, O(n log n) time, O(n) space.
+Ties are broken to keep the earliest possible tail, which preserves reconstruction.
+"""
+
+from __future__ import annotations
+
+import bisect
+from typing import Sequence, TypeVar
+
+T = TypeVar("T")
+
+
+def _to_list(seq) -> list:
+    try:
+        return list(seq)
+    except TypeError as e:
+        raise TypeError("seq must be iterable") from e
+
+
+def lis_length(seq) -> int:
+    """Length of a longest strictly increasing subsequence of seq."""
+    a = _to_list(seq)
+    tails = []
+    for x in a:
+        i = bisect.bisect_left(tails, x)
+        if i == len(tails):
+            tails.append(x)
+        else:
+            tails[i] = x
+    return len(tails)
+
+
+def lis(seq) -> list:
+    """Return one longest strictly increasing subsequence of seq (as a list)."""
+    a = _to_list(seq)
+    n = len(a)
+    if n == 0:
+        return []
+    # tails[i] = smallest tail value of an increasing subsequence of length i+1
+    tails = []
+    # parent[k] = index in a of predecessor of a[k] in the chosen subsequence
+    parent = [-1] * n
+    # tail_idx[i] = index in a of the element currently at tails[i]
+    tail_idx = []
+    for k, x in enumerate(a):
+        i = bisect.bisect_left(tails, x)
+        if i == len(tails):
+            tails.append(x)
+            tail_idx.append(k)
+        else:
+            tails[i] = x
+            tail_idx[i] = k
+        if i > 0:
+            parent[k] = tail_idx[i - 1]
+    # reconstruct from the last tail
+    last = tail_idx[-1]
+    out = []
+    while last != -1:
+        out.append(a[last])
+        last = parent[last]
+    out.reverse()
+    return out
+
+
+def lis_indices(seq) -> list[int]:
+    """Return indices into seq of one longest strictly increasing subsequence."""
+    a = _to_list(seq)
+    n = len(a)
+    if n == 0:
+        return []
+    tails = []
+    tail_idx = []
+    parent = [-1] * n
+    for k, x in enumerate(a):
+        i = bisect.bisect_left(tails, x)
+        if i == len(tails):
+            tails.append(x)
+            tail_idx.append(k)
+        else:
+            tails[i] = x
+            tail_idx[i] = k
+        if i > 0:
+            parent[k] = tail_idx[i - 1]
+    last = tail_idx[-1]
+    out = []
+    while last != -1:
+        out.append(last)
+        last = parent[last]
+    out.reverse()
+    return out
+
+
+def is_increasing(sub) -> bool:
+    s = _to_list(sub)
+    return all(s[i] < s[i + 1] for i in range(len(s) - 1))
+
+
+if __name__ == "__main__":
+    import random
+    # brute-force check
+    def brute(seq):
+        n = len(seq)
+        best = 0
+        # DP
+        dp = [1] * n
+        for i in range(n):
+            for j in range(i):
+                if seq[j] < seq[i]:
+                    dp[i] = max(dp[i], dp[j] + 1)
+        return max(dp) if n else 0
+    for _ in range(2000):
+        L = random.randint(0, 12)
+        seq = [random.randint(-5, 5) for _ in range(L)]
+        assert lis_length(seq) == brute(seq), (seq, lis_length(seq), brute(seq))
+        sub = lis(seq)
+        assert len(sub) == brute(seq)
+        assert is_increasing(sub)
+        # indices valid
+        idx = lis_indices(seq)
+        assert [seq[i] for i in idx] == sub
+    print("self-check ok")

wizard-work/beam/test_lis.py

@@ -0,0 +1,71 @@
+"""Tests for lis.py."""
+import unittest
+import random
+from lis import lis_length, lis, lis_indices, is_increasing
+
+
+class TestLIS(unittest.TestCase):
+    def test_empty(self):
+        self.assertEqual(lis_length([]), 0)
+        self.assertEqual(lis([]), [])
+        self.assertEqual(lis_indices([]), [])
+
+    def test_single(self):
+        self.assertEqual(lis_length([5]), 1)
+        self.assertEqual(lis([5]), [5])
+
+    def test_sorted(self):
+        self.assertEqual(lis_length([1,2,3,4]), 4)
+        self.assertEqual(lis([1,2,3,4]), [1,2,3,4])
+
+    def test_reverse(self):
+        self.assertEqual(lis_length([4,3,2,1]), 1)
+
+    def test_duplicates_strict(self):
+        self.assertEqual(lis_length([1,1,1]), 1)
+        self.assertEqual(lis([2,2,1,1,3]), [1,3])
+
+    def test_known(self):
+        seq = [10, 9, 2, 5, 3, 7, 101, 18]
+        self.assertEqual(lis_length(seq), 4)
+        sub = lis(seq)
+        self.assertEqual(len(sub), 4)
+        self.assertTrue(is_increasing(sub))
+
+    def test_indices(self):
+        seq = [3, 1, 2]
+        idx = lis_indices(seq)
+        self.assertEqual([seq[i] for i in idx], lis(seq))
+        self.assertEqual(lis_length(seq), 2)
+
+    def test_brute_random(self):
+        def brute(seq):
+            n = len(seq)
+            if n == 0:
+                return 0
+            dp = [1]*n
+            for i in range(n):
+                for j in range(i):
+                    if seq[j] < seq[i]:
+                        dp[i] = max(dp[i], dp[j]+1)
+            return max(dp)
+        rng = random.Random(42)
+        for _ in range(300):
+            L = rng.randint(0, 14)
+            seq = [rng.randint(-6, 6) for _ in range(L)]
+            self.assertEqual(lis_length(seq), brute(seq))
+            sub = lis(seq)
+            self.assertEqual(len(sub), brute(seq))
+            self.assertTrue(is_increasing(sub))
+            idx = lis_indices(seq)
+            self.assertEqual([seq[i] for i in idx], sub)
+
+    def test_tuple_input(self):
+        self.assertEqual(lis_length((3,1,4,1,5)), 3)
+        sub = lis((3,1,4,1,5))
+        self.assertEqual(len(sub), 3)
+        self.assertTrue(is_increasing(sub))
+
+
+if __name__ == "__main__":
+    unittest.main()