关于旅行商问题的求助

最近有一个作业是用遗传算法求解旅行商问题,在网上找了博主的代码,但结果方面还是存在不小的误差,希望有大神能解答一下

主函数:

clc
clear
 load("berlin52.mat"); % 加在数据
 location = load("berlin52.mat"); % 开始对数据进行可视化
x = berlin52(:,1);
y = berlin52(:,2);
plot(x,y,'bo','MarkerSize', 4);
xlabel('城市横坐标');
ylabel('城市纵坐标');
grid on;
 
NIND = 300;         %种群大小
MAXGEN = 1200;       %最大迭代次数
Pc = 0.9;           %交叉概率,相当于基因遗传的时候染色体交叉
Pm = 0.05;          %染色体变异
GGAP = 0.9;         %这个是代沟,通过遗传方式得到的子代数为父代数*GGAP
D = Distance(berlin52);    %通过这个函数可以计算i,j两点之间的距离
N = size(D,1);      %计算有多少个坐标点

Chrom = InitPop(NIND,N);   %%初始化种群

gen = 0;
trace = zeros(1,MAXGEN);
ObjV = PathLength(D,Chrom);%计算当前路线长度,即上面随机产生的那些个体路径
preObjV = min(ObjV);%当前最优解
while gen<MAXGEN
    %%计算适应度
    ObjV = PathLength(D,Chrom);     %计算路线长度
    pause(0.0001);
    preObjV = min(ObjV);
    trace(1,gen+1)=preObjV;
    FitnV = Fitness(ObjV);
    %选择
    SelCh = Select(Chrom,FitnV,GGAP);
    SelCh = Reverse(SelCh,D);
    %交叉操作
    SelCh = Recombin(SelCh,Pc);
    SelCh = Reverse(SelCh,D);
    %变异
    SelCh = Mutate(SelCh,D,Pm);
    Chrom = Reins(Chrom,SelCh,ObjV);
    gen = gen + 1;
    SelCh = Reverse(SelCh,D);
end

%画出最优解的路线图
ObjV = PathLength(D,Chrom);     %计算路线长度
[minObjV,minInd] = min(ObjV);
DrawPath(Chrom(minInd(1),:),berlin52)
figure();
plot([1:MAXGEN],trace(1,:));
title('优化过程');
xlabel('迭代次数');
ylabel('当前最优解');
%%输出最优解的距离
disp(['旅行商走过的总距离:',num2str(ObjV(minInd(1)))]);
disp('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~')

Distance函数:

function D = Distance(x)
row = size(x,1);
D = zeros(row,row);
for i = 1:row
    for j = i+1:row
        D(i,j) = ((x(i,1)-x(j,1))^2+(x(i,2)-x(j,2))^2)^0.5;
        D(j,i) = D(i,j);
    end
end

DrawPath函数:

function DrawPath(Chrom,X)
 
R =  [Chrom(1,:) Chrom(1,1)]; %第一个随机解(个体)【Chrom(1,:)取第一行数据】,一共有n个城市,但是这里R有n+1个值,因为后面又补了一个Chrom(1,1),“是为了让路径最后再回到起点”
figure;
hold on
plot(X(:,1),X(:,2),'bo','MarkerSize', 4);%X(:,1),X(:,2)分别代表的X轴坐标和Y轴坐标
plot(X(Chrom(1,1),1),X(Chrom(1,1),2),'rv','MarkerSize',9)%标记起点
A = X(R,:);                         %A是将之前的坐标顺序用R打乱后重新存入A中
row = size(A,1);                    %row为坐标数+1
for i = 2:row
    % 直接使用 plot 函数连接两个点
        plot(A(i-1:i,1), A(i-1:i,2), 'g-', 'LineWidth', 1);
    end
hold off
xlabel('横坐标x')
ylabel('纵坐标y')
title('旅行商轨迹图')
box on

初始化函数:

function Chrom = InitPop(NIND,N)
Chrom = zeros(NIND,N); % 定义存储种群的变量
for i = 1:NIND
    Chrom(i,:) = randperm(N);
end

Fitness函数

    function FitnV = Fitness(len) %适应度函数
    %输入:
    %len 个体的长度(TSP的距离)
    %输出:
    %FitnV 个体的适应度值
    FitnV = 1./len;

Intercross:

function [a,b] = intercross(a,b)
%输入:
%a和b为两个待交叉的个体
%输出:
%a和b为交叉后得到的两个个体
L = length(a);
%随机产生交叉区段
r1 = randsrc(1,1,[1:L]); % 这里先随机选出两个位置,位置不能超过
r2 = randsrc(1,1,[1:L]); % 
if r1~=r2
    a0 = a;
    b0 = b;
    s = min([r1,r2]);
    e = max([r1,r2]);
    for i = s:e
        a1 = a;
        b1 = b;
        %第一次互换
        a(i) = b0(i);
        b(i) = a0(i);
        %寻找相同的城市
        x = find(a==a(i)); % 如果有相同的城市,x会得到包含i的两个值,y同理
        y = find(b==b(i));
        %第二次互换产生新的解
        i1 = x(x~=i);      % 将位置不是i但相同的城市标记出来
        i2 = y(y~=i);
        if ~isempty(i1)
            a(i1)=a1(i);
        end
        if ~isempty(i2)
            b(i2)=b1(i);   % 注意,原文这里有误,应该是b(i2)
        end
    end
end
% 这里增加一段代码,r1=r2时,两个个体只在一点交叉
if r1 == r2
    a0 = a;
    b0 = b;
    a(r1) = b0(r1);
    b(r1) = a0(r1);
    x = find(a==a(r1));
    y = find(b==b(r1));
    i1 = x(x~=r1);
    i2 = y(y~=r1);
    if ~isempty(i1)
        a(i1) = a0(r1);
    end
    if ~isempty(i2)
        b(i2) = b0(r1);
    end
end

mutate

function SelCh = Mutate(SelCh,D,Pm)
%变异操作
%SelCh  被选择的个体
%Pm  变异概率
%输出:
%SelCh  变异后的个体
index = SelCh;
col = size(SelCh,2); %返回SelCh的列数
lenSelCh = PathLength(D,SelCh);
[NSel,L] = size(SelCh);
for i = 1:NSel
    if Pm >= rand
        R = randperm(L);
        index(i,R(1:2)) = index(i,R(2:-1:1)); % 将个体i中R(1)和R(2)这两个位置的城市互换
        p = [index(i,:) index(i,1)];
        i1 = p(1:end-1);
        i2 = p(2:end);
        DIndexi = sum(D((i1-1)*col+i2));              % 计算出变异后个体的路径距离
        if DIndexi < lenSelCh(i)                 % 如果变异后路径距离比变异前更小,则保留变异
            SelCh(i,:) = index(i,:);
        end
    end
end

Pathlength

function len = PathLength(D,Chrom)
%计算所有个体的路线长度
%输入
%D两两城市之间的距离
%Chrom个体的轨迹
 [~,col] = size(D); %返回D的列数 
NIND = size(Chrom,1);%NIND等于初始种群
len = zeros(NIND,1);%初始化一个大小等于NIND的len来记录长度
for i = 1:NIND
    p = [Chrom(i,:) Chrom(i,1)];%构造p矩阵保存路线图 将第一行路线提出 再加上第一个构成回路
    i1 = p(1:end-1);%i1从第一个开始遍历到倒数第二个
    i2 = p(2:end);%i2从第二个开始遍历到倒数第一个
    len(i,1) = sum(D((i1-1)*col+i2));%计算出每种路线(种群的个体)的长度,这里相当于把D矩阵沿行展开
end

recombin:

function SelCh = Recombin(SelCh,Pc)
%交叉操作
%输入:
%SelCh 被选择的个体
%Pc    交叉概率
%输出:
%SelCh 交叉后的个体
 
NSel = size(SelCh,1);
for i = 1:2:NSel - mod(NSel,2) % 若不减去余数且NSel如果是奇数,最后一个i会等于NSel,i+1报错
    if Pc>=rand %交叉概率PC
        [SelCh(i,:),SelCh(i+1,:)] = intercross(SelCh(i,:),SelCh(i+1,:));
    end
end

Reins

function Chrom = Reins(Chrom,SelCh,ObjV)
%重插入子代的种群
%输入:
%Chrom      父代的种群
%SelCh      子代的种群
%ObjV       父代适应度
%输出:
%Chrom      组合父代与子代后得到的新种群
NIND = size(Chrom,1);
NSel = size(SelCh,1);
[TobjV,index] = sort(ObjV);
Chrom =  [Chrom(index(1:NIND-NSel),:);SelCh];

Reverse:

function SelCh = Reverse(SelCh,D)
%进化逆转函数
%输入:
%SelCh  被选择的个体
%D  各城市的距离矩阵
%输出:
%SelCh  进化逆转后的个体
 
