Branch-and-cut and Concorde

Read first

• Post 5: Branch and bound. A search tree splits tours into groups, a lower bound on each group lets you prune it, and the bound’s quality decides how much of the tree you explore.
• Post 6: Integer programming formulations. The TSP as an integer program. Post 6 wrote it with one variable per directed leg; here I use its symmetric version: one 0/1 variable per edge, a degree constraint saying every stop has exactly two tour edges, and subtour elimination constraints (the Dantzig-Fulkerson-Johnson, or DFJ, form) saying every proper subset S of stops must have at least two tour edges leaving it. The LP relaxation lets the variables take any value between 0 and 1, and its optimum is a lower bound.

The question this post answers: the DFJ formulation has exponentially many subtour constraints, so how can anyone solve it?

Intuition

The DFJ model is correct but enormous. For n stops, every subset S with at least two stops on each side gives a subtour constraint. Each constraint and its mirror image (S and the stops outside S) say the same thing, so the count of distinct constraints is about 2ⁿ⁻¹.

Writing down 10¹⁵ constraints is out of the question. The idea of this post is that you don’t have to. Almost all of those constraints are already satisfied by the solution you would get without them. So start with a small model, solve it, look for a constraint the solution breaks, add only that one, and solve again. This is the cutting-plane loop. An added constraint is a cut, because it cuts the current solution out of the feasible region.

I will run it by hand on six of our stops: A, B, D, E, G and H. The coordinates are in kilometres and the distances in units of 100 m, rounded, as in every post so far.

Round 1: degree constraints only. Ask for the cheapest set of edges that gives every stop exactly two of them, with no subtour constraints at all. The answer is two triangles: A-B-G-A (41 + 45 + 36 = 122) and D-E-H-D (36 + 50 + 28 = 114), total 236. Every stop has degree 2, but it is two loops, not one tour.

Separation. Now find a subtour constraint that this solution breaks. That search has a name, the separation problem, and here it is easy: the triangle {D, E, H} has zero edges leaving it, and its subtour constraint demands at least two. Add that single constraint.

Round 2: degree constraints plus one cut. Solve again. The cheapest answer now is the tour A-B-D-H-E-G-A, cost 259. It is a single loop, so no subtour constraint is broken, and the separation step finds nothing. A solution that satisfies every constraint of the full model, found by solving a relaxation of it, is optimal for the full model. One cut out of 25 was enough.

Step through the diagram below to watch the loop. The drawn lengths are unrounded kilometres, so they differ from the matrix totals by a few tenths.

The matrix below prices the final tour. Each outlined cell is one leg, and the strip adds them in driving order.

On all eight stops the loop is even shorter. Round 1 already returns the optimal tour A-B-C-D-H-E-F-G-A at 315, so it needs no cuts at all. That is luck of the geometry, not a rule. On berlin52 it took one round of cuts, as the table further down shows.

The formal idea

The model. For a symmetric instance with stops V, write E for the set of edges, cₑ for the cost of edge e, and xₑ for its variable. For a subset S of V, δ(S) is the set of edges with exactly one end in S, and x(δ(S)) is the sum of their variables. The DFJ formulation from Post 6 is:

minimize    Σ  cₑ · xₑ
           e∈E

subject to  x(δ({v})) = 2          for every stop v           (degree)
            x(δ(S))   ≥ 2          for every S ⊂ V, 2 ≤ |S| ≤ n − 2   (subtour)
            xₑ ∈ {0, 1}            for every edge e

Cutting planes. Replace the full constraint list with a small working list L, which starts empty, and relax xₑ ∈ {0, 1} to 0 ≤ xₑ ≤ 1. Then repeat:

1. Solve the LP with the degree constraints and the cuts in L. Call its solution x*.
2. Separation: find a set S with x*(δ(S)) < 2. If you find one, add it to L and go back to step 1.
3. If no such S exists and x* is a 0/1 tour, it is optimal.
4. If no such S exists but x* is fractional, you are stuck with this family of cuts. Either look for another family, or branch.

Separating subtour constraints. If x* is a 0/1 solution, its edges form disjoint cycles, and each cycle that is not the whole tour is a violated set S. If x* is fractional, treat x*ₑ as the capacity of edge e. A violated subtour constraint exists exactly when the minimum cut of that weighted graph is below 2. Minimum cuts can be found in polynomial time, for example with the Stoer-Wagner algorithm (1997), so the exponential size of the constraint family is no obstacle: you only ever look at the cuts that matter.

