diff --git a/CLAUDE.md b/CLAUDE.md
index 1bab6b4..6636a85 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -36,7 +36,7 @@ Optimal-solver pattern databases are baked in-app (Settings › Solving): the 7-
 
 Two layers with a hard boundary:
 
-- **CubeKit/** — local Swift package, zero UI dependencies, all cube logic. `CubeState` is the cubie-level source of truth (corner/edge permutation + orientation, centers fixed). `FaceletCube` is the 54-sticker view with legality validation. Solvers: `KociembaSolver` (two-phase, coordinate tables cached on disk, time-budgeted search) and `BeginnerSolver` (layer-by-layer, returns a `StagedSolution`). Everything here is unit-tested; UI code never mutates cube state directly.
+- **CubeKit/** — local Swift package, zero UI dependencies, all cube logic. `CubeState` is the cubie-level source of truth (corner/edge permutation + orientation, centers fixed). `FaceletCube` is the 54-sticker view with legality validation. Solvers: `KociembaSolver` (two-phase, coordinate tables cached on disk, time-budgeted search), `ThistlethwaiteSolver` (four-phase group reduction, greedy descent over BFS distance tables, each phase provably shortest), and `BeginnerSolver` (layer-by-layer). The staged solvers return a `StagedSolution<StageKind>` (generic over the stage enum). Everything here is unit-tested; UI code never mutates cube state directly.
 - **cube/** — the app target (folder-synchronized Xcode group; new files under `cube/` join the target automatically). `AppModel` owns the logical state; `CubeSceneController` owns the RealityKit scene. The scene *renders and animates*; every completed turn flows back through `onMoveCommitted` into `AppModel.commit`, which is the only place `cubeState` advances.
 
 Key mechanisms to understand before touching play/scene code:
diff --git a/CubeKit/Sources/CubeKit/BeginnerSolver.swift b/CubeKit/Sources/CubeKit/BeginnerSolver.swift
index a056db9..ac759a8 100644
--- a/CubeKit/Sources/CubeKit/BeginnerSolver.swift
+++ b/CubeKit/Sources/CubeKit/BeginnerSolver.swift
@@ -1,5 +1,5 @@
 /// The stages of the classic layer-by-layer (beginner) method, in order.
-public enum BeginnerStage: String, CaseIterable, Sendable {
+public enum BeginnerStage: String, CaseIterable, Sendable, SolverStage {
     case bottomCross
     case bottomCorners
     case middleEdges
@@ -21,24 +21,13 @@ public enum BeginnerStage: String, CaseIterable, Sendable {
     }
 }
 
-/// A solution annotated with the method stage each move belongs to.
-public struct StagedSolution: Sendable {
-    public struct Stage: Sendable {
-        public let stage: BeginnerStage
-        public let moves: [Move]
-    }
-
-    public let stages: [Stage]
-    public var moves: [Move] { stages.flatMap(\.moves) }
-}
-
 /// Solves the way people are taught to: layer by layer, with short
 /// memorable algorithms. Solutions are long (~150–250 moves) but every
 /// stage has a understandable goal — the point is to follow along.
 public struct BeginnerSolver: Sendable {
     public init() {}
 
-    public func solve(_ state: CubeState) -> StagedSolution? {
+    public func solve(_ state: CubeState) -> StagedSolution<BeginnerStage>? {
         guard state.isLegal else { return nil }
         var worker = Worker(state: state)
         do {
@@ -59,7 +48,7 @@ private enum SolverFailure: Error {
 
 private struct Worker {
     var state: CubeState
-    var finishedStages: [StagedSolution.Stage] = []
+    var finishedStages: [StagedSolution<BeginnerStage>.Stage] = []
     private var currentMoves: [Move] = []
 
     init(state: CubeState) {
@@ -82,8 +71,7 @@ private struct Worker {
     ) rethrows {
         currentMoves = []
         try body(&self)
-        finishedStages.append(
-            StagedSolution.Stage(stage: stage, moves: Self.merged(currentMoves)))
+        finishedStages.append(.init(stage: stage, moves: Self.merged(currentMoves)))
         currentMoves = []
     }
 
diff --git a/CubeKit/Sources/CubeKit/SolverTables.swift b/CubeKit/Sources/CubeKit/SolverTables.swift
index 4bd6b89..485d193 100644
--- a/CubeKit/Sources/CubeKit/SolverTables.swift
+++ b/CubeKit/Sources/CubeKit/SolverTables.swift
@@ -35,7 +35,7 @@ public final class SolverTables: Sendable {
         // basic moves, flattened for the generators.
         let moves = CubeState.basicMoves
 
-        twistMove = Self.buildOrientationTable(
+        twistMove = TableBuilder.orientationTable(
             count: Coordinates.twistCount, pieceCount: 8, modulus: 3,
             decode: Coordinates.cornerOrientations(forTwist:),
             encode: { orientations in
@@ -47,7 +47,7 @@ public final class SolverTables: Sendable {
             moveOrientation: { moves[$0].cornerOrientation }
         )
 
-        flipMove = Self.buildOrientationTable(
+        flipMove = TableBuilder.orientationTable(
             count: Coordinates.flipCount, pieceCount: 12, modulus: 2,
             decode: Coordinates.edgeOrientations(forFlip:),
             encode: { orientations in
@@ -61,50 +61,29 @@ public final class SolverTables: Sendable {
 
         // Slice table: representative places slice edges at the ranked
         // positions; the filler arrangement is irrelevant to the result.
-        var sliceTable = [UInt16](repeating: 0, count: Coordinates.sliceCount * 18)
-        for coordinate in 0..<Coordinates.sliceCount {
-            let positions = Set(Coordinates.slicePositions(forSlice: coordinate))
-            var permutation = [Int](repeating: 0, count: 12)
-            var nextSlice = 8
-            var nextOther = 0
-            for slot in 0..<12 {
-                if positions.contains(slot) {
-                    permutation[slot] = nextSlice
-                    nextSlice += 1
-                } else {
-                    permutation[slot] = nextOther
-                    nextOther += 1
-                }
-            }
-            for move in 0..<18 {
-                let moved = Self.permutationApplying(
-                    permutation, moves[move / 3].edgePermutation, times: move % 3 + 1)
-                var rank = 0
-                var found = 0
-                for slot in 0..<12 where moved[slot] >= 8 {
-                    found += 1
-                    rank += Coordinates.binomial[slot][found]
-                }
-                sliceTable[coordinate * 18 + move] = UInt16(rank)
-            }
-        }
-        sliceMove = sliceTable
+        sliceMove = TableBuilder.occupancyTable(
+            slotCount: 12, markerCount: 4, totalSlots: 12, moves: Array(0..<18),
+            movePermutation: { moves[$0].edgePermutation }
+        )
 
-        cornerPermutationMove = Self.buildPermutationTable(
+        cornerPermutationMove = TableBuilder.permutationTable(
             count: Coordinates.cornerPermutationCount, pieceCount: 8,
+            moves: Self.phase2Moves,
             movePermutation: { moves[$0].cornerPermutation },
             project: { $0 }, embed: { $0 }
         )
 
-        udEdgePermutationMove = Self.buildPermutationTable(
+        udEdgePermutationMove = TableBuilder.permutationTable(
             count: Coordinates.udEdgePermutationCount, pieceCount: 8,
+            moves: Self.phase2Moves,
             movePermutation: { moves[$0].edgePermutation },
             project: { Array($0[0..<8]) },
             embed: { $0 + [8, 9, 10, 11] }
         )
 
-        sliceEdgePermutationMove = Self.buildPermutationTable(
+        sliceEdgePermutationMove = TableBuilder.permutationTable(
             count: Coordinates.sliceEdgePermutationCount, pieceCount: 4,
+            moves: Self.phase2Moves,
             movePermutation: { moves[$0].edgePermutation },
             project: { $0[8..<12].map { $0 - 8 } },
             embed: { Array(0..<8) + $0.map { $0 + 8 } }
@@ -118,10 +97,10 @@ public final class SolverTables: Sendable {
         let udEdgePermutationMove = self.udEdgePermutationMove
         let sliceEdgePermutationMove = self.sliceEdgePermutationMove
 
-        prune1SliceFlip = Self.breadthFirstDistances(
+        prune1SliceFlip = TableBuilder.breadthFirstDistances(
             stateCount: Coordinates.sliceCount * Coordinates.flipCount,
             moveCount: 18,
-            start: Coordinates.solvedSlice * Coordinates.flipCount
+            starts: [Coordinates.solvedSlice * Coordinates.flipCount]
         ) { state, move in
             let slice = state / Coordinates.flipCount
             let flip = state % Coordinates.flipCount
@@ -129,10 +108,10 @@ public final class SolverTables: Sendable {
                 + Int(flipMove[flip * 18 + move])
         }
 
-        prune1SliceTwist = Self.breadthFirstDistances(
+        prune1SliceTwist = TableBuilder.breadthFirstDistances(
             stateCount: Coordinates.sliceCount * Coordinates.twistCount,
             moveCount: 18,
-            start: Coordinates.solvedSlice * Coordinates.twistCount
+            starts: [Coordinates.solvedSlice * Coordinates.twistCount]
         ) { state, move in
             let slice = state / Coordinates.twistCount
             let twist = state % Coordinates.twistCount
@@ -140,10 +119,10 @@ public final class SolverTables: Sendable {
                 + Int(twistMove[twist * 18 + move])
         }
 
-        prune2SliceCorner = Self.breadthFirstDistances(
+        prune2SliceCorner = TableBuilder.breadthFirstDistances(
             stateCount: Coordinates.sliceEdgePermutationCount * Coordinates.cornerPermutationCount,
             moveCount: 10,
-            start: 0
+            starts: [0]
         ) { state, move in
             let slicePerm = state / Coordinates.cornerPermutationCount
             let cornerPerm = state % Coordinates.cornerPermutationCount
@@ -152,10 +131,10 @@ public final class SolverTables: Sendable {
                 + Int(cornerPermutationMove[cornerPerm * 10 + move])
         }
 
-        prune2SliceEdge = Self.breadthFirstDistances(
+        prune2SliceEdge = TableBuilder.breadthFirstDistances(
             stateCount: Coordinates.sliceEdgePermutationCount * Coordinates.udEdgePermutationCount,
             moveCount: 10,
-            start: 0
+            starts: [0]
         ) { state, move in
             let slicePerm = state / Coordinates.udEdgePermutationCount
             let edgePerm = state % Coordinates.udEdgePermutationCount
@@ -184,93 +163,6 @@ public final class SolverTables: Sendable {
         self.prune2SliceEdge = prune2SliceEdge
     }
 
-    // MARK: Generation helpers
-
-    private static func permutationApplying(
-        _ permutation: [Int], _ movePermutation: [Int], times: Int
-    ) -> [Int] {
-        var result = permutation
-        for _ in 0..<times {
-            var next = result
-            for i in 0..<result.count { next[i] = result[movePermutation[i]] }
-            result = next
-        }
-        return result
-    }
-
-    private static func buildOrientationTable(
-        count: Int, pieceCount: Int, modulus: Int,
-        decode: (Int) -> [Int], encode: ([Int]) -> Int,
-        movePermutation: (Int) -> [Int], moveOrientation: (Int) -> [Int]
-    ) -> [UInt16] {
-        var table = [UInt16](repeating: 0, count: count * 18)
-        for coordinate in 0..<count {
-            let orientations = decode(coordinate)
-            for move in 0..<18 {
-                let face = move / 3
-                let permutation = movePermutation(face)
-                let orientation = moveOrientation(face)
-                var current = orientations
-                for _ in 0..<(move % 3 + 1) {
-                    var next = current
-                    for i in 0..<pieceCount {
-                        next[i] = (current[permutation[i]] + orientation[i]) % modulus
-                    }
-                    current = next
-                }
-                table[coordinate * 18 + move] = UInt16(encode(current))
-            }
-        }
-        return table
-    }
-
-    /// Builds a phase-2 permutation move table. `project` extracts the
-    /// ranked sub-permutation from a full arrangement; `embed` rebuilds a
-    /// representative full arrangement from it.
-    private static func buildPermutationTable(
-        count: Int, pieceCount: Int,
-        movePermutation: (Int) -> [Int],
-        project: ([Int]) -> [Int], embed: ([Int]) -> [Int]
-    ) -> [UInt16] {
-        var table = [UInt16](repeating: 0, count: count * 10)
-        for coordinate in 0..<count {
-            let full = embed(Coordinates.unrankPermutation(coordinate, count: pieceCount))
-            for (index, moveValue) in phase2Moves.enumerated() {
-                let moved = permutationApplying(
-                    full, movePermutation(moveValue / 3), times: moveValue % 3 + 1)
-                table[coordinate * 10 + index] =
-                    UInt16(Coordinates.rankPermutation(project(moved)))
-            }
-        }
-        return table
-    }
-
-    private static func breadthFirstDistances(
-        stateCount: Int, moveCount: Int, start: Int,
-        neighbor: (Int, Int) -> Int
-    ) -> [Int8] {
-        var distances = [Int8](repeating: -1, count: stateCount)
-        var queue = [Int32]()
-        queue.reserveCapacity(stateCount)
-        distances[start] = 0
-        queue.append(Int32(start))
-        var head = 0
-        while head < queue.count {
-            let state = Int(queue[head])
-            head += 1
-            let next = distances[state] + 1
-            for move in 0..<moveCount {
-                let target = neighbor(state, move)
-                if distances[target] < 0 {
-                    distances[target] = next
-                    queue.append(Int32(target))
-                }
-            }
-        }
-        assert(!distances.contains(-1), "pruning table has unreachable states")
-        return distances
-    }
-
     // MARK: Disk cache
 
     private static let cacheVersion: UInt32 = 1
diff --git a/CubeKit/Sources/CubeKit/StagedSolution.swift b/CubeKit/Sources/CubeKit/StagedSolution.swift
new file mode 100644
index 0000000..a103c4c
--- /dev/null
+++ b/CubeKit/Sources/CubeKit/StagedSolution.swift
@@ -0,0 +1,25 @@
+/// A named phase of a multi-stage solving method.
+public protocol SolverStage: Hashable, Sendable {
+    var displayName: String { get }
+}
+
+/// A solution annotated with the method stage each move belongs to.
+public struct StagedSolution<StageKind: SolverStage>: Sendable {
+    public struct Stage: Sendable {
+        public let stage: StageKind
+        public let moves: [Move]
+
+        public init(stage: StageKind, moves: [Move]) {
+            self.stage = stage
+            self.moves = moves
+        }
+    }
+
+    public let stages: [Stage]
+
+    public init(stages: [Stage]) {
+        self.stages = stages
+    }
+
+    public var moves: [Move] { stages.flatMap(\.moves) }
+}
diff --git a/CubeKit/Sources/CubeKit/TableBuilder.swift b/CubeKit/Sources/CubeKit/TableBuilder.swift
new file mode 100644
index 0000000..31ff9e0
--- /dev/null
+++ b/CubeKit/Sources/CubeKit/TableBuilder.swift
@@ -0,0 +1,156 @@
+import Foundation
+
+/// Shared builders for coordinate move tables and breadth-first distance
+/// tables, used by both the Kociemba and Thistlethwaite table sets.
+enum TableBuilder {
+    static func permutationApplying(
+        _ permutation: [Int], _ movePermutation: [Int], times: Int
+    ) -> [Int] {
+        var result = permutation
+        for _ in 0..<times {
+            var next = result
+            for i in 0..<result.count { next[i] = result[movePermutation[i]] }
+            result = next
+        }
+        return result
+    }
+
+    /// Move table for an orientation coordinate, indexed
+    /// `[coordinate * 18 + move]` over all eighteen moves.
+    static func orientationTable(
+        count: Int, pieceCount: Int, modulus: Int,
+        decode: (Int) -> [Int], encode: ([Int]) -> Int,
+        movePermutation: (Int) -> [Int], moveOrientation: (Int) -> [Int]
+    ) -> [UInt16] {
+        var table = [UInt16](repeating: 0, count: count * 18)
+        for coordinate in 0..<count {
+            let orientations = decode(coordinate)
+            for move in 0..<18 {
+                let face = move / 3
+                let permutation = movePermutation(face)
+                let orientation = moveOrientation(face)
+                var current = orientations
+                for _ in 0..<(move % 3 + 1) {
+                    var next = current
+                    for i in 0..<pieceCount {
+                        next[i] = (current[permutation[i]] + orientation[i]) % modulus
+                    }
+                    current = next
+                }
+                table[coordinate * 18 + move] = UInt16(encode(current))
+            }
+        }
+        return table
+    }
+
+    /// Move table for a (sub-)permutation coordinate, indexed
+    /// `[coordinate * moves.count + moveIndex]` over the given `Move` raw
+    /// values. `project` extracts the ranked sub-permutation from a full
+    /// arrangement; `embed` rebuilds a representative full arrangement.
+    static func permutationTable(
+        count: Int, pieceCount: Int, moves: [Int],
+        movePermutation: (Int) -> [Int],
+        project: ([Int]) -> [Int], embed: ([Int]) -> [Int]
+    ) -> [UInt16] {
+        var table = [UInt16](repeating: 0, count: count * moves.count)
+        for coordinate in 0..<count {
+            let full = embed(Coordinates.unrankPermutation(coordinate, count: pieceCount))
+            for (index, moveValue) in moves.enumerated() {
+                let moved = permutationApplying(
+                    full, movePermutation(moveValue / 3), times: moveValue % 3 + 1)
+                table[coordinate * moves.count + index] =
+                    UInt16(Coordinates.rankPermutation(project(moved)))
+            }
+        }
+        return table
+    }
+
+    /// Move table for an occupancy coordinate: the colex rank of which of
+    /// the first `slotCount` slots (of a `totalSlots`-piece permutation)
+    /// hold the `markerCount` tracked pieces. Representatives place marker
+    /// values (≥ `totalSlots - markerCount`) at the occupied slots; every
+    /// allowed move must keep markers within the first `slotCount` slots.
+    static func occupancyTable(
+        slotCount: Int, markerCount: Int, totalSlots: Int, moves: [Int],
+        movePermutation: (Int) -> [Int]
+    ) -> [UInt16] {
+        let stateCount = Coordinates.binomial[slotCount][markerCount]
+        let markerFloor = totalSlots - markerCount
+        var table = [UInt16](repeating: 0, count: stateCount * moves.count)
+        for coordinate in 0..<stateCount {
+            let positions = Set(occupiedSlots(
+                forRank: coordinate, slots: slotCount, count: markerCount))
+            var permutation = [Int](repeating: 0, count: totalSlots)
+            var nextMarker = markerFloor
+            var nextFiller = 0
+            for slot in 0..<totalSlots {
+                if positions.contains(slot) {
+                    permutation[slot] = nextMarker
+                    nextMarker += 1
+                } else {
+                    permutation[slot] = nextFiller
+                    nextFiller += 1
+                }
+            }
+            for (index, moveValue) in moves.enumerated() {
+                let moved = permutationApplying(
+                    permutation, movePermutation(moveValue / 3), times: moveValue % 3 + 1)
+                var rank = 0
+                var found = 0
+                for slot in 0..<slotCount where moved[slot] >= markerFloor {
+                    found += 1
+                    rank += Coordinates.binomial[slot][found]
+                }
+                table[coordinate * moves.count + index] = UInt16(rank)
+            }
+        }
+        return table
+    }
+
+    /// The occupied slot list (ascending) for a colex occupancy rank.
+    static func occupiedSlots(forRank rank: Int, slots: Int, count: Int) -> [Int] {
+        var positions: [Int] = []
+        var r = rank
+        for i in stride(from: count, through: 1, by: -1) {
+            var p = slots - 1
+            while Coordinates.binomial[p][i] > r { p -= 1 }
+            positions.append(p)
+            r -= Coordinates.binomial[p][i]
+        }
+        return positions.reversed()
+    }
+
+    /// Exact distances from the `starts` set to every reachable state.
+    /// Unreached states stay -1; pass `requireFullCoverage: false` for
+    /// spaces that are legitimately only partially reachable.
+    static func breadthFirstDistances(
+        stateCount: Int, moveCount: Int, starts: [Int],
+        requireFullCoverage: Bool = true,
+        neighbor: (Int, Int) -> Int
+    ) -> [Int8] {
+        var distances = [Int8](repeating: -1, count: stateCount)
+        var queue = [Int32]()
+        queue.reserveCapacity(stateCount)
+        for start in starts where distances[start] < 0 {
+            distances[start] = 0
+            queue.append(Int32(start))
+        }
+        var head = 0
+        while head < queue.count {
+            let state = Int(queue[head])
+            head += 1
+            let next = distances[state] + 1
+            for move in 0..<moveCount {
+                let target = neighbor(state, move)
+                if distances[target] < 0 {
+                    distances[target] = next
+                    queue.append(Int32(target))
+                }
+            }
+        }
+        assert(
+            !requireFullCoverage || !distances.contains(-1),
+            "pruning table has unreachable states")
+        return distances
+    }
+}
diff --git a/CubeKit/Sources/CubeKit/ThistlethwaiteSolver.swift b/CubeKit/Sources/CubeKit/ThistlethwaiteSolver.swift
new file mode 100644
index 0000000..41c6673
--- /dev/null
+++ b/CubeKit/Sources/CubeKit/ThistlethwaiteSolver.swift
@@ -0,0 +1,73 @@
+/// The four phases of Thistlethwaite's algorithm, in order. Each phase
+/// restricts the move set further and lands the cube in a smaller
+/// subgroup, until only half turns remain.
+public enum ThistlethwaiteStage: CaseIterable, Sendable, SolverStage {
+    case orientEdges
+    case orientCornersAndSlice
+    case separatePieces
+    case finish
+
+    public var displayName: String {
+        switch self {
+        case .orientEdges: "Orient edges"
+        case .orientCornersAndSlice: "Orient corners + slice"
+        case .separatePieces: "Separate pieces"
+        case .finish: "Half-turn finish"
+        }
+    }
+}
+
+/// Thistlethwaite's 1981 four-phase group reduction: a fixed nesting of
+/// move-set restrictions, each phase solved provably shortest from its
+/// distance table. Solutions run ~30–45 moves and every phase has a
+/// clear goal, which makes the method pleasant to follow.
+public struct ThistlethwaiteSolver: Sendable {
+    private let tables: ThistlethwaiteTables
+
+    public init(tables: ThistlethwaiteTables) {
+        self.tables = tables
+    }
+
+    public func solve(_ state: CubeState) -> StagedSolution<ThistlethwaiteStage>? {
+        guard state.isLegal else { return nil }
+        var current = state
+        var stages: [StagedSolution<ThistlethwaiteStage>.Stage] = []
+        for (phase, stage) in ThistlethwaiteStage.allCases.enumerated() {
+            guard let moves = descend(&current, phase: phase) else {
+                assertionFailure("thistlethwaite stuck in phase \(phase + 1)")
+                return nil
+            }
+            stages.append(.init(stage: stage, moves: moves))
+        }
+        guard current.isSolved else {
+            assertionFailure("thistlethwaite end state not solved")
+            return nil
+        }
+        return StagedSolution(stages: stages)
+    }
+
+    /// Greedy descent over one phase's distance table: from distance d,
+    /// some allowed move always reaches d − 1 (the move sets are closed
+    /// under inverses, so a BFS parent edge runs backwards from every
+    /// state), and any such move lies on a shortest path.
+    private func descend(_ state: inout CubeState, phase: Int) -> [Move]? {
+        let allowed = ThistlethwaiteTables.phaseMoves[phase].map { Move(rawValue: $0)! }
+        guard var distance = tables.distance(of: state, phase: phase) else { return nil }
+        var moves: [Move] = []
+        while distance > 0 {
+            var stepped = false
+            for move in allowed {
+                let next = state.applying(move)
+                if tables.distance(of: next, phase: phase) == distance - 1 {
+                    state = next
+                    moves.append(move)
+                    distance -= 1
+                    stepped = true
+                    break
+                }
+            }
+            guard stepped else { return nil }
+        }
+        return moves
+    }
+}
diff --git a/CubeKit/Sources/CubeKit/ThistlethwaiteTables.swift b/CubeKit/Sources/CubeKit/ThistlethwaiteTables.swift
new file mode 100644
index 0000000..218dcda
--- /dev/null
+++ b/CubeKit/Sources/CubeKit/ThistlethwaiteTables.swift
@@ -0,0 +1,346 @@
+import Foundation
+
+/// Exact BFS distance tables for the four Thistlethwaite phases (~5 MB).
+/// Each table measures, in the move set allowed for that phase, the
+/// distance to the next subgroup in the chain G0 ⊃ G1 ⊃ G2 ⊃ G3 ⊃ G4,
+/// so a solve is a greedy descent — every phase comes out
+/// phase-optimal. Generation takes about a second in release builds;
+/// use `cached(in:)` to persist across launches.
+public final class ThistlethwaiteTables: Sendable {
+    /// Moves allowed in each phase, as `Move` raw values. Every set is
+    /// closed under inverses — that is what guarantees the greedy
+    /// descent always finds a distance-reducing move (the move graph is
+    /// undirected, so a BFS parent edge runs backwards from any state).
+    static let phaseMoves: [[Int]] = [
+        Array(0..<18),                                      // G0: everything
+        [0, 1, 2, 9, 10, 11, 12, 13, 14, 3, 4, 5, 7, 16],   // G1: U* D* L* R* F2 B2
+        SolverTables.phase2Moves,                           // G2: U* D* R2 L2 F2 B2
+        [1, 10, 4, 13, 7, 16],                              // G3: half turns only
+    ]
+
+    /// The M-slice edge pieces (UF, UB, DF, DB) live in the U/D layers,
+    /// slots 0–7, once the cube is in G2.
+    static let mEdges: Set<Int> = [1, 3, 5, 7]
+    static let separationCount = 70  // C(8,4)
+
+    static let corner96Count = 96
+    static let sliceEdgeArrangements = 13824  // 24³
+    static let phase4Count = 96 * 13824
+
+    // Distance tables. Phase 4 legitimately covers only half its index
+    // space (half turns are even permutations, so odd combined edge
+    // parity is unreachable); unreached entries stay -1.
+    let phase1: [Int8]  // [flip]
+    let phase2: [Int8]  // [twist * 495 + slice]
+    let phase3: [Int8]  // [cornerPermutation * 70 + mSeparation]
+    let phase4: [Int8]  // [corner96 * 13824 + ((mPerm * 24) + sPerm) * 24 + ePerm]
+
+    /// Lehmer ranks of the 96 corner permutations reachable by half
+    /// turns, mapped to dense indices for the phase-4 coordinate.
+    let corner96Index: [Int: Int]
+
+    // MARK: Generation
+
+    public static func generate() -> ThistlethwaiteTables {
+        ThistlethwaiteTables()
+    }
+
+    private init() {
+        let moves = CubeState.basicMoves
+        let corner96 = Self.corner96Indices()
+        corner96Index = corner96
+
+        // Phase 1 — orient every edge: the flip coordinate alone.
+        let flipMove = TableBuilder.orientationTable(
+            count: Coordinates.flipCount, pieceCount: 12, modulus: 2,
+            decode: Coordinates.edgeOrientations(forFlip:),
+            encode: { orientations in
+                var value = 0
+                for i in 0..<11 { value = value * 2 + orientations[i] }
+                return value
+            },
+            movePermutation: { moves[$0].edgePermutation },
+            moveOrientation: { moves[$0].edgeOrientation }
+        )
+        phase1 = TableBuilder.breadthFirstDistances(
+            stateCount: Coordinates.flipCount, moveCount: 18, starts: [0]
+        ) { state, move in
+            Int(flipMove[state * 18 + move])
+        }
+
+        // Phase 2 — orient corners and bring the E-slice edges into the
+        // E slice, using G1 moves only.
+        let twistMove = TableBuilder.orientationTable(
+            count: Coordinates.twistCount, pieceCount: 8, modulus: 3,
+            decode: Coordinates.cornerOrientations(forTwist:),
+            encode: { orientations in
+                var value = 0
+                for i in 0..<7 { value = value * 3 + orientations[i] }
+                return value
+            },
+            movePermutation: { moves[$0].cornerPermutation },
+            moveOrientation: { moves[$0].cornerOrientation }
+        )
+        let sliceMove = TableBuilder.occupancyTable(
+            slotCount: 12, markerCount: 4, totalSlots: 12, moves: Array(0..<18),
+            movePermutation: { moves[$0].edgePermutation }
+        )
+        let g1 = Self.phaseMoves[1]
+        phase2 = TableBuilder.breadthFirstDistances(
+            stateCount: Coordinates.twistCount * Coordinates.sliceCount,
+            moveCount: g1.count,
+            starts: [Coordinates.solvedSlice]  // twist 0, slice solved
+        ) { state, index in
+            let move = g1[index]
+            let twist = state / Coordinates.sliceCount
+            let slice = state % Coordinates.sliceCount
+            return Int(twistMove[twist * 18 + move]) * Coordinates.sliceCount
+                + Int(sliceMove[slice * 18 + move])
+        }
+
+        // Phase 3 — into the square group: corner permutation and
+        // M-edge separation, using G2 moves only. The goal is the whole
+        // coset (every half-turn-reachable corner permutation with the
+        // M edges home), seeded as a multi-source BFS — that sidesteps
+        // hand-deriving Thistlethwaite's tetrad-twist condition.
+        let g2 = Self.phaseMoves[2]
+        let cornerPermutationMove = TableBuilder.permutationTable(
+            count: Coordinates.cornerPermutationCount, pieceCount: 8, moves: g2,
+            movePermutation: { moves[$0].cornerPermutation },
+            project: { $0 }, embed: { $0 }
+        )
+        let separationMove = TableBuilder.occupancyTable(
+            slotCount: 8, markerCount: 4, totalSlots: 12, moves: g2,
+            movePermutation: { moves[$0].edgePermutation }
+        )
+        let solvedSeparation = Self.mSeparation(of: .solved)
+        phase3 = TableBuilder.breadthFirstDistances(
+            stateCount: Coordinates.cornerPermutationCount * Self.separationCount,
+            moveCount: g2.count,
+            starts: corner96.keys.map {
+                $0 * Self.separationCount + solvedSeparation
+            }
+        ) { state, index in
+            let corner = state / Self.separationCount
+            let separation = state % Self.separationCount
+            return Int(cornerPermutationMove[corner * g2.count + index])
+                * Self.separationCount
+                + Int(separationMove[separation * g2.count + index])
+        }
+
+        // Phase 4 — finish inside the square group with half turns.
+        var corner96Move = [UInt16](repeating: 0, count: Self.corner96Count * 6)
+        for (rank, index) in corner96 {
+            let permutation = Coordinates.unrankPermutation(rank, count: 8)
+            for face in 0..<6 {
+                let moved = TableBuilder.permutationApplying(
+                    permutation, moves[face].cornerPermutation, times: 2)
+                corner96Move[index * 6 + face] =
+                    UInt16(corner96[Coordinates.rankPermutation(moved)]!)
+            }
+        }
+        let mMove = Self.sliceClassMove(slots: [1, 3, 5, 7])
+        let sMove = Self.sliceClassMove(slots: [0, 2, 4, 6])
+        let eMove = Self.sliceClassMove(slots: [8, 9, 10, 11])
+        phase4 = TableBuilder.breadthFirstDistances(
+            stateCount: Self.phase4Count, moveCount: 6,
+            starts: [corner96[0]! * Self.sliceEdgeArrangements],
+            requireFullCoverage: false
+        ) { state, face in
+            let e = state % 24
+            let s = (state / 24) % 24
+            let m = (state / 576) % 24
+            let corner = state / Self.sliceEdgeArrangements
+            return Int(corner96Move[corner * 6 + face]) * Self.sliceEdgeArrangements
+                + Int(mMove[m * 6 + face]) * 576
+                + Int(sMove[s * 6 + face]) * 24
+                + Int(eMove[e * 6 + face])
+        }
+        assert(
+            phase4.lazy.filter { $0 >= 0 }.count == Self.phase4Count / 2,
+            "square group should fill exactly half the phase-4 space")
+    }
+
+    private init(phase1: [Int8], phase2: [Int8], phase3: [Int8], phase4: [Int8]) {
+        self.phase1 = phase1
+        self.phase2 = phase2
+        self.phase3 = phase3
+        self.phase4 = phase4
+        corner96Index = Self.corner96Indices()
+    }
+
+    // MARK: Coordinates
+
+    /// Colex rank of which of the eight U/D slots hold the M-slice
+    /// edges. Defined once the cube is in G2.
+    static func mSeparation(of state: CubeState) -> Int {
+        var rank = 0
+        var found = 0
+        for slot in 0..<8 where mEdges.contains(state.edgePermutation[slot]) {
+            found += 1
+            rank += Coordinates.binomial[slot][found]
+        }
+        return rank
+    }
+
+    /// The distance left in `phase` (0–3) for a state that has completed
+    /// the phases before it, or nil if the state is outside the phase's
+    /// domain (which means an earlier phase failed).
+    func distance(of state: CubeState, phase: Int) -> Int? {
+        switch phase {
+        case 0:
+            return Int(phase1[Coordinates.flip(of: state)])
+        case 1:
+            return Int(phase2[
+                Coordinates.twist(of: state) * Coordinates.sliceCount
+                    + Coordinates.slice(of: state)])
+        case 2:
+            return Int(phase3[
+                Coordinates.rankPermutation(state.cornerPermutation) * Self.separationCount
+                    + Self.mSeparation(of: state)])
+        default:
+            guard let corner = corner96Index[
+                Coordinates.rankPermutation(state.cornerPermutation)],
+                let m = Self.classPermutation(of: state, slots: [1, 3, 5, 7]),
+                let s = Self.classPermutation(of: state, slots: [0, 2, 4, 6]),
+                let e = Self.classPermutation(of: state, slots: [8, 9, 10, 11])
+            else { return nil }
+            let index = corner * Self.sliceEdgeArrangements
+                + (Coordinates.rankPermutation(m) * 24 + Coordinates.rankPermutation(s)) * 24
+                + Coordinates.rankPermutation(e)
+            let value = phase4[index]
+            return value < 0 ? nil : Int(value)
+        }
+    }
+
+    /// The within-class permutation of one edge slice (indices into
+    /// `slots`), or nil if a foreign piece sits in the class.
+    private static func classPermutation(of state: CubeState, slots: [Int]) -> [Int]? {
+        var permutation: [Int] = []
+        for slot in slots {
+            guard let index = slots.firstIndex(of: state.edgePermutation[slot])
+            else { return nil }
+            permutation.append(index)
+        }
+        return permutation
+    }
+
+    /// The 96 corner permutations reachable by half turns, as Lehmer
+    /// rank → dense index. Sub-millisecond, so recomputed on every init
+    /// rather than persisted.
+    private static func corner96Indices() -> [Int: Int] {
+        var visited: Set<Int> = [0]
+        var queue: [[Int]] = [Array(0..<8)]
+        var head = 0
+        while head < queue.count {
+            let permutation = queue[head]
+            head += 1
+            for face in 0..<6 {
+                let moved = TableBuilder.permutationApplying(
+                    permutation, CubeState.basicMoves[face].cornerPermutation, times: 2)
+                if visited.insert(Coordinates.rankPermutation(moved)).inserted {
+                    queue.append(moved)
+                }
+            }
+        }
+        assert(visited.count == corner96Count, "half-turn corner closure should be 96")
+        var indices: [Int: Int] = [:]
+        for (index, rank) in visited.sorted().enumerated() {
+            indices[rank] = index
+        }
+        return indices
+    }
+
+    /// Half-turn move table for the within-class permutation of one edge
+    /// slice (each half turn maps every slice class to itself).
+    private static func sliceClassMove(slots: [Int]) -> [UInt16] {
+        var table = [UInt16](repeating: 0, count: 24 * 6)
+        for rank in 0..<24 {
+            let classPermutation = Coordinates.unrankPermutation(rank, count: 4)
+            var full = Array(0..<12)
+            for (i, slot) in slots.enumerated() { full[slot] = slots[classPermutation[i]] }
+            for face in 0..<6 {
+                let moved = TableBuilder.permutationApplying(
+                    full, CubeState.basicMoves[face].edgePermutation, times: 2)
+                let projected = slots.map { slot in slots.firstIndex(of: moved[slot])! }
+                table[rank * 6 + face] = UInt16(Coordinates.rankPermutation(projected))
+            }
+        }
+        return table
+    }
+
+    // MARK: Disk cache
+
+    private static let cacheVersion: UInt32 = 1
+    private static let cacheMagic: UInt32 = 0x434B_3157  // "CK1W"
+
+    /// Loads tables from `directory`, generating and saving them on a
+    /// cache miss. The file is device-local (native byte order).
+    public static func cached(in directory: URL) throws -> ThistlethwaiteTables {
+        let url = directory.appendingPathComponent("thistlethwaite-tables.bin")
+        if let data = try? Data(contentsOf: url), let tables = decode(data) {
+            return tables
+        }
+        let tables = generate()
+        try? FileManager.default.createDirectory(
+            at: directory, withIntermediateDirectories: true)
+        try tables.encoded().write(to: url, options: .atomic)
+        return tables
+    }
+
+    private func encoded() -> Data {
+        var data = Data()
+        func append(_ value: UInt32) {
+            withUnsafeBytes(of: value) { data.append(contentsOf: $0) }
+        }
+        func append8(_ table: [Int8]) {
+            append(UInt32(table.count))
+            table.withUnsafeBytes { data.append(contentsOf: $0) }
+        }
+        append(Self.cacheMagic)
+        append(Self.cacheVersion)
+        append8(phase1)
+        append8(phase2)
+        append8(phase3)
+        append8(phase4)
+        return data
+    }
+
+    private static func decode(_ data: Data) -> ThistlethwaiteTables? {
+        var offset = 0
+        func readUInt32() -> UInt32? {
+            guard offset + 4 <= data.count else { return nil }
+            let value = data.withUnsafeBytes {
+                $0.loadUnaligned(fromByteOffset: offset, as: UInt32.self)
+            }
+            offset += 4
+            return value
+        }
+        func read8(expectedCount: Int) -> [Int8]? {
+            guard let count = readUInt32(), Int(count) == expectedCount,
+                  offset + expectedCount <= data.count else { return nil }
+            let table = [Int8](unsafeUninitializedCapacity: expectedCount) { buffer, written in
+                data.withUnsafeBytes { raw in
+                    let source = raw.baseAddress!.advanced(by: offset)
+                    memcpy(buffer.baseAddress!, source, expectedCount)
+                }
+                written = expectedCount
+            }
+            offset += expectedCount
+            return table
+        }
+
+        guard readUInt32() == cacheMagic, readUInt32() == cacheVersion,
+              let phase1 = read8(expectedCount: Coordinates.flipCount),
+              let phase2 = read8(
+                expectedCount: Coordinates.twistCount * Coordinates.sliceCount),
+              let phase3 = read8(
+                expectedCount: Coordinates.cornerPermutationCount * separationCount),
+              let phase4 = read8(expectedCount: phase4Count),
+              offset == data.count
+        else { return nil }
+
+        return ThistlethwaiteTables(
+            phase1: phase1, phase2: phase2, phase3: phase3, phase4: phase4)
+    }
+}
diff --git a/CubeKit/Tests/CubeKitTests/Support/TestTables.swift b/CubeKit/Tests/CubeKitTests/Support/TestTables.swift
index 73b24ea..fcf4dae 100644
--- a/CubeKit/Tests/CubeKitTests/Support/TestTables.swift
+++ b/CubeKit/Tests/CubeKitTests/Support/TestTables.swift
@@ -11,4 +11,11 @@ enum TestTables {
     }()
 
     static let solver = KociembaSolver(tables: shared)
+
+    static let thistlethwaite: ThistlethwaiteTables = {
+        let directory = FileManager.default.temporaryDirectory
+            .appendingPathComponent("CubeKitTestThistlethwaiteTables")
+        return (try? ThistlethwaiteTables.cached(in: directory))
+            ?? ThistlethwaiteTables.generate()
+    }()
 }
diff --git a/CubeKit/Tests/CubeKitTests/ThistlethwaiteSolverTests.swift b/CubeKit/Tests/CubeKitTests/ThistlethwaiteSolverTests.swift
new file mode 100644
index 0000000..16d8ff0
--- /dev/null
+++ b/CubeKit/Tests/CubeKitTests/ThistlethwaiteSolverTests.swift
@@ -0,0 +1,124 @@
+import Foundation
+import Testing
+@testable import CubeKit
+
+@Suite struct ThistlethwaiteSolverTests {
+    private var solver: ThistlethwaiteSolver {
+        ThistlethwaiteSolver(tables: TestTables.thistlethwaite)
+    }
+
+    /// Moves each phase may use, written out independently of the
+    /// solver's own tables.
+    private static let phaseMoves: [Set<Move>] = [
+        Set(Move.allCases),
+        Set([Move](notation: "U U2 U' D D2 D' L L2 L' R R2 R' F2 B2")!),
+        Set([Move](notation: "U U2 U' D D2 D' R2 L2 F2 B2")!),
+        Set([Move](notation: "U2 D2 R2 L2 F2 B2")!),
+    ]
+
+    /// The corner permutations reachable by half turns alone, built here
+    /// by direct closure so the solver's own 96-set isn't trusted.
+    private static let squareGroupCornerPermutations: Set<[Int]> = {
+        let generators = [Move.u2, .d2, .r2, .l2, .f2, .b2].map {
+            CubeState.solved.applying($0).cornerPermutation
+        }
+        var visited: Set<[Int]> = [Array(0..<8)]
+        var queue: [[Int]] = [Array(0..<8)]
+        var head = 0
+        while head < queue.count {
+            let permutation = queue[head]
+            head += 1
+            for generator in generators {
+                let next = (0..<8).map { permutation[generator[$0]] }
+                if visited.insert(next).inserted { queue.append(next) }
+            }
+        }
+        return visited
+    }()
+
+    @Test func squareGroupClosureHas96CornerPermutations() {
+        #expect(Self.squareGroupCornerPermutations.count == 96)
+    }
+
+    @Test func solvesSolvedCubeWithFourEmptyStages() throws {
+        let solution = try #require(solver.solve(.solved))
+        #expect(solution.stages.count == 4)
+        #expect(solution.stages.allSatisfy { $0.moves.isEmpty })
+    }
+
+    @Test func rejectsIllegalState() {
+        var state = CubeState.solved
+        state.cornerOrientation[0] = 1
+        #expect(solver.solve(state) == nil)
+    }
+
+    /// The known per-phase worst cases: 7, 10, 13, 15 (total 45). These
+    /// pin the coordinate conventions — any drift in the tables shows up
+    /// here first.
+    @Test func tableDepthsMatchKnownBounds() {
+        let tables = TestTables.thistlethwaite
+        #expect(tables.phase1.max() == 7)
+        #expect(tables.phase2.max() == 10)
+        #expect(tables.phase3.max() == 13)
+        #expect(tables.phase4.max() == 15)
+    }
+
+    @Test func solvesRandomStatesWithValidStages() throws {
+        var rng = SeededRandom(seed: 1981)
+        for iteration in 0..<200 {
+            let start = Scrambler.randomState(using: &rng)
+            let solution = try #require(
+                solver.solve(start), "iteration \(iteration) returned nil")
+            #expect(solution.stages.count == 4)
+            #expect(solution.moves.count <= 52, "iteration \(iteration) too long")
+
+            var state = start
+            for (phase, stage) in solution.stages.enumerated() {
+                #expect(
+                    stage.moves.allSatisfy { Self.phaseMoves[phase].contains($0) },
+                    "phase \(phase + 1) used a forbidden move at iteration \(iteration)")
+                state = state.applying(stage.moves)
+                checkInvariant(phase: phase, state: state, iteration: iteration)
+            }
+            #expect(state.isSolved, "iteration \(iteration) not solved")
+        }
+    }
+
+    private func checkInvariant(
+        phase: Int, state: CubeState, iteration: Int,
+        sourceLocation: SourceLocation = #_sourceLocation
+    ) {
+        var holds = state.edgeOrientation.allSatisfy { $0 == 0 }
+        if phase >= 1 {
+            holds = holds && state.cornerOrientation.allSatisfy { $0 == 0 }
+                && (8..<12).allSatisfy { state.edgePermutation[$0] >= 8 }
+        }
+        if phase >= 2 {
+            holds = holds && [1, 3, 5, 7].allSatisfy { [1, 3, 5, 7].contains(state.edgePermutation[$0]) }
+                && Self.squareGroupCornerPermutations.contains(state.cornerPermutation)
+        }
+        if phase >= 3 {
+            holds = holds && state.isSolved
+        }
+        #expect(
+            holds, "phase \(phase + 1) invariant failed at iteration \(iteration)",
+            sourceLocation: sourceLocation)
+    }
+
+    @Test func cachedTablesMatchGenerated() throws {
+        let directory = FileManager.default.temporaryDirectory
+            .appendingPathComponent("CubeKitThistlethwaiteRoundTrip-\(UUID().uuidString)")
+        defer { try? FileManager.default.removeItem(at: directory) }
+        let generated = TestTables.thistlethwaite
+        // First call writes the cache, second one must load it.
+        let written = try ThistlethwaiteTables.cached(in: directory)
+        let reloaded = try ThistlethwaiteTables.cached(in: directory)
+        for tables in [written, reloaded] {
+            #expect(tables.phase1 == generated.phase1)
+            #expect(tables.phase2 == generated.phase2)
+            #expect(tables.phase3 == generated.phase3)
+            #expect(tables.phase4 == generated.phase4)
+            #expect(tables.corner96Index == generated.corner96Index)
+        }
+    }
+}
diff --git a/cube/AppModel.swift b/cube/AppModel.swift
index b378d4b..a34eb29 100644
--- a/cube/AppModel.swift
+++ b/cube/AppModel.swift
@@ -246,6 +246,9 @@ final class AppModel {
 
     private(set) var solveSession: SolveSession?
     private(set) var isComputingSolution = false
+    /// Built lazily on the first Thistlethwaite solve (the tables take
+    /// about a second to generate, then load from disk), kept for reuse.
+    private var thistlethwaiteSolver: ThistlethwaiteSolver?
 
     func startSolve() {
         guard canSolve, let solver else { return }
@@ -256,33 +259,58 @@ final class AppModel {
         isComputingSolution = true
         let state = cubeState
         let method = settings.solvingMethod
+        let prebuiltThistlethwaite = thistlethwaiteSolver
+        let directory = Self.supportDirectory
         Task.detached(priority: .userInitiated) {
             let solution: [Move]?
             var markers: [SolveSession.StageMarker]?
+            var builtThistlethwaite: ThistlethwaiteSolver?
             switch method {
             case .fast, .optimal:
                 solution = solver.solve(state, timeBudget: .milliseconds(300))
+            case .thistlethwaite:
+                let thistlethwaite = prebuiltThistlethwaite
+                    ?? (try? ThistlethwaiteTables.cached(in: directory))
+                        .map(ThistlethwaiteSolver.init)
+                builtThistlethwaite = thistlethwaite
+                if let staged = thistlethwaite?.solve(state) {
+                    solution = staged.moves
+                    markers = Self.stageMarkers(of: staged)
+                } else {
+                    solution = nil
+                }
             case .beginner:
                 if let staged = BeginnerSolver().solve(state) {
                     solution = staged.moves
-                    var start = 0
-                    markers = staged.stages.compactMap { stage in
-                        guard !stage.moves.isEmpty else { return nil }
-                        defer { start += stage.moves.count }
-                        return SolveSession.StageMarker(
-                            name: stage.stage.displayName,
-                            range: start..<(start + stage.moves.count))
-                    }
+                    markers = Self.stageMarkers(of: staged)
                 } else {
                     solution = nil
                 }
             }
-            await MainActor.run { [weak self, markers] in
+            await MainActor.run { [weak self, markers, builtThistlethwaite] in
+                if let builtThistlethwaite {
+                    self?.thistlethwaiteSolver = builtThistlethwaite
+                }
                 self?.beginSolveSession(solution, stages: markers)
             }
         }
     }
 
+    /// Stage chips for the solution panel: one marker per non-empty
+    /// stage of a staged solution.
+    private nonisolated static func stageMarkers<Stage>(
+        of staged: StagedSolution<Stage>
+    ) -> [SolveSession.StageMarker] {
+        var start = 0
+        return staged.stages.compactMap { stage in
+            guard !stage.moves.isEmpty else { return nil }
+            defer { start += stage.moves.count }
+            return SolveSession.StageMarker(
+                name: stage.stage.displayName,
+                range: start..<(start + stage.moves.count))
+        }
+    }
+
     private func beginSolveSession(
         _ solution: [Move]?, stages: [SolveSession.StageMarker]?, isOptimal: Bool = false
     ) {
diff --git a/cube/Stores/SettingsStore.swift b/cube/Stores/SettingsStore.swift
index 754a9dd..eebaca9 100644
--- a/cube/Stores/SettingsStore.swift
+++ b/cube/Stores/SettingsStore.swift
@@ -104,6 +104,9 @@ final class SettingsStore {
         /// Korf IDA* over pattern databases: provably shortest, takes
         /// minutes and requires prepared tables.
         case optimal
+        /// Thistlethwaite's four-phase group reduction: each phase
+        /// provably shortest within its move set.
+        case thistlethwaite
         /// Classic layer-by-layer, the way people learn: long but
         /// followable, with named stages.
         case beginner
@@ -113,6 +116,7 @@ final class SettingsStore {
             switch self {
             case .fast: "Fast"
             case .optimal: "Optimal"
+            case .thistlethwaite: "Thistlethwaite"
             case .beginner: "Beginner"
             }
         }
@@ -121,6 +125,7 @@ final class SettingsStore {
             switch self {
             case .fast: "~20 moves, found instantly"
             case .optimal: "Proven shortest — takes minutes, needs tables"
+            case .thistlethwaite: "Four phases, each provably shortest (~30–45 moves)"
             case .beginner: "Step by step, the way people learn (~200 moves)"
             }
         }
@@ -138,7 +143,11 @@ final class SettingsStore {
         var orbitSensitivity: Double?
         var rotationLock: Bool?
         var lockMode: LockMode?
-        var solvingMethod: SolvingMethod?
+        /// Stored as a raw string so a file written by a newer version
+        /// (with methods this build doesn't know) still decodes — an
+        /// unknown method falls back to `.fast` instead of throwing and
+        /// silently resetting every setting.
+        var solvingMethod: String?
         var twoFingerOrbit: Bool?
         var faceColors: [StoredColor]?
     }
@@ -199,7 +208,7 @@ final class SettingsStore {
         orbitSensitivity = snapshot.orbitSensitivity ?? 1.0
         lockMode = snapshot.lockMode
             ?? (snapshot.rotationLock == true ? .layersOnly : .free)
-        solvingMethod = snapshot.solvingMethod ?? .fast
+        solvingMethod = snapshot.solvingMethod.flatMap(SolvingMethod.init(rawValue:)) ?? .fast
         twoFingerOrbit = snapshot.twoFingerOrbit ?? false
         if let colors = snapshot.faceColors, colors.count == 6 {
             faceColors = colors
@@ -211,7 +220,7 @@ final class SettingsStore {
         let snapshot = Snapshot(
             soundEnabled: soundEnabled, hapticsEnabled: hapticsEnabled,
             turnSpeed: turnSpeed, orbitSensitivity: orbitSensitivity,
-            rotationLock: nil, lockMode: lockMode, solvingMethod: solvingMethod,
+            rotationLock: nil, lockMode: lockMode, solvingMethod: solvingMethod.rawValue,
             twoFingerOrbit: twoFingerOrbit, faceColors: faceColors)
         do {
             try FileManager.default.createDirectory(
diff --git a/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone-seeded.xcscheme b/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone-seeded.xcscheme
index 70d3250..8f6662a 100644
--- a/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone-seeded.xcscheme
+++ b/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone-seeded.xcscheme
@@ -1,7 +1,7 @@
 <?xml version="1.0" encoding="UTF-8"?>
 <Scheme
    LastUpgradeVersion = "2700"
-   version = "1.7">
+   version = "1.3">
    <BuildAction
       parallelizeBuildables = "YES"
       buildImplicitDependencies = "YES">
@@ -27,6 +27,8 @@
       selectedDebuggerIdentifier = "Xcode.DebuggerFoundation.Debugger.LLDB"
       selectedLauncherIdentifier = "Xcode.DebuggerFoundation.Launcher.LLDB"
       shouldUseLaunchSchemeArgsEnv = "YES">
+      <Testables>
+      </Testables>
    </TestAction>
    <LaunchAction
       buildConfiguration = "Release-Seeded"
diff --git a/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone.xcscheme b/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone.xcscheme
index d0e46e2..97873f1 100644
--- a/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone.xcscheme
+++ b/cubeone.xcodeproj/xcshareddata/xcschemes/cubeone.xcscheme
@@ -1,7 +1,7 @@
 <?xml version="1.0" encoding="UTF-8"?>
 <Scheme
    LastUpgradeVersion = "2700"
-   version = "1.7">
+   version = "1.3">
    <BuildAction
       parallelizeBuildables = "YES"
       buildImplicitDependencies = "YES">
@@ -27,6 +27,8 @@
       selectedDebuggerIdentifier = "Xcode.DebuggerFoundation.Debugger.LLDB"
       selectedLauncherIdentifier = "Xcode.DebuggerFoundation.Launcher.LLDB"
       shouldUseLaunchSchemeArgsEnv = "YES">
+      <Testables>
+      </Testables>
    </TestAction>
    <LaunchAction
       buildConfiguration = "Debug"
diff --git a/docs/ARCHITECTURE.md b/docs/ARCHITECTURE.md
index 36b6c0a..7ade881 100644
--- a/docs/ARCHITECTURE.md
+++ b/docs/ARCHITECTURE.md
@@ -40,6 +40,16 @@ The classic layer-by-layer method, seven stages, each emitting its moves into a
 
 Algorithm conventions are *self-calibrating* where possible (e.g. the A-perm's fixed corner is computed from the cube algebra at startup, not assumed). Every loop has an iteration guard that throws instead of hanging. The test sweep solves 200 seeded random states and asserts the stage invariant after **every** stage, not just the end state.
 
+### ThistlethwaiteSolver (four-phase group reduction)
+
+Thistlethwaite's 1981 algorithm: restrict the move set in four steps, each landing the cube in a smaller subgroup — G0 ⊃ G1 = ⟨U, D, L, R, F2, B2⟩ ⊃ G2 = ⟨U, D, R2, L2, F2, B2⟩ ⊃ G3 = ⟨U2, D2, R2, L2, F2, B2⟩ ⊃ solved. Each phase has its own exact BFS distance table over a small coordinate, so solving is a **greedy descent**: from distance d, apply any allowed move whose table entry reads d − 1. Every phase move set is closed under inverses, so the move graph is undirected and such a move always exists — each phase comes out provably shortest for its goal. The known worst cases are 7 + 10 + 13 + 15 = 45 moves (the test suite pins all four maxima); typical solves run ~30–45.
+
+Per-phase coordinates: ① edge flip (2048); ② corner twist × slice-edge locations (2187 × 495); ③ corner permutation × M-edge separation (40320 × C(8,4)); ④ within the square group, the 96 half-turn-reachable corner permutations × the three within-slice edge permutations (96 × 24³, of which exactly half is reachable — half turns are even permutations).
+
+Phase 3's goal is the subtle one (Thistlethwaite's "tetrad twist" condition). Instead of deriving it by hand, the table seeds a **multi-source BFS from the whole goal coset**: every half-turn-reachable corner permutation with the M edges home is distance 0. That is provably the exact G3 membership test for a G2 state — legality ties edge parity to corner parity, and 96 × 6912 matches the square group's order exactly.
+
+All four tables (~5 MB) generate in about a second (release) and cache as `thistlethwaite-tables.bin` in Application Support, built lazily on the first Thistlethwaite solve. Playback shows the four named stages just like the beginner method (`StagedSolution` is generic over the stage kind: `StagedSolution<ThistlethwaiteStage>`).
+
 ### OptimalSolver (proven-shortest solutions)
 
 Korf-style IDA* over **pattern databases**: nibble-packed tables holding the exact solve distance of three projections — all corner configurations (88M entries), and two overlapping edge subsets (7-edge tier: 2×511M entries ≈ 0.5 GB; 8-edge tier: 2×5.1B entries ≈ 4.8 GB). The heuristic is the max of the three lookups (admissible), so the first solution found by per-bound exhaustive deepening is provably optimal.