Intractability/Approximation Algorithms

Lesson 12.45,680 words

Approximation Algorithms

When a problem is NP\mathsf{NP}-hard we can still ask for a solution provably close to optimal. This lesson makes the approximation ratio ρ\rho precise, separates absolute from relative guarantees, and proves the ratios of four classic algorithms: greedy set cover (HnlnnH_n \approx \ln n), the MST-doubling 22-approximation for metric TSP, load balancing, and the knapsack FPTAS.

╌╌╌╌

The previous lesson, Coping with NP-Hardness, surveyed four honest responses to a problem we cannot solve exactly and quickly: approximate, use heuristics, pay for a smart exponential search, or exploit special structure. It worked the first response through a single example, the -approximation for vertex cover. This lesson develops that response in detail. We define the approximation ratio with care, separate two flavors of guarantee, and then prove ratios for a sequence of canonical algorithms, each illustrating a different proof technique. We finish by classifying which problems can be approximated arbitrarily well, and which cannot be approximated at all.

The recurring move, established last time, is worth restating because every proof below is a variation on it: we bound our output not against the optimum , which we cannot compute, but against a surrogate quantity, a lower bound for a minimization (or an upper bound for a maximization) that we can reason about. Most of approximation analysis is finding the right surrogate.

What approximate means precisely

Fix an optimization problem. On an instance , let be the cost of an optimal solution and the cost of the solution our algorithm returns. There are two distinct ways to say is close to optimal, and conflating them causes endless confusion.

Absolute guarantees are rare: most -hard problems scale, so an additive error that holds on small instances cannot hold on large ones. (If you could always come within of the optimal graph coloring, you could take disjoint copies of an instance and divide out the additive slack to solve it exactly.) So the relative ratio is the standard guarantee, and we pin it down exactly.

Two conventions are in use. Some authors quote the maximization ratio as a fraction (a -approximation returns at least half of optimal); others keep throughout. They say the same thing; we will use for minimization and the fraction for maximization, and state which we mean each time.

A ratio may be a constant ( for vertex cover), may grow with the input ( for set cover, as we are about to prove), or may be tunable to any value above at the cost of running time (an approximation scheme). The last section organizes these into a hierarchy.

Greedy set cover: a logarithmic ratio

Our first proof gives a ratio that is not constant — it grows like the logarithm of the input — and it introduces the charging argument, where we account for the algorithm's cost by distributing it onto the elements it serves.

Set cover is -hard, and it is the abstract form of many resource-selection problems: pick the fewest cell towers covering every neighborhood, the fewest tests catching every fault, the fewest committees spanning every skill. The natural greedy rule is simple: repeatedly take the set that covers the most still-uncovered elements.1

Algorithm 1:Greedy-Set-Cover(U,F)\textsc{Greedy-Set-Cover}(U, \mathcal{F}) — repeatedly grab the largest uncovered set
  1. 1
    C\mathcal{C} \gets \emptyset
  2. 2
    RUR \gets U
    elements still uncovered
  3. 3
    while RR \neq \emptyset do
  4. 4
    pick SFS \in \mathcal{F} maximizing SR|S \cap R|
    most new coverage
  5. 5
    CC{S}\mathcal{C} \gets \mathcal{C} \cup \{S\}
  6. 6
    RRSR \gets R \setminus S
    remove newly covered elements
  7. 7
    return C\mathcal{C}

Each iteration covers at least one new element, so the loop runs at most times; finding the best set is polynomial, so the whole algorithm is. The figure shows the first greedy choice: among the available sets, the one covering the largest part of the uncovered universe wins.

One step of : among candidate sets, pick the one covering the most still-uncovered elements (here with vs. and ).
Greedy run over three rounds: each round grabs the set covering the most still-uncovered elements (green = covered so far), until the universe is full

Greedy can use more sets than necessary, and the worst case is genuinely logarithmic. The classic bad instance is a universe of elements that the optimum covers with two sets, but greedy peels off in steps, each grabbing barely more than half of what remains. The precise bound uses the harmonic number .

The logarithmic factor is not an artifact of the analysis: greedy is essentially optimal here. A deep result shows that, unless , no polynomial-time algorithm approximates set cover within .2 Greedy is the best we can hope for, up to lower-order terms.

greedy_set_cover.pypython
from collections.abc import Hashable
from math import log
from typing import Optional, TypeVar

Element = TypeVar("Element", bound=Hashable)

