MATLAB实验报告,遗传算法解最短路径以及函数最小值问题

硕士生考查课程考试试卷

考试科目：MATLAB教程

考生姓名：考生学号：

学院：专业：

考生成绩：

任课老师(签名)

考试日期：20 年月日午时至时

《MATLAB 教程》试题：

A 、利用MATLA

B 设计遗传算法程序，寻找下图11个端点的最短路径，其中没有连接的端点表示没有路径。要求设计遗传算法对该问题求解。

a

d e

h

k

B 、设计遗传算法求解f (x)极小值，具体表达式如下：

3

21231(,,)5.12 5.12,1,2,3

i i i f x x x x x i =?=???-≤≤=?

∑ 要求必须使用m 函数方式设计程序。

C 、利用MATLAB 编程实现：三名商人各带一个随从乘船渡河，一只小船只能容纳二人，由他们自己划行，随从们密约，在河的任一岸，一旦随从的人数比商人多，就杀人越货，但是如何乘船渡河的大权掌握在商人手中，商人们怎样才能安全渡河？

D 、结合自己的研究方向选择合适的问题，利用MATLAB 进行实验。

以上四题任选一题进行实验，并写出实验报告。

选择题目： A 一、问题分析（10分）

1

4

10

11

如图如示，将节点编号，依次为 1.2.3.4.5.6.7.8.9.10.11，由图论知识，则可写出其带权邻接矩阵为：

0 2 8 1 500 500 500 500 500 500 500 2 0 6 500 1 500 500 500 500 500 500 8 6 0 7 500 1 500 500 500 500 500 1 500 7 0 500 500 9 500 500 500 500 500 1 500 500 0 3 500 2 500 500 500 500 500 1 500 3 0 4 500 6 500 500 500 500 500 9 500 4 0 500 500 1 500 500 500 500 500 2 500 500 0 7 500 9 500 500 500 500 500 6 500 7 0 1 2 500 500 500 500 500 500 1 500 1 0 4 500 500 500 500 500 500 500 9 2 4 0 注：为避免计算时无穷大数吃掉小数，此处为令inf=500。

问题要求求出任意两点间的最短路径，Floyd 算法采用的是在两点间尝试插入顶点，比较距离长短的方法。我思考后认为，用遗传算法很难找到一个可以统一表示最短路径的函数，但是可以对每一对点分别计算，然后加入for 循环，可将相互之间的所有情况解出。观察本题可发现，所有节点都是可双向行走，则可只计算i 到j 的路径与距离，然后将矩阵按主对角线翻折即可得到全部数据。 二、实验原理与数学模型（20分）

实现原理为遗传算法原理：

按所选择的适应度函数并通过遗传中的复制、交叉及变异对个体进行筛选，使得适应度高的个体被保留下来，组成新的群体，新的群体既继承了上一代的信息，又优于上一代。这样周而复始，群体中个体适应度不断提高，直到满足一定的条件。

数学模型如下： 设图G 由非空点集合12{,...}n V V V V = 和边集合12{,...}m E e e e = 组成，其中121221(,)e ,P ,)(P ,P ),

i i i i i i i i e P P E P =∈≠且若（则G 为一个有向图； 又设i e 的值为i a ，12{,...},m A a a a = 故G 可表示为一个三元组{,,}G P E A = 则求最短路径的数学模型可以描述为：

1

min *..

n

i i

i i i A E s t A A E E

=?

???∈∈?∑ 实验具体：

第一：编码与初始化

因采用自然编码，且产生的编码不能重复，于是我采用了randperm 函数产生不重复的随机自然数。因解题方法是使用的是计算每一对点，则我们编码时将第一个节点单独放入，合并成完整编码。

因为节点有11个，可采用一个1行11列的矩阵储存数据，同时，由于编号为数字，可直接使用数字编码表示路径的染色体。具体如下：

采用等长可变染色体的方式，例如由2到9的路径，染色体编码可能为（2，5，1，8，4，6，9，3，10，7，11），超过9之后的编码，用来进行算子的运算，不具备实际意义。

第二：计算适应度，因取最短路径值，即最小值，常用方法为C-F(x)或C/F(x)（C 为一常数），此处采用前一种方式。于是，可进一步计算相对适应度。

