# Local Search And Neighborhoods > When the user wants to design or implement neighborhood-based local search - choosing moves (swap, insertion, 2-opt, Or-opt, exchange), writing O(1)/O(n) delta evaluation, first vs best improvement, scan order, and move data structures. Also use when the user mentions "local search," "2-opt," "neighborhood," "delta evaluation," "hill climbing," "first improvement," or when a heuristic recomputes… - **Type:** Skill - **Install:** `agentstack add skill-hajibabaie-combinatorial-optimization-skills-local-search-and-neighborhoods` - **Verified:** Yes — security-reviewed for prompt injection and unsafe behavior - **Seller:** [hajibabaie](https://agentstack.voostack.com/s/hajibabaie) - **Installs:** 0 - **Category:** [Search](https://agentstack.voostack.com/c/search) - **Latest version:** 0.1.0 - **License:** MIT - **Upstream author:** [hajibabaie](https://github.com/hajibabaie) - **Source:** https://github.com/hajibabaie/combinatorial-optimization-skills/tree/main/skills/local-search-and-neighborhoods ## Install ```sh agentstack add skill-hajibabaie-combinatorial-optimization-skills-local-search-and-neighborhoods ``` Requires the [AgentStack CLI](https://agentstack.voostack.com/docs/cli). Works with Claude Code, Cursor, and any MCP-compatible agent. ## About # Local Search and Neighborhoods You are an expert in neighborhood-based search for combinatorial optimization. This skill covers neighborhood design (swap, insertion, 2-opt, Or-opt, exchange), constant- and linear-time delta evaluation, first- vs best-improvement pivoting, neighborhood scanning order, move data structures, and the limits of plain hill climbing. Use the framework below to build descent engines that are correct (deltas audited against full recomputation) and fast (no full objective evaluation inside the move loop), and that drop in unchanged as the inner loop of simulated annealing, tabu search, and iterated local search. ## Initial Assessment Establish these facts before proposing a neighborhood or writing any move code: - **Representation.** Permutation (tour, schedule), binary vector (selection), integer assignment array, or partition into sets (routes, machine loads)? The move catalog and the delta formulas depend entirely on this. See **solution-encodings**-style criteria: every move must map feasible representations to feasible representations, or you must plan a penalty scheme. - **Objective structure.** Additive over solution elements (sum of edge lengths, sum of assignment costs, linear penalties)? Additive objectives almost always admit O(1) or O(n) deltas. Max-type / critical-path objectives (flow-shop or job-shop makespan) do not; plan for an acceleration scheme (Taillard heads/tails, critical-path filtering) instead of a true delta. - **Hard vs soft constraints.** Decide per constraint: keep moves feasibility-preserving, or relax the constraint into the objective and evaluate the penalty change inside the delta. Mixing the two without a plan is the most common source of wrong deltas. - **Instance size.** n ≤ 500: any O(n²) neighborhood with full sweeps is fine in plain Python. n in 10³–10⁴: you need candidate lists, don't-look bits, or vectorized move scoring. n above 10⁴: you also need O(1) position queries and careful apply costs (segment reversal is O(n) per move). - **Evaluation cost.** Time one full objective evaluation and one delta evaluation early. If the delta is not at least ~n times cheaper than the full evaluation, the neighborhood implementation is wrong or the objective does not decompose. - **Role of the local search.** Standalone descent, multistart, or inner loop of a metaheuristic? The wrapper changes the acceptance logic; the move/delta code should not change at all. Keep them in separate functions from day one. - **Quality target.** 2-opt local optima on random Euclidean TSP sit roughly 5% above optimal from good starts (Johnson & McGeoch 1997, "The traveling salesman problem: a case study in local optimization"). If the target is tighter, plan a metaheuristic wrapper or a richer neighborhood, not more sweeps. - **Data format.** Full distance matrix (O(n²) memory) vs coordinates plus k-nearest neighbors via a KD-tree? At n = 20,000, a float64 matrix is 3.2 GB — candidate lists are not optional there. - **Determinism.** Seed every random scan order and tie-break (`np.random.default_rng`) so descents are reproducible run to run. - **Validation plan.** An independent objective recompute and a delta audit (compare each claimed delta against a from-scratch evaluation on random moves) must exist before any performance tuning. ## Neighborhood Anatomy A neighborhood is a function $N: S \to 2^S$ mapping each solution to the set of solutions reachable by one move. Local search descends until it reaches a local optimum: $$ s^\* \text{ is a local optimum of } f \text{ w.r.t. } N \iff f(s^\*) \le f(s) \quad \forall s \in N(s^\*). $$ A move $m$ applied to $s$ is written $s \oplus m$, and the entire craft of fast local search is computing the delta $$ \Delta(s, m) = f(s \oplus m) - f(s) $$ without ever computing $f(s \oplus m)$ from scratch. Local optimality is always *relative to N*: a 2-opt optimum is usually not an Or-opt optimum, which is the basis of compound neighborhoods and variable neighborhood descent. ### Move catalog | Move | Acts on | Neighborhood size | Delta cost | Notes | |---|---|---|---|---| | Swap (exchange two elements) | permutation, assignment array | n(n−1)/2 | O(1) for additive costs; O(n) for QAP | Weak on tours (touches 4 edges but rarely improving); strong on assignments | | Insertion (shift one element) | permutation | (n−1)² | O(1) additive; O(nm) per job in flow shop via Taillard | The dominant scheduling neighborhood | | 2-opt (reverse a segment) | tours | n(n−1)/2 | O(1) symmetric distances; **no O(1) delta if asymmetric** | Apply costs O(n) (reversal); Croes (1958) | | Or-opt (relocate segment of 1–3) | tours | ≈ 3n² | O(1) | Orientation-preserving, so valid for asymmetric instances; Or (1976) | | k-flip | binary vectors | C(n,k) | O(degree) per flip | Knapsack, max-cut, SAT-style problems | | Inter-set exchange / relocate | partitions, routes | O(n²) | O(1)–O(route length) | CVRP inter-route moves; feasibility (load, time) checked in the delta | ### Delta evaluation For the symmetric TSP, reversing positions $i{+}1 \dots j$ of tour $\pi$ removes two edges and adds two edges; every interior edge keeps its (symmetric) length: $$ \Delta_{\text{2opt}}(i, j) = d_{\pi_i \pi_j} + d_{\pi_{i+1} \pi_{j+1}} - d_{\pi_i \pi_{i+1}} - d_{\pi_j \pi_{j+1}}. $$ This O(1) formula is why 2-opt scans millions of moves per second. On *asymmetric* instances the reversed segment changes direction, so the formula is invalid — use orientation-preserving moves (Or-opt, 3-opt segment reinsertion) instead. The same discipline applies everywhere: write the delta as (cost of edges/terms created) − (cost of edges/terms destroyed), and prove to yourself that nothing else changes. ### First vs best improvement - **Best improvement** scans the whole neighborhood, applies the steepest move. Fewer, larger steps; deterministic given the scan; pairs naturally with vectorized scoring of all moves at once. - **First improvement** applies the first improving move found. Much cheaper sweeps early in the descent when improving moves are dense; final quality is usually on par with best improvement, sometimes better, especially from random starts (Hansen & Mladenović 2006, "First vs. best improvement: an empirical study"). - Practical default: first improvement with a *randomized scan order per sweep* (removes lexicographic bias), switching to best improvement only when moves are scored in bulk with numpy. - The termination certificate is one full improving-move-free sweep. That last sweep costs |N| delta evaluations no matter what — pruning rules and candidate lists are the only way to shrink it. ### Scanning order and move data structures Scan order matters because first improvement commits to whatever it sees first. Lexicographic order biases the search toward low indices; a seeded random permutation per sweep is the cheap fix. Queue-driven orders (re-examine elements whose surroundings just changed) are the systematic version and lead directly to don't-look bits (Bentley 1992, "Fast algorithms for geometric traveling salesman problems"). Keep these structures current under every apply: - `pos` array with `pos[element] = position`, the inverse of the permutation — O(1) move legality and delta lookups, updated in O(segment) on reversal. - Cached per-element quantities used by deltas (route load, completion times, gains), each with an O(1)/O(k) incremental update rule. - For very large tours, array reversal is the bottleneck; doubly linked or two-level list representations bound the apply cost (Fredman, Johnson, McGeoch & Ostheimer 1995, "Data structures for traveling salesmen"). ### Hill climbing and its limits Plain hill climbing terminates at the first local optimum of $N$ — typically in O(n) to O(n log n) accepted moves empirically, and there it stops. Larger neighborhoods give better local optima at higher sweep cost; no fixed polynomial neighborhood removes local optima for NP-hard problems. Escaping is the wrapper's job: probabilistic acceptance (simulated annealing), memory (tabu search), perturbation and restarts (iterated local search). Build the descent so those wrappers reuse the identical move and delta code. ## Generic Descent Engine The reusable skeleton, independent of problem and move type: ```text LOCAL-SEARCH(s, N, rule): f_s tuple[float, dict[str, int]]: """Generic descent. enumerate_moves() yields (move, delta) pairs for the CURRENT solution; apply_move mutates the solution in place.""" stats = {"sweeps": 0, "moves": 0, "delta_evals": 0} while stats["moves"] None: """Hill climbing for weighted max-cut: flip one vertex side, O(n) delta.""" rng = np.random.default_rng(seed) w = rng.uniform(0.0, 1.0, (n, n)) w = np.triu(w, 1) + np.triu(w, 1).T # symmetric, zero diagonal side = rng.integers(0, 2, n) def cut_value(s: np.ndarray) -> float: return float(w[s[:, None] != s[None, :]].sum() / 2.0) def enumerate_moves() -> Iterable[tuple[int, float]]: same = side[:, None] == side[None, :] gain = (w * same).sum(axis=1) - (w * ~same).sum(axis=1) for i in np.argsort(-gain): # promising flips first yield int(i), -float(gain[i]) # we minimize -cut def apply_move(i: int) -> None: side[i] ^= 1 obj, stats = local_search(-cut_value(side), enumerate_moves, apply_move, improvement="first") assert abs(obj + cut_value(side)) 20 rarely pays off | | Or-opt segment lengths | 1–3 | Longer segments: bigger neighborhood, diminishing returns | | Max sweeps / move budget | bound by time, not count, in experiments | A budget cap protects wrappers (ILS, SA) from runaway descents | | Restarts (standalone use) | 5–50 seeded starts, keep best | Cheap variance reduction; superseded by ILS-style perturbation | | Delta audit frequency | every run in tests; every ~10⁴ moves in production | Catches drift and wrong delta formulas at negligible cost | ## Worked Example: TSP 2-opt with O(1) Delta Evaluation Tour as a numpy position array `tour` (city at each position). The move `(i, j)` reverses `tour[i+1 : j+1]`; the delta is the four-edge formula above. First, the scalar first-improvement descent with a randomized scan order: ```python import numpy as np def tour_length(tour: np.ndarray, dist: np.ndarray) -> float: """Length of the closed tour (array of city indices by position).""" return float(dist[tour, np.roll(tour, -1)].sum()) def two_opt_delta(tour: np.ndarray, dist: np.ndarray, i: int, j: int) -> float: """O(1) length change of reversing tour[i+1 : j+1], for positions i tuple[np.ndarray, float]: """First-improvement 2-opt descent with a randomized scan order per sweep.""" rng = np.random.default_rng(seed) tour = tour.copy() n = len(tour) length = tour_length(tour, dist) improved = True while improved: improved = False for i in rng.permutation(n - 2): i = int(i) for j in range(i + 2, n - (i == 0)): # skip the wrap pair (0, n-1) delta = two_opt_delta(tour, dist, i, j) if delta tuple[np.ndarray, np.ndarray]: """n points on the unit circle; the optimal tour visits them in angular order.""" rng = np.random.default_rng(seed) theta = np.sort(rng.uniform(0.0, 2.0 * np.pi, n)) pts = np.column_stack([np.cos(theta), np.sin(theta)]) dist = np.linalg.norm(pts[:, None] - pts[None, :], axis=2) return pts, dist _, dist = circle_instance(12, seed=42) start = np.random.default_rng(7).permutation(12) tour, length = two_opt_first(start, dist) print(f"start={tour_length(start, dist):.4f} final={length:.4f} " f"optimal={tour_length(np.arange(12), dist):.4f}") # Expected: start=14.6762 final=5.9567 optimal=5.9567. For points in convex # position every 2-opt local optimum is crossing-free, and the only # crossing-free tour through points in convex position is the hull tour. ``` Best improvement becomes attractive once the whole delta matrix is computed in one shot. With `succ = np.roll(tour, -1)` and `edge[i]` the cost of the edge leaving position i, the n×n delta matrix is four broadcasts — one (n, n) temporary per sweep (n = 5,000 → ~200 MB, so this variant is for n up to a few thousand): ```python import numpy as np def two_opt_best_vectorized(tour: np.ndarray, dist: np.ndarray, tol: float = 1e-9) -> tuple[np.ndarray, float]: """Best-improvement 2-opt; each sweep scores all O(n^2) moves at once.""" tour = tour.copy() n = len(tour) length = float(dist[tour, np.roll(tour, -1)].sum()) iu = np.triu_indices(n, k=2) # positions with j >= i+2 while True: succ = np.roll(tour, -1) edge = dist[tour, succ] # cost of edge at each position delta = (dist[tour[:, None], tour[None, :]] + dist[succ[:, None], succ[None, :]] - edge[:, None] - edge[None, :]) k = int(np.argmin(delta[iu])) i, j = int(iu[0][k]), int(iu[1][k]) if delta[i, j] >= -tol: return tour, length tour[i + 1 : j + 1] = tour[i + 1 : j + 1][::-1] length += float(delta[i, j]) rng = np.random.default_rng(3) pts = rng.random((100, 2)) dist = np.linalg.norm(pts[:, None] - pts[None, :], axis=2) start = rng.permutation(100) tour, length = two_opt_best_vectorized(start, dist) print(f"start={float(dist[start, np.roll(start, -1)].sum()):.3f} final={length:.3f}") # Expected: start=46.152 final=8.733 -- about 81% shorter; a best-improvement # 2-opt local optimum for 100 uniform points (optimal is typically ~7.7-8.0, # so this sits 5-12% above it, normal for 2-opt from a random start). ``` The wrap pair (i=0, j=n−1) stays in the index set because its delta is exactly zero (reversing everything after position 0 re-traverses the same closed tour), so it can never be selected as improving. Every delta formula gets an audit before it gets trusted. This is non-negotiable: a wrong delta produces a descent that silently reports objectives that are not the objective. ```python import numpy as np def audit_two_opt_deltas(n: int = 60, trials: int = 500, seed: int = 1) -> float: """Largest |claimed delta - true length change| over random 2-opt moves.""" rng = np.random.default_rng(seed) pts = rng.random((n, 2)) dist = np.linalg.norm(pts[:, None] - pts[None, :], axis=2) tour = rng.permutation(n) base = float(dist[tour, np.roll(tour, -1)].sum()) worst = 0.0 for _ in range(trials): i = int(rng.integers(0, n - 3)) j = int(rng.integers(i + 2, n - (i == 0))) a, b, c, d = tour[i], tour[i + 1], tour[j], tour[(j + 1) % n] delta = float(dist[a, c] + dist[b, d] - dist[a, b] - dist[c, d]) trial = tour.copy() trial[i + 1 : j + 1] = trial[i + 1 : j + 1][::-1] true_change = float(dist[trial, np.roll(trial, -1)].sum()) - base worst = max(worst, abs(true_change - delta)) return worst print(f"max |delta - true change| = {audit_two_opt_deltas():.2e}") # Expected: below 1e-12 -- pure floating-point rounding; the formula is exact. ``` 2-opt local optima still contain misplaced *chains* of cities that no segment reversal can fix. Or-opt relocates segments of length 1–3 (optionally reversed) elsewhere in the tour — an O(1) delta with three edges destroyed and three created, and since the segment keeps its orientation, the formula is also valid for asymmetric instances: ```python import numpy as np def closed_tour_length(tour: list[int], dist: np.ndarray) -> float: """Length of a closed tour stored as a Python list of city indices.""" n = len(tour) return float(sum(dist[tour[k], tour[(k + 1) % n]] for k in range(n))) def first_or_opt_move(tour: list[int], d … ## Source & license This open-source skill is cataloged on AgentStack and links to its original source — we do not rehost the code. - **Author:** [hajibabaie](https://github.com/hajibabaie) - **Source:** [hajibabaie/combinatorial-optimization-skills](https://github.com/hajibabaie/combinatorial-optimization-skills) - **License:** MIT Install and usage instructions live in the source repository linked above. ## Pricing - **Free** — Free ## Security capabilities Automated source analysis of v0.1.0 — what this tool can access: - **Network access:** no - **Filesystem access:** no - **Shell / process execution:** no - **Environment & secrets:** no - **Dynamic code execution:** no *"Yes" means the capability is present in the source — more access means more to trust, not that it is unsafe.* ## Versions - **0.1.0** — security scan: passed — Imported from the upstream source. ## Links - Listing page: https://agentstack.voostack.com/l/skill-hajibabaie-combinatorial-optimization-skills-local-search-and-neighborhoods - Seller: https://agentstack.voostack.com/s/hajibabaie - Browse the marketplace: https://agentstack.voostack.com/browse --- Listed on AgentStack — the marketplace for AI agent skills and MCP servers. Every listing is security-reviewed. Creators keep 70%.