def greedy_set_cover(
  universe: set[Element],
  family: list[set[Element]],
) -> list[int]:
  """
    A subfamily of `family` whose union is `universe`, chosen greedily.\n
    Returns the indices (into `family`) of the chosen sets, in pick order.\n
    Raises ValueError when `family` cannot cover `universe`.\n
  """
  # reject up front if the whole family still leaves some element uncovered.
  coverable: set[Element] = set()
  coverable.update(*family)
  if not universe <= coverable:
    raise ValueError("family does not cover the universe")

  remaining: set[Element] = set(universe)
  chosen: list[int] = []

  while remaining:

    # the set biting off the largest chunk of what is still uncovered.
    best_index: Optional[int] = None
    best_gain: int = 0
    for index, candidate in enumerate(family):
      gain: int = len(candidate & remaining)
      if gain > best_gain:
        best_gain = gain
        best_index = index

    # no set adds new coverage, yet elements remain: impossible to finish.
    if best_index is None:
      raise ValueError("family does not cover the universe")

    # commit the winner and drop everything it covers from the leftovers.
    chosen.append(best_index)
    remaining -= family[best_index]

  return chosen

def harmonic(count: int) -> float:
  """
    The harmonic number H_n = 1 + 1/2 + ... + 1/n (zero for n <= 0).\n
    This is the greedy set-cover approximation ratio.\n
  """
  return sum(1.0 / term for term in range(1, count + 1))

def cover_ratio_bound(universe_size: int) -> float:
  """
    The greedy guarantee H_n <= ln n + 1 expressed as the looser ln bound.\n
  """
  if universe_size <= 0:
    return 0.0
  return log(universe_size) + 1.0

Metric TSP: the MST-doubling 2-approximation

The traveling-salesman problem asks for a minimum-cost cycle visiting every vertex once. In full generality it admits no constant-factor approximation (established at the end of this lesson). But almost every TSP that arises from geography or physical distance satisfies the triangle inequality, — going direct is never longer than detouring through a third point. This metric case admits a -approximation built on the minimum spanning tree.3

The coping lesson sketched this construction; here we develop it in full and prove the ratio. The whole idea fits in one line: an MST is a cheap connected backbone, and a tour is just a connected backbone with the extra demand of being a single cycle, so the tour cannot cost much more than the tree.

Algorithm 2:MST-TSP(G,c)\textsc{MST-TSP}(G, c) — metric TSP tour within 2×2\times optimal
  1. 1
    TT \gets minimum spanning tree of GG
    e.g. Prim or Kruskal
  2. 2
    pick any vertex rr as root
  3. 3
    LL \gets vertices in the order first visited by a preorder walk of TT from rr
  4. 4
    return the cycle that visits the vertices in the order LL, then closes back to rr

The algorithm never literally doubles edges; the preorder walk is the clean way to realize the double then shortcut intuition. Doubling each tree edge makes every degree even, so an Euler tour exists that traverses every (doubled) edge once and has cost exactly . Walking that Euler tour and skipping vertices already seen yields precisely the preorder sequence . The figure shows both halves: the MST backbone, then the shortcut Hamiltonian cycle.

Metric-TSP -approximation. Left: MST (cost ). Right: shortcut the doubled Euler walk into a Hamiltonian cycle of cost .
Doubling each MST edge makes every degree even, so an Euler walk traverses all edges at cost ; shortcutting past repeats yields the tour

Both inequalities are tight in isolation, so is the honest guarantee for this particular algorithm. It is not the best ratio known: Christofides' algorithm replaces the wasteful doubling with a cleverer parity fix. Instead of duplicating every tree edge, it adds a minimum-weight perfect matching on just the odd-degree vertices of — the only vertices whose degree blocks an Euler tour. That matching costs at most , and the same shortcutting argument then yields a tour of cost at most .4 For decades the resulting was the best constant known for metric TSP.

mst_tsp.pypython
from collections.abc import Hashable
from typing import Callable, Generic, TypeVar

Point = TypeVar("Point", bound=Hashable)

# a symmetric distance function on the points (must obey the triangle ineq.).
DistanceFn = Callable[[Point, Point], float]

class MSTNode(Generic[Point]):
  """
    One vertex in the spanning tree: its point and its tree children.\n
    The children are filled as Prim attaches each new vertex to its parent.\n
  """

  def __init__(self, point: Point) -> None:
    self.point: Point = point
    self.children: list[MSTNode[Point]] = []

  def __repr__(self) -> str:
    return f"MSTNode({self.point!r})"

