Lipton-Tarjan separator algorithm in Python
Finds a planar separator in linear time using a BFS tree and fundamental cycles.
- File
lipton_tarjan_separator.py, 389 lines- Level
- I0 not checked yet
- Toolchain
- Python 3.14, pytest 9
- Test framework
- pytest
- Run the tests
python3 -m pytest
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/lipton-tarjan-separator/. The authors are listed in prov/YCV3MB.jsonl.
from collections import deque
Rotation = dict[int, list[int]]
def _components(vertices: list[int], adj: Rotation, removed: set[int]) -> list[list[int]]:
"""Connected components of the graph minus removed, by iterative search."""
seen = set(removed)
comps: list[list[int]] = []
for s in vertices:
if s in seen:
continue
seen.add(s)
comp = [s]
stack = [s]
while stack:
u = stack.pop()
for w in adj[u]:
if w not in seen:
seen.add(w)
comp.append(w)
stack.append(w)
comps.append(comp)
return comps
def _check_planar(rot: Rotation) -> None:
"""Raise ValueError unless rot is a planar rotation system of a simple graph."""
pos: dict[tuple[int, int], int] = {}
m2 = 0
for v, around in rot.items():
if len(set(around)) != len(around) or v in around:
raise ValueError("the graph must be simple")
for i, w in enumerate(around):
pos[(v, w)] = i
m2 += 1
for v, w in pos:
if (w, v) not in pos:
raise ValueError("the rotation system is not symmetric")
seen: set[tuple[int, int]] = set()
faces = 0
for (u, v) in pos:
if (u, v) in seen:
continue
faces += 1
d = (u, v)
while d not in seen:
seen.add(d)
a, b = d
ring = rot[b]
d = (b, ring[pos[(b, a)] - 1])
n = len(rot)
isolated = sum(1 for v in rot if not rot[v])
c = len(_components(list(rot), rot, set()))
if n - m2 // 2 + faces + isolated != 2 * c:
raise ValueError("the rotation system is not planar")
def _group(pieces: list[list[int]], n: int) -> tuple[list[int], list[int]]:
"""Split pieces of at most 2n/3 vertices each into two sides of at most 2n/3."""
a: list[int] = []
b: list[int] = []
big = [p for p in pieces if 3 * len(p) >= n]
if big:
a = list(big[0])
for p in pieces:
if p is not big[0]:
b.extend(p)
return a, b
for p in pieces:
if 3 * len(a) < n:
a.extend(p)
else:
b.extend(p)
return a, b
def _within(cost: int, budget: int) -> bool:
"""True if cost <= 2 * sqrt(budget), in exact integer arithmetic."""
return cost * cost <= 4 * budget
def _cycle_separator(comp: list[int], rot: Rotation, level: dict[int, int], parent: dict[int, int],
l0: int, l2: int, root_q: int) -> list[int]:
"""Lines L7 to L12: the vertices of M on a balanced fundamental cycle."""
in_m = {v for v in comp if l0 < level[v] < l2} # @line L7: H is M plus the levels <= l0 contracted into x
in_r = {v for v in comp if level[v] <= l0}
ids: dict[int, int] = {}
names: list[int | None] = []
if l0 >= 0:
ids[-1] = 0 # the contracted root x has no name in G
names.append(None)
for v in comp:
if v in in_m:
ids[v] = len(names)
names.append(v)
hrot: list[list[int]] = [[] for _ in names]
if l0 >= 0:
x = 0
kept: dict[int, int] = {} # vertex of M -> vertex of R whose edge is kept
stack = [(root_q, 0, len(rot[root_q]))] # Euler tour of the BFS tree inside R
while stack:
v, start, count = stack.pop()
ring = rot[v]
while count > 0:
w = ring[start % len(ring)]
start += 1
count -= 1
if w in in_r and parent.get(w) == v:
stack.append((v, start, count))
i = ring_index(rot, w, v)
stack.append((w, i + 1, len(rot[w]) - 1))
break
if w in in_m and w not in kept:
kept[w] = v
hrot[x].append(ids[w])
for v in comp:
if v in in_m:
out = hrot[ids[v]]
for w in rot[v]:
if w in in_m:
out.append(ids[w])
elif w in in_r and kept.get(v) == w:
out.append(x)
root = x
else:
for v in comp:
if v in in_m:
hrot[ids[v]] = [ids[w] for w in rot[v] if w in in_m]
root = ids[root_q]
h = len(names)
# Darts of H: dart 2e and 2e + 1 are the two directions of edge e.
tail: list[int] = []
head: list[int] = []
dart_of: dict[tuple[int, int], int] = {}
for a in range(h):
for b in hrot[a]:
if (a, b) not in dart_of:
dart_of[(a, b)] = len(tail)
tail.append(a)
head.append(b)
dart_of[(b, a)] = len(tail)
tail.append(b)
head.append(a)
drot = [[dart_of[(a, b)] for b in hrot[a]] for a in range(h)]
nh = len(tail) # @line L8: put a new vertex in every face of H, joined to each corner
pos = [0] * nh
for a in range(h):
for i, d in enumerate(drot[a]):
pos[d] = i
hface = [-1] * nh
face_darts: list[list[int]] = []
for d0 in range(nh):
if hface[d0] >= 0:
continue
walk = []
d = d0
while hface[d] < 0:
hface[d] = len(face_darts)
walk.append(d)
b = head[d]
ring = drot[b]
d = ring[pos[d ^ 1] - 1]
face_darts.append(walk)
total = h + len(face_darts)
spoke: list[int] = [0] * nh # the dart from tail(d) to the new vertex of face(d)
new_rot: list[list[int]] = [[] for _ in range(total)]
for d in range(nh):
c = h + hface[d]
spoke[d] = len(tail)
tail.append(tail[d])
head.append(c)
tail.append(c)
head.append(tail[d])
for a in range(h):
for d in drot[a]:
new_rot[a].append(d)
new_rot[a].append(spoke[d])
for f, walk in enumerate(face_darts):
new_rot[h + f] = [spoke[d] ^ 1 for d in walk]
nd = len(tail)
tparent = [-1] * total # @line L9: the tree T, inherited from the search; weights 1 on M
tdart = [-1] * total # the dart from a vertex to its parent
for v in comp:
if v in in_m:
p = parent[v]
pv = 0 if (l0 >= 0 and level[v] == l0 + 1) else ids[p]
tparent[ids[v]] = pv
tdart[ids[v]] = dart_of[(ids[v], pv)]
for f, walk in enumerate(face_darts):
c = h + f
tparent[c] = tail[walk[0]]
tdart[c] = spoke[walk[0]] ^ 1
weight = [0] * total
for a in range(h):
if names[a] is not None:
weight[a] = 1
wtotal = sum(weight)
in_tree = [False] * nd
for v in range(total):
if tdart[v] >= 0:
in_tree[tdart[v]] = True
in_tree[tdart[v] ^ 1] = True
order = [root] # vertices of T in breadth-first order
children: list[list[int]] = [[] for _ in range(total)]
for v in range(total):
if tparent[v] >= 0:
children[tparent[v]].append(v)
for v in order:
order.extend(children[v])
depth = [0] * total
pre = [0] * total # weight of the tree path from v up to the root
for v in order:
if v != root:
depth[v] = depth[tparent[v]] + 1
pre[v] = pre[tparent[v]] + weight[v]
else:
pre[v] = weight[v]
up = [[tparent[v] if tparent[v] >= 0 else v for v in range(total)]]
while (1 << len(up)) <= max(depth) + 1:
last = up[-1]
up.append([last[last[v]] for v in range(total)])
def lca(u: int, v: int) -> int:
if depth[u] < depth[v]:
u, v = v, u
diff = depth[u] - depth[v]
j = 0
while diff:
if diff & 1:
u = up[j][u]
diff >>= 1
j += 1
if u == v:
return u
for j in range(len(up) - 1, -1, -1):
if up[j][u] != up[j][v]:
u, v = up[j][u], up[j][v]
return tparent[u]
# Faces of the triangulation: every face is a triangle.
tpos = [0] * nd
for a in range(total):
for i, d in enumerate(new_rot[a]):
tpos[d] = i
face = [-1] * nd
nf = 0
for d0 in range(nd):
if face[d0] >= 0:
continue
d = d0
while face[d] < 0:
face[d] = nf
ring = new_rot[head[d]]
d = ring[tpos[d ^ 1] - 1]
nf += 1
dual: list[list[int]] = [[] for _ in range(nf)]
for d in range(0, nd, 2):
if not in_tree[d]:
dual[face[d]].append(d)
dual[face[d ^ 1]].append(d ^ 1)
owner_w = [0] * nf # @line L10: the dual tree of non-tree edges, with subtree sums
for v in range(total):
od = tdart[v] if tdart[v] >= 0 else new_rot[v][0]
owner_w[face[od]] += weight[v]
fparent_dart = [-1] * nf # the dart whose left face is the child face
tin = [0] * nf
tout = [0] * nf
sub_w = owner_w[:]
sub_n = [1] * nf
clock = 0
visited = [False] * nf
visited[0] = True
stack2: list[tuple[int, int]] = [(0, 0)]
tin[0] = clock
clock += 1
while stack2:
f, i = stack2.pop()
if i < len(dual[f]):
stack2.append((f, i + 1))
d = dual[f][i]
g = face[d ^ 1]
if not visited[g]:
visited[g] = True
fparent_dart[g] = d ^ 1
tin[g] = clock
clock += 1
stack2.append((g, 0))
else:
tout[f] = clock
if f != 0:
p = face[fparent_dart[f] ^ 1]
sub_w[p] += sub_w[f]
sub_n[p] += sub_n[f]
best = -1 # @line L11: weights inside and on each cycle; the candidate with fewest faces
best_size = nf + 1
best_inside = 0
for g in range(1, nf):
d = fparent_dart[g] # the inside face g lies on the left of d = (u, v)
u, v = tail[d], head[d]
z = lca(u, v)
on_cycle = pre[u] + pre[v] - 2 * pre[z] + weight[z]
oz = tdart[z] if tdart[z] >= 0 else new_rot[z][0]
fz = face[oz]
z_inside = tin[g] <= tin[fz] < tout[g]
owners_on_cycle = pre[v] - pre[z] + (weight[z] if z_inside else 0)
inside = sub_w[g] - owners_on_cycle
if 3 * (inside + on_cycle) >= wtotal and sub_n[g] < best_size:
best, best_size, best_inside = g, sub_n[g], inside
assert 3 * best_inside <= 2 * wtotal # the lemma of the separator theorem
d = fparent_dart[best]
u, v = tail[d], head[d]
z = lca(u, v)
cycle = [z] # @line L12: the vertices of G on the cycle
for s in (u, v):
while s != z:
cycle.append(s)
s = tparent[s]
return [names[a] for a in cycle if a < h and names[a] is not None]
def ring_index(rot: Rotation, v: int, w: int) -> int:
return rot[v].index(w)
def lipton_tarjan_separator(rot: Rotation) -> tuple[list[int], list[int], list[int]]:
"""Return (A, B, S) for a planar graph given by a planar rotation system.
rot[v] lists the neighbours of v in counterclockwise order. The result is
a partition of the vertices with no edge between A and B, len(A) and
len(B) at most 2n/3, and len(S) at most 2 * sqrt(2) * sqrt(n).
Raises ValueError if rot is not a planar rotation system of a simple graph.
"""
_check_planar(rot)
vertices = list(rot)
n = len(vertices)
if n == 0:
return [], [], []
comps = _components(vertices, rot, set())
comp = max(comps, key=len)
sep: list[int] = []
if 3 * len(comp) > 2 * n: # @line L1
nk = len(comp)
q = comp[0]
level = {q: 0} # @line L2
parent: dict[int, int] = {}
levels: list[list[int]] = [[q]]
queue = deque([q])
while queue:
u = queue.popleft()
for w in rot[u]:
if w not in level:
level[w] = level[u] + 1
parent[w] = u
if level[w] == len(levels):
levels.append([])
levels[level[w]].append(w)
queue.append(w)
r = len(levels) - 1
size = [len(lv) for lv in levels]
def lsize(i: int) -> int:
return size[i] if 0 <= i <= r else 0
count = 0
l1 = 0
while 2 * (count + size[l1]) < nk:
count += size[l1]
l1 += 1
k = count + size[l1]
if size[l1] * size[l1] <= 8 * nk: # @line L3
sep = list(levels[l1])
else:
l0 = l1 # @line L4
while not _within(lsize(l0) + 2 * (l1 - l0), k):
l0 -= 1
l2 = l1 + 1 # @line L5
while not _within(lsize(l2) + 2 * (l2 - l1 - 1), nk - k):
l2 += 1
sep = [v for i in (l0, l2) if 0 <= i <= r for v in levels[i]] # @line L6
m_size = sum(lsize(i) for i in range(l0 + 1, l2))
if 3 * m_size > 2 * nk:
sep += _cycle_separator(comp, rot, level, parent, l0, l2, q)
pieces = _components(vertices, rot, set(sep)) # @line L13
a, b = _group(pieces, n)
return sorted(a), sorted(b), sorted(sep) # @line L14
test_lipton_tarjan_separator.py 182 lines
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/lipton-tarjan-separator/. The authors are listed in prov/YCV3MB.jsonl.
import itertools
import random
import pytest
from lipton_tarjan_separator import lipton_tarjan_separator
Rotation = dict[int, list[int]]
def check(rot: Rotation) -> tuple[list[int], list[int], list[int]]:
"""Run the algorithm and verify the output independently."""
a, b, s = lipton_tarjan_separator(rot)
n = len(rot)
assert sorted(a + b + s) == sorted(rot)
side_b = set(b)
for u in a:
for w in rot[u]:
assert w not in side_b
assert 3 * len(a) <= 2 * n and 3 * len(b) <= 2 * n
assert len(s) * len(s) <= 8 * n
return a, b, s
def fan(n: int) -> Rotation:
rot: Rotation = {0: list(range(n - 1, 0, -1))}
for i in range(1, n):
rot[i] = [j for j in (i + 1, i - 1) if 1 <= j < n] + [0]
return rot
def grid(k: int, first: int = 0) -> Rotation:
"""The k x k grid; neighbours counterclockwise: right, up, left, down."""
rot: Rotation = {}
cells = [(i, j) for i in range(k) for j in range(k)]
cells.sort(key=lambda c: c != divmod(first, k))
for i, j in cells:
around = []
for di, dj in ((0, 1), (1, 0), (0, -1), (-1, 0)):
if 0 <= i + di < k and 0 <= j + dj < k:
around.append((i + di) * k + j + dj)
rot[i * k + j] = around
return rot
def stacked(rng: random.Random, n: int) -> Rotation:
"""A random stacked triangulation with n >= 3 vertices."""
rot: Rotation = {0: [1, 2], 1: [2, 0], 2: [0, 1]}
faces = [(0, 1, 2), (0, 2, 1)]
for v in range(3, n):
i = rng.randrange(len(faces))
x, y, z = faces[i]
faces[i] = (x, y, v)
faces += [(y, z, v), (z, x, v)]
for p, q in ((x, y), (y, z), (z, x)):
rot[p].insert(rot[p].index(q) + 1, v)
rot[v] = [x, y, z]
return rot
def random_planar(rng: random.Random, n: int) -> Rotation:
rot = stacked(rng, n) if n >= 3 else {v: [] for v in range(n)}
if n == 2 and rng.random() < 0.5:
rot = {0: [1], 1: [0]}
for _ in range(rng.randint(0, 2 * n)):
v = rng.randrange(n)
if rot[v]:
w = rng.choice(rot[v])
rot[v].remove(w)
rot[w].remove(v)
perm = list(range(n))
rng.shuffle(perm)
return {perm[v]: [perm[w] for w in rot[v]] for v in rng.sample(range(n), n)}
def least_separator(rot: Rotation) -> int:
"""Brute force: the least size of a separator with sides of at most 2n/3."""
vs = list(rot)
n = len(vs)
for size in range(n + 1):
for s in itertools.combinations(vs, size):
seen = set(s)
parts = []
for v in vs:
if v in seen:
continue
seen.add(v)
stack, count = [v], 0
while stack:
u = stack.pop()
count += 1
for w in rot[u]:
if w not in seen:
seen.add(w)
stack.append(w)
parts.append(count)
sums = {0}
for c in parts:
sums |= {x + c for x in sums}
total = sum(parts)
if any(3 * x <= 2 * n and 3 * (total - x) <= 2 * n for x in sums):
return size
return n
def test_page_example() -> None:
rot = fan(11)
assert rot[0] == [10, 9, 8, 7, 6, 5, 4, 3, 2, 1] and rot[5] == [6, 4, 0] and rot[1] == [2, 0]
assert lipton_tarjan_separator(rot) == ([4, 5, 6, 7, 8, 9], [1, 2], [0, 3, 10])
def test_small_cases() -> None:
assert lipton_tarjan_separator({}) == ([], [], [])
assert lipton_tarjan_separator({7: []}) == ([], [], [7])
assert lipton_tarjan_separator({0: [1], 1: [0]}) == ([1], [], [0])
k4 = {0: [1, 2, 3], 1: [2, 0, 3], 2: [0, 1, 3], 3: [0, 2, 1]}
a, b, s = check(k4)
assert s == [1, 2, 3] and a == [0] # the level of the neighbours of 0
star = {0: [1, 2, 3, 4, 5], 1: [0], 2: [0], 3: [0], 4: [0], 5: [0]}
assert check(star)[2] == [1, 2, 3, 4, 5] # the median level, not the least separator
path = {i: [j for j in (i - 1, i + 1) if 0 <= j < 9] for i in range(9)}
check(path)
isolated = {v: [] for v in range(5)}
assert check(isolated)[2] == []
def test_grid_from_center() -> None:
a, b, s = check(grid(3, first=4))
assert s == [1, 3, 5, 7] and a == [0, 2, 4] and b == [6, 8]
for k in range(2, 30):
check(grid(k))
def test_k5_and_k33_are_rejected() -> None:
k5 = {v: [w for w in range(5) if w != v] for v in range(5)}
with pytest.raises(ValueError):
lipton_tarjan_separator(k5)
k33 = {v: ([3, 4, 5] if v < 3 else [0, 1, 2]) for v in range(6)}
with pytest.raises(ValueError):
lipton_tarjan_separator(k33)
toroidal_k4 = {0: [1, 2, 3], 1: [0, 2, 3], 2: [0, 1, 3], 3: [0, 1, 2]}
with pytest.raises(ValueError):
lipton_tarjan_separator(toroidal_k4)
def test_fans_use_the_cycle_step() -> None:
for n in range(11, 200, 7):
a, b, s = check(fan(n))
assert 0 in s and len(s) == 3
def test_random_against_oracle() -> None:
rng = random.Random(2024)
for trial in range(600):
n = rng.randint(0, 9) if trial < 200 else rng.randint(9, 160)
rot = random_planar(rng, n)
a, b, s = check(rot)
if n <= 9:
assert least_separator(rot) <= len(s)
def test_disjoint_unions() -> None:
rng = random.Random(7)
for _ in range(100):
parts = [stacked(rng, rng.randint(3, 40)) for _ in range(rng.randint(1, 4))]
rot: Rotation = {}
for part in parts:
base = len(rot)
for v, around in part.items():
rot[base + v] = [base + w for w in around]
check(rot)
def test_large_inputs() -> None:
rng = random.Random(1)
a, b, s = check(stacked(rng, 20000))
assert len(s) <= 400
check(fan(5000))
check(grid(100))
What the next level needs
- the tests did not run in this build