愿你搜索到人生的最优解
文案：挽月
排版：随心390
hello，大家好。各位可点击左下方阅读原文，访问公众号官方店铺。谨防上当受骗，感谢各位支持！
在上周的运筹学课中，听老师吹了一波差分进化算法，小编阅读了很多文献发现，差分进化算法在求解连续型函数问题上有了许多应用，但是在离散优化上鲜有资料。于是小编动手实现了基于差分进化算法的TSP求解，效果与遗传算法一样但是速度显著提高，代码量显著减少。
本期推文目录：
01 差分进化算法基本思路
02 差分进化算法求解TSP问题MATLAB代码
01 | 差分进化算法基本思路
首先来看看差分进化算法背景，差分进化算法首先由K.V. Price于1997年提出，算法的背景基于随机游走算法。随机游走算法基础思想如下：当找到一个函数值优于原来函数值时则会替代原有函数值，否则不变。
02 | 差分进化算法求解TSP问题MATLAB代码
MAIN函数
%% Main函数
F0=0.4; %变异因子
CR=0.1; %交叉概率
MaxGens=1000;
x_high=500;
x_low=-500;
X =[16.47,96.10
16.47,94.44
20.09,92.54
22.39,93.37
25.23,97.24
22.00,96.05
20.47,97.02
17.20,96.29
16.30,97.38
14.05,98.12
16.53,97.38
21.52,95.59
19.41,97.13
20.09,92.55];
DM=Distanse(X);
D=length(X);
NP=8*D;
pop=rand(NP,D)*(x_high-x_low)+x_low; %初始化种群
g=1; %迭代标记
fit=PathLength(DM,pop);
trace(1)=min(fit);
v=zeros(NP,D); %竞争矩阵
for gen=1:MaxGens
% 变异操作
for i=1:NP
r1=randi([1,NP],1,1);
pop_s = fit(r1);
while(pop_s<fit(i))
r1 =randi([1,NP],1,1);
pop_s = fit(r1);
end
%产生r2,r3
r2=randi([1,NP],1,1);
while(r2==r1)||(r2==i)
r2=randi([1,NP],1,1);
end
r3=randi([1,NP],1,1);
while(r3==i)||(r3==r2)||(r3==r1)
r3=randi([1,NP],1,1);
end
%交叉
Z=randi([1,D],1,1);
r=rand;
for j = 1:D
if (r<=CR) || (j==Z)
v(i,j)=(pop(r1,j)+pop(i,j))/2+F0*(pop(r1,j)-pop(i,j)+pop(r2,j)-pop(r3,j));
else
v(i,j)=pop(i,j);
end
end
if PathLength(DM,v(i,:))<PathLength(DM,pop(i,:));
pop(i,:)=v(i,:);
end
end
fit=PathLength(DM,pop);
[trace(gen+1),minId]=min(fit);
tt=min(fit);
end
[~,pop(minId,:)]=sort(pop(minId,:),2,'ascend');
DrawPath(pop(minId,:),X)
figure(2);
title(['差分进化算法(DE)', '最小值: ', num2str(tt)]);
xlabel('迭代次数');
ylabel('目标函数值');
plot(trace);
PathLength.m函数
%% PathLength.m
%% 计算各个体的路径长度
% 输入：
% D 两两城市之间的距离
% Chrom 个体的轨迹
function len=PathLength(D,Chrom)
[B,Chrom]=sort(Chrom,2,'ascend');
[row,col]=size(D);
NIND=size(Chrom,1);
len=zeros(NIND,1);
for i=1:NIND
p=[Chrom(i,:) Chrom(i,1)];
i1=p(1:end-1);
i2=p(2:end);
len(i,1)=sum(D((i1-1)*col+i2));
end
Distanse.m函数
%% Distanse.m
%% 计算两两城市之间的距离
%输入 a 各城市的位置坐标
%输出 D 两两城市之间的距离
function D=Distanse(a)
row=size(a,1);
D=zeros(row,row);
for i=1:row
for j=i+1:row
D(i,j)=((a(i,1)-a(j,1))^2+(a(i,2)-a(j,2))^2)^0.5;
D(j,i)=D(i,j);
end
end
dsxy2figxy.m函数
%% dsxy2figxy.m
function varargout = dsxy2figxy(varargin)
if length(varargin{1}) == 1 && ishandle(varargin{1}) ...
&& strcmp(get(varargin{1},'type'),'axes')
hAx = varargin{1};
varargin = varargin(2:end);
else
hAx = gca;
end;
if length(varargin) == 1
pos = varargin{1};
else
[x,y] = deal(varargin{:});
end
axun = get(hAx,'Units');
set(hAx,'Units','normalized');
axpos = get(hAx,'Position');
axlim = axis(hAx);
axwidth = diff(axlim(1:2));
axheight = diff(axlim(3:4));
if exist('x','var')
varargout{1} = (x - axlim(1)) * axpos(3) / axwidth + axpos(1);
varargout{2} = (y - axlim(3)) * axpos(4) / axheight + axpos(2);
else
pos(1) = (pos(1) - axlim(1)) / axwidth * axpos(3) + axpos(1);
pos(2) = (pos(2) - axlim(3)) / axheight * axpos(4) + axpos(2);
pos(3) = pos(3) * axpos(3) / axwidth;
pos(4) = pos(4) * axpos(4 )/ axheight;
varargout{1} = pos;
end
set(hAx,'Units',axun)
DrawPath.m函数
%% DrawPath.m
%% 画路径函数
%输入
% Chrom 待画路径
% X 各城市坐标位置
function DrawPath(Chrom,X)
R=[Chrom(1,:) Chrom(1,1)];
figure;
hold on
plot(X(:,1),X(:,2),'o','color',[0.5,0.5,0.5])
plot(X(Chrom(1,1),1),X(Chrom(1,1),2),'rv','MarkerSize',20)
for i=1:size(X,1)
text(X(i,1)+0.05,X(i,2)+0.05,num2str(i),'color',[1,0,0]);
end
A=X(R,:);
row=size(A,1);
for i=2:row
[arrowx,arrowy] = dsxy2figxy(gca,A(i-1:i,1),A(i-1:i,2)); %坐标转换
annotation('textarrow',arrowx,arrowy,'HeadWidth',8,'color',[0,0,1]);
end
hold off
xlabel('横坐标')
ylabel('纵坐标')
title('轨迹图')
box on
求解结果如下图所示：
各位小伙伴可在留言板留言，未来我们讲解的具体内容由你做主。如果可以的话，可以把希望讲解的文献也在留言板上写出来。