Branch-and-cut. Step 4 is where Post 5 comes back. When the LP solution is fractional and no more cuts can be found, branch on a fractional edge: one child sets xₑ = 0, the other xₑ = 1. Each child runs its own cutting-plane loop, and the LP value at each node is the lower bound for pruning. Branch and bound with a cutting-plane loop at every node is branch-and-cut. Padberg and Rinaldi named it in 1987, when they used it to solve a 532-city instance. The same year they solved a 2,392-city instance with it, published in 1991.

The diagram below shows how the pieces fit together.

flowchart TD
  S["Start: degree constraints only"] --> LP["Solve the LP"]
  LP --> SEP{"Separation: a violated cut?"}
  SEP -- yes --> ADD["Add the cut"]
  ADD --> LP
  SEP -- no --> INT{"Is the LP solution a tour?"}
  INT -- yes --> DONE["Optimal at this node"]
  INT -- no --> BR["Branch on a fractional edge"]
  BR --> LP

Beyond subtours: combs. Subtour constraints are not enough on their own. Some fractional solutions satisfy every one of them and are still not tours. The standard example has six stops in two triangles, {U, V, W} and {X, Y, Z}. Put xₑ = 1/2 on the six triangle edges and xₑ = 1 on the three edges U-X, V-Y and W-Z. Every stop has degree 1/2 + 1/2 + 1 = 2, and every one of the 62 possible cuts carries at least 2, yet the solution is half of two different tours. The diagram below shows it.

A comb inequality catches it. A comb has a handle H, a set of stops, and an odd number k ≥ 3 of teeth T₁, …, Tₖ, each a set that overlaps the handle and also reaches outside it, with the teeth disjoint from each other. Every tour satisfies

x(δ(H)) + Σᵢ x(δ(Tᵢ)) ≥ 3k + 1

With handle {U, V, W} and teeth {U, X}, {V, Y}, {W, Z}, the left side is 3 + 2 + 2 + 2 = 9, and the right side is 3 · 3 + 1 = 10. The fractional solution breaks the comb, so adding it as a cut moves the LP away. Combs were introduced by Chvátal (1973) and generalized by Grötschel and Padberg (1979). They are one of several cut families that make the LP bound tighter than subtour constraints alone.

How Concorde does it. Concorde, written by David Applegate, Robert Bixby, Vašek Chvátal and William Cook, is the branch-and-cut code that holds the records for exact TSP solving. Its main pieces are below. Each one fixes a reason the simple loop above would not scale.

flowchart TD
  T["Good tour first: chained Lin-Kernighan"] --> UB["Upper bound"]
  E["Core LP: a small set of likely edges"] --> LP["Solve LP"]
  LP --> C["Separation: subtours, combs, more families"]
  C --> LP
  LP --> P["Pricing: add any missing edge that could help"]
  P --> LP
  LP --> B{"Fractional and no more cuts?"}
  B -- yes --> BR["Branch"]
  BR --> LP
  B -- no --> OPT["Optimal tour, with proof"]
  UB --> B

With these pieces Concorde has solved all 110 instances in TSPLIB, the largest being pla85900 with 85,900 cities, in a computation run in 2005 and 2006. The 81,998-stop Korea tour from Post 2 was proven optimal in 2025 the same way: the LKH heuristic found the tour, and Concorde’s cutting-plane code proved no shorter one exists.

The milestones below show how far the idea has carried, from the original paper to today.

Dantzig, Fulkerson and Johnson took one city in each of the 48 states plus Washington, D.C. They removed seven East Coast cities, solved the 42-city problem that was left, and noted that its optimal tour used the Washington-Boston edge, whose shortest road passes through all seven. That gave the 49-city tour too.

Common misconceptions

“An exponential number of constraints makes the model useless.” The model is never written out in full. The cutting-plane loop only generates the constraints that the current LP solution breaks, and the separation problem for subtour constraints is a minimum cut, which is fast. On berlin52, 4 cuts out of about 2.25 × 10¹⁵ were enough.

“Cutting planes replace branching.” For most instances they shrink the search tree enormously, but they rarely remove it. Once the LP solution is fractional and the separation routines find nothing, the only way forward is to branch. Branch-and-cut is branch and bound with a much better bound, not an alternative to it.

