"""
combinatorial_explosion.py -- Lighthill's core argument made numerical.
Demonstrates the combinatorial explosion argument from the 1973 Lighthill
Report on Artificial Intelligence, and shows how heuristic search (A*) and
iterative-deepening A* (IDA*) reduce the growth rate in practice.
The domain is the 15-puzzle (4x4 sliding tile puzzle), a classic benchmark
for heuristic search with known optimal solution depths.
Run: python combinatorial_explosion.py
Reference:
Lighthill, J. (1973). Artificial Intelligence: A General Survey.
Science Research Council of Great Britain.
Hart, P.E., Nilsson, N.J., Raphael, B. (1968). A Formal Basis for the
Heuristic Determination of Minimum Cost Paths. IEEE Trans. SSC, 4(2).
Korf, R.E. (1985). Depth-first iterative-deepening: An optimal admissible
tree search. Artificial Intelligence, 27(1), 97-109.
"""
from __future__ import annotations
import heapq
from collections import deque
from typing import Optional
# ---------------------------------------------------------------------------
# 15-puzzle representation
# ---------------------------------------------------------------------------
# State: tuple of 16 ints, 0 = blank tile. Index = row*4 + col.
State = tuple[int, ...]
GOAL_STATE: State = tuple(range(16)) # 0,1,2,...,15; 0 is blank in top-left
MOVES = { # blank position -> list of neighbour positions
i: [
n
for n in [i - 4, i + 4, i - 1 if i % 4 > 0 else -1, i + 1 if i % 4 < 3 else 16]
if 0 <= n < 16
]
for i in range(16)
}
def successors(state: State) -> list[State]:
"""Generate all states reachable from state in one move."""
blank = state.index(0)
result = []
for nbr in MOVES[blank]:
lst = list(state)
lst[blank], lst[nbr] = lst[nbr], lst[blank]
result.append(tuple(lst))
return result
def manhattan(state: State) -> int:
"""
Sum of Manhattan distances of each tile from its goal position.
Admissible heuristic for A* on the 15-puzzle.
h(state) <= true distance to goal for all states.
"""
dist = 0
for idx, tile in enumerate(state):
if tile == 0:
continue
goal_row, goal_col = tile // 4, tile % 4
cur_row, cur_col = idx // 4, idx % 4
dist += abs(cur_row - goal_row) + abs(cur_col - goal_col)
return dist
def make_puzzle(depth: int, seed: int = 42) -> State:
"""
Create a puzzle exactly `depth` moves from the goal
by doing `depth` random (non-reversing) moves.
"""
import random
rng = random.Random(seed)
state = GOAL_STATE
prev_blank = -1
for _ in range(depth):
blank = state.index(0)
nbrs = [n for n in MOVES[blank] if n != prev_blank]
nbr = rng.choice(nbrs)
lst = list(state)
lst[blank], lst[nbr] = lst[nbr], lst[blank]
prev_blank = blank
state = tuple(lst)
return state
# ---------------------------------------------------------------------------
# Breadth-first search (no heuristic -- bare exponential growth)
# ---------------------------------------------------------------------------
def bfs(
start: State, goal: State, max_nodes: int = 200_000
) -> tuple[int, Optional[int]]:
"""
BFS on the 15-puzzle.
Returns (nodes_expanded, solution_depth) or (max_nodes, None) if limit hit.
"""
if start == goal:
return 0, 0
queue: deque[tuple[State, int]] = deque([(start, 0)])
visited: set[State] = {start}
count = 0
while queue:
state, depth = queue.popleft()
count += 1
if count >= max_nodes:
return count, None
for nxt in successors(state):
if nxt == goal:
return count, depth + 1
if nxt not in visited:
visited.add(nxt)
queue.append((nxt, depth + 1))
return count, None
# ---------------------------------------------------------------------------
# A* with Manhattan heuristic
# ---------------------------------------------------------------------------
def astar(
start: State, goal: State, max_nodes: int = 200_000
) -> tuple[int, Optional[int]]:
"""
A* on the 15-puzzle using Manhattan distance.
Returns (nodes_expanded, solution_depth) or (max_nodes, None).
"""
h0 = manhattan(start)
# heap entries: (f, g, state)
heap: list[tuple[int, int, State]] = [(h0, 0, start)]
g_cost: dict[State, int] = {start: 0}
count = 0
while heap:
f, g, state = heapq.heappop(heap)
count += 1
if count >= max_nodes:
return count, None
if state == goal:
return count, g
if g > g_cost.get(state, float("inf")):
continue # stale entry
for nxt in successors(state):
ng = g + 1
if ng < g_cost.get(nxt, float("inf")):
g_cost[nxt] = ng
heapq.heappush(heap, (ng + manhattan(nxt), ng, nxt))
return count, None
# ---------------------------------------------------------------------------
# Iterative-deepening A* (Korf 1985)
# ---------------------------------------------------------------------------
def idastar(
start: State, goal: State, max_nodes: int = 500_000
) -> tuple[int, Optional[int]]:
"""
IDA* on the 15-puzzle.
Memory-efficient: stores only the current path.
Returns (nodes_expanded, solution_depth) or (max_nodes, None).
"""
threshold = manhattan(start)
counter = [0]
def search(path: list[State], g: int, bound: int) -> int | float:
state = path[-1]
f = g + manhattan(state)
if f > bound:
return f
if state == goal:
return -1 # found
minimum = float("inf")
for nxt in successors(state):
if nxt in path:
continue
counter[0] += 1
if counter[0] >= max_nodes:
return float("inf")
path.append(nxt)
t = search(path, g + 1, bound)
if t == -1:
return -1
if t < minimum:
minimum = t
path.pop()
return minimum
while True:
counter[0] = 0
path: list[State] = [start]
result = search(path, 0, threshold)
if result == -1:
return counter[0], threshold
if result == float("inf") or counter[0] >= max_nodes:
return counter[0], None
threshold = int(result)
# ---------------------------------------------------------------------------
# Theoretical growth rates
# ---------------------------------------------------------------------------
def theoretical_bfs_nodes(depth: int, branching: float = 2.13) -> float:
"""
Expected BFS nodes at depth d for a tree with average branching factor b.
BFS expands all nodes up to depth d: sum_{i=0}^{d} b^i = (b^{d+1}-1)/(b-1)
Branching factor for 15-puzzle is approximately 2.13 (corners have 2,
edges have 3, interior has 4; blank is uniform over all positions).
"""
return (branching ** (depth + 1) - 1) / (branching - 1)
# ---------------------------------------------------------------------------
# Main: comparison table
# ---------------------------------------------------------------------------
if __name__ == "__main__":
print("=" * 70)
print("Combinatorial Explosion: Lighthill (1973) argument, made numerical")
print("15-puzzle search comparison: BFS vs A* vs IDA*")
print("=" * 70)
print("""
Lighthill's core claim: for a problem with n objects, the state space
has order n! or 2^n elements. Any uninformed search expands this space
exhaustively. Heuristic search reduces expansion but cannot escape the
worst case. The table below shows nodes expanded for increasing puzzle
depths (= solution lengths), demonstrating exponential growth for BFS
and dramatically reduced growth for A* and IDA*.
""")
depths = [4, 6, 8, 10, 12, 14, 16, 18, 20]
print(
f" {'Depth':>5} {'Theoretical BFS':>16} {'Actual BFS':>11} "
f"{'A* nodes':>9} {'IDA* nodes':>11} {'A*/BFS':>7}"
)
print(
f" {'-----':>5} {'---------------':>16} {'-----------':>11} "
f"{'--------':>9} {'-----------':>11} {'------':>7}"
)
for d in depths:
puzzle = make_puzzle(d)
theory = int(theoretical_bfs_nodes(d))
# BFS (capped at 200K for deep puzzles)
bfs_count, bfs_sol = bfs(puzzle, GOAL_STATE, max_nodes=200_000)
bfs_str = f"{bfs_count:>10,}" if bfs_sol is not None else ">200,000"
# A*
as_count, as_sol = astar(puzzle, GOAL_STATE, max_nodes=200_000)
as_str = f"{as_count:>8,}" if as_sol is not None else ">200,000"
# IDA*
ida_count, ida_sol = idastar(puzzle, GOAL_STATE, max_nodes=500_000)
ida_str = f"{ida_count:>10,}" if ida_sol is not None else ">500,000"
ratio = (
(as_count / bfs_count)
if bfs_count > 0 and bfs_sol is not None
else float("nan")
)
ratio_str = f"{ratio:.3f}" if not (ratio != ratio) else " ---"
theory_str = f"{theory:>15,}" if theory < 10**13 else f" ~{theory:.2e}"
print(
f" {d:>5} {theory_str} {bfs_str:>11} {as_str:>9} "
f"{ida_str:>11} {ratio_str:>7}"
)
print("""
Interpretation:
- BFS grows roughly as (2.13)^depth -- pure exponential
- A* nodes grow far more slowly due to the Manhattan distance heuristic
pruning branches that cannot improve on the current best path
- IDA* achieves comparable node counts with O(d) memory vs O(b^d) for A*
- Even with heuristics, very deep puzzles (depth > 50) exceed any
computer's capacity: this is Lighthill's point. The 15-puzzle has
solutions up to depth 80; complete optimal search requires IDA*
with pattern database heuristics (Korf & Felner 2002) -- far beyond
the techniques available in 1973.
""")
# Show growth rates explicitly for the first few depths
print("=" * 70)
print("Raw node counts side by side (first 6 depths):")
print("=" * 70)
for d in depths[:6]:
puzzle = make_puzzle(d, seed=7)
bfs_n, _ = bfs(puzzle, GOAL_STATE, 50_000)
as_n, _ = astar(puzzle, GOAL_STATE, 50_000)
ida_n, _ = idastar(puzzle, GOAL_STATE, 50_000)
bfs_bar = "#" * min(40, bfs_n // 500)
as_bar = "#" * min(40, as_n // 500)
ida_bar = "#" * min(40, ida_n // 500)
print(f" depth={d:>2}")
print(f" BFS [{bfs_bar:<40}] {bfs_n:,}")
print(f" A* [{as_bar:<40}] {as_n:,}")
print(f" IDA* [{ida_bar:<40}] {ida_n:,}")