def _minimum_spanning_tree(
  points: list[Point],
  distance: DistanceFn[Point],
) -> MSTNode[Point]:
  """
    A rooted MST of the complete graph on `points` built by Prim's rule.\n
    Returns the root node; tree edges hang as parent -> children links.\n
  """
  root = MSTNode(points[0])
  nodes: dict[Point, MSTNode[Point]] = {points[0]: root}

  # closest tree vertex (its node) to each outside point, with that distance.
  closest_node: dict[Point, MSTNode[Point]] = {
    point: root for point in points[1:]
  }
  closest_distance: dict[Point, float] = {
    point: distance(points[0], point) for point in points[1:]
  }
  outside: set[Point] = set(points[1:])

  while outside:

    # pull in the outside point nearest to the current tree.
    nearest: Point = min(outside, key=lambda candidate: closest_distance[candidate])
    parent: MSTNode[Point] = closest_node[nearest]
    attached = MSTNode(nearest)
    parent.children.append(attached)
    nodes[nearest] = attached
    outside.discard(nearest)

    # relax every remaining outside point against the freshly added vertex.
    for point in outside:
      candidate_distance: float = distance(nearest, point)
      if candidate_distance < closest_distance[point]:
        closest_distance[point] = candidate_distance
        closest_node[point] = attached

  return root

def _preorder(root: MSTNode[Point]) -> list[Point]:
  """
    The points of the tree in preorder (root before its subtrees).\n
    This realizes the doubled-then-shortcut Euler walk.\n
  """
  # iterative DFS from the root, emitting each node as it is popped.
  order: list[Point] = []
  stack: list[MSTNode[Point]] = [root]
  while stack:
    node: MSTNode[Point] = stack.pop()
    order.append(node.point)

    # push children reversed so they pop in their listed order.
    stack.extend(reversed(node.children))
  return order

def mst_tsp(points: list[Point], distance: DistanceFn[Point]) -> list[Point]:
  """
    A Hamiltonian tour over `points` of cost at most 2 * OPT on a metric.\n
    Returns the visiting order (each point once); close it back to the first\n
    point to form the cycle. `distance` must be symmetric and obey the\n
    triangle inequality for the 2-approximation guarantee to hold.\n
  """
  if not points:
    return []
  root: MSTNode[Point] = _minimum_spanning_tree(points, distance)
  return _preorder(root)

def tour_cost(tour: list[Point], distance: DistanceFn[Point]) -> float:
  """
    The total length of the closed cycle visiting `tour` in order.\n
    Empty and single-vertex tours have cost zero.\n
  """
  if len(tour) < 2:
    return 0.0

  # sum each leg's length, wrapping the last point back to the first.
  total: float = 0.0
  for index in range(len(tour)):
    current: Point = tour[index]
    following: Point = tour[(index + 1) % len(tour)]
    total += distance(current, following)
  return total

def euclidean_distance(
  left: tuple[float, float],
  right: tuple[float, float],
) -> float:
  """
    Plane distance between two points — a metric that fits the guarantee.\n
  """
  return ((left[0] - right[0]) ** 2 + (left[1] - right[1]) ** 2) ** 0.5

Load balancing: a second 2-approximation, by averaging

To show a different surrogate at work, consider makespan scheduling: assign jobs with processing times to identical machines so as to minimize the makespan, the finishing time of the busiest machine. Deciding the optimum is -hard. The list-scheduling greedy assigns each job, in turn, to whichever machine is currently least loaded.

Algorithm 3:List-Schedule(p1pn,m)\textsc{List-Schedule}(p_1 \dots p_n, m) — greedy makespan within 2×2\times optimal
  1. 1
    load[1m]0\text{load}[1 \dots m] \gets 0
    current load of each machine
  2. 2
    for j1j \gets 1 to nn do
  3. 3
    iargminkload[k]i \gets \arg\min_{k} \text{load}[k]
    least-loaded machine
  4. 4
    assign job jj to machine ii
  5. 5
    load[i]load[i]+pj\text{load}[i] \gets \text{load}[i] + p_j
  6. 6
    return the assignment
List-scheduling on machines: each job (a bar of height ) lands on the currently shortest stack; the tallest stack is the makespan

This greedy needs only two facts that lower-bound any schedule, including the optimum.

The same two surrogates, sharpened, show that sorting jobs longest-first before list-scheduling improves the ratio to — a reminder that the order in which a greedy commits often matters more than the greedy rule itself.5

list_scheduling.pypython
import heapq
from typing import NamedTuple