“If the LP solution satisfies every subtour constraint, it is a tour.” Only if it is also a 0/1 solution. The two-triangle example above satisfies all 62 cut conditions and every degree constraint with half-valued edges. That is why Concorde separates combs and other families too.

“Concorde is just a very fast brute force.” Concorde never enumerates tours. Its power comes from the LP bound: on many instances the LP with good cuts is within a fraction of a percent of the optimum, so the search tree that remains is small.

Try it in Python

This script runs the cutting-plane loop from the Intuition section on the six stops. It uses no LP solver. Instead, the “solve” step is a brute-force search for the cheapest set of disjoint loops that gives every stop two edges and respects every cut found so far. On an instance this small that gives the same answers as the LP. The separation step is the real thing: a loop that is not the whole tour is a violated subtour constraint. The last part of the script checks the two-triangle example against every subtour constraint and against the comb.

import math
from itertools import combinations, permutations

STOPS = {"A": (2, 2), "B": (6, 1), "D": (8, 7), "E": (5, 9), "G": (4, 5), "H": (10, 9)}  # km

def d(a, b):  # straight-line distance in units of 100 m, rounded
    (x1, y1), (x2, y2) = STOPS[a], STOPS[b]
    return int(10 * math.hypot(x1 - x2, y1 - y2) + 0.5)

def cycles_through(v, others):
    """Every cycle that starts at v and uses 2 or more of `others`."""
    for k in range(2, len(others) + 1):
        for chosen in combinations(others, k):
            for order in permutations(chosen):
                if order[0] < order[-1]:             # skip the same cycle driven backwards
                    yield (v, *order)

def edges(cycle):
    return {frozenset((cycle[i], cycle[(i + 1) % len(cycle)])) for i in range(len(cycle))}

def crossing(edge_set, S):  # how many chosen edges have exactly one end in S
    return sum(len(e & S) == 1 for e in edge_set)

def cheapest_cover(cuts):
    """Cheapest set of disjoint cycles covering every stop (degree 2 everywhere)
    that respects every cut found so far. Brute force stands in for the LP/IP solver."""
    best = (math.inf, None)
    def grow(left, chosen, cost):
        nonlocal best
        if cost >= best[0]:
            return
        if not left:
            es = set().union(*(edges(c) for c in chosen))
            if all(crossing(es, S) >= 2 for S in cuts):
                best = (cost, chosen)
            return
        v, rest = min(left), sorted(left - {min(left)})
        for c in cycles_through(v, rest):
            grow(left - set(c), chosen + [c], cost + sum(d(a, b) for a, b in zip(c, c[1:] + c[:1])))
    grow(set(STOPS), [], 0)
    return best

# The cutting-plane loop: solve, look for a violated subtour constraint, add it, repeat
cuts = []
while True:
    cost, cycles = cheapest_cover(cuts)
    shown = " + ".join("-".join(c + c[:1]) for c in cycles)
    if len(cycles) == 1:
        print(f"round {len(cuts) + 1}: {cost:>4}  {shown}  -> one tour, optimal")
        break
    S = frozenset(cycles[-1])                      # separation: a subtour is a violated cut
    print(f"round {len(cuts) + 1}: {cost:>4}  {shown}  -> add cut: at least 2 edges leave {{{', '.join(sorted(S))}}}")
    cuts.append(S)

# A fractional point that passes every subtour test but fails a comb inequality
x = {frozenset(e): 0.5 for e in ["UV", "VW", "UW", "XY", "YZ", "XZ"]}
x.update({frozenset(e): 1.0 for e in ["UX", "VY", "WZ"]})
cut = lambda S: sum(v for e, v in x.items() if len(e & S) == 1)
nodes = "UVWXYZ"
weakest = min(cut(frozenset(S)) for k in range(1, 6) for S in combinations(nodes, k))
handle, teeth = frozenset("UVW"), [frozenset("UX"), frozenset("VY"), frozenset("WZ")]
print(f"smallest cut over all {2 ** 6 - 2} subsets: {weakest}  (subtour constraints need >= 2)")
print(f"comb: {cut(handle)} + {' + '.join(str(cut(t)) for t in teeth)} = "
      f"{cut(handle) + sum(cut(t) for t in teeth)}  (comb inequality needs >= {3 * len(teeth) + 1})")

It prints:

Change STOPS to all eight stops and round 1 prints the 315 tour straight away, with no cuts.

I also ran the same loop on berlin52 with an LP solver and a minimum-cut separation routine. The LP code isn’t shown, because this part of the series sticks to plain Python, but the numbers are worth seeing next to the lower bounds from Post 5.

The degree-only LP split berlin52 into 4 separate loops, 6 of its edges at fractional values. Adding one subtour cut for each loop, 4 cuts in all, and solving again gave a 0/1 solution: a single tour of length 7,542, TSPLIB’s optimum, proven optimal without any branching. 1,326 edge variables, 52 degree constraints and 4 cuts did what 2.25 × 10¹⁵ constraints describe.

Exercise

Round 1 on the six stops produced two triangles, and the loop added the cut for {D, E, H}. Show that the cut for {A, B, G} is exactly the same constraint. Then use that fact to confirm the table in the Intuition section: six stops have 25 distinct subtour constraints, and eight stops have 119.

Hint: which edges are in δ(S), and which are in δ of the stops outside S? For the count, take subsets of size 2 to n − 2 and pair each one with its complement.

What this method can’t do yet

Branch-and-cut is at its best on the pure TSP, because decades of research have mapped which cuts the TSP polytope needs and how to separate them fast. Add the rules real routes carry and that knowledge runs out.

Take the time windows from Post 3: reach B by minute 60 and F by minute 140. An integer program needs an arrival time tᵢ at each stop and a constraint like tⱼ ≥ tᵢ + d(i, j) − M(1 − xᵢⱼ) for every pair, where M is a large constant that switches the constraint off when the edge isn’t used. In the LP relaxation, xᵢⱼ can be 0.5, and with M = 1,000 the constraint becomes tⱼ ≥ tᵢ + d(i, j) − 500, which says nothing. The LP bound collapses, and there is no comb-like cut family that repairs it in general. Driver breaks, skills, vehicle-order compatibility and multi-trip reloads all add more of the same. The model is hard to write and its relaxation is weak. Post 8, “Constraint programming for routing”, takes a different approach: state the rules directly and let propagation remove impossible choices, instead of hoping an LP bound will catch them.

References

• Dantzig, G., Fulkerson, R., and Johnson, S. (1954). “Solution of a Large-Scale Traveling-Salesman Problem.” Journal of the Operations Research Society of America 2(4): 393–410.
• Padberg, M., and Rinaldi, G. (1987). “Optimization of a 532-City Symmetric Traveling Salesman Problem by Branch and Cut.” Operations Research Letters 6(1): 1–7.
• Padberg, M., and Rinaldi, G. (1991). “A Branch-and-Cut Algorithm for the Resolution of Large-Scale Symmetric Traveling Salesman Problems.” SIAM Review 33(1): 60–100.
• Grötschel, M., and Holland, O. (1991). “Solution of Large-Scale Symmetric Travelling Salesman Problems.” Mathematical Programming 51: 141–202.
• Chvátal, V. (1973). “Edmonds Polytopes and Weakly Hamiltonian Graphs.” Mathematical Programming 5: 29–40.
• Grötschel, M., and Padberg, M. W. (1979). “On the Symmetric Travelling Salesman Problem I: Inequalities.” Mathematical Programming 16: 265–280.
• Stoer, M., and Wagner, F. (1997). “A Simple Min-Cut Algorithm.” Journal of the ACM 44(4): 585–591.
• Applegate, D., Bixby, R., Chvátal, V., and Cook, W. (2003). “Implementing the Dantzig-Fulkerson-Johnson Algorithm for Large Traveling Salesman Problems.” Mathematical Programming 97: 91–153.
• Applegate, D. L., Bixby, R. E., Chvátal, V., and Cook, W. J. (2006). The Traveling Salesman Problem: A Computational Study. Princeton University Press.
• Concorde TSP Solver. University of Waterloo. https://www.math.uwaterloo.ca/tsp/concorde.html
• Cook, W., Espinoza, D., Goycoolea, M., and Helsgaun, K. (2025). “Korea TSPs.” University of Waterloo. https://www.math.uwaterloo.ca/tsp/korea/
• Reinelt, G. (1991). “TSPLIB—A Traveling Salesman Problem Library.” ORSA Journal on Computing 3(4): 376–384; berlin52 files and list of optimal solutions.

See the whole series · ← All writing