matlab实现遗传算法解决TSP问题
文章基于 遗传算法--旅行商问题(TSP问题)-Matlab 上进行,优化了计算的时间成本问题
1、旅行商问题(TSP问题)
假设有一个旅行商人要拜访全国31个省会城市,它需要选择所要走的路径,路径的限制是每个城市只能拜访一次,而且最后要回到原来出发的城市。对路径选择的要求是:所选路径的路成为所有路径之中的最小值。
全国31个省会城市的坐标为
[1304 2312;3639 1315;4177 2244;3712 1399;3488 1535;3326 1556;3238 1229;4196 1044;
4312 790;4386 570;3007 1970;2562 1756;2788 1491;2381 1676;1332 695;3715 1678;
3918 2179;4061 2370;3780 2212;3676 2578;4029 2838;4263 2931;3429 1908;3507 2376;
3394 2643;3439 3201;2935 3240;3140 3550;2545 2357;2778 2826;2370 2975]
2、仿真过程
原思路:
(1)初始化种群数目NP=200,染色体基因维数为N=31,最大进化代数G=1000.
(2)产生初始种群,计算个体适应度值,即路径长度:采用基于概率的方式选择进行操作的个体;对选中的成对个体,随机交叉所选中的成对城市坐标,以确保交叉后路径每个城市只到访一次;对选中的单个个体,随机交换其一对城市坐标作为变异操作,产生新的种群,进行下一次遗传操作。
(3)判断是否满足终止条件:若满足,则结束搜索过程,输出优化值,若不满足,则继续进迭代优化。
(4)为迭代过程预设一个文件,用于存储保存每个迭代的图像帧(仅当当前图与上次不同),制作迭代过程图
clear all; % 清除所有变量
close all; % 清图
clc; % 清屏
C = [1304 2312; 3639 1315; 4177 2244; 3712 1399; 3488 1535; 3326 1556; ...
3238 1229; 4196 1044; 4312 790; 4386 570; 3007 1970; 2562 1756; ...
2788 1491; 2381 1676; 1332 695; 3715 1678; 3918 2179; 4061 2370; ...
3780 2212; 3676 2578; 4029 2838; 4263 2931; 3429 1908; 3507 2376; ...
3394 2643; 3439 3201; 2935 3240; 3140 3550; 2545 2357; 2778 2826; ...
2370 2975]; %31个省会城市坐标
N = size(C,1); % TSP问题的规模, 即城市数目
D = zeros(N); % 任意两个城市距离间隔矩阵
% 计算任意两个城市的距离
for i = 1:N
for j = 1:N
D(i,j) = sqrt((C(i,1) - C(j,1))^2 + (C(i,2) - C(j,2))^2);
end
end
NP = 200; % 种群规模
G = 2500; % 最大遗传代数
f = zeros(NP, N); % 用于存储种群
F = []; % 种群更新中间存储
for i = 1:NP
f(i,:) = randperm(N); % 随机生成初始种群
end
R = f(1,:); % 存储最优种群
len = zeros(NP, 1); % 存储路径长度
fitness = zeros(NP, 1); % 存储归一化适应值
gen = 0;
% 为动画预设一个文件
video = VideoWriter('TSP_GA_Animation', 'MPEG-4');
open(video);
% 上一次保存的路径图
last_frame = [];
% 遗传算法循环
while gen < G
% 计算路径长度
for i = 1:NP
len(i,1) = D(f(i,N), f(i,1));
for j = 1:(N-1)
len(i,1) = len(i,1) + D(f(i,j), f(i,j+1));
end
end
maxlen = max(len); % 最长路径
minlen = min(len); % 最短路径
% 更新最短路径
rr = find(len == minlen);
R = f(rr(1,1),:);
% 计算归一化适应值
for i = 1:length(len)
fitness(i,1) = (1 - ((len(i,1) - minlen) / (maxlen - minlen + 0.001)));
end
% 选择操作
nn = 0;
for i = 1:NP
if fitness(i,1) >= rand
nn = nn + 1;
F(nn,:) = f(i,:);
end
end
[aa, bb] = size(F);
while aa < NP
nnper = randperm(nn);
A = F(nnper(1), :);
B = F(nnper(2), :);
% 交叉操作
W = ceil(N / 10); % 交叉点个数
p = unidrnd(N - W + 1); % 随机选择交叉范围,从p到p+W
for i = 1:W
x = find(A == B(p+i-1));
y = find(B == A(p+i-1));
temp = A(p+i-1);
A(p+i-1) = B(p+i-1);
B(p+i-1) = temp;
temp = A(x);
A(x) = B(y);
B(y) = temp;
end
% 变异操作
p1 = floor(1 + N * rand());
p2 = floor(1 + N * rand());
while p1 == p2
p1 = floor(1 + N * rand());
p2 = floor(1 + N * rand());
end
tmp = A(p1);
A(p1) = A(p2);
A(p2) = tmp;
tmp = B(p1);
B(p1) = B(p2);
B(p2) = tmp;
F = [F; A; B];
[aa, bb] = size(F);
end
if aa > NP
F = F(1:NP,:); %保持种群规模为n
end
f = F; % 更新种群
f(1,:) = R; % 保留每代最优个体
clear F;
gen = gen + 1;
% 记录最优路径长度
Rlength(gen) = minlen;
% 绘制当前的路径,并添加迭代次数
figure(1);
clf;
hold on;
for i = 1:N-1
plot([C(R(i),1), C(R(i+1),1)], [C(R(i),2), C(R(i+1),2)], 'bo-');
end
plot([C(R(N),1), C(R(1),1)], [C(R(N),2), C(R(1),2)], 'ro-');
title(['优化最短距离: ', num2str(minlen)]);
% 添加当前迭代次数到图形上
text(0.5, 0.05, ['Iteration: ', num2str(gen)], 'Units', 'normalized', 'HorizontalAlignment', 'center', 'FontSize', 12);
% 保存每个迭代的图像帧(仅当当前图与上次不同)
frame = getframe(gcf);
if ~isequal(last_frame, frame)
writeVideo(video, frame);
last_frame = frame; % 更新上次帧
end
end
% 关闭视频文件
close(video);
% 绘制适应度进化曲线
figure(2);
plot(Rlength);
xlabel('迭代次数');
ylabel('目标函数值');
title('适应度进化曲线');
TSP_GA_Animation-1

