数模论坛

 找回密码
 注-册-帐-号
搜索
热搜: 活动 交友 discuz
查看: 12282|回复: 15

TSP问题的遗传算法

[复制链接]
发表于 2003-8-4 03:16:58 | 显示全部楼层 |阅读模式
我这几天做了一个队员选择问题,其中一个问题我是用遗传算法做的,现在我把它整理成解决TSP问题的遗传算法:旅行商问题(Traveling Saleman Problem,简称TSP):
    已知n个城市之间的相互距离,现有一个推销员必须遍访这n个城市,并且每个城市只能访问一次,最后又必须返回出发城市。如何安排他对这些城市的访问次序,可使其旅行路线的总长度最短?
    用图论的术语来说,假设有一个图G=(V,E),其中V是顶点集,E是边集,设D=(dij)是由顶点i和顶点j之间的距离所组成的距离矩阵,旅行商问题就是求出一条通过所有顶点且每个顶点只通过一次的具有最短距离的回路。
    这个问题可分为对称旅行商问题(dij=dji,,任意i,j=1,2,3,…,n)和非对称旅行商问题(dij≠dji,,任意i,j=1,2,3,…,n)。
    若对于城市V={v1,v2,v3,…,vn}的一个访问顺序为T=(t1,t2,t3,…,ti,…,tn),其中ti∈V(i=1,2,3,…,n),且记tn+1= t1,则旅行商问题的数学模型为:
     min     L=Σd(t(i),t(i+1))  (i=1,…,n)
    旅行商问题是一个典型的组合优化问题,并且是一个NP难问题,其可能的路径数目与城市数目n是成指数型增长的,所以一般很难精确地求出其最优解,本文采用遗传算法求其近似解。
    遗传算法:
初始化过程:用v1,v2,v3,…,vn代表所选n个城市。定义整数pop-size作为染色体的个数,并且随机产生pop-size个初始染色体,每个染色体为1到18的整数组成的随机序列。
适应度f的计算:对种群中的每个染色体Vi,计算其适应度,f=Σd(t(i),t(i+1)).
评价函数eval(Vi):用来对种群中的每个染色体Vi设定一个概率,以使该染色体被选中的可能性与其种群中其它染色体的适应性成比例,既通过轮盘赌,适应性强的染色体被选择产生后台的机会要大,设alpha∈(0,1),本文定义基于序的评价函数为eval(Vi)=alpha*(1-alpha).^(i-1) 。[随机规划与模糊规划]
选择过程:选择过程是以旋转赌轮pop-size次为基础,每次旋转都为新的种群选择一个染色体。赌轮是按每个染色体的适应度进行选择染色体的。
   step1 、对每个染色体Vi,计算累计概率qi,q0=0;qi=Σeval(Vj)   j=1,…,i;i=1,…pop-size.
   step2、从区间(0,pop-size)中产生一个随机数r;
   step3、若qi-1<r<qi,则选择第i个染色体 ;
   step4、重复step2和step3共pop-size次,这样可以得到pop-size个复制的染色体。
Grefenstette编码:由于常规的交叉运算和变异运算会使种群中产生一些无实际意义的染色体,本文采用Grefenstette编码《遗传算法原理及应用》可以避免这种情况的出现。所谓的Grefenstette编码就是用所选队员在未选(不含淘汰)队员中的位置,如:
          8        15        2        16        10        7        4        3        11        14        6        12        9        5        18        13        17        1
          对应:
          8        14        2        13        8        6        3        2        5        7        3        4        3        2        4        2        2        1。
交叉过程:本文采用常规单点交叉。为确定交叉操作的父代,从 到pop-size重复以下过程:从[0,1]中产生一个随机数r,如果r<pc ,则选择Vi作为一个父代。
           将所选的父代两两组队,随机产生一个位置进行交叉,如:
          8        14        2        13        8        6        3        2        5        7        3        4        3        2        4        2        2        1
          6        12        3        5        6        8        5        6        3        1        8        5        6        3        3        2        1        1
交叉后为:
         8        14        2        13        8        6        3        2        5        1        8        5        6        3        3        2        1        1
         6        12        3        5        6        8        5        6        3        7        3        4        3        2        4        2        2        1
