基于kmeans聚类和遗传优化算法的TSP问题求解matlab仿真

目录

1. K-means聚类

2.遗传算法求解TSP

3.MATLAB程序

4.仿真结果

---

       旅行商问题（Traveling Salesman Problem, TSP）是一个典型的组合优化问题，旨在寻找遍历所有城市并返回原点的最短路径。

1. K-means聚类

       K-means聚类算法是一种广泛应用的无监督学习方法，主要用于数据分类，通过迭代过程将数据集划分为预设数量的聚类（或类别）。

      首先，我们引入K-means聚类算法对城市进行分组，减少问题的规模。设有一组城市坐标C={c1​,c2​,...,cn​}，其中ci​∈R2。K-means的目标是将城市划分为K 个簇，使得每个簇内部的城市尽可能接近。K-means的主要步骤如下：

初始化：随机选择K 个中心点M={m1​,m2​,...,mK​}。

迭代：重复分配和更新中心操作，直至簇中心不再显著变化或达到预设的最大迭代次数。分配：对于每个城市ci​，计算其与各个簇中心的距离，将其分配到最近的簇。更新中心：每个簇的新中心为该簇内所有城市的坐标平均值。

2.遗传算法求解TSP

       接下来，我们将使用遗传算法来求解简化后的TSP问题，这里的简化在于我们先用K-means将城市聚类，然后在每个簇内部建立子TSP问题。遗传算法（Genetic Algorithm, GA）是一种模拟自然进化和生物遗传过程的优化算法，它通过模拟自然选择、遗传重组、突变等生物学现象，搜索问题的潜在解空间，寻求最优解或近似最优解。

编码方案：

设有一个染色体P={p1​,p2​,...,pK​}，其中pk​ 是第k 个簇中城市的排列。

对于每个簇，我们可以用一个排列（顺序）来表示城市访问顺序，例如，对于K 个簇，我们可以得到K 个排列，组成一个染色体。

适应度函数：

对于每个子TSP问题，计算其哈密尔顿回路的总距离，作为适应度函数值。总适应度可以通过将各个子问题的适应度相加得到。

其中D(pk​) 是第k 个簇中城市按照排列pk​ 访问的总距离。

遗传操作：

选择：根据适应度函数值，选择优秀的染色体进入下一代。
交叉：对两个父代染色体执行部分匹配交叉（Order Crossover, OX）或其它适合的交叉操作，产生新的子代染色体。
变异：随机改变某些子TSP排列中的城市访问顺序，保持多样性。
精英保留：确保每一代中最优秀的一些个体直接传递到下一代。
终止条件：

达到预设的最大迭代次数。
连续若干代没有显著的适应度改善。

        综合以上两部分，首先应用K-means对城市进行聚类，然后用遗传算法针对每个簇内的子TSP问题寻找近似最优解。最后，将各个簇内部的最优路径拼接成全局路径，得到最终的旅行商问题解决方案。

3.MATLAB程序

...........................................................................

if jcpoint1<jcpoint2
for k = 1:jcpoint1
    Route_best2(k,:)=Route_best(k,:);
end
 if jcpoint3 < jcpoint4%没问题
    for kk = jcpoint3:-1:1
    Route_best2((jcpoint1+jcpoint3-kk+1),:)=date3(optRoute3(kk),:);
    end
    for  kkk = n3:-1:jcpoint4
    Route_best2((jcpoint1+jcpoint3+n3-kkk+1),:)=date3(optRoute3(kkk),:);
    end
    for  kkkkk = jcpoint2 : (n1+n2)
    Route_best2((n3+kkkkk),:)=Route_best(kkkkk,:);
    end
 else     %没问题
    for k = jcpoint3 : n3
    Route_best2((jcpoint1+k-jcpoint3+1),:)=date3(optRoute3(k),:);
    end
     for k = 1 : jcpoint4
    Route_best2((jcpoint1+n3-jcpoint4+k),:)=date3(optRoute3(k),:);
     end
    for k = jcpoint2 : n1+n2
    Route_best2((n3+k),:)=Route_best(k,:);
    end
 end
