免费代码 | 头脑风暴优化(BSO)算法求解旅行商问题(TSP)MATLAB代码

原创
08/07 21:05
阅读数 177


各位小伙伴可在闲鱼搜索 优化算法交流地,即可搜索到官方闲鱼账号, 谨防上当受骗
hello,大家立秋快乐!今天为各位免费分享 头脑风暴优化(BSO)算法求解旅行商问题(TSP)(附MATLAB代码) 这篇推文的源代码
近期有很多粉丝向我们反应闲鱼上有商家高价售卖我们公众号的代码,对于这种情况,我们真诚地希望各位小伙伴 多点击右下角 ,让更多的小伙伴能够看到我们,防止上当受骗。你们的支持是我们继续推出原创推文的动力。

这份源代码包含7个函数,分别如下:

01 | 主函数
主函数的输入是文本文件pr226.txt 第1列是序号,第2列是x坐标,第3列是y坐标 ),输出是最优路线 文本文件可根据自己需要进行替换,只要保持3列的这种形式即可,文本文件格式如下

序号

x坐标

y坐标


%%      @作者:随心390%      @微信公众号:优化算法交流地%ticclearclc%% 导入数据pr226=importdata('pr226.txt');N=size(pr226,1);                %城市数目vertexs=pr226(:,2:3);           %各点xy坐标x=vertexs(:,1);                 %x坐标y=vertexs(:,2);                 %y坐标h=pdist(vertexs);dist=squareform(h);             %距离矩阵%% 参数初始化MAXGEN=1000;                    %最大迭代次数NIND=50;                        %种群数目cluster_num=5;                  %聚类数目p_replace=0.1;                  %用随机解替换一个聚类中心的概率p_one=0.5;                      %选择1个聚类的概率p_two=1-p_one;                  %选择2个聚类的概率,p_two=1-p_onep_one_center=0.3;               %选择1个聚类中聚类中心的概率p_two_center=0.2;               %选择2个聚类中聚类中心的概率%% 种群初始化Population=InitPop(NIND,N);%% 主循环gen=1;                          %计数器初始化BestPop=zeros(MAXGEN,N);        %记录每次迭代过程中全局最优个体BestL=zeros(MAXGEN,1);          %记录每次迭代过程中全局最优个体的总距离while gen<=MAXGEN    %% 计算目标函数值    Obj=ObjFunction(Population,dist);    %% K-means聚类    [Idx,C,sumD,D]=kmeans(Obj,cluster_num,'Distance','cityblock','Replicates',2);    cluster=cell(cluster_num,2);                    %将解储存在每一个聚类中    order_cluster=cell(cluster_num,2);              %将储存在每一个聚类中的解按照目标函数值排序    for i=1:cluster_num        cluster{i,1}=Population(Idx==i,:);          %将个体按照所处的聚类编号储存到对应的聚类中        cluster_row(i)=size(cluster{i,1},1);           %计算当前聚类中个体数目        for j=1:cluster_row(i)            Individual=cluster{i,1}(j,:);           %当前聚类中第j个个体            cluster{i,2}(j,1)=ObjFunction(Individual,dist); %计算当前聚类中第j个个体的目标函数值        end        [order_cluster{i,2},order_index]=sort(cluster{i,2}) ; %将当前聚类中的所有个体按照目标函数值从小到大的顺序进行排序        order_cluster{i,1}=cluster{i,1}(order_index,:); %将当前聚类中的所有个体按照排序结果重新排列        order_index=0;                                  %重置排序序号    end    cluster_fit=cell2mat(order_cluster);                %将聚类的元胞数组转换为矩阵,最后一列为个体的目标函数值    %% 以一定的概率随机从m个聚类中心中选择出一个聚类中心,并用一个新产生的随机解更新这个被选中的聚类中心    R1=rand(1,1);    if R1<=p_replace        %随机选择一个聚类中心        repalce_cluster_num=randi([1,cluster_num],1,1);        %随机产生一个解        replace_solution=randperm(N);        %并用这个新产生的随机解更新这个被选中的聚类中心        order_cluster{repalce_cluster_num,1}(1,:)=replace_solution;        %计算新解的目标函数值        replace_solution_fitness=ObjFunction(replace_solution,dist);        %将新解的目标函数值储存到order_cluster中        order_cluster{repalce_cluster_num,2}(1,:)=replace_solution_fitness;    end    %% 更新这n个个体    for i=1:NIND        %如果随机数小于选择1个聚类的概率,则随机选择1个聚类        if rand()<p_one            select_one_cluster=randi([1,cluster_num],1,1); %选择选择一个聚类            %如果随机数小于选择1个聚类中聚类中心的概率,或当前聚类中只有一个个体            if rand()<p_one_center||cluster_row(select_one_cluster)==1                select_ind=order_cluster{select_one_cluster,1}(1,:);        %选择当前解的聚类中心(只有一个个体时,就选择该个体)            else                r_1= randi([2,cluster_row(select_one_cluster)],1,1);        %随机选择当前聚类中除聚类中心外的其它个体序号                select_ind=order_cluster{select_one_cluster,1}(r_1,:);      %随机选择当前聚类中除聚类中心外的其它个体            end            indi_temp=Swap(select_ind);                                     %将该个体进行交换操作        else %如果随机数不小于选择1个聚类的概率,则随机选择2个聚类            cluster_two=[0,0];                                              %随机产生两个聚类的序号            while cluster_two(1,1)==cluster_two(1,2)                cluster_two=randi([1,cluster_num],1,2);                     %这两个聚类序号不能相同            end            %如果随机数小于选择2个聚类中聚类中心的概率,或当前两个聚类中都只有一个个体            if (rand()<p_two_center)||(cluster_row(cluster_two(1,1))==1&&cluster_row(cluster_two(1,2))==1)                %选择这2个聚类中聚类中心                select_ind1=order_cluster{cluster_two(1,1),1}(1,:);      %第一个被选择的聚类中心                select_ind2=order_cluster{cluster_two(1,2),1}(1,:);      %第二个被选择的聚类中心            else                %如果第1个选择的聚类中只有一个个体                if cluster_row(cluster_two(1,1))==1                    r_2=randi([2,cluster_row(cluster_two(1,2))],1,1);        %选择第2个聚类中除聚类中心外的其它个体                    select_ind1=order_cluster{cluster_two(1,1),1}(1,:);      %第一个聚类中被选择的个体                    select_ind2=order_cluster{cluster_two(1,2),1}(r_2,:);    %第二个聚类中被选择的个体                elseif cluster_row(cluster_two(1,2))==1  %如果第2个选择的聚类中只有一个个体                    r_3=randi([2,cluster_row(cluster_two(1,1))],1, 1);        %选择第1个聚类中除聚类中心外的其它个体                    select_ind1=order_cluster{cluster_two(1,1),1}(r_3,:);     %第一个聚类中被选择的个体                    select_ind2=order_cluster{cluster_two(1,2),1}(1,:);       %第二个聚类中被选择的个体                elseif cluster_row(cluster_two(1,1))>1&&cluster_row(cluster_two(1,2))>1                    %%选择这2个聚类中除聚类中心外的其它个体                    r_4=randi([2,cluster_row(cluster_two(1,1))],1,1);                    r_5=randi([2,cluster_row(cluster_two(1,2))],1,1);                    select_ind1=order_cluster{cluster_two(1,1),1}(r_4,:);       %第一个聚类中被选择的个体                    select_ind2=order_cluster{cluster_two(1,2),1}(r_5,:);       %第二个聚类中被选择的个体                end            end            [child1,child2,min_index,start_c]=heuristic_crossover(select_ind1,select_ind2,dist);  %启发式交叉算子            %如果子代个体1目标函数值更小,则将indi_temp赋值为child1,否则赋值为child2            if min_index==1                indi_temp=child1;            else                indi_temp=child2;            end        end        fit_indi_temp=ObjFunction(indi_temp,dist);        %如果fit_indi_temp比原来位置上的个体目标函数值更小,则更新该位置上的个体        if fit_indi_temp<order_cluster{Idx(i),2}(1,:)             Population(i,:) =indi_temp(1,:);        end    end    %% 计算目标函数值    Obj=ObjFunction(Population,dist);    [min_len,min_index]=min(Obj);               %当前种群中最优个体以及所对应的序号    BestL(gen,1)=min_len;                       %记录当前迭代中最优个体       BestPop(gen,:)=Population(min_index,:);     %记录当前迭代中最优个体的总距离    %% 打印各代最优解    disp(['第',num2str(gen),'代最优个体的总距离为',num2str(min_len)]);    %% 计数器加1    gen=gen+1;end%% 计算最终种群的目标函数值Obj=ObjFunction(Population,dist);[best_len,best_index]=min(Obj);         bestR=Population(best_index,:);         %全局最优个体PlotRoute(bestR,x,y);%% 打印每次迭代的全局最优个体的总距离变化趋势图figure;plot(BestL,'LineWidth',1);title('优化过程')xlabel('迭代次数');ylabel('总距离');
toc