变异过程:本文采用均匀多点变异。类似交叉操作中选择父代的过程,在r<pm 的标准下选择多个染色体Vi作为父代。对每一个选择的父代,随机选择多个位置,使其在每位置按均匀变异(该变异点xk的取值范围为[Ukmin,Ukmax],产生一个[0,1]中随机数r,该点变异为x'k=Ukmin+r(Ukmax-Ukmin))操作。如:
         8        14        2        13        8        6        3        2        5        7        3        4        3        2        4        2        2        1
      变异后:
        8        14        2        13        10        6        3        2        2        7        3        4        5        2        4        1        2        1
反Grefenstette编码:交叉和变异都是在Grefenstette编码之后进行的,为了循环操作和返回最终结果,必须逆Grefenstette编码过程,将编码恢复到自然编码。
循环操作:判断是否满足设定的带数xzome,否,则跳入适应度f的计算;是,结束遗传操作,跳出。



function [bestpop,trace]=ga(D,termOps,num,pc,cxOps,pm,alpha)
%
%————————————————————————
%[bestpop,trace]=ga(D,termOps,num,pc,cxOps,pm,alpha)
%D:距离矩阵
%termOps:种群带数
%num:每带染色体的个数
%pc:交叉概率
%cxOps:由于本程序采用单点交叉,交叉点的设置在本程序中没有很好的解决,所以本文了采用定点,即第cxOps,可以随机产生。
%pm:变异概率
%alpha:评价函数eval(Vi)=alpha*(1-alpha).^(i-1).
%bestpop:返回的最优种群
%trace:进化轨迹
%------------------------------------------------
%####@@@##版权所有!欢迎广大网友改正,改进!##@@@####
%E-mail:tobysidney33@sohu.com
%####################################################
%
citynum=size(D,2);
n=nargin;
if n<2
    disp('缺少变量!!')
    disp('^_^开个玩笑^_^')
end
if n<2
    termOps=500;
    num=50;
    pc=0.25;
    cxOps=3;
    pm=0.30;
    alpha=0.10;
end
if n<3
    num=50;
    pc=0.25;
    cxOps=3;
    pm=0.30;
    alpha=0.10;
end
if n<4
    pc=0.25;
    cxOps=3;
    pm=0.30;
    alpha=0.10;
end
if n<5
    cxOps=3;
    pm=0.30;
    alpha=0.10;
end
if n<6
    pm=0.30;
    alpha=0.10;
end
if n<7
    alpha=0.10;
end
if isempty(cxOps)
    cxOps=3;
end

[T]=initializega(num,citynum);
for i=1:termOps
[L]=f(D,T);
[x,y]=find(L==max(L));
trace(i)=-L(y(1));
bestpop=T(y(1),;
[T]=select(T,L,alpha);
[G]=grefenstette(T);
[G1]=crossover(G,pc,cxOps);
[G]=mutation(G1,pm);  %均匀变异
[T]=congrefenstette(G);
end

---------------------------------------------------------
function [T]=initializega(num,citynum)
for i=1:num
    T(i,=randperm(citynum);
end
-----------------------------------------------------------
function [L]=f(D,T)
[m,n]=size(T);
for k=1:m
    for i=1:n-1
      l(k,i)=D(T(k,i),T(k,i+1));
    end
      l(k,n)=D(T(k,n),T(k,1));
      L(k)=-sum(l(k,);
end
-----------------------------------------------------------
function [T]=select(T,L,alpha)
[m,n]=size(L);
T1=T;
[beforesort,aftersort1]=sort(L,2);%fsort from l to u
for i=1:n
    aftersort(i)=aftersort1(n+1-i);      %change
end
for k=1:n;
    T(k,:)=T1(aftersort(k),:);
    L1(k)=L(aftersort(k));
end
T1=T;
L=L1;
for i=1:size(aftersort,2)
    evalv(i)=alpha*(1-alpha).^(i-1);
end
m=size(T,1);
q=cumsum(evalv);
qmax=max(q);
for k=1:m
    r=qmax*rand(1);
    for j=1:m
        if j==1&r<=q(1)
            T(k,:)=T1(1,:);
        elseif j~=1&r>q(j-1)&r<=q(j)
            T(k,:)=T1(j,:);
        end
    end
end
--------------------------------------------------
function [G]=grefenstette(T)
[m,n]=size(T);
for k=1:m
    T0=1:n;
   for i=1:n
       for j=1:length(T0)
           if T(k,i)==T0(j)
              G(k,i)=j;
              T0(j)=[];
              break
           end
       end
    end
end
-------------------------------------------
function [G]=crossover(G,pc,cxOps)
[m,n]=size(G);
ran=rand(1,m);
r=cxOps;
[x,ru]=find(ran<pc);
if ru>=2
    for k=1:2:length(ru)-1
       G1(ru(k),:)=[G(ru(k),[1:r]),G(ru(k+1),[(r+1):n])];
       G(ru(k+1),:)=[G(ru(k+1),[1:r]),G(ru(k),[(r+1):n])];
       G(ru(k),:)=G1(ru(k),:);
    end
end
--------------------------------------------
function [G]=mutation(G,pm)    %均匀变异
[m,n]=size(G);
ran=rand(1,m);
r=rand(1,3);      %dai gai jin
rr=floor(n*rand(1,3)+1);
[x,mu]=find(ran<pm);
for k=1:length(mu)
    for i=1:length(r)
        umax(i)=n+1-rr(i);
        umin(i)=1;
        G(mu(k),rr(i))=umin(i)+floor((umax(i)-umin(i))*r(i));
    end
end
---------------------------------------------------
function [T]=congrefenstette(G)
[m,n]=size(G);
for k=1:m
    T0=1:n;
   for i=1:n
      T(k,i)=T0(G(k,i));
      T0(G(k,i))=[];
   end
end
-------------------------------------------------
请多指教!

[em06][em06]
发表于 2003-8-4 03:45:03 | 显示全部楼层
恩,不错。楼主编程功底很好。
偶也做过TSP的遗传算法,但是效果不够满意,是用C++写的。
你有兴趣可以挑战世界记录,偶做过71009个城市的TSP,和世界记录差了13%,太丢人了。
记录是美国的普林斯顿大学的,他们用的混合整数规划来做,用C语言写的,技术之高令人震惊。记录在http://www.math.princeton.edu/tsp/world/chlog.html。
 楼主| 发表于 2003-8-4 03:52:19 | 显示全部楼层
谢谢鼓励!
我想程序还有很多地方需要改进,请多指教!
发表于 2003-8-4 03:57:35 | 显示全部楼层
当时princeton大学的结果:

偶的结果:


距离差了13%。
 楼主| 发表于 2003-8-4 04:08:39 | 显示全部楼层
看了!
你已经很厉害了,我现在大受鼓舞!
一方面发现自己的知识太少了,另一方面我希望更多的年轻人挑战极限(包括我)
我现在有点震惊!非常感谢!
发表于 2003-8-4 04:18:55 | 显示全部楼层
呵呵,去挑战极限吧,现在在很多国外的网站上举行过有很多数学问题的记录比赛,几乎没有几个中国人保持过记录,说明我国的学术水平还不高啊,见过的一个是杭州某大学的一个本科生参与过。还有就是上面的图中竟然没有台湾,可见美国人对中国的蔑视,争取做一下超过他们,并附上完整的中国地图也是偶的心愿,可是自己还是太忙了,几乎被逼得没有时间去选择自己想研究的问题。楼上一起努力啊。
 楼主| 发表于 2003-8-4 04:31:02 | 显示全部楼层
一定!大学还有两年,我想在只要有激情和大虾的指导(引),会出现新的局面
还有台湾一定会加上去的!
每一个中国人都不希望中国的地图少任何一块地方,包括台湾
发表于 2003-8-4 05:00:51 | 显示全部楼层
呵呵,真是高深呀......精彩精彩!
发表于 2003-8-5 08:44:33 | 显示全部楼层
所谓的grefenstette编码就是用所选队员在未选(不含淘汰)队员中的位置,如:
          8 15 2 16 10 7 4 3 11 14 6 12 9 5 18 13 17 1
          对应:
          8 14 2 13 8 6 3 2 5 7 3 4 3 2 4 2 2 1。

其他的都看了,就这个过程不太懂。
请教!
 楼主| 发表于 2003-8-6 03:28:34 | 显示全部楼层
比如说:
如果把城市A,B,C,D的顺序记为1,2,3,4的话
那么ACDB的标号为1342;
每个城市在未到过的所有城市列表中的位置为
A在ABCD中是1
C 在BCD的位置为2
D在DC的位置为1
C在C的位置为1
故ACDB的编码为1211
您需要登录后才可以回帖 登录 | 注-册-帐-号

本版积分规则

小黑屋|手机版|Archiver|数学建模网 ( 湘ICP备11011602号 )

GMT+8, 2024-11-27 06:20 , Processed in 0.077804 second(s), 18 queries .

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表