第三：选择与复制

采用轮盘赌算法，产生一个随机值，比较它与累计相对适应度的关系，从而选择出优良个体进入下一代。

第四：交叉。

因编码是不重复的数字，所以采用传统的交叉方法，即上一行与下一行对位交叉，会产生无效路径，于是，采用了不同的交叉方法，具体如下：

（1）在表示路径的染色体Tx 和Ty 中，随机选取两个基因座（不能为起点基因座）i 和j, 即将i 个基因座和第j 个基因座之间的各个基因座定义为交叉域，并将交叉的内容分别记忆为temp1和temp2。

（2）根据交叉区域中的映射关系，在个体Tx 中找出所有与temp2相同的元素，在个体Ty 中找出所有与temp1相同的元素，全部置为0。

（3）将个体Tx 、Ty 进行循环左移，遇到0就删除，直到编码串中交叉区域的左端不再有0：然后将所有空位集中到交叉区域，而将交叉区域内原有的基因依次向后移动。因0元素可能较多，在程序实现时，我是将非零元素提出，后面再合成。

（4）将temp2插入到Tx 的交叉区域，temp1插入到Ty 的交叉区域。形成新的染色体[1]。

第五：变异

染色体编码为从1到11的无重复编码，所以不能采用一般的生成一个随机数替代的办法。此处采用交换变异法。即随机产生两个数，交换两个节点的顺序。例：[1,2,3,4,5,6,7,8,9,10,11],13,28p K K === 则新染色体编码为：

[1,2,8,4,5,6,7,3,9,10,11]p =

三、实验过程记录（含基本步骤、程序代码及异常情况记录等）（60分）

首先，写程序，修复Bug 。

然后，调试种群数量，遗传代数，交叉概率，变异概率等，不断运行程序，以达到较理想的状态。

有一次异常情况：算出来的最短距离均为0，最短路径没有终点出现，经过分析发现，因为交叉处的代码较复杂，弄错了一点，导致新基因有部分为非法基因。最后采用提出非零数值组成向量，再合成新基因的方式解决。

Matlab程序代码如下：

clc;clear;

%初始化参数

%注：popsize=200,MaxGeneration=100,约跑2分钟。若不要求太精确，可减少循环次数。

pointnumber=11; %节点个数

Popsize=200; %种群规模，只能取偶数（因67行的循环）MaxGeneration=100; %最大代数

Pc=0.8;Pm=0.3; %交叉概率和变异概率

A=[0 2 8 1 50 50 50 50 50 50 50

2 0 6 50 1 50 50 50 50 50 50

8 6 0 7 50 1 50 50 50 50 50

1 50 7 0 50 50 9 50 50 50 50

50 1 50 50 0 3 50 2 50 50 50

50 50 1 50 3 0 4 50 6 50 50

50 50 50 9 50 4 0 50 50 1 50

50 50 50 50 2 50 50 0 7 50 9

50 50 50 50 50 6 50 7 0 1 2

50 50 50 50 50 50 1 50 1 0 4

50 50 50 50 50 50 50 9 2 4 0]; %带权邻接矩阵。

A(A==50)=500; %取值50过小而修正为500；

Bestindividual=zeros(MaxGeneration,1);

outdistance=zeros(11,11);

outpath=cell(11,11); %用于存放11个点相互之间的最短路径

%****** 生成初始种群******

for a=1:pointnumber %起点的编号

%a=1;

tempvary=[1 2 3 4 5 6 7 8 9 10 11];

tempvary(a)=[]; %暂时剔除起点

tempmatrix=a*ones(Popsize,1); %将起点单独放一矩阵

path=zeros(Popsize,pointnumber-1); %声明矩阵大小，避免减慢速度

for i=1:Popsize

temprand=randperm(pointnumber-1);

path(i,:)=tempvary(temprand(1:end)); %生成一系列剔除起点的随机路线

end

path=[tempmatrix path]; %合成包括起点的完整路线

[row,col]=size(path);

for b=a:pointnumber %终点的编号

%b=10;

for k=1:1:MaxGeneration

