《MATLAB智能优化算法》学习笔记三
一、TSP概述
见上个学习笔记。
二、算法简介
虽然VNS在搜索过程中使用若干个不同的邻域以扩大搜索范围,但是当TSP中城市数目增多时,VNS搜索范围有限这一缺点就会暴露出来。为了能够进一步扩大搜索范围,本文使用大规模邻域搜索算法(LNS)求解TSP。LNS的思想是先“破坏”解,然后将破坏后的解进行“修复”,最终获得更高质量的解。
还是以上个学习笔记为例,假设TSP中城市数目为4,初始路线为1234,现在分别使用交换操作,逆转操作和插入操作求出初始解1234对应的邻域,可以得到10种不同邻域解。因此,除去初始解1234外还有13种可能解没有通过上述3种操作获得。由此可见,VNS在使用这3种邻域操作时搜索范围的局限性。
(1)在LNS求解问题过程中,如何“破坏”解,以及如何“修复”解?
“破坏”解其实就是将若干城市从当前路线中移除。“修复”解就是将被移除的城市再重新插入到“破坏”解中去。
比如当前路线为1234,我要移除城市2、3、4,此时“破坏”后的解就是1。
“修复”的过程可以看下图:
每次插入的时候都有不同的位置可以选择。
(2)LNS的大规模是如何体现的呢?
上述“修复”解的顺序为先将3插回,然后将4插回,最后将2插回,即插回顺序为342。实际上也可以按照234、243、324、423和432这5种插回顺序“修复”解。每种插回顺序同样也能得到这24个不同的“修复解”。因为当 TSP 中城市数目为4时,一共有$A_4^4=24$种排序方式,所以采用任何一种插回顺序都能得到全部排序方式,这其实就体现了大规模的思想,即能搜索到更多不同类型的解。
综上所述,LNS求解TSP的流程如下图所示,其中$f(S)$表示解$S$的行走距离。
三、求解策略
(1)构造初始解
与VNS构造TSP初始解相同,即采用贪婪算法构造TSP的初始解。
(2)“破坏”解
假设TSP中城市数目为$N$,当前解为$route$,且要从当前解中最多移除相连接的$L_{max}$个城市,最少移除相连接的$L_{min}$个城市,则具体的移除步骤如下。
STEP1:从$L_{min}\sim L_{max}$中随机取一个整数作为要移除的城市数目。
STEP2:从当前解$route$中随机选择一个城市$visit$作为即将移除的相连接$L$个城市的参考城市。
STEP3:计算$visit$在$route$中的位置$findv$,以起始点$route(1)$为界限计算在$route$中$visit$左侧的城市数目$v_{LN}=findv-1$,同理也计算出在$route$中$visit$右侧的城市数目$v_{RN}=N-findv$。
STEP4:如果$v_{LN}\le v_{RN}$,则转至STEP5,否则转至STEP6。
STEP5:分以下3种情况计算$visit$右侧要移除城市的数目$n_R$和左侧要移除城市的数目$n_L$。
(1)如果$v_{LN}<L-1$且$v_{RN}<L-1$,则$n_R=L-1-v_{LN}+(0\sim (v_{RN}-L+1+v_{RN}))$之间的随机整数,$n_L=L-1-n_R$。
(2)如果$v_{LN}>L-1$且$v_{RN}>L-1$,则$n_R=(0\sim v_{LN})$之间的随机整数,$n_L=L-1-n_R$。
(3)如果$v_{LN}\le L-1$且$v_{RN}\ge L-1$,则$n_R=L-1-v_{LN}+(0\sim v_{LN})$之间的随机整数,$n_L=L-1-N_R$。
SREP6:分以下3种情况计算$visit$右侧要移除城市的数目$n_R$和左侧要移除城市的数目$n_L$。
(1)如果$v_{LN}<L-1$且$v_{RN}<L-1$,则$n_L=L-1-v_{RN}+(0\sim (v_{LN}-L+1+v_{RN}))$之间的随机整数,$n_R=L-1-n_L$。
(2)如果$v_{LN}>L-1$且$v_{RN}>L-1$,则$n_L=(0\sim v_{RN})$之间的随机整数,$n_R=L-1-n_L$。
(3)如果$v_{LN}\le L-1$且$v_{RN}\ge L-1$,则$n_L=L-1-v_{RN}+(0\sim v_{RN})$之间的随机整数,$n_R=L-1-N_L$。
STEP7:从$route$中提取被移除的城市集合$removed=route(findv-n_L:findv+n_R)$,并将removed中所有城市从$route$中删除,得到“破坏”后的解$S_{destroy}$。
(3)“修复”解
在讲解“修复”解的步骤前,先阐述“插入成本”这一概念。如果当前“破坏”后的解为$S_{destroy}$,那么在将removed中的一个城市插回到$S_{destroy}$中的某个插入位置后,此时“修复”后的路线总距离减去$S_{destroy}$的总距离即为该城市插入该位置的“插入成本”。
在介绍“插入成本”之后,我们再来了解一下“遗憾值”这一概念。如果$S_{destroy}$中城市的数目为$lr$,那么将在removed中的一个城市插回到$S_{destroy}$时共有$lr+1$个可能的插入位置,这$lr+1$个插入位置对应$lr+1$个“插入成本”。
接下来将这$lr+1$个“插入成本”从小到大进行排序,若排序后的“插入成本”为up_delta,则该城市插回到$S_{repair}$的“遗憾值”即为排序后排在第二位的“插入成本”减去排在第一位的“插入成本”,即up_delta(2)-up_delta(1)。
“修复”的步骤如下:
STEP1:初始为“修复后”的解$S_{repair}$赋值,即$S_{repair}=S_{destroy}$。
STEP2:如果removed非空,转至STEP3,否则转至STEP6。
STEP3:计算当前removed中的城市数目$nr$,计算将removed中各个城市插回到$S_{repair}$的“遗憾值”regret,即regret是$nr$行1列的矩阵。
STEP4:找出regret中最大“遗憾值”对应的序号max_index,确定出即将被插回的城市renisert_city=removed(max_index),将renisert_city插回到$S_{repair}$中“插入成本”最小的位置。
STEP5:更新removed(max_index)=[ ],转至STEP2。
STEP6:“修复”结束,返回“修复”后 的解$S_{repair}$。
四、MATLAB程序实现
(1)构造初始解函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
|
function [init_route,init_len]=construct_route(dist) N=size(dist,1);
for i=1:N for j=1:N if i==j dist(i,j)=inf; end end end
unvisited=1:N; visited=[];
min_dist=min(min(dist)); [row,col]=find(dist==min_dist); first=row(1);
unvisited(unvisited==first)=[]; visited=[visited,first]; pre_point=first; while ~isempty(unvisited) pre_dist=dist(pre_point,:); pre_dist(visited)=inf; [~,pre_point]=min(pre_dist); unvisited(unvisited==pre_point)=[]; visited=[visited,pre_point]; end init_route=visited; init_len=route_length(init_route,dist); end
|
(2)路线总距离计算函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14
|
function len=route_length(route,dist) N=numel(route); route=[route route(1)]; len=0; for k=1:N i=route(k); j=route(k+1); len=len+dist(i,j); end end
|
(3)“破坏”函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
|
function [Sdestroy,removed]=destroy(route) N=numel(route); Lmin=1; Lmax=min(ceil(N/2),25); visit=ceil(rand*N); L=Lmin+ceil((Lmax-Lmin)*rand);
findv=find(route==visit,1,'first'); vLN=findv-1; vRN=N-findv;
if vLN<=vRN if (vRN<L-1)&&(vLN<L-1) nR=L-1-vLN+round(rand*(vRN-L+1+vLN)); nL=L-1-nR; elseif (vRN>L-1)&&(vLN>L-1) nR=round(rand*vLN); nL=L-1-nR; else nR=L-1-vLN+round(rand*vLN); nL=L-1-nR; end else if (vLN<L-1)&&(vRN<L-1) nL=L-1-vRN+round(rand*(vLN-L+1+vRN)); nR=L-1-nL; elseif (vLN>L-1)&&(vRN>L-1) nL=round(rand*vRN); nR=L-1-nL; else nL=L-1-vRN+round(rand*vRN); nR=L-1-nL; end end removed=route(findv-nL:findv+nR);
Sdestroy=route; for i=1:L Sdestroy(Sdestroy==removed(i))=[]; end end
|
(4)“插入成本”计算函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
|
function [new_route,up_delta]=ins_route(visit,dist,route) lr=numel(route); rc0=zeros(lr+1,lr+1); delta0=zeros(lr+1,1); for i=1:lr+1 if i==lr+1 rc=[route visit]; elseif i==1 rc=[visit route]; else rc=[route(1:i-1) visit route(i:end)]; end rc0(i,:)=rc; dif=route_length(rc,dist)-route_length(route,dist); delta0(i,1)=dif; end up_delta=sort(delta0); [~,ind]=min(delta0); new_route=rc0(ind,:); end
|
(5)“修复”函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
|
function [Srepair,repair_length]=repair(removed,Sdestroy,dist) Srepair=Sdestroy;
while ~isempty(removed) nr=numel(removed); regret=zeros(nr,1); for i=1:nr visit=removed(i); [~,up_delta]=ins_route(visit,dist,Srepair); de12=up_delta(2)-up_delta(1); regret(i)=de12; end [~,max_index]=max(regret); reinsert_city=removed(max_index); Srepair=ins_route(reinsert_city,dist,Srepair); removed(max_index)=[]; end repair_length=route_length(Srepair,dist); end
|
(6)TSP路线图函数
1 2 3 4 5 6 7 8 9 10
|
function plot_route(route,x,y) figure route=[route route(1)]; plot(x(route),y(route),'k-o','MarkerSize',10,'MarkerFaceColor','w','LineWidth',1.5); xlabel('x'); ylabel('y'); end
|
(7)主函数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55
| tic clear clc
dataset=importdata('input.txt'); x=dataset(:,2); y=dataset(:,3); vertexes=dataset(:,2:3); h=pdist(vertexes); dist=squareform(h);
MAXGEN=300;
[Sinit,init_len]=construct_route(dist); init_length=route_length(Sinit,dist); str1=['初始总路线长度 = ' num2str(init_length)]; disp(str1)
Scurr=Sinit; curr_length=init_length; Sbest=Sinit; best_length=init_length;
gen=1; BestL=zeros(MAXGEN,1); while gen<=MAXGEN [Sdestroy,removed]=destroy(Scurr); [Srepair,repair_length]=repair(removed,Sdestroy,dist); if repair_length<curr_length Scurr=Srepair; curr_length=repair_length; end if curr_length<best_length Sbest=Scurr; best_length=curr_length; end disp(['第',num2str(gen),'代最优路线总长度 = ' num2str(best_length)]) BestL(gen,1)=best_length; gen=gen+1; end str2=['搜索完成! 最优路线总长度 = ' num2str(best_length)]; disp(str2)
figure; plot(BestL,'LineWidth',1); title('优化过程') xlabel('迭代次数'); ylabel('总距离');
plot_route(Sbest,x,y); toc
|
五、实验结果展示
1. 《MATLAB智能优化算法》(曹旺著) ↩
感谢鼓励