自重启伪遗传改良算法解决TSP问题(matlab实现)
说明:解决旅行商问题。所有代码和数据上传到github,可自行下载和修改。
地址:代码和数据
1.实验数据说明
实验数据为ch130.mat,格式如下:
序号 | x坐标 | y坐标 |
---|---|---|
1 | 100.312 | 140.213 |
2 | 122.343 | 124.343 |
3 | 563.378 | 123.676 |
... | ... | ... |
表示130个点的序号,x坐标和y坐标。
2.算法概述
本算法的灵感来自改良圈算法,可参考该博客中对改良圈的描述。改良圈算法描述
改良圈算法运用了一个巧妙的思路,将初始随机路径快速地改良为一个较为良好的路径,然后经过遗传算法,继续收敛。但是我们在实验中发现,由于改良圈算法得到的路径已经很优秀,遗传算法很难收敛下去。于是有了将改良圈算法与遗传算法结合的想法。之所以在算法的名称中加入“自启”和“伪”两个词,是因为针对传统遗传算法做了修改?!白云簟敝杆惴ㄔ诳赡苁樟驳骄植孔钚≈?,无法再收敛下去的时候,自行重启整个流程。“伪”指这个算法中将交叉率和变异率都设为了1,且做了较大的改动。
算法流程图:
3.算法优势
快。
4.实验结果
ch130在网站上给出的最优解为6110.8,目前的程序能在较短时间内快速收敛到6110.7。用其他数据集测试的时候效果也比较良好。
代码附录:
(1)改良圈
function [ A ] = circle_modification( initial_population, rows, L, distance_matrix )
% 改良圈改良一个初始种群initial_population,返回改良后的种群A
A = zeros([rows, L]);
for i=1:rows
c = initial_population(i, :);
flag = 1;
while flag > 0
flag = 0;
for m = 1:L-3
for n = m+2:L-1
if distance_matrix(c(m), c(n))+distance_matrix(c(m+1), c(n+1)) < distance_matrix(c(m), c(m+1))+distance_matrix(c(n), c(n+1))
flag=1;
c(m+1:n) = c(n:-1:m+1);
end
end
end
end
A(i,c) = 1:L;
end
end
(2)归一化
function [ A ] = normalization( A, L )
% 对A进行归一化
A = A / L;
A(:, 1) = 0;
A(:, L) = 1;
end
(3)产生初始种群
function [ initial_population ] = generate_population( n, L )
% 产生指定数量n的种群
initial_population = zeros([n, L]);
for i=1:n
c = randperm(L - 2); % 生成1---(L-2)的乱序排列, 本题中为1-129
current_path = [1, c+1, L]; % 生成真正的个体,最后一个数为c的第一个数,c1范围为1-131
initial_population(i, :) = current_path;
end
end
(4)获取距离矩阵
function [ distance_matrix ] = get_distance_matrix( point_info )
% 获取距离矩阵
number_of_point = length(point_info);
D = pdist(point_info,'euclidean'); % 每两点之间的欧式距离,矩阵大小为(1,n)
distance_matrix = zeros(number_of_point + 1); % 每次都要回到起点,因此,最后一列为到起点的距离
sum = 0;
for i = 1:number_of_point
if i < number_of_point
sum = sum + number_of_point - i;
distance_matrix(i, i+1:number_of_point) = D(1, sum-(number_of_point-i-1):sum); % 将已知的D距离改为矩阵形式
end
end
distance_matrix = distance_matrix + distance_matrix';
distance_matrix(1:number_of_point, number_of_point+1) = distance_matrix(1:number_of_point, 1); % 对最后一列进行填值
distance_matrix(number_of_point+1, :) = distance_matrix(1, :); % 对最后一行进行填值,完整的距离矩阵
end
(5)变异
function [ C ] = mutation( A, w, L, distance_matrix )
% A变异产生子代C
C = A;
% 对变异行的元素重新排列
for i=1:w
mutation_location = 2 + floor((L-2)*rand(1,3));
mutation_location = sort(mutation_location);
C(i,:)=C(i,[1:mutation_location(1)-1, mutation_location(2)+1:mutation_location(3), mutation_location(1):mutation_location(2), mutation_location(3)+1:L]); % 对该个体重新排列组合
end
[~, mutation_path_table] = sort(C, 2); % 对每一行的元素进行排序,获取所有的路径
C = normalization(circle_modification(mutation_path_table, w, L, distance_matrix), L);
end
(6)交叉
function [ B ] = cross( A, w, L, distance_matrix )
% A交叉产生子代B
B = A(randperm(w), :);
for i = 1:2:w
cross_position = 2 + floor((L-2) * rand(1, 2)); % 2<F<L-1,为整数[2, L-1],首尾的点不去改变
sort(cross_position);
cp1 = cross_position(1);
cp2 = cross_position(2);
temp = B(i, cp1:cp2); % 交换第i个和第i+1个子代的F:L列
B(i, cp1:cp2) = B(i+1, cp1:cp2);
B(i+1, cp1:cp2) = temp;
end % 交换完毕,产生w个新个体
[~, cross_path_table] = sort(B, 2); % 对每一行的元素进行排序,获取所有的路径
B = normalization(circle_modification(cross_path_table, w, L, distance_matrix), L);
end
(7)选取子代
function [ A_next, optimal_path, optimal_path_length ] = select_next_generation( A, B, C, w, L, distance_matrix )
% 选择遗传的下一代
fp = 0.3; % 外来人口比例
np = 0.2; % 非最优比例
G = [A;B;C]; % A为原始的,B为交叉后的,C为变异后的
total_length = size(G, 1); % 父代和子代的总长度
%在父代和子代中选择优良品种作为新的父代
[~, path_table] = sort(G,2); % 对每一行的元素进行排序,获取所有的路径
all_path_length = [];
all_path_length(1:total_length) = 0; % 初始化所有路径的长度为0
for path_index = 1:total_length % 计算每条路径的长度
for i = 1:L-1
all_path_length(path_index) = all_path_length(path_index) + distance_matrix(path_table(path_index, i), path_table(path_index, i+1));
end
end
[sorted_length, sorted_index] = sort(all_path_length);
A_next = G(sorted_index(1:w), :); %选出距离最短的w个子代
optimal_path = path_table(sorted_index(1), :); % 此时的最短路径为第1个个体
optimal_path_length = sorted_length(1); % Dz为该最短路径的实际距离
% 人人能娶老婆
up_n = w*(1-np-fp);
for i = up_n+1:w*(1-fp)
r = up_n + unidrnd(w - up_n);
A_next(i) = G(r);
end
% 加入外来人口,但是前提是要对外来人口进行优化
A_next(int16(w*(1-fp))+1:w, :) = normalization(circle_modification(generate_population(int16(w*fp), L), int16(w*fp), L, distance_matrix), L);
(8)获取路径长度
function [ path_length ] = get_path_length( path, L, distance_matrix )
% 获取路径长度
path_length = 0;
for i = 1:L-1
path_length = path_length + distance_matrix(path(i), path(i+1));
end
end
(9)画出路径
function [ ] = plot_path( point_position_x_and_y, path, length )
% ?????·???·
x_=point_position_x_and_y(path,1);
y_=point_position_x_and_y(path,2);
plot(x_,y_,'-o');
title(length);
drawnow;
end
(10)主函数
tic
clc,clear
rng('shuffle'); %改变随机数的初始状态
% -----------------参数------------------
w = 500; % 种群规模
restart_times = 200; % 重启次数
iterations = 300; % 迭代次数
repeat_time_threshold = 100; % 重复次数阈值
% ---------------------------------------
% 从文件中读取信息
load ch130.mat % 载入数据集
point_info = ch130(:, 2:3);
point_position_x_and_y = [point_info; point_info(1,:)];
distance_matrix = get_distance_matrix(point_info);
L = length(ch130) + 1; % 为了保证最终能回到起点,实际的个体长度设为L,L的最后一个数和第一个数相同,保证回到起点
optimal_path = zeros([1, L]); % 记录全局最佳路径
optimal_path_length = 999999; % 记录全局最佳路径的长度
for r_index = 1:restart_times
current_optimal_path = zeros([1, L]); % 记录每次重启的最佳路径
current_optimal_path_length = 999999; % 记录每次重启的最佳路径长度
last_optimal_path_length = current_optimal_path_length; % 记录单次遗传算法中上一次最佳路径长度,用于判断是否陷入局部最优解
same_time = 0; % 单次遗传算法陷入局部最优解的次数
% 产生初始种群
initial_population = generate_population(w, L);
% 改良圈改良初始种群
A = circle_modification(initial_population, w, L, distance_matrix);
% 归一化
A = normalization(A, L);
% 以下为遗传算法实现过程
for k=1:iterations
% 交叉产生子代 B
B = cross(A, w, L, distance_matrix);
% 变异产生子代 C
C = mutation(A, w, L, distance_matrix);
% 选择下一代
[A, current_optimal_path, current_optimal_path_length] = select_next_generation(A, B, C, w, L, distance_matrix);
% 更新全局最优路径长度
if current_optimal_path_length < optimal_path_length
optimal_path_length = current_optimal_path_length;
optimal_path = current_optimal_path;
end
% 绘制最优路径图
plot_path(point_position_x_and_y, current_optimal_path, current_optimal_path_length)
% 输出每次迭代的信息
fprintf('第%0d次重启,迭代次数%04d,最优路径长度%.5f,全局最优路径长度%.5f\n' , r_index, k, current_optimal_path_length, optimal_path_length);
% 记录重复次数(重复次数过高表示陷入局部最优解)
if current_optimal_path_length == last_optimal_path_length
same_time = same_time + 1;
if same_time >= repeat_time_threshold
break;
end
else
same_time = 0; % 重置same_time
end
% 更新last_optimal_path_length
last_optimal_path_length = current_optimal_path_length;
end
end
toc % 输出程序运行的时间