end


if jcpoint1>jcpoint2
for k = 1:jcpoint2
    Route_best2(k,:)=Route_best(k,:);
end
 if jcpoint3 < jcpoint4
    for k = jcpoint4 : n3
    Route_best2((jcpoint2+k-jcpoint4+1),:)=date3(optRoute3(k),:);
    end
     for  k = 1 : jcpoint3
    Route_best2((jcpoint2+n3+k-jcpoint3-jcpoint4+1),:)=date3(optRoute3(k),:);
     end
     for  k = jcpoint1 :(n1+n2)
    Route_best2((n3+k),:)=date1(optRoute1(k),:);
    end
  else
    for k = jcpoint4:-1 : 1
        Route_best2((jcpoint2+jcpoint4-k+1),:)=date3(optRoute3(k),:);
    end
     for k = n3 : -1 :jcpoint3
        Route_best2((jcpoint2+jcpoint4+n3-k+1),:)=date3(optRoute3(k),:);
     end
    for k = jcpoint1 : n1+n2
    Route_best2((n3+k),:)=Route_best2(k,:);
    end
 end
end

d = sqrt((Route_best2(1,1)-(Route_best2(R,1)))^2+(Route_best2(1,2)-(Route_best2(R,2)))^2);
        for k = 2:R
            d = d + sqrt((Route_best2(k-1,1)-(Route_best2(k,1)))^2+(Route_best2(k-1,2)-(Route_best2(k,2)))^2);  %交叉算子、变异算子
        end

          figure(pfig);
          rte = 1:R;
          rte(R+1)=1;
          plot(Route_best2(rte,1),Route_best2(rte,2),'r.-');
          title(sprintf('Total Distance = %1.4f, Iteration = %d',d,iter));
   
% 画图
if showResult
    % Plots the GA Results
    figure('Name','TSP_GA | Results','Numbertitle','off');
    % subplot(2,2,1);
    figure(1),
    pclr = ~get(0,'DefaultAxesColor');
    if dims > 2, 
        plot3(date1(:,1),date1(:,2),date1(:,3),'.','Color',pclr);
    else
        plot(date1(:,1),date1(:,2),'.','Color',pclr)
        hold on
        plot(date3(:,1),date3(:,2),'.','Color',pclr)
        hold on
        plot(date3(:,1),date3(:,2),'.','Color',pclr)
    end
    title('城市位置');
    grid on
%     subplot(2,2,2);
    figure(3),
    rte1 = optRoute1([1:n1 1]);
    rte2 = optRoute2([1:n2 1]);
    rte3 = optRoute3([1:n3 1]);
    if dims > 2, plot3(date1(rte,1),date1(rte,2),date1(rte,3),'r.-');
    else
        plot(date1(rte1,1),date1(rte1,2),'r.-')
        hold on
        plot(date2(rte2,1),date2(rte2,2),'g.-')
        hold on
        plot(date3(rte3,1),date3(rte3,2),'b.-')
        hold on
    for i=1:R 
    if pos2(i,3)==1   
         plot(pos2(i,1),pos2(i,2),'ro'); 
    elseif pos2(i,3)==2
         plot(pos2(i,1),pos2(i,2),'go'); 
    else 
         plot(pos2(i,1),pos2(i,2),'bo'); 
    end
    end
    end
    title(sprintf('最短距离 = %1.4f',minDist1,minDist2,minDist3));
    grid on
%     subplot(2,2,4);
    for k=1:Iter
        distHistory(k)=distHistory1(k)+distHistory2(k)+distHistory3(k);
    end
    figure(4),
    plot(distHistory,'b','LineWidth',2);
    hold on
    title('Best fitness curve');
    grid on
  
end
up4064
 

4.仿真结果