02 | 种群初始化函数
种群初始化函数InitPop的输入是种群大小NIND、城市数目N,输出是随机生成的初始种群Population。
%%      @作者:随心390%      @微信公众号:优化算法交流地%%% 种群初始化%输入NIND:种群大小%输入N:城市数目%输出Population:随机生成的初始种群function Population=InitPop(NIND,N)Population=zeros(NIND,N);                %种群初始化为NIND行N列的零矩阵for i=1:NIND    Population(i,:)=randperm(N);         %每个个体为1~N的随机排列endend


03 | 目标函数
目标函数ObjFunction的输入是种群Population、距离矩阵dist,输出是每个个体的目标函数值Obj,即每个个体的总距离。
%%      @作者:随心390%      @微信公众号:优化算法交流地%%% 计算种群目标函数值,即每个个体的总距离%输入Population:种群%输入dist:距离矩阵%输出Obj:每个个体的目标函数值,即每个个体的总距离function Obj=ObjFunction(Population,dist)NIND=size(Population,1);                %种群大小Obj=zeros(NIND,1);                      %目标函数初始化为0for i=1:NIND    route=Population(i,:);              %当前个体    Obj(i,1)=RouteLength(route,dist);   %计算当前个体的总距离endend


04 | 路线长度计算函数
路线长度计算函数RouteLength的输入是路线route、距离矩阵dist,输出是该路线总距离L。
%%      @作者:随心390%      @微信公众号:优化算法交流地%%输入route:               路线%输入dist:                距离矩阵%输出L:                   该路线总距离function L=RouteLength(route,dist)    n=length(route);    route=[route route(1)];    L=0;    for k=1:n         i=route(k);        j=route(k+1);         L=L+dist(i,j);     endend