一个在2000次迭代后优化到1700以下的实例
代码改进思路
通过引入更复杂的交叉和变异操作,改进适应度函数,并增强种群多样性,可以大幅提高遗传算法解决TSP问题的效率和质量。此外,通过并行计算加速路径计算,可以显著减少算法的运行时间,尤其是在大规模问题中。
交叉和变异操作的改进: 当前的交叉操作是随机选择一段基因进行交换,可以考虑使用更复杂的交叉策略,如部分映射交叉(PMX)或顺序交叉(OX)。 变异操作可以增加更多的变异策略,如逆转变异或插入变异。
- 部分映射交叉(PMX):这种方法通过在父代中选择一个子串,保留父代的一部分,并且调整剩余部分以保持合法的城市顺序。
- 顺序交叉(OX):在这个方法中,首先选择父代中的一部分基因,然后在另一个父代中填充剩余部分,保持原有顺序。
% 逆转变异 function mutated = reverse_mutation(route) N = length(route); p1 = randi([1, N]); p2 = randi([p1, N]); mutated = route; mutated(p1:p2) = route(p2:-1:p1); end % 插入变异 function mutated = insertion_mutation(route) N = length(route); p1 = randi([1, N]); p2 = randi([1, N]); while p1 == p2 p2 = randi([1, N]); end mutated = route; gene = mutated(p1); mutated(p1) = mutated(p2); mutated(p2) = gene; end
适应度函数的改进: 当前的适应度函数是简单的归一化处理,可以考虑使用指数函数或其他非线性函数来增强选择压力。这样做能够使得适应度较高的个体拥有更大的选择概率。
- 计算归一化适应值,使用指数函数增强选择压力
% 计算归一化适应值,使用指数函数增强选择压力 function fitness = compute_fitness(len, minlen, maxlen) fitness = exp(-((len - minlen) / (maxlen - minlen + 0.001))); end
种群多样性维护: 在进化过程中,种群可能会过早收敛。可以通过引入多样性维护机制,如小生境技术或移民策略,来保持种群的多样性。 并行计算: 遗传算法的种群评估可以并行化,以加速计算过程。
- 小生境技术:通过将种群分成多个子种群,每个子种群内部进行进化,定期将优秀个体从一个子种群迁移到另一个子种群,保持种群多样性。
- 移民策略:每隔一段时间将一部分个体从外部加入到种群中,防止种群过早收敛。
% 种群多样性维护:移民策略 function new_population = immigration_strategy(population, best_individuals) immigration_rate = 0.1; % 迁移率 num_immigrants = round(immigration_rate * size(population, 1)); new_population = population; for i = 1:num_immigrants immigrant = best_individuals(randi([1, length(best_individuals)]), :); new_population(randi([1, size(new_population, 1)]), :) = immigrant; end end并行计算:加速种群评估
可以使用 MATLAB 的并行计算工具箱来加速路径计算,尤其是当种群较大时。你可以将适应度评估部分并行化。
% 并行计算种群适应度 parfor i = 1:NP len(i, 1) = D(f(i, N), f(i, 1)); for j = 1:(N - 1) len(i, 1) = len(i, 1) + D(f(i, j), f(i, j + 1)); end end整合后代码
clear all; % 清除所有变量 close all; % 清图 clc; % 清屏 C = [1304 2312; 3639 1315; 4177 2244; 3712 1399; 3488 1535; 3326 1556; ... 3238 1229; 4196 1044; 4312 790; 4386 570; 3007 1970; 2562 1756; ... 2788 1491; 2381 1676; 1332 695; 3715 1678; 3918 2179; 4061 2370; ... 3780 2212; 3676 2578; 4029 2838; 4263 2931; 3429 1908; 3507 2376; ... 3394 2643; 3439 3201; 2935 3240; 3140 3550; 2545 2357; 2778 2826; ... 2370 2975]; % 31个省会城市坐标 N = size(C,1); % TSP问题的规模, 即城市数目 D = zeros(N); % 任意两个城市距离间隔矩阵 % 计算任意两个城市的距离 for i = 1:N for j = 1:N D(i,j) = sqrt((C(i,1) - C(j,1))^2 + (C(i,2) - C(j,2))^2); % 使用欧几里得距离计算城市间距离 end end NP = 200; % 种群规模 G = 2500; % 最大遗传代数 f = zeros(NP, N); % 用于存储种群 F = []; % 种群更新中间存储 for i = 1:NP f(i,:) = randperm(N); % 随机生成初始种群 end R = f(1,:); % 存储最优种群 len = zeros(NP, 1); % 存储路径长度 fitness = zeros(NP, 1); % 存储归一化适应值 gen = 0; % 为动画预设一个文件 video = VideoWriter('TSP_GA_Animation', 'MPEG-4'); open(video); % 上一次保存的路径图 last_frame = []; % 遗传算法循环 while gen < G % 计算路径长度 parfor i = 1:NP len(i,1) = D(f(i,N), f(i,1)); % 计算从末尾城市到起始城市的路径长度 for j = 1:(N-1) len(i,1) = len(i,1) + D(f(i,j), f(i,j+1)); % 计算路径中其他城市之间的总距离 end end maxlen = max(len); % 最长路径 minlen = min(len); % 最短路径 % 更新最短路径 rr = find(len == minlen); % 找到最短路径 R = f(rr(1,1),:); % 存储最优路径 % 计算归一化适应值,使用指数函数增强选择压力 fitness = exp(-((len - minlen) / (maxlen - minlen + 0.001))); % 将路径长度映射为适应度 % 选择操作 nn = 0; for i = 1:NP if fitness(i,1) >= rand % 适应度大于随机值则选择 nn = nn + 1; F(nn,:) = f(i,:); % 将选择的个体加入中间种群 end end [aa, bb] = size(F); while aa < NP % 保证新种群的大小为NP nnper = randperm(nn); % 随机排列 A = F(nnper(1), :); B = F(nnper(2), :); % 部分映射交叉 (PMX) [A, B] = PMX_crossover(A, B); % 交叉操作 % 变异操作:逆转变异 A = reverse_mutation(A); % 对个体A进行逆转变异 B = reverse_mutation(B); % 对个体B进行逆转变异 % 移民策略(增强种群多样性) best_individuals = f(1:round(0.1*NP), :); % 选择前10%最优个体 F = immigration_strategy(F, best_individuals); % 增加种群的多样性 F = [F; A; B]; % 将交叉后的两个个体加入到新种群中 [aa, bb] = size(F); end if aa > NP F = F(1:NP,:); % 保持种群规模为NP end f = F; % 更新种群 f(1,:) = R; % 保留每代最优个体 clear F; gen = gen + 1; % 记录最优路径长度 Rlength(gen) = minlen; % 绘制当前的路径,并添加迭代次数 figure(1); clf; hold on; for i = 1:N-1 plot([C(R(i),1), C(R(i+1),1)], [C(R(i),2), C(R(i+1),2)], 'bo-'); % 绘制路径 end plot([C(R(N),1), C(R(1),1)], [C(R(N),2), C(R(1),2)], 'ro-'); % 连接最后一个城市与第一个城市 title(['优化最短距离: ', num2str(minlen)]); % 显示当前的最短路径长度 % 添加当前迭代次数到图形上 text(0.5, 0.05, ['Iteration: ', num2str(gen)], 'Units', 'normalized', 'HorizontalAlignment', 'center', 'FontSize', 12); % 保存每个迭代的图像帧(仅当当前图与上次不同) frame = getframe(gcf); if ~isequal(last_frame, frame) writeVideo(video, frame); % 写入视频文件 last_frame = frame; % 更新上次帧 end end % 关闭视频文件 close(video); % 绘制适应度进化曲线 figure(2); plot(Rlength); xlabel('迭代次数'); ylabel('目标函数值'); title('适应度进化曲线'); % 部分映射交叉 (PMX) function [offspring1, offspring2] = PMX_crossover(parent1, parent2) N = length(parent1); point1 = randi([1, N-1]); % 随机选择交叉点 point2 = randi([point1+1, N]); offspring1 = parent1; offspring2 = parent2; % 交换选定区间的部分 for i = point1:point2 temp = offspring1(i); offspring1(i) = offspring2(i); offspring2(i) = temp; end % 修复交叉后的冲突 for i = 1:N if ismember(offspring1(i), offspring1(point1:point2)) for j = 1:N if ~ismember(parent2(j), offspring1) offspring1(i) = parent2(j); break; end end end end for i = 1:N if ismember(offspring2(i), offspring2(point1:point2)) for j = 1:N if ~ismember(parent1(j), offspring2) offspring2(i) = parent1(j); break; end end end end end % 逆转变异 function mutated = reverse_mutation(route) N = length(route); p1 = randi([1, N]); % 随机选择变异区间 p2 = randi([p1, N]); mutated = route; mutated(p1:p2) = route(p2:-1:p1); % 将选择的部分反转 end % 移民策略(增强种群多样性) function new_population = immigration_strategy(population, best_individuals) immigration_rate = 0.1; % 迁移率 num_immigrants = round(immigration_rate * size(population, 1)); % 计算迁移个体数量 % 防止超出索引 num_best = size(best_individuals, 1); % best_individuals 的个体数 new_population = population; for i = 1:num_immigrants immigrant_idx = randi([1, min(num_best, length(best_individuals))]); immigrant = best_individuals(immigrant_idx, :); % 选择迁移的个体 new_population(randi([1, size(new_population, 1)]), :) = immigrant; % 将其加入新种群中 end end

- 一些在600次以内优化到17000一下的实例
TSP_GA_after_Animation-2
结论:遗传算法解决TSP问题难免陷入局部最优解的境地,往往与理想中的精确解有很大一段距离。
魔乐社区(Modelers.cn) 是一个中立、公益的人工智能社区,提供人工智能工具、模型、数据的托管、展示与应用协同服务,为人工智能开发及爱好者搭建开放的学习交流平台。社区通过理事会方式运作,由全产业链共同建设、共同运营、共同享有,推动国产AI生态繁荣发展。
更多推荐


所有评论(0)