class Schedule(NamedTuple):
  """
    The result of placing jobs: the per-machine job lists and the makespan.\n
    `assignment[i]` holds the processing times routed to machine i.\n
  """
  assignment: list[list[float]]
  makespan: float

def list_schedule(
  processing_times: list[float],
  machines: int,
) -> Schedule:
  """
    Greedily assign each job in input order to the least-loaded machine.\n
    Returns the assignment and resulting makespan, which is <= 2 * OPT.\n
    Raises ValueError when `machines` is not positive.\n
  """
  return _schedule(processing_times, machines, sort_longest_first=False)

def longest_processing_time(
  processing_times: list[float],
  machines: int,
) -> Schedule:
  """
    List scheduling after sorting jobs longest-first (the LPT rule).\n
    Committing the big jobs early sharpens the guarantee to <= 4/3 * OPT.\n
  """
  return _schedule(processing_times, machines, sort_longest_first=True)

def _schedule(
  processing_times: list[float],
  machines: int,
  sort_longest_first: bool,
) -> Schedule:
  """
    Shared placement loop for both the plain and LPT orderings.\n
  """
  if machines <= 0:
    raise ValueError("need at least one machine")

  # one bucket per machine; under LPT the big jobs are committed first.
  assignment: list[list[float]] = [[] for _ in range(machines)]
  order: list[float] = list(processing_times)
  if sort_longest_first:
    order.sort(reverse=True)

  # min-heap of (current load, machine index): the top is least loaded.
  loads: list[tuple[float, int]] = [
    (0.0, index) for index in range(machines)
  ]
  heapq.heapify(loads)

  # route each job to the least-loaded machine, then re-heap with its new load.
  for job in order:
    load, machine_index = heapq.heappop(loads)
    assignment[machine_index].append(job)
    heapq.heappush(loads, (load + job, machine_index))

  # the makespan is the finishing time of the busiest machine.
  makespan: float = max(
    (sum(jobs) for jobs in assignment),
    default=0.0,
  )
  return Schedule(assignment=assignment, makespan=makespan)

def optimal_makespan(
  processing_times: list[float],
  machines: int,
) -> float:
  """
    The exact minimum makespan by exhaustive assignment — the reference.\n
    Only for small inputs: it tries every machine for every job.\n
  """
  if machines <= 0:
    raise ValueError("need at least one machine")
  if not processing_times:
    return 0.0

  # best[0] is the best complete-schedule makespan found; loads tracks the dive.
  best: list[float] = [float("inf")]
  loads: list[float] = [0.0 for _ in range(machines)]

  def place(index: int) -> None:

    # every job placed: record this schedule's makespan if it beats the best.
    if index == len(processing_times):
      best[0] = min(best[0], max(loads))
      return

    # prune branches that already exceed the best full schedule seen.
    if max(loads) >= best[0]:
      return

    job: float = processing_times[index]
    seen_empty: bool = False
    for machine_index in range(machines):

      # symmetric empty machines are interchangeable: try just the first.
      if loads[machine_index] == 0.0:
        if seen_empty:
          continue
        seen_empty = True

      # try this machine, recurse, then undo for the next choice.
      loads[machine_index] += job
      place(index + 1)
      loads[machine_index] -= job

  place(0)
  return best[0]

A knapsack FPTAS: approximation you can dial in

The algorithms so far hit a fixed ratio. For some problems, any accuracy is achievable: a - (or -) approximation, paying in running time. The 0/1 knapsack problem is the textbook case, and its scheme is fully polynomial — polynomial in both the input size and .6

Recall knapsack: items with weights and profits , a capacity ; choose a subset of weight at most maximizing total profit. There is an exact dynamic program indexed by profit running in time , where . That is fast when profits are small but blows up when they are large — it is pseudo-polynomial, polynomial in the value , not in its bit-length.

The FPTAS exploits this dependence. Large profits are what inflate the profit-indexed table, so we scale them down and round, shrinking the table to a controllable size while distorting the answer only slightly.

Algorithm 4:Knapsack-FPTAS(w,p,W,ε)\textsc{Knapsack-FPTAS}(w, p, W, \varepsilon)(1ε)(1-\varepsilon)-optimal profit
  1. 1
    PmaxipiP \gets \max_i p_i
    largest single profit
  2. 2
    μεPn\mu \gets \dfrac{\varepsilon P}{n}
    scaling unit
  3. 3
    for i1i \gets 1 to nn do
  4. 4
    p^ipi/μ\hat p_i \gets \big\lfloor p_i / \mu \big\rfloor
    round profits down
  5. 5
    solve knapsack exactly on profits p^\hat p by the profit-indexed DP
  6. 6
    return that chosen subset (valued by the true pip_i)