05 | 交换函数
交换函数Swap的输入是路线1 route1,输出是经过交换结构变换后的路线2 route2
%%      @作者:随心390%      @微信公众号:优化算法交流地%
%比如说有6个城市,当前解为123456,我们随机选择两个位置,然后将这两个位置上的元素进行交换。%比如说,交换2和5两个位置上的元素,则交换后的解为153426。
%输入route1: 路线1%输出route2: 经过交换结构变换后的路线2function route2=Swap(route1)n=length(route1);seq=randperm(n);I=seq(1:2);i1=I(1);i2=I(2);route2=route1;route2([i1 i2])=route1([i2 i1]);end


06 | 启发式交叉函数
启发式交叉函数heuristic_crossover的输入是父代个体1 parent1、父代个体2 parent2、距离矩阵 dist,输出是子代个体1 child1、子代个体2 child2、目标函数值更小的个体序号min_index、初始选择一个城市作为子代个体的起点start_c
%%      @作者:随心390%      @微信公众号:优化算法交流地%%% 启发式交叉算子%输入parent1:                 父代个体1%输入parent2:                 父代个体2%输入dist:                    距离矩阵%输出child1:                  子代个体1%输出child2:                  子代个体2%输出min_index:               目标函数值更小的个体序号%输出start_c:                 初始选择一个城市作为子代个体的起点function [child1,child2,min_index,start_c]=heuristic_crossover(parent1,parent2,dist)n=numel(parent1);               %城市数目child1=zeros(1,n);              %初始化子代个体1child2=zeros(1,n);              %初始化子代个体2start=randperm(n,1);            %随机选择一个城市作为子代个体的起点start_c=start;                  %复制startchild1(1)=start;                %更新子代个体1的起点为startchild2(1)=start;                %更新子代个体2的起点为start
%% 终止条件是parent1_c1只剩一个点时parent1_c1=parent1; %复制parent1parent2_c1=parent2; %复制parent2pos=2; %记录child中的位置while numel(parent1_c1)~=1 start1=find(parent1_c1==start,1,'first'); %找出start在parent1的位置 start2=find(parent2_c1==start,1,'first'); %找出start在parent2的位置 %% 找出start在parent1_c1的右紧邻点 if start1==numel(parent1_c1) right1=parent1_c1(1); %如果start1是parent1最后一个点,则start1右边紧邻的点是parent1第一个点 else right1=parent1_c1(start1+1); %start1右边紧邻的点 end %% 找出start在parent2_c1的右紧邻点 if start2==numel(parent2_c1) right2=parent2_c1(1); %如果start2是parent2最后一个点,则start2右边紧邻的点是parent2第一个点 else right2=parent2_c1(start2+1); %start2右边紧邻的点 end
%% 比较dist(start,right1)和dist(start,right2),目的是更新chlid1 if dist(start,right1)<=dist(start,right2) %如果dist(start,right1)<=dist(start,right2),则说明start离right1更近 %则将child1的下一个位置更新为right1 child1(pos)=right1; else %如果dist(start,right1)>dist(start,right2),则说明start离right2更近 %则将child1的下一个位置更新为right2 child1(pos)=right2; end %% 更新parent1_c1和parent2_c1 parent1_c1(parent1_c1==start)=[]; parent2_c1(parent2_c1==start)=[]; %% 更新start start=child1(pos); %% 更新pos pos=pos+1;end
%% 终止条件是parent1_c2只剩一个点时parent1_c2=parent1; %复制parent1parent2_c2=parent2; %复制parent2start=start_c; %更新start为初始startpos=2; %记录child中的位置while numel(parent1_c2)~=1 start1=find(parent1_c2==start,1,'first'); %找出start在parent1的位置 start2=find(parent2_c2==start,1,'first'); %找出start在parent2的位置 %% 找出start在parent1_c1的左紧邻点 if start1==1 left1=parent1_c2(end); %如果start1是parent1第一个点,则start1左边紧邻的点是parent1最后一个点 else left1=parent1_c2(start1-1); %start1左边紧邻的点 end %% 找出start在parent2_c1的左紧邻点 if start2==1 left2=parent2_c2(end); %如果start2是parent2第一个点,则start2左边紧邻的点是parent2最后一个点 else left2=parent2_c2(start2-1); %start2左边紧邻的点 end %% 比较dist(left1,start)和dist(left2.start),目的是更新chlid2 if dist(left1,start)<=dist(left2,start) %如果dist(left1,start)<=dist(left2.start),则说明start离left1更近 %则将child2的下一个位置更新为left1 child2(pos)=left1; else %如果dist(left1,start)>dist(left2.start),则说明start离left2更近 %则将child2的下一个位置更新为left2 child2(pos)=left2; end %% 更新parent1_c2和parent2_c2 parent1_c2(parent1_c2==start)=[]; parent2_c2(parent2_c2==start)=[]; %% 更新start start=child2(pos); %% 更新pos pos=pos+1;endlen1=ObjFunction(child1,dist);len2=ObjFunction(child2,dist);[~,min_index]=min([len1,len2]);end


07 | 画图函数
画图函数PlotRoute的输入是路线 route、横纵坐标x和y,输出是最优路线图。
%%      @作者:随心390%      @微信公众号:优化算法交流地%%输入route:           路线%输入x,y:             x,y坐标function PlotRoute(route,x,y)route=[route route(1)];plot(x(route),y(route),'k-o','MarkerSize',10,'MarkerFaceColor','w','LineWidth',1.5);xlabel('x');ylabel('y');end

最后,各位小伙伴有什么疑问可点击 写留言 小程序为我们留言,我们会尽可能为各位解答。

点击下方小程序写留言
写留言

往期精选





知乎 | bilibili:随心390

长按识别二维码关注我们
觉得内容还不错的话,给我点个“在看”呗


本文分享自微信公众号 - 优化算法交流地(ROSECW123)。
如有侵权,请联系 support@oschina.cn 删除。
本文参与“OSC源创计划”,欢迎正在阅读的你也加入,一起分享。

展开阅读全文
打赏
0
0 收藏
分享
加载中
更多评论
打赏
0 评论
0 收藏
0
分享
在线直播报名
返回顶部
顶部