[row,col] = size(SelCh);
ObjV = PathLength(D,SelCh);
SelCh1 = SelCh;
for i = 1:row
    r1 = randsrc(1,1,[1:col]);
    r2 = randsrc(1,1,[1:col]);
    mininverse = min([r1 r2]);
    maxinverse = max([r1 r2]);
    SelCh1(i,mininverse:maxinverse) = SelCh1(i,maxinverse:-1:mininverse);
end
ObjV1 = PathLength(D,SelCh1);%计算路线长度
index = ObjV1<ObjV;
SelCh(index,:)=SelCh1(index,:);

select:

function SelCh = Select(Chrom,FitnV,GGAP)
%输入:
%Chrom 种群
%FitnV 适应度值
%GGAP 选择概率
%输出:
%SelCh 被选择的个体
NIND = size(Chrom,1);%种群个体总数
NSel = max(floor(NIND * GGAP+.5),2);%确定下一代种群的个体数,如果不足二个就计为二个
ChrIx = Sus(FitnV,NSel);
SelCh = Chrom(ChrIx,:);

Sus

function NewChrIx = Sus(FitnV,Nsel)
%输入:
%FitnV 个体的是适应度值
%Nsel 被选个体的数目
%输出:
%NewChrIx 被选择个体的索引号
[Nind,ans] = size(FitnV);%Nind为FitnV的行数也就是个体数 ans为列数1
cumfit = cumsum(FitnV);%对适应度累加 例如 1 2 3 累加后 1 3 6 用来计算累积概率
trials = cumfit(Nind)/Nsel * (rand + (0:Nsel-1)');%cumfit(Nind)表示的是矩阵cumfit的第Nind个元素 A.'是一般转置,A'是共轭转置 rand返回一个在区间 (0,1) 内均匀分布的随机数
%cumfit(Nind)/Nsel平均适应度 * (rand +(0:Nsel-1)')从0开始到89的转置矩阵(行矩阵变列矩阵)加上每一项加上一个0-1的随机数
%生成随机数矩阵 用来筛选
Mf = cumfit(:,ones(1,Nsel));%将生成的累积概率 
Mt = trials(:,ones(1,Nind))';
[NewChrIx,ans] = find(Mt<Mf & [zeros(1,Nsel);Mf(1:Nind-1,:)]<= Mt);%寻找满足条件的个体 返回返回数组 X 中每个元素的行和列下标
[ans,shuf] = sort(rand(Nsel,1));%生成Nsel*1的随机数矩阵  按升序对 A 的元素进行排序 返回选择的shuf 随机打乱个体顺序 保证后续遗传算子操作的随机性
NewChrIx = NewChrIx(shuf);%返回shuf索引的矩阵元素

数据集为berlin52,

NAME: berlin52
TYPE: TSP
COMMENT: 52 locations in Berlin (Groetschel)
DIMENSION: 52
EDGE_WEIGHT_TYPE: EUC_2D
NODE_COORD_SECTION
1 565.0 575.0
2 25.0 185.0
3 345.0 750.0
4 945.0 685.0
5 845.0 655.0
6 880.0 660.0
7 25.0 230.0
8 525.0 1000.0
9 580.0 1175.0
10 650.0 1130.0
11 1605.0 620.0 
12 1220.0 580.0
13 1465.0 200.0
14 1530.0 5.0
15 845.0 680.0
16 725.0 370.0
17 145.0 665.0
18 415.0 635.0
19 510.0 875.0  
20 560.0 365.0
21 300.0 465.0
22 520.0 585.0
23 480.0 415.0
24 835.0 625.0
25 975.0 580.0
26 1215.0 245.0
27 1320.0 315.0
28 1250.0 400.0
29 660.0 180.0
30 410.0 250.0
31 420.0 555.0
32 575.0 665.0
33 1150.0 1160.0
34 700.0 580.0
35 685.0 595.0
36 685.0 610.0
37 770.0 610.0
38 795.0 645.0
39 720.0 635.0
40 760.0 650.0
41 475.0 960.0
42 95.0 260.0
43 875.0 920.0
44 700.0 500.0
45 555.0 815.0
46 830.0 485.0
47 1170.0 65.0
48 830.0 610.0
49 605.0 625.0
50 595.0 360.0
51 1340.0 725.0
52 1740.0 245.0
EOF

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值