Rounding the profits to multiples of caps the scaled profits at , so the profit-indexed DP now runs in time — polynomial in and in . It remains to bound the lost profit.

The two error sources are now visible. Rounding down loses at most per chosen item, hence at most in total; the scaling choice makes that total exactly . Smaller means finer , a larger table, and a slower-but-sharper answer — the running time grows as , never as an exponential.

knapsack_fptas.pypython
from typing import NamedTuple, Sequence

class KnapsackItem(NamedTuple):
  """
    One knapsack item: the profit it earns and the weight it costs.\n
  """
  profit: float
  weight: int

def _profit_indexed_knapsack(
  items: Sequence[KnapsackItem],
  scaled_profits: list[int],
  capacity: int,
) -> list[int]:
  """
    Exact 0/1 knapsack indexed by integer profit, returning chosen indices.\n
    `min_weight[value]` is the least weight achieving exactly that scaled\n
    profit; we sweep items and pick the largest value that still fits.\n
  """
  total_profit: int = sum(scaled_profits)

  # least weight to reach each scaled-profit value; infinity if unreachable.
  min_weight: list[float] = [float("inf") for _ in range(total_profit + 1)]
  min_weight[0] = 0

  # remember which item realized each (value, achievable) state to backtrack.
  chosen_at: list[list[int]] = [[] for _ in range(total_profit + 1)]

  for index, item in enumerate(items):
    item_profit: int = scaled_profits[index]
    if item_profit <= 0:
      continue

    # high-to-low so each item is used at most once.
    for value in range(total_profit, item_profit - 1, -1):
      candidate_weight: float = min_weight[value - item_profit] + item.weight
      if candidate_weight < min_weight[value]:
        min_weight[value] = candidate_weight
        chosen_at[value] = chosen_at[value - item_profit] + [index]

  # the best attainable scaled profit whose realizing subset fits.
  best_value: int = 0
  for value in range(total_profit + 1):
    if min_weight[value] <= capacity:
      best_value = value
  return chosen_at[best_value]

