算法-TSP启发式算法-04
旅行商问题 (Traveling Salesperson Problem, TSP)，是计算机科学和运筹学领域一个非常经典且著名的组合优化问题。简单来说，它的核心任务是：
给定一系列城市和每对城市之间的距离，找到一条最短的哈密顿回路，使得旅行商从一个城市出发，经过所有其他城市恰好一次，最后返回到起点。
问题拆解成几个关键点：
-
城市（Cities）: 问题的节点，你可以把它们看作图中的顶点。
-
路径/距离（Paths/Distances）: 连接城市之间的边，每条边都有一个权重，通常代表距离、成本或时间。
-
哈密顿回路（Hamiltonian Cycle）: 一条路径，它从一个城市出发，经过所有其他城市恰好一次，然后回到起点。
-
目标（Objective）: 在所有可能的哈密顿回路中，找到总距离最短的那一条。
TSP 的目标就是比较路径的总距离，然后选择较短的那一条。
但当城市的数量增加时，问题的难度会呈指数级增长。
假设有 n 个城市。从一个城市出发，你可以在剩下的 n−1 个城市中选择一个，然后从剩下的 n−2 个城市中选择一个，以此类推。因此，总的可能路径数量是 (n−1)!，这是一个NP-hard 问题，这意味着目前没有已知的算法可以在多项式时间内（即时间复杂度为 n 的多项式，如 n2,n3）找到最优解。
常见的解法
小规模问题
-
暴力穷举法（Brute-force）: 这是最直观的方法，通过生成所有可能的路径，然后计算每条路径的总距离，最后找出最短的那条。这种方法能保证找到最优解，但只适用于城市数量很少的情况（通常小于10）。
-
动态规划（Dynamic Programming）: 动态规划是一种更高效的精确解法，其中最著名的是Held-Karp算法。它的时间复杂度大约为 O(n2⋅2n)，空间复杂度为 O(n⋅2^n)。它利用子问题的最优解来构建整个问题的最优解，可以求解城市数量在20个左右的问题
当城市数量超过20个时，精确解法变得不可行。这时，我们会转向近似算法（Approximation Algorithms）或启发式算法（Heuristic Algorithms）。这些算法不能保证找到最优解，但可以在合理的时间内找到一个足够好的近似解。
最近邻算法
最近邻算法（Nearest Neighbor Algorithm）: 从一个随机城市出发，每次都前往离当前城市最近的未访问城市，直到所有城市都被访问。这种方法简单快速，但结果可能不是最优的。这种方法可以用于快速的生成初始路径。
遗传算法
遗传算法（Genetic Algorithm）: 一种基于生物进化原理的搜索算法。它通过模拟“基因突变”和“自然选择”来不断优化解，适用于大规模复杂问题。
遗传算法灵感来自生物进化。把“路径（tour）”看成个体（染色体），把路径长度看作“适应度”（fitness，通常适应度和代价成反比）。通过“选择、交叉、变异、替代（世代更新）”这些算子，不断产生新一代，期望得到更好的个体（更短的路径）。
核心步骤：
tour的表示：用一个整数列表表示顺序（例如 [0,3,1,2,4] 表示访问城市顺序）。这种表示天然适合 TSP。
初始种群：随机生成若干条合法路径（或用启发式如最近邻生成一些良好起点）。
选择（Selection）：按“优胜劣汰”选择父代，常见方法：锦标赛选择（Tournament）、轮盘赌（Roulette）。锦标赛简单且稳定：随机抽 k 个个体，选最好的作为父本。
交叉（Crossover）：把两个父本合并产生后代。TSP 常用的有效交叉有 Order Crossover (OX)、Partially Mapped Crossover (PMX)。它们保证后代仍为合法排列（无重复城市）。
变异（Mutation）：引入随机扰动（swap、inversion、scramble 等），避免早熟收敛。
替代 / 精英策略（Replacement / Elitism）：保留最优个体直接进入下一代（防止好的解丢失）。
终止条件：达到固定代数、长时间无改进、或运行时间到达上限。
易实现、对复杂目标灵活、并行友好。参数较多（种群、交叉率、变异率等）、收敛速度与质量依赖细致调参。
常把 GA 与局部搜索（如 2-opt/3-opt）结合成“混合算法”，能显著提升解的质量。
代码实现：
-
距离矩阵计算与路径代价函数
-
初始种群（随机或用最近邻生成部分）
-
锦标赛选择、OX（Order Crossover）、PMX（可选）
-
变异：swap / inversion / scramble
-
可选 2-opt 局部改进（混合策略）
-
精英保留、历史记录与示例运行与可视化
"""
ga_tsp.py — 遗传算法 (GA) 求解 TSP 的示例实现
主要入口：genetic_algorithm_tsp(dist, ...)
返回 best_tour, best_cost, history（每代的最优成本）
作者：ChatGPT（示例代码，供学习与改造）
"""
import random
import math
import time
from typing import List, Tuple, Optional
import copy
# ---------- 基础函数 ----------
def compute_distance_matrix(points: List[Tuple[float, float]]) -> List[List[float]]:
n = len(points)
dist = [[0.0] * n for _ in range(n)]
for i in range(n):
xi, yi = points[i]
for j in range(i + 1, n):
xj, yj = points[j]
d = math.hypot(xi - xj, yi - yj)
dist[i][j] = d
dist[j][i] = d
return dist
def tour_length(tour: List[int], dist: List[List[float]]) -> float:
n = len(tour)
if n == 0:
return 0.0
total = 0.0
for i in range(n):
total += dist[tour[i]][tour[(i + 1) % n]]
return total
# --------- 选择（锦标赛） ----------
def tournament_selection(population: List[List[int]], costs: List[float], k: int = 3) -> List[int]:
# 随机抽 k 个（不重复）进行锦标赛，返回最好者（复制）
n = len(population)
if k <= 1:
# just pick best randomly among two
idx = random.randrange(n)
return population[idx][:]
contestants = random.sample(range(n), min(k, n))
best = min(contestants, key=lambda i: costs[i])
return population[best][:]
# --------- 交叉算子 ----------
def order_crossover(p1: List[int], p2: List[int]) -> List[int]:
# OX (Order Crossover)
n = len(p1)
child = [None] * n
a = random.randint(0, n - 2)
b = random.randint(a + 1, n - 1)
# copy segment from p1
child[a:b + 1] = p1[a:b + 1]
# fill remaining from p2 in order
pos = (b + 1) % n
for gene in p2:
if gene not in child:
child[pos] = gene
pos = (pos + 1) % n
return child
def pmx_crossover(p1: List[int], p2: List[int]) -> List[int]:
# PMX (Partially Mapped Crossover)
n = len(p1)
a = random.randint(0, n - 2)
b = random.randint(a + 1, n - 1)
child = [-1] * n
# copy slice
for i in range(a, b + 1):
child[i] = p1[i]
# mapping function
for i in range(a, b + 1):
if p2[i] not in child:
pos = i
val = p2[i]
while True:
# find position of val in p1
pos = p1.index(val)
if child[pos] == -1:
child[pos] = p2[i]
break
else:
val = p2[pos]
# fill remaining
for i in range(n):
if child[i] == -1:
child[i] = p2[i]
return child
# --------- 变异 ----------
def swap_mutation(ind: List[int], mutation_rate: float):
n = len(ind)
if random.random() < mutation_rate:
i, j = random.sample(range(n), 2)
ind[i], ind[j] = ind[j], ind[i]
def inversion_mutation(ind: List[int], mutation_rate: float):
n = len(ind)
if random.random() < mutation_rate:
a = random.randint(0, n - 2)
b = random.randint(a + 1, n - 1)
ind[a:b + 1] = ind[a:b + 1][::-1]
def scramble_mutation(ind: List[int], mutation_rate: float):
n = len(ind)
if random.random() < mutation_rate:
a = random.randint(0, n - 2)
b = random.randint(a + 1, n - 1)
segment = ind[a:b + 1][:]
random.shuffle(segment)
ind[a:b + 1] = segment
# --------- 可选：2-opt 局部改进（混合策略） ----------
def two_opt_delta(tour: List[int], i: int, j: int, dist: List[List[float]]) -> float:
n = len(tour)
a = tour[i - 1]
b = tour[i]
c = tour[j]
d = tour[(j + 1) % n]
return dist[a][c] + dist[b][d] - dist[a][b] - dist[c][d]
def two_opt_local_search(tour: List[int], dist: List[List[float]], max_iter: int = 50) -> List[int]:
# 简单的 2-opt 改进（贪心寻找首个改进）
n = len(tour)
if n < 4:
return tour
improved = True
it = 0
while improved and it < max_iter:
improved = False
it += 1
for i in range(1, n - 1):
for j in range(i + 1, n):
delta = two_opt_delta(tour, i, j, dist)
if delta < -1e-12:
# apply 2-opt swap
tour[i:j + 1] = reversed(tour[i:j + 1])
improved = True
break
if improved:
break
return tour
# --------- GA 主过程 ----------
def genetic_algorithm_tsp(
dist: List[List[float]],
pop_size: int = 100,
generations: int = 500,
crossover_rate: float = 0.9,
mutation_rate: float = 0.2,
elite_size: int = 1,
selection_k: int = 3,
crossover_type: str = "ox", # "ox" or "pmx"
mutation_type: str = "swap", # "swap", "inversion", "scramble"
initial_method: str = "random", # "random" or "nearest"
use_local_search: bool = False,
ls_max_iter: int = 20,
seed: Optional[int] = None,
verbose: bool = False
):
"""
返回 best_tour, best_cost, history（每代最优cost）
"""
if seed is not None:
random.seed(seed)
n = len(dist)
if n == 0:
return [], 0.0, []
# 初始化种群
population: List[List[int]] = []
if initial_method == "nearest":
# 产生若干最近邻起点的解（不同起点）
for start in range(min(pop_size, n)):
tour = nearest_neighbor_solution(dist, start=start)
population.append(tour)
# 若不够，再随机补齐
while len(population) < pop_size:
tour = list(range(n))
random.shuffle(tour)
population.append(tour)
# 计算初代代价
costs = [tour_length(ind, dist) for ind in population]
best_idx = min(range(pop_size), key=lambda i: costs[i])
best_tour = population[best_idx][:]
best_cost = costs[best_idx]
history = [best_cost]
# 迭代进化
start_time = time.time()
for gen in range(1, generations + 1):
# 精英保留：排序并保留 elite_size 个
pop_with_cost = list(zip(population, costs))
pop_with_cost.sort(key=lambda x: x[1])
new_population = [copy.deepcopy(ind) for ind, _c in pop_with_cost[:elite_size]]
# 生成其余个体
while len(new_population) < pop_size:
parent1 = tournament_selection(population, costs, k=selection_k)
parent2 = tournament_selection(population, costs, k=selection_k)
# crossover
if random.random() < crossover_rate:
if crossover_type == "pmx":
child1 = pmx_crossover(parent1, parent2)
child2 = pmx_crossover(parent2, parent1)
else:
child1 = order_crossover(parent1, parent2)
child2 = order_crossover(parent2, parent1)
else:
child1 = parent1[:]
child2 = parent2[:]
# mutation
if mutation_type == "swap":
swap_mutation(child1, mutation_rate)
swap_mutation(child2, mutation_rate)
elif mutation_type == "inversion":
inversion_mutation(child1, mutation_rate)
inversion_mutation(child2, mutation_rate)
else:
scramble_mutation(child1, mutation_rate)
scramble_mutation(child2, mutation_rate)
# optional local search
if use_local_search:
child1 = two_opt_local_search(child1, dist, max_iter=ls_max_iter)
child2 = two_opt_local_search(child2, dist, max_iter=ls_max_iter)
new_population.append(child1)
if len(new_population) < pop_size:
new_population.append(child2)
# 更新种群与代价
population = new_population
costs = [tour_length(ind, dist) for ind in population]
gen_best_idx = min(range(pop_size), key=lambda i: costs[i])
gen_best_cost = costs[gen_best_idx]
gen_best_tour = population[gen_best_idx][:]
if gen_best_cost < best_cost:
best_cost = gen_best_cost
best_tour = gen_best_tour[:]
history.append(best_cost)
if verbose and (gen % max(1, generations // 10) == 0 or gen <= 5):
elapsed = time.time() - start_time
print(f"Gen {gen}/{generations} best_cost={best_cost:.4f} elapsed={elapsed:.1f}s")
return best_tour, best_cost, history
# --------- 辅助：最近邻初解 ----------
def nearest_neighbor_solution(dist: List[List[float]], start: Optional[int] = None) -> List[int]:
n = len(dist)
if n == 0:
return []
if start is None:
start = random.randrange(n)
unvisited = set(range(n))
tour = [start]
unvisited.remove(start)
current = start
while unvisited:
next_node = min(unvisited, key=lambda x: dist[current][x])
tour.append(next_node)
unvisited.remove(next_node)
current = next_node
return tour
# ---------------- Demo / example ----------------
if __name__ == "__main__":
# 生成随机点并运行 GA
random.seed(1)
N = 50 # 城市数量，可修改尝试（例如 20, 50, 100）
points = [(random.random() * 100, random.random() * 100) for _ in range(N)]
dist = compute_distance_matrix(points)
best_tour, best_cost, history = genetic_algorithm_tsp(
dist,
pop_size=120,
generations=400,
crossover_rate=0.9,
mutation_rate=0.2,
elite_size=2,
selection_k=3,
crossover_type="ox",
mutation_type="inversion",
initial_method="nearest",
use_local_search=True, # 启用 2-opt 局部改进会显著提高解质量但更慢
ls_max_iter=10,
seed=1,
verbose=True
)
print("Best cost:", best_cost)
# 可视化（如果安装 matplotlib）
try:
import matplotlib.pyplot as plt
xs = [points[i][0] for i in best_tour] + [points[best_tour[0]][0]]
ys = [points[i][1] for i in best_tour] + [points[best_tour[0]][1]]
plt.figure(figsize=(6, 6))
plt.plot(xs, ys, '-o')
plt.title(f"GA TSP best (cost={best_cost:.2f})")
plt.show()
plt.figure()
plt.plot(history)
plt.xlabel("Generation")
plt.ylabel("Best cost so far")
plt.title("Convergence")
plt.show()
except Exception:
pass
其他经验：
种群大小（pop_size）：通常 50–300；种群大收敛更稳定但更慢。
代数（generations）：视问题规模，几十到几千代不等。
交叉率（crossover_rate）：常设 0.7–0.95。
变异率（mutation_rate）：一般较小 0.01–0.3（对排列编码，30%也常见），如果多代无改进可稍微调高。
选择方式与 k（selection_k）：锦标赛 k=2~5，k 越大选择压力越强（更快收敛但易陷入局部）。
精英数（elite_size）：1~5，保证好解不会丢失。
是否加局部搜索：把 GA 与 2-opt/3-opt 混合通常大幅提升解质量（代价是慢）。在后代上做一次轻量 2-opt 是常用实践。
多次运行取最好：GA 有随机性，重复多次取最好结果简单且有效。
为什么GA+局部搜索 混合方式较好？
GA 的全局搜索，与局部搜索的快速局部改进，优势互补，
GA 的全局搜索可以交叉产生多样性，局部搜索可以进行快速局部改进。前者探索结构化的解空间，后者快速消除局部不优。
tips：实际强算法通常是“元启发式 + 局部搜索”的组合（memetic algorithm）。
模拟退火
模拟退火算法（Simulated Annealing）: 借鉴了物理学中固体退火的过程，通过接受一些次优的解来跳出局部最优，从而寻找全局最优解。
目标：把 TSP 中的一条路径（tour）看成“系统状态”，把路径长度看成“能量（cost）”，我们希望使能量尽可能低。
邻域（如何从一个解跳到另一个解）：常用的邻域操作有 swap（交换两个城市）和 2-opt（选取 i<j，然后把 i..j 段翻转）。2-opt 在 TSP 中通常很有效。
温度（Temperature）：一个控制随机接受劣解概率的参数。温度高时更容易接受变差的解（允许跳出局部最优）；温度低时趋向只接受更好的解。
接受准则：若新解更好（cost 减少），直接接受；若更差，按概率 exp(-Δ / T) 接受（Δ = 新cost - 旧cost）。
降温（Cooling schedule）：常用指数冷却：T = alpha * T，alpha 在 [0.90, 0.9999] 中，越接近 1 冷却越慢。
停止条件：温度降到某个阈值或达到最大迭代次数或长时间没有改进。
2-opt（选取 i<j，然后把 i..j 段翻转）翻转一段能消除两条交叉边，常能显著改善路径，而且用解析式可以快速计算成本变化（不用重算整条路径），因此速度与效果都不错。
import math
import random
import time
from typing import List, Tuple, Optional
def compute_distance_matrix(points: List[Tuple[float,float]]) -> List[List[float]]:
"""Compute symmetric Euclidean distance matrix."""
n = len(points)
dist = [[0.0]*n for _ in range(n)]
for i in range(n):
xi, yi = points[i]
for j in range(i+1, n):
xj, yj = points[j]
d = math.hypot(xi-xj, yi-yj)
dist[i][j] = d
dist[j][i] = d
return dist
def tour_length(tour: List[int], dist: List[List[float]]) -> float:
n = len(tour)
total = 0.0
for i in range(n):
total += dist[tour[i]][tour[(i+1)%n]]
return total
def nearest_neighbor_solution(dist: List[List[float]], start: Optional[int]=None) -> List[int]:
"""Greedy nearest-neighbor initial solution (fast, usually reasonable)."""
n = len(dist)
if n == 0:
return []
if start is None:
start = random.randrange(n)
unvisited = set(range(n))
tour = [start]
unvisited.remove(start)
current = start
while unvisited:
next_node = min(unvisited, key=lambda x: dist[current][x])
tour.append(next_node)
unvisited.remove(next_node)
current = next_node
return tour
def two_opt_delta(tour: List[int], i: int, j: int, dist: List[List[float]]) -> float:
"""
Compute cost difference if we reverse tour[i:j+1] (2-opt).
delta = new_cost - old_cost. Positive -> worse.
Uses only 4 edge lookups (constant time).
"""
n = len(tour)
a = tour[i-1]
b = tour[i]
c = tour[j]
d = tour[(j+1) % n]
# Remove edges (a-b) and (c-d), add (a-c) and (b-d)
return dist[a][c] + dist[b][d] - dist[a][b] - dist[c][d]
def simulated_annealing_tsp(
dist: List[List[float]],
initial_solution: str = "nearest", # "nearest" or "random"
initial_temp: Optional[float] = None,
stopping_temp: float = 1e-3,
alpha: float = 0.995,
iter_per_temp: Optional[int] = None,
max_iterations: int = 1_000_000,
seed: Optional[int] = None,
verbose: bool = False
):
"""
Simulated annealing for symmetric TSP using 2-opt neighbor.
Returns best_tour, best_cost, history (list of best_costs over time).
"""
if seed is not None:
random.seed(seed)
n = len(dist)
if n == 0:
return [], 0.0, []
# initial solution
if initial_solution == "nearest":
current_tour = nearest_neighbor_solution(dist, start=0)
else:
current_tour = list(range(n))
random.shuffle(current_tour)
current_cost = tour_length(current_tour, dist)
best_tour = current_tour[:]
best_cost = current_cost
# automatic temperature initialization based on average edge length
if initial_temp is None:
sum_edges = 0.0
for i in range(n):
for j in range(i+1, n):
sum_edges += dist[i][j]
avg_edge = (2 * sum_edges) / (n * (n - 1)) if n > 1 else 0.0
initial_temp = avg_edge * n if avg_edge > 0 else 1.0
if iter_per_temp is None:
iter_per_temp = 100 * n # heuristic: try ~100*n attempts per temperature
T = initial_temp
iteration = 0
history = [best_cost]
start_time = time.time()
while T > stopping_temp and iteration < max_iterations:
for _ in range(iter_per_temp):
# choose i < j for 2-opt
i = random.randint(0, n-2)
j = random.randint(i+1, n-1)
# skip reversing the entire tour (gives same tour)
if i == 0 and j == n-1:
continue
delta = two_opt_delta(current_tour, i, j, dist)
# accept if better or by probability
if delta < 0 or random.random() < math.exp(-delta / T):
# apply 2-opt: reverse the segment i..j
current_tour = current_tour[:i] + current_tour[i:j+1][::-1] + current_tour[j+1:]
current_cost += delta
# update best
if current_cost < best_cost:
best_cost = current_cost
best_tour = current_tour[:]
history.append(best_cost)
iteration += 1
if iteration >= max_iterations:
break
T *= alpha
if verbose and iteration % (10 * iter_per_temp) == 0:
elapsed = time.time() - start_time
print(f"iter={iteration}, T={T:.5f}, best_cost={best_cost:.4f}, elapsed={elapsed:.1f}s")
return best_tour, best_cost, history
# ---------------- Example usage ----------------
if __name__ == "__main__":
import random
random.seed(0)
N = 50
# generate random points in 0..100 square
pts = [(random.random()*100, random.random()*100) for _ in range(N)]
d = compute_distance_matrix(pts)
# run SA
best, cost, hist = simulated_annealing_tsp(
d,
initial_solution="nearest",
alpha=0.995,
iter_per_temp=200*len(d), # more attempts per T -> better but slower
seed=0,
verbose=True
)
print("best cost:", cost)
# optional plotting (requires matplotlib)
try:
import matplotlib.pyplot as plt
def plot_tour(points, tour, title=None):
xs = [points[i][0] for i in tour] + [points[tour[0]][0]]
ys = [points[i][1] for i in tour] + [points[tour[0]][1]]
plt.figure(figsize=(6,6))
plt.plot(xs, ys, '-o')
if title:
plt.title(title)
plt.show()
plot_tour(pts, best, f"SA best (cost={cost:.2f})")
except Exception:
pass
蚁群优化算法
蚁群优化算法（Ant Colony Optimization）: 灵感来源于蚂蚁寻找食物的路径行为。蚂蚁在寻找食物时会在路径上留下信息素（pheromone），其他蚂蚁更倾向于沿着信息素浓度高的路径走。好的路径（短）会被更多蚂蚁强化，从而逐渐成为主流解。
一个“蚂蚁”构造一条 TSP 路径（tour），路径长度为目标（越短越好）。
两种信息：
-
信息素（pheromone）τ(i,j)：基于历史经验，某条边的吸引力（学习到的好/坏）。
-
启发式信息（visibility）η(i,j)：瞬时启发式，通常取
1 / distance(i,j)，也就是更倾向走短边。
转移概率：蚂蚁从城市 i 选择下一个城市 j 的概率与 (τ(i,j)^α) * (η(i,j)^β) 成正比，α 控制信息素重要性，β 控制启发式重要性。
信息素更新：
-
挥发（evaporation）：所有边的 τ 乘上
(1 - ρ)，避免过早收敛。 -
沉积（deposit）：每只（或最优）蚂蚁在其经过的边上增加
Δτ = Q / L（L 为该蚂蚁路径长度），优良路径贡献更大。
优缺点：
适合组合优化、容易并行化、能自然平衡探索/利用。缺点：参数较敏感、对大规模问题单纯 ACO 收敛慢。
把 ACO 和局部搜索（如 2-opt/3-opt）结合，或采用 Max-Min ACO/Ant Colony System 的变种以改进收敛。
"""
aco_tsp.py — 教学版 ACO 求解 TSP
返回值：
best_tour, best_cost, history
主要参数说明见函数注释（ant_colony_tsp）。
"""
import random
import math
import time
from typing import List, Tuple, Optional
# ---------- 基础工具 ----------
def compute_distance_matrix(points: List[Tuple[float, float]]) -> List[List[float]]:
n = len(points)
dist = [[0.0]*n for _ in range(n)]
for i in range(n):
xi, yi = points[i]
for j in range(i+1, n):
xj, yj = points[j]
d = math.hypot(xi - xj, yi - yj)
dist[i][j] = d
dist[j][i] = d
return dist
def tour_length(tour: List[int], dist: List[List[float]]) -> float:
n = len(tour)
if n == 0:
return 0.0
total = 0.0
for i in range(n):
total += dist[tour[i]][tour[(i+1) % n]]
return total
def nearest_neighbor_solution(dist: List[List[float]], start: int = 0) -> List[int]:
n = len(dist)
unvisited = set(range(n))
tour = [start]
unvisited.remove(start)
current = start
while unvisited:
next_node = min(unvisited, key=lambda x: dist[current][x])
tour.append(next_node)
unvisited.remove(next_node)
current = next_node
return tour
# ---------- 可选的 2-opt 局部改进（提升解质量） ----------
def two_opt_delta(tour: List[int], i: int, j: int, dist: List[List[float]]) -> float:
n = len(tour)
a = tour[i-1]
b = tour[i]
c = tour[j]
d = tour[(j+1) % n]
return dist[a][c] + dist[b][d] - dist[a][b] - dist[c][d]
def two_opt_local_search(tour: List[int], dist: List[List[float]], max_iter: int = 50) -> List[int]:
n = len(tour)
if n < 4:
return tour
improved = True
it = 0
while improved and it < max_iter:
improved = False
it += 1
for i in range(1, n-1):
for j in range(i+1, n):
delta = two_opt_delta(tour, i, j, dist)
if delta < -1e-12:
tour[i:j+1] = reversed(tour[i:j+1])
improved = True
break
if improved:
break
return tour
# ---------- ACO 主算法 ----------
def ant_colony_tsp(
dist: List[List[float]],
n_ants: Optional[int] = None,
n_iterations: int = 200,
alpha: float = 1.0, # 信息素重要度
beta: float = 2.0, # 启发函数(visibility)重要度
rho: float = 0.1, # 挥发率
q: float = 1.0, # 信息素强度因子（沉积时 Q / L）
deposit: str = "iteration_best", # "all" / "iteration_best" / "global_best"
initial_pheromone: Optional[float] = None,
use_2opt: bool = False,
ls_max_iter: int = 10,
seed: Optional[int] = None,
verbose: bool = False
):
"""
dist: 对称距离矩阵
n_ants: 蚂蚁数，默认等于城市数
n_iterations: 迭代次数
alpha, beta: 控制转移概率中信息素和启发信息的权重
rho: 信息素挥发率 (0 < rho < 1)
q: 信息素沉积常数（通常与问题规模和距离尺度相关）
deposit: 沉积策略，"all" 所有蚂蚁沉积, "iteration_best" 本代最佳沉积, "global_best" 全局最佳沉积
initial_pheromone: 如果 None，会根据最近邻解自动估计一个合理值
use_2opt: 是否对每个蚂蚁（或最佳）做 2-opt 局部改进
ls_max_iter: 2-opt 的最大迭代次数（每次局部改进的限制）
seed: 随机种子（便于复现）
"""
if seed is not None:
random.seed(seed)
n = len(dist)
if n == 0:
return [], 0.0, []
if n_ants is None:
n_ants = n
# 计算启发信息（visibility = 1 / distance）
eps = 1e-12
visibility = [[0.0]*n for _ in range(n)]
for i in range(n):
for j in range(n):
if i != j:
visibility[i][j] = 1.0 / (dist[i][j] + eps)
# 初始信息素估计（用最近邻长度做一个启发）
if initial_pheromone is None:
try:
nn = nearest_neighbor_solution(dist, start=0)
L_nn = tour_length(nn, dist)
initial_pheromone = 1.0 / (n * L_nn) if L_nn > 0 else 1.0
except Exception:
initial_pheromone = 1.0
# 信息素矩阵（对称）
tau = [[initial_pheromone]*n for _ in range(n)]
for i in range(n):
tau[i][i] = 0.0
best_tour = None
best_cost = float('inf')
history = []
start_time = time.time()
for it in range(1, n_iterations+1):
all_tours = []
all_costs = []
# 每只蚂蚁构造解
for a in range(n_ants):
start_city = random.randrange(n)
tour = [start_city]
visited = set(tour)
current = start_city
while len(tour) < n:
# 计算可选城市的选择概率（未访问城市）
probs = []
denom = 0.0
for j in range(n):
if j in visited:
probs.append(0.0)
else:
val = (tau[current][j] ** alpha) * (visibility[current][j] ** beta)
probs.append(val)
denom += val
if denom == 0.0:
# 若所有值都为 0（极端情况），随机选择一个未访问城市
choices = [c for c in range(n) if c not in visited]
next_city = random.choice(choices)
else:
# 轮盘赌选择
r = random.random()
cumulative = 0.0
next_city = None
for j, p in enumerate(probs):
if j in visited:
continue
cumulative += p / denom
if r <= cumulative:
next_city = j
break
if next_city is None:
# 数值误差导致未选中，退化为随机选择
choices = [c for c in range(n) if c not in visited]
next_city = random.choice(choices)
tour.append(next_city)
visited.add(next_city)
current = next_city
# 可选局部改进
if use_2opt:
tour = two_opt_local_search(tour, dist, max_iter=ls_max_iter)
cost = tour_length(tour, dist)
all_tours.append(tour)
all_costs.append(cost)
# 找出本代最优解
iter_best_idx = min(range(len(all_costs)), key=lambda i: all_costs[i])
iter_best_cost = all_costs[iter_best_idx]
iter_best_tour = all_tours[iter_best_idx]
if iter_best_cost < best_cost:
best_cost = iter_best_cost
best_tour = iter_best_tour[:]
history.append(best_cost)
# 信息素挥发
for i in range(n):
for j in range(n):
tau[i][j] *= (1.0 - rho)
if tau[i][j] < 1e-16:
tau[i][j] = 1e-16 # 避免数值归零
# 信息素沉积
if deposit == "all":
# 所有蚂蚁沉积
for k, tour in enumerate(all_tours):
Lk = all_costs[k]
delta = q / (Lk + 1e-12)
for idx in range(n):
i = tour[idx]
j = tour[(idx+1) % n]
tau[i][j] += delta
tau[j][i] += delta
elif deposit == "iteration_best":
# 本代最好蚂蚁沉积
delta = q / (iter_best_cost + 1e-12)
for idx in range(n):
i = iter_best_tour[idx]
j = iter_best_tour[(idx+1) % n]
tau[i][j] += delta
tau[j][i] += delta
elif deposit == "global_best":
# 全局最好蚂蚁沉积
delta = q / (best_cost + 1e-12)
for idx in range(n):
i = best_tour[idx]
j = best_tour[(idx+1) % n]
tau[i][j] += delta
tau[j][i] += delta
else:
# 默认与 iteration_best 相同
delta = q / (iter_best_cost + 1e-12)
for idx in range(n):
i = iter_best_tour[idx]
j = iter_best_tour[(idx+1) % n]
tau[i][j] += delta
tau[j][i] += delta
if verbose and (it % max(1, n_iterations//10) == 0 or it <= 5):
elapsed = time.time() - start_time
print(f"Iter {it}/{n_iterations} best_cost_so_far={best_cost:.4f} iter_best={iter_best_cost:.4f} elapsed={elapsed:.1f}s")
return best_tour, best_cost, history
# ------------------ 示例与可视化 ------------------
if __name__ == "__main__":
# 生成随机点并运行示例
random.seed(0)
N = 40
points = [(random.random()*100, random.random()*100) for _ in range(N)]
dist = compute_distance_matrix(points)
best_tour, best_cost, history = ant_colony_tsp(
dist,
n_ants=30,
n_iterations=200,
alpha=1.0,
beta=3.0,
rho=0.1,
q=1.0,
deposit="iteration_best",
use_2opt=True,
ls_max_iter=8,
seed=0,
verbose=True
)
print("Best cost:", best_cost)
# 可视化（需要 matplotlib）
try:
import matplotlib.pyplot as plt
xs = [points[i][0] for i in best_tour] + [points[best_tour[0]][0]]
ys = [points[i][1] for i in best_tour] + [points[best_tour[0]][1]]
plt.figure(figsize=(6,6))
plt.plot(xs, ys, '-o')
plt.title(f"ACO best (cost={best_cost:.2f})")
plt.show()
plt.figure()
plt.plot(history)
plt.xlabel("Iteration")
plt.ylabel("Best cost so far")
plt.title("Convergence")
plt.show()
except Exception:
pass
<!DOCTYPE html>
<html lang="zh-CN">
<head>
<meta charset="UTF-8">
<meta name="viewport" content="width=device-width, initial-scale=1.0">
<title>TSP 蚁群算法动画演示</title>
<style>
body {
margin: 0;
padding: 20px;
font-family: 'Arial', sans-serif;
background: linear-gradient(135deg, #667eea 0%, #764ba2 100%);
min-height: 100vh;
color: white;
}
.container {
max-width: 1200px;
margin: 0 auto;
background: rgba(255, 255, 255, 0.1);
backdrop-filter: blur(10px);
border-radius: 20px;
padding: 30px;
box-shadow: 0 8px 32px rgba(0, 0, 0, 0.2);
}
h1 {
text-align: center;
margin-bottom: 30px;
font-size: 2.5em;
text-shadow: 2px 2px 4px rgba(0, 0, 0, 0.3);
}
.controls {
display: flex;
justify-content: center;
gap: 15px;
margin-bottom: 30px;
flex-wrap: wrap;
}
.control-group {
display: flex;
flex-direction: column;
align-items: center;
gap: 5px;
}
button {
padding: 12px 24px;
border: none;
border-radius: 25px;
background: linear-gradient(45deg, #4CAF50, #45a049);
color: white;
font-size: 16px;
cursor: pointer;
transition: all 0.3s ease;
box-shadow: 0 4px 15px rgba(0, 0, 0, 0.2);
}
button:hover {
transform: translateY(-2px);
box-shadow: 0 6px 20px rgba(0, 0, 0, 0.3);
}
button:disabled {
background: linear-gradient(45deg, #666, #555);
cursor: not-allowed;
transform: none;
}
.reset-btn {
background: linear-gradient(45deg, #f44336, #d32f2f);
}
input[type="range"] {
width: 150px;
margin: 0 10px;
}
label {
font-size: 14px;
color: rgba(255, 255, 255, 0.9);
}
#canvas {
border: 2px solid rgba(255, 255, 255, 0.3);
border-radius: 15px;
background: rgba(0, 0, 0, 0.1);
display: block;
margin: 0 auto;
box-shadow: 0 8px 32px rgba(0, 0, 0, 0.3);
}
.stats {
display: flex;
justify-content: space-around;
margin-top: 20px;
gap: 20px;
flex-wrap: wrap;
}
.stat-card {
background: rgba(255, 255, 255, 0.1);
padding: 15px 25px;
border-radius: 15px;
text-align: center;
min-width: 120px;
backdrop-filter: blur(5px);
}
.stat-title {
font-size: 14px;
opacity: 0.8;
margin-bottom: 5px;
}
.stat-value {
font-size: 24px;
font-weight: bold;
color: #4CAF50;
}
.info {
margin-top: 20px;
padding: 15px;
background: rgba(255, 255, 255, 0.1);
border-radius: 10px;
font-size: 14px;
line-height: 1.6;
}
.legend {
display: flex;
justify-content: center;
gap: 30px;
margin-top: 15px;
flex-wrap: wrap;
}
.legend-item {
display: flex;
align-items: center;
gap: 8px;
}
.legend-color {
width: 20px;
height: 20px;
border-radius: 50%;
border: 2px solid white;
}
</style>
</head>
<body>
<div class="container">
<h1>🐜 TSP 蚁群算法动画演示</h1>
<div class="controls">
<div class="control-group">
<label>城市数量</label>
<input type="range" id="cityCount" min="5" max="20" value="10">
<span id="cityCountValue">10</span>
</div>
<div class="control-group">
<label>蚂蚁数量</label>
<input type="range" id="antCount" min="10" max="50" value="20">
<span id="antCountValue">20</span>
</div>
<div class="control-group">
<label>动画速度</label>
<input type="range" id="speed" min="1" max="10" value="5">
<span id="speedValue">5</span>
</div>
</div>
<div class="controls">
<button id="generateBtn">🏙️ 生成新城市</button>
<button id="startBtn">🚀 开始算法</button>
<button id="pauseBtn" disabled>⏸️ 暂停</button>
<button id="resetBtn" class="reset-btn">🔄 重置</button>
</div>
<canvas id="canvas" width="800" height="600"></canvas>
<div class="stats">
<div class="stat-card">
<div class="stat-title">当前代数</div>
<div class="stat-value" id="generation">0</div>
</div>
<div class="stat-card">
<div class="stat-title">最佳距离</div>
<div class="stat-value" id="bestDistance">∞</div>
</div>
<div class="stat-card">
<div class="stat-title">活跃蚂蚁</div>
<div class="stat-value" id="activeAnts">0</div>
</div>
<div class="stat-card">
<div class="stat-title">信息素强度</div>
<div class="stat-value" id="pheromoneLevel">1.0</div>
</div>
</div>
<div class="legend">
<div class="legend-item">
<div class="legend-color" style="background: #FFD700;"></div>
<span>城市</span>
</div>
<div class="legend-item">
<div class="legend-color" style="background: #FF6B6B;"></div>
<span>蚂蚁</span>
</div>
<div class="legend-item">
<div class="legend-color" style="background: #4ECDC4;"></div>
<span>最佳路径</span>
</div>
<div class="legend-item">
<div class="legend-color" style="background: rgba(100, 200, 100, 0.3);"></div>
<span>信息素</span>
</div>
</div>
<div class="info">
<strong>🧠 算法原理：</strong>
蚁群算法模拟真实蚂蚁寻找食物的行为。蚂蚁在路径上留下信息素，越短的路径信息素浓度越高，
后续蚂蚁更倾向于选择信息素浓度高的路径。通过多代迭代，最终找到近似最优解。
<br><br>
<strong>📊 参数说明：</strong>
信息素重要性(α)=1，启发式重要性(β)=2，信息素挥发率(ρ)=0.1，信息素增强系数(Q)=100
</div>
</div>
<script>
class TSPAntColony {
constructor() {
this.canvas = document.getElementById('canvas');
this.ctx = this.canvas.getContext('2d');
this.width = this.canvas.width;
this.height = this.canvas.height;
// 算法参数
this.alpha = 1; // 信息素重要性
this.beta = 2; // 启发式重要性
this.rho = 0.1; // 信息素挥发率
this.Q = 100; // 信息素增强系数
// 状态变量
this.cities = [];
this.ants = [];
this.pheromones = [];
this.distances = [];
this.bestPath = [];
this.bestDistance = Infinity;
this.generation = 0;
this.isRunning = false;
this.animationId = null;
this.initializeEventListeners();
this.generateRandomCities();
}
initializeEventListeners() {
// 控制按钮
document.getElementById('generateBtn').onclick = () => this.generateRandomCities();
document.getElementById('startBtn').onclick = () => this.startAlgorithm();
document.getElementById('pauseBtn').onclick = () => this.pauseAlgorithm();
document.getElementById('resetBtn').onclick = () => this.resetAlgorithm();
// 滑动条更新
['cityCount', 'antCount', 'speed'].forEach(id => {
const slider = document.getElementById(id);
const valueSpan = document.getElementById(id + 'Value');
slider.oninput = () => {
valueSpan.textContent = slider.value;
if (id === 'cityCount' && !this.isRunning) {
this.generateRandomCities();
}
};
});
}
generateRandomCities() {
const count = parseInt(document.getElementById('cityCount').value);
this.cities = [];
// 生成随机城市位置
for (let i = 0; i < count; i++) {
this.cities.push({
x: Math.random() * (this.width - 100) + 50,
y: Math.random() * (this.height - 100) + 50,
id: i
});
}
this.calculateDistances();
this.initializePheromones();
this.resetStats();
this.draw();
}
calculateDistances() {
const n = this.cities.length;
this.distances = Array(n).fill().map(() => Array(n).fill(0));
for (let i = 0; i < n; i++) {
for (let j = 0; j < n; j++) {
if (i !== j) {
const dx = this.cities[i].x - this.cities[j].x;
const dy = this.cities[i].y - this.cities[j].y;
this.distances[i][j] = Math.sqrt(dx * dx + dy * dy);
}
}
}
}
initializePheromones() {
const n = this.cities.length;
this.pheromones = Array(n).fill().map(() => Array(n).fill(1));
}
startAlgorithm() {
if (this.cities.length < 3) return;
this.isRunning = true;
document.getElementById('startBtn').disabled = true;
document.getElementById('pauseBtn').disabled = false;
document.getElementById('generateBtn').disabled = true;
this.runGeneration();
}
pauseAlgorithm() {
this.isRunning = false;
if (this.animationId) {
cancelAnimationFrame(this.animationId);
}
document.getElementById('startBtn').disabled = false;
document.getElementById('pauseBtn').disabled = true;
document.getElementById('generateBtn').disabled = false;
}
resetAlgorithm() {
this.pauseAlgorithm();
this.resetStats();
this.initializePheromones();
this.draw();
}
resetStats() {
this.bestPath = [];
this.bestDistance = Infinity;
this.generation = 0;
this.ants = [];
this.updateStats();
}
runGeneration() {
if (!this.isRunning) return;
// 创建蚂蚁
const antCount = parseInt(document.getElementById('antCount').value);
this.ants = [];
for (let i = 0; i < antCount; i++) {
this.ants.push({
currentCity: Math.floor(Math.random() * this.cities.length),
visitedCities: [],
path: [],
totalDistance: 0,
isComplete: false
});
}
this.generation++;
this.animateAnts();
}
animateAnts() {
const speed = parseInt(document.getElementById('speed').value);
let antIndex = 0;
let stepCount = 0;
const animate = () => {
if (!this.isRunning) return;
// 移动蚂蚁
let allComplete = true;
for (let ant of this.ants) {
if (!ant.isComplete) {
allComplete = false;
this.moveAnt(ant);
}
}
this.draw();
this.updateStats();
if (allComplete) {
// 更新信息素
this.updatePheromones();
// 继续下一代
setTimeout(() => {
if (this.isRunning) {
this.runGeneration();
}
}, 1000 / speed);
} else {
this.animationId = requestAnimationFrame(() => {
setTimeout(animate, 100 / speed);
});
}
};
animate();
}
moveAnt(ant) {
if (ant.path.length === 0) {
ant.path.push(ant.currentCity);
ant.visitedCities.push(ant.currentCity);
}
if (ant.visitedCities.length === this.cities.length) {
// 回到起始城市
const startCity = ant.path[0];
ant.path.push(startCity);
ant.totalDistance += this.distances[ant.currentCity][startCity];
ant.isComplete = true;
// 更新最佳路径
if (ant.totalDistance < this.bestDistance) {
this.bestDistance = ant.totalDistance;
this.bestPath = [...ant.path];
}
return;
}
// 选择下一个城市
const nextCity = this.selectNextCity(ant);
if (nextCity !== -1) {
ant.totalDistance += this.distances[ant.currentCity][nextCity];
ant.currentCity = nextCity;
ant.path.push(nextCity);
ant.visitedCities.push(nextCity);
}
}
selectNextCity(ant) {
const unvisited = [];
for (let i = 0; i < this.cities.length; i++) {
if (!ant.visitedCities.includes(i)) {
unvisited.push(i);
}
}
if (unvisited.length === 0) return -1;
// 计算概率
const probabilities = [];
let totalProbability = 0;
for (let city of unvisited) {
const pheromone = Math.pow(this.pheromones[ant.currentCity][city], this.alpha);
const heuristic = Math.pow(1.0 / this.distances[ant.currentCity][city], this.beta);
const probability = pheromone * heuristic;
probabilities.push(probability);
totalProbability += probability;
}
// 轮盘赌选择
let random = Math.random() * totalProbability;
for (let i = 0; i < unvisited.length; i++) {
random -= probabilities[i];
if (random <= 0) {
return unvisited[i];
}
}
return unvisited[unvisited.length - 1];
}
updatePheromones() {
const n = this.cities.length;
// 信息素挥发
for (let i = 0; i < n; i++) {
for (let j = 0; j < n; j++) {
this.pheromones[i][j] *= (1 - this.rho);
}
}
// 蚂蚁留下信息素
for (let ant of this.ants) {
if (ant.isComplete) {
const delta = this.Q / ant.totalDistance;
for (let i = 0; i < ant.path.length - 1; i++) {
const from = ant.path[i];
const to = ant.path[i + 1];
this.pheromones[from][to] += delta;
this.pheromones[to][from] += delta;
}
}
}
}
draw() {
this.ctx.clearRect(0, 0, this.width, this.height);
// 绘制信息素
this.drawPheromones();
// 绘制最佳路径
if (this.bestPath.length > 0) {
this.drawPath(this.bestPath, '#4ECDC4', 4);
}
// 绘制城市
this.drawCities();
// 绘制蚂蚁
this.drawAnts();
}
drawPheromones() {
const maxPheromone = Math.max(...this.pheromones.flat());
for (let i = 0; i < this.cities.length; i++) {
for (let j = i + 1; j < this.cities.length; j++) {
const intensity = this.pheromones[i][j] / maxPheromone;
if (intensity > 0.1) {
this.ctx.strokeStyle = `rgba(100, 200, 100, ${intensity * 0.5})`;
this.ctx.lineWidth = intensity * 3;
this.ctx.beginPath();
this.ctx.moveTo(this.cities[i].x, this.cities[i].y);
this.ctx.lineTo(this.cities[j].x, this.cities[j].y);
this.ctx.stroke();
}
}
}
}
drawCities() {
this.ctx.fillStyle = '#FFD700';
this.ctx.strokeStyle = '#FFA500';
this.ctx.lineWidth = 2;
for (let i = 0; i < this.cities.length; i++) {
const city = this.cities[i];
this.ctx.beginPath();
this.ctx.arc(city.x, city.y, 8, 0, 2 * Math.PI);
this.ctx.fill();
this.ctx.stroke();
// 城市标签
this.ctx.fillStyle = 'white';
this.ctx.font = '12px Arial';
this.ctx.textAlign = 'center';
this.ctx.fillText(i.toString(), city.x, city.y + 4);
this.ctx.fillStyle = '#FFD700';
}
}
drawAnts() {
for (let ant of this.ants) {
if (!ant.isComplete && ant.path.length > 0) {
const city = this.cities[ant.currentCity];
// 蚂蚁身体
this.ctx.fillStyle = '#FF6B6B';
this.ctx.beginPath();
this.ctx.arc(city.x, city.y, 4, 0, 2 * Math.PI);
this.ctx.fill();
// 蚂蚁路径
if (ant.path.length > 1) {
this.ctx.strokeStyle = 'rgba(255, 107, 107, 0.3)';
this.ctx.lineWidth = 1;
this.ctx.beginPath();
for (let i = 0; i < ant.path.length - 1; i++) {
const from = this.cities[ant.path[i]];
const to = this.cities[ant.path[i + 1]];
if (i === 0) {
this.ctx.moveTo(from.x, from.y);
}
this.ctx.lineTo(to.x, to.y);
}
this.ctx.stroke();
}
}
}
}
drawPath(path, color, width) {
if (path.length < 2) return;
this.ctx.strokeStyle = color;
this.ctx.lineWidth = width;
this.ctx.setLineDash([]);
this.ctx.beginPath();
const startCity = this.cities[path[0]];
this.ctx.moveTo(startCity.x, startCity.y);
for (let i = 1; i < path.length; i++) {
const city = this.cities[path[i]];
this.ctx.lineTo(city.x, city.y);
}
this.ctx.stroke();
}
updateStats() {
document.getElementById('generation').textContent = this.generation;
document.getElementById('bestDistance').textContent =
this.bestDistance === Infinity ? '∞' : Math.round(this.bestDistance);
document.getElementById('activeAnts').textContent =
this.ants.filter(ant => !ant.isComplete).length;
if (this.pheromones.length > 0) {
const avgPheromone = this.pheromones.flat().reduce((a, b) => a + b, 0) /
(this.pheromones.length * this.pheromones.length);
document.getElementById('pheromoneLevel').textContent = avgPheromone.toFixed(2);
}
}
}
// 初始化应用
window.addEventListener('load', () => {
new TSPAntColony();
});
</script>
</body>
</html>
https://www.ceas3.uc.edu/ret/archive/2016/ret/docs/abstract/reading33.pdf
浙公网安备 33010602011771号