for i=1:row

position2=find(path(i,:)==b); %找出终点在路线中的位置

pathlong(i)=0;

for j=1:position2-1

pathlong(i)=pathlong(i)+A(path(i,j),path(i,j+1));

end

end

%计算适应度

Fitness=length(A)*max(max(A))-pathlong; %因要求最小值，采且常数减函数值构造适应度Fitness=Fitness./sum(Fitness);

%****** Step 1 : 选择最优个体******

Bestindividual(k)=min(pathlong);

[Orderfi,Indexfi]=sort(Fitness); %按照适应度大小排序

Bestfi=Orderfi(Popsize); %Oderfi中最后一个即是最大的适应度

BestS=path(Indexfi(Popsize),:); %记录每一代中最优个体的路线

%****** Step 2 : 选择与复制操作******

temppath=path;

roulette=cumsum(Fitness);

for i=1:Popsize

tempP=rand(1);

for j=1:length(roulette)

if tempP

break;

end

end

path(i,:)=temppath(j,:);

end

%************ Step 3 : 交叉操作************

temppath2=path;

for i=1:2:row

tempP2=rand(1);

if(tempP2

temPm2=fix((rand(1)+0.2)*10); %因起点基因不能改变

temPm3=fix((rand(1)+0.2)*10); %随机取出两个位置为2到11基因座

temPm4=min(temPm2,temPm3);

temPm5=max(temPm2,temPm3);

temp1=path(i,temPm4:temPm5); %将两点之间的基因储存，方便交叉

temp2=path(i+1,temPm4:temPm5);

[c d]=find(ismember(path(i,:),temp2));

path(i,d)=0; %找出i行在i+1行取出区域中的数，置为0

[e f]=find(ismember(path(i+1,:),temp1));

path(i+1,f)=0; %找出i+1行在i行取出区域中的数，置为0

[g h]=find(path(i,:)~=0);

v1=path(i,h); %取出i行的非零元素，成一向量

[l m]=find(path(i+1,:)~=0);

v2=path(i+1,m); %取出i+1行的非零元素，成一向量

path(i,:)=[v1(1:temPm4-1) temp2 v1(temPm4-1+size(temp1):end)];

path(i+1,:)=[v2(1:temPm4-1) temp1 v2(temPm4-1+size(temp2):end)]; %基因交叉end

end

path(Popsize,:)=BestS;

%************ Step 4: 变异操作**************

for i=1:Popsize

tempPm=rand(1);

if(tempPm

temPm6=fix((rand(1)+0.2)*10);

temPm7=fix((rand(1)+0.2)*10); %产生两个用于交换的随机数

tempvessel=path(i,temPm6); %交换前用一临时容器存放数据

path(i,temPm6)=path(i,temPm7);

path(i,temPm7)=tempvessel; %变异交换

end

end

path(Popsize,:)=BestS;

end

[aa bb]=find(BestS==b); %找出终点

Bestpath=BestS(1:bb); %剔除后面无用的点，留下实际路线

outdistance(a,b)=Bestindividual(k); %将最短距离写入矩阵

outpath{a,b}=Bestpath; %写入路径，因数据类型为矩阵，所以采用元胞数组储存

end

end

for i=1:pointnumber

for j=1:i

outdistance(i,j)=outdistance(j,i); %实现距离的对称

outpath{i,j}=fliplr(outpath{j,i}); %实现路径的对称与翻转

end

end

%*************** 结果输出*****************

outdistance

celldisp(outpath)

%xlswrite('tempdata.xls', outpath) %存入excel中进行操作

四、实验结果与总结（10分）

距离矩阵：a(i,j) i表示起点，j表示终点。

outdistance =

0 2 7 1 3 6 10 5 12 11 14

2 0 5

3 1

4 8 3 10 9 12

7 5 0 7 4 1 5 6 7 6 9

1 3 7 0 4 8 9 6 11 10 13

3 1

4 4 0 3 7 2 9 8 11

6 4 1 8 3 0 4 5 6 5 8

10 8 5 9 7 4 0 9 2 1 4

5 3

6 6 2 5 9 0

7

8 9

12 10 7 11 9 6 2 7 0 1 2

11 9 6 10 8 5 1 8 1 0 3

14 12 9 13 11 8 4 9 2 3 0

路径：b(i,j) i表示起点，j表示终点。

outpath:

此程序运算速度有待提高，程序的收敛速度不是很快。可能的原因如下：

（1）在变异操作时，可能将本来很好的解弃掉，换来更差的染色体，导致收敛速度不佳。解决办法：可以在变异操作时，增加个体求优的自学习过程。

即在某位基因变异后，计算新染色体的适应函数值，若适应值变大，即路径更短，则保留；否则，保持原来的染色体不变。

（2）算法的进一步改进，例如可加入Floyd算法的思想，在父代产生子代的过程中，不是单纯的交叉，可以考虑随机加入顶点是否路径变短。

参考文献：

[1]康晓军，王茂才.基于遗传算法的最短路径问题的求解.计算机工程与应用[J],2008,44(23)

第二题代码：

clc;clear;

%Rosenbrock函数的极大值0-1编码的GA算法

%初始参数

tic;

Size=80;

G=100;

CodeL=10;

umax=5.12;

umin=-5.12;

E=round(rand(Size,3*CodeL)); %生成初始种群

%主程序

for k=1:1:G

time(k)=k;

for s=1:1:Size

m=E(s,:);

y1=0;y2=0; y3=0;

%解码方法

m1=m(1:1:CodeL);

for i=1:1:CodeL

y1=y1+m1(i)*2^(i-1);

end

x1=(umax-umin)*y1/1023+umin;

m2=m(CodeL+1:1:2*CodeL);

for i=1:1:CodeL

y2=y2+m2(i)*2^(i-1);

end

x2=(umax-umin)*y2/1023+umin;

m3=m(2*CodeL+1:1:end);

for i=1:1:CodeL

y3=y3+m3(i)*2^(i-1);

end

x3=(umax-umin)*y3/1023+umin;

F(s)=x1^2+x2^2+x3^3;

end

%****** Step 1 : 选择最优个体******

BestJ(k)=min(F); %记录每一代中最大个体的函数值

fi=F; %适应度函数

[Oderfi,Indexfi]=sort(fi); %按照适应度大小排序

Bestfi=Oderfi(1); %Oderfi中最后一个即是最大的适应度BestS=E(Indexfi(1),:); %记录每一代中最优个体的0-1编码

bfi(k)=Bestfi; %记录每一代中最优个体的适应度

%****** Step 2 : 选择与复制操作******

fi_sum=sum(fi);

fi_Size=(Oderfi/fi_sum)*Size; %计算个体相对适应度

fi_S=floor(fi_Size); %对80个个体依据相对适应度进行划分等级

kk=1;

for i=1:1:Size

for j=1:1:fi_S(i) %选择等级高的个体，等级越高被选次数越多

TempE(kk,:)=E(Indexfi(i),:);

kk=kk+1; %选择进入下一代个体的个数，显然不够80个个体end

end

%************ Step 3 : 交叉操作************

pc=0.60;

n=ceil(20*rand);

for i=1:2:(Size-1)

temp=rand;

if pc>temp %交叉条件

TempE(i,n:end)=E(i+1,n:end);

TempE(i+1,n:end)=E(i,n:end);

end

end

TempE(Size,:)=BestS; %最优个体保留

E=TempE; %种群替换

%************ Step 4: 变异操作**************

%pm=0.001;

%pm=0.001-[1:1:Size]*(0.001)/Size; %自适应变异概率

%pm=0.0; %没有变异

pm=0.1; %较大的变异概率

for i=1:1:Size

for j=1:1:2*CodeL

temp=rand;

if pm>temp %变异条件

if TempE(i,j)==0

TempE(i,j)=1;

else

TempE(i,j)=0;

end

end

end

end

TempE(Size,:)=BestS; %保留当代最优个体

E=TempE; %种群替换end

%*************** 结果输出***************** Max_V alue=Bestfi

BestS

x1

x2

figure(1);

plot(time,BestJ);

xlabel('Times');ylabel('Best J');

figure(2);

plot(time,bfi);

xlabel('times');ylabel('Best F');

toc;