def knapsack_fptas(
  items: Sequence[KnapsackItem],
  capacity: int,
  epsilon: float,
) -> list[int]:
  """
    Indices of a feasible subset with profit at least (1-epsilon)*OPT.\n
    Items heavier than `capacity` are discarded first (they fit no subset),\n
    profits are scaled and floored, then solved exactly. `epsilon` must lie\n
    in (0, 1].\n
  """
  if not 0.0 < epsilon <= 1.0:
    raise ValueError("epsilon must be in (0, 1]")

  # only items that fit on their own can ever be chosen.
  feasible: list[int] = [
    index for index, item in enumerate(items) if item.weight <= capacity
  ]
  if not feasible:
    return []

  # no positive profit anywhere means the empty set is optimal.
  feasible_items: list[KnapsackItem] = [items[index] for index in feasible]
  largest_profit: float = max(item.profit for item in feasible_items)
  if largest_profit <= 0:
    return []

  # scale by mu = eps*P/n and floor, then solve the small instance exactly.
  scaling_unit: float = epsilon * largest_profit / len(feasible_items)
  scaled_profits: list[int] = [
    int(item.profit // scaling_unit) for item in feasible_items
  ]
  local_choice: list[int] = _profit_indexed_knapsack(
    feasible_items, scaled_profits, capacity
  )

  # translate indices back to positions in the original `items` sequence.
  return [feasible[index] for index in local_choice]

def subset_profit(items: Sequence[KnapsackItem], chosen: list[int]) -> float:
  """
    The total true profit of the items at indices `chosen`.\n
  """
  return sum(items[index].profit for index in chosen)

def subset_weight(items: Sequence[KnapsackItem], chosen: list[int]) -> int:
  """
    The total weight of the items at indices `chosen`.\n
  """
  return sum(items[index].weight for index in chosen)

The hierarchy, and where approximation runs out

Naming the qualitative differences we have seen organizes the whole subject.

An FPTAS is the strongest guarantee: accuracy becomes a parameter. The guarantees form a hierarchy, each level strictly stronger than the one below.

Approximability ladder, strongest to weakest. FPTAS (knapsack) gives any ratio in time poly in and ; PTAS (Euclidean TSP) is poly in only; constant ratio (vertex cover , metric TSP ); growing ratio (set cover, ); inapproximable (general TSP) admits no constant ratio unless .

The bottom level is real: some problems admit no constant-factor approximation at all, and the proofs are themselves reductions. The canonical example is the general traveling-salesman problem, with no triangle inequality.7

The gap-creating trick — pricing non-edges so high that any decent approximation is forced to avoid them — is the prototype for hardness-of-approximation. Its modern descendant, the PCP theorem, constructs such gaps for a vast range of problems, settling the exact approximation threshold of MAX-3SAT, set cover, and many others, and explaining why the ratios proved above so often turn out to be the best possible.

tsp_gap_reduction.pypython
from __future__ import annotations

from collections.abc import Hashable
from itertools import permutations
from typing import Callable, Optional, TypeVar

from graph import Graph

Label = TypeVar("Label", bound=Hashable)

# a tour solver: given labels and a symmetric cost(a, b), return a visit order.
TourSolver = Callable[[list[Label], Callable[[Label, Label], float]], list[Label]]

def gap_instance(
  graph: Graph[Label],
  ratio: float,
) -> tuple[list[Label], Callable[[Label, Label], float]]:
  """
    The complete TSP instance encoding `graph` for a target ratio `ratio`.\n
    Returns the vertex labels and a cost function: 1 on edges of `graph`,\n
    and ratio*n + 1 on non-edges, so any tour below ratio*n uses only real\n
    edges. `graph` is treated as undirected.\n
  """
  labels: list[Label] = graph.labels
  vertex_count: int = len(labels)
  expensive: float = ratio * vertex_count + 1.0

  # adjacency lookup over the original graph's edges (both directions).
  adjacency: dict[Label, set[Label]] = {label: set() for label in labels}
  for edge in graph.edges():
    adjacency[edge.source.label].add(edge.target.label)
    adjacency[edge.target.label].add(edge.source.label)

  def cost(first: Label, second: Label) -> float:
    if second in adjacency[first]:
      return 1.0
    return expensive

  return labels, cost

def decides_hamiltonicity_via_tsp(
  graph: Graph[Label],
  ratio: float,
  solver: TourSolver[Label],
) -> bool:
  """
    Whether `graph` is Hamiltonian, decided through a rho-approximate TSP.\n
    Builds the gap instance, runs `solver` (assumed a ratio-approximation),\n
    and reports True exactly when the returned tour stays below ratio*n.\n
  """
  labels, cost = gap_instance(graph, ratio)
  if len(labels) < 2:
    return True

  tour: list[Label] = solver(labels, cost)
  length: float = _tour_length(tour, cost)
  threshold: float = ratio * len(labels)
  return length < threshold

def has_hamiltonian_cycle(graph: Graph[Label]) -> bool:
  """
    Whether `graph` has a cycle visiting every vertex exactly once.\n
    Brute force over vertex orderings — the reference for small graphs.\n
  """
  # empty and single-vertex graphs are vacuously Hamiltonian.
  labels: list[Label] = graph.labels
  vertex_count: int = len(labels)
  if vertex_count <= 1:
    return True

  # undirected adjacency over the original edges.
  adjacency: dict[Label, set[Label]] = {label: set() for label in labels}
  for edge in graph.edges():
    adjacency[edge.source.label].add(edge.target.label)
    adjacency[edge.target.label].add(edge.source.label)

  # fix the first vertex; permute the rest to break the cyclic symmetry.
  first: Label = labels[0]
  for tail in permutations(labels[1:]):
    cycle: tuple[Label, ...] = (first,) + tail
    if _is_cycle(cycle, adjacency):
      return True
  return False

def _is_cycle(
  order: tuple[Label, ...],
  adjacency: dict[Label, set[Label]],
) -> bool:
  """
    Whether consecutive vertices in `order` (closing back) are all adjacent.\n
  """
  # walk each vertex to its successor (wrapping at the end); all must be edges.
  count: int = len(order)
  for index in range(count):
    current: Label = order[index]
    following: Label = order[(index + 1) % count]
    if following not in adjacency[current]:
      return False
  return True

def brute_force_tsp(
  labels: list[Label],
  cost: Callable[[Label, Label], float],
) -> list[Label]:
  """
    The exact minimum-cost tour by trying every ordering — a TSP solver.\n
    Being optimal, it is in particular a rho-approximation for every rho,\n
    so it slots into `decides_hamiltonicity_via_tsp`.\n
  """
  if len(labels) < 2:
    return list(labels)

  # fix the first vertex, permute the rest, and keep the cheapest closed tour.
  first: Label = labels[0]
  best_order: Optional[list[Label]] = None
  best_length: float = float("inf")
  for tail in permutations(labels[1:]):
    order: list[Label] = [first, *tail]
    length: float = _tour_length(order, cost)
    if length < best_length:
      best_length = length
      best_order = order

  assert best_order is not None
  return best_order

def _tour_length(
  tour: list[Label],
  cost: Callable[[Label, Label], float],
) -> float:
  """
    The total cost of the closed cycle visiting `tour` in order.\n
  """
  if len(tour) < 2:
    return 0.0

  # sum each step's cost, wrapping the last vertex back to the first.
  total: float = 0.0
  for index in range(len(tour)):
    current: Label = tour[index]
    following: Label = tour[(index + 1) % len(tour)]
    total += cost(current, following)
  return total
graph.pypython
from collections.abc import Hashable, Iterator
from typing import Generic, Optional, TypeVar


Label = TypeVar("Label", bound=Hashable)


class Edge(Generic[Label]):
  """
    A directed connection from `source` to `target`, carrying a weight.\n
  """

  def __init__(
    self,
    source: Vertex[Label],
    target: Vertex[Label],
    weight: float = 1.0,
  ) -> None:
    self.source: Vertex[Label] = source
    self.target: Vertex[Label] = target
    self.weight: float = weight

  def __repr__(self) -> str:
    return f"Edge({self.source.label!r} -> {self.target.label!r}, w={self.weight})"


class Vertex(Generic[Label]):
  """
    A graph vertex: a label plus the list of edges leaving it.\n
  """

  def __init__(self, label: Label) -> None:
    self.label: Label = label
    self.outgoing: list[Edge[Label]] = []

  def neighbors(self) -> list[Vertex[Label]]:
    """
      The vertices reachable from this one by a single edge.\n
    """
    return [edge.target for edge in self.outgoing]

  def edge_to(self, label: Label) -> Optional[Edge[Label]]:
    """
      The outgoing edge to the vertex with `label`, or None.\n
    """
    for edge in self.outgoing:
      if edge.target.label == label:
        return edge
    return None

  def __repr__(self) -> str:
    return f"Vertex({self.label!r})"


class Graph(Generic[Label]):
  """
    A graph of Vertex objects linked by Edge objects.\n
    Pass `directed=True` for a digraph; otherwise each `add_edge` inserts\n
    the reverse edge too.\n
  """

  def __init__(self, directed: bool = False) -> None:
    self.directed: bool = directed
    self._vertices: dict[Label, Vertex[Label]] = {}

  def add_vertex(self, label: Label) -> Vertex[Label]:
    """
      Return the vertex for `label`, creating it if it is absent.\n
    """
    # reuse the existing vertex, or mint and register a fresh one.
    vertex = self._vertices.get(label)
    if vertex is None:
      vertex = Vertex(label)
      self._vertices[label] = vertex
    return vertex

  def add_edge(
    self,
    source_label: Label,
    target_label: Label,
    weight: float = 1.0,
  ) -> None:
    """
      Connect two labels (creating either vertex as needed).\n
      Adds the reverse edge as well when the graph is undirected.\n
    """
    source = self.add_vertex(source_label)
    target = self.add_vertex(target_label)

    # link source to target, and mirror it back when undirected.
    source.outgoing.append(Edge(source, target, weight))
    if not self.directed:
      target.outgoing.append(Edge(target, source, weight))

  def vertex(self, label: Label) -> Vertex[Label]:
    """
      The vertex carrying `label` (raises KeyError if absent).\n
    """
    return self._vertices[label]

  @property
  def vertices(self) -> list[Vertex[Label]]:
    """
      Every vertex, in insertion order.\n
    """
    return list(self._vertices.values())

  @property
  def labels(self) -> list[Label]:
    """
      Every vertex label, in insertion order.\n
    """
    return list(self._vertices)

  def edges(self) -> Iterator[Edge[Label]]:
    """
      Each edge once — an undirected edge is yielded a single time.\n
    """
    # track undirected endpoint pairs so each is emitted only once.
    seen: set[frozenset[Label]] = set()

    for vertex in self._vertices.values():
      for edge in vertex.outgoing:
        # skip an undirected edge already yielded from the other endpoint.
        if not self.directed:
          endpoints = frozenset((edge.source.label, edge.target.label))
          if endpoints in seen:
            continue
          seen.add(endpoints)

        yield edge

  def __contains__(self, label: Label) -> bool:
    return label in self._vertices

  def __iter__(self) -> Iterator[Vertex[Label]]:
    return iter(self._vertices.values())

  def __len__(self) -> int:
    return len(self._vertices)

Tight thresholds and the approximability frontier

The gap reduction above shows some ratio is unachievable; the modern result is that the guarantees this lesson proved are often exactly optimal — the algorithm and the hardness bound meet at the same constant. Three landmark results draw that frontier:

  • Set cover is , tight. The greedy ratio from earlier in this lesson cannot be improved: Feige (1998) proved that a -approximation for set cover would imply , essentially .8 No polynomial algorithm for set cover can do better than greedy.
  • MAX-3SAT is , tight. A trivial random assignment satisfies of the clauses of a MAX-3SAT instance in expectation. Håstad (2001) proved via the PCP machinery that beating is -hard — so the simplest possible algorithm is already optimal.9
  • Dinur's combinatorial PCP. The original PCP theorem's proof was famously intricate; Dinur (2007) gave a different gap-amplification proof that builds the required hardness gap incrementally, making inapproximability accessible without the heavy algebra.10
Algorithm and hardness meet: for set cover and MAX-3SAT the proved approximation ratio equals the proved inapproximability threshold.

One frontier remains genuinely open. The Unique Games Conjecture (Khot, 2002) posits the hardness of a specific constraint problem; if true, it would fix the exact approximation threshold of a wide range of problems at once, including vertex cover (making our factor- optimal) and MAX-CUT.11 Whether UGC holds is one of the central open questions of complexity theory today.

Takeaways

  • An approximation algorithm runs in polynomial time with a provable bound on how far from optimal it can land. The relative ratio (cost for minimization, for maximization) is the standard guarantee; absolute additive guarantees are rare because hard problems scale.
  • Every proof bounds the output against a surrogate for the unknown optimum: a charging argument for greedy set cover (, and no polynomial algorithm beats ), the MST for metric TSP ( gives a -approximation; Christofides sharpens it to ), and two universal lower bounds (average load, biggest job) for list-scheduling makespan (-approximation).
  • A PTAS achieves any ratio in time polynomial in ; an FPTAS is also polynomial in . The knapsack FPTAS scales and rounds profits to shrink the pseudo-polynomial DP, losing at most .
  • Some problems are inapproximable: unless , general TSP has no constant-factor approximation, shown by pricing non-edges to create a gap that would otherwise decide Hamiltonicity. Hardness of approximation is proved by reductions, just like -hardness itself.

Footnotes

  1. CLRS, Ch. 35 — Approximation Algorithms (§35.3): the greedy set-cover heuristic and its approximation bound via the per-element charging argument.
  2. CLRS, Ch. 35 — Approximation Algorithms: greedy set cover is essentially optimal; under no polynomial algorithm achieves ratio .
  3. CLRS, Ch. 35 — Approximation Algorithms (§35.2): the MST-doubling -approximation for metric (triangle-inequality) TSP, with and shortcutting bounded by the triangle inequality.
  4. Skiena, §11.10 — Approximation Algorithms: Christofides' -approximation for metric TSP, fixing odd-degree parity with a minimum-weight perfect matching instead of doubling the whole tree.
  5. CLRS, Ch. 35 — Approximation Algorithms (Problem 35-5): list scheduling for makespan is a -approximation; the longest-processing-time rule improves the bound to .
  6. CLRS, Ch. 35 — Approximation Algorithms (§35.5): the fully polynomial-time approximation scheme for the subset-sum / knapsack problem by trimming and scaling, giving in time polynomial in and .
  7. Erickson, Ch. 12 — NP-Hardness: gap reductions and inapproximability, including the proof that general TSP admits no constant-factor approximation unless .
  8. Uriel Feige, A Threshold of for Approximating Set Cover, Journal of the ACM 45(4), 1998 — the matching hardness lower bound making greedy optimal for set cover.
  9. Johan Håstad, Some Optimal Inapproximability Results, Journal of the ACM 48(4), 2001 — the tight threshold for MAX-3SAT and other optimal PCP-based inapproximability results.
  10. Irit Dinur, The PCP Theorem by Gap Amplification, Journal of the ACM 54(3), 2007 — a combinatorial gap-amplification proof of the PCP theorem.
  11. Subhash Khot, On the Power of Unique 2-Prover 1-Round Games, STOC 2002 — the Unique Games Conjecture, which if true fixes the optimal approximation ratio for vertex cover, MAX-CUT, and many others.
Practice

╌╌ END ╌╌