自重启伪遗传改良算法解决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                                                                   % 输出程序运行的时间