1.算法描述

在PSO中，群中的每个粒子表示为向量。在投资组合优化的背景下，这是一个权重向量，表示每个资产的分配资本。矢量转换为多维搜索空间中的位置。每个粒子也会记住它最好的历史位置。对于PSO的每次迭代，找到全局最优位置。这是群体中最好的最优位置。一旦找到全局最优位置，每个粒子都会更接近其局部最优位置和全局最优位置。当在多次迭代中执行时，该过程产生一个解决该问题的良好解决方案，因为粒子会聚在近似最优解上。

Ray等人通过将PSO算法和Pareto排序机制想结合起来。采用Pareto排序法来选择一组精英解，全局最优粒子的选择则是采用轮盘赌方式从中选择。实际运行时，只有少量的个体选择概率大，种群多样性保持不好。Coello等在PSO算法中选择群体最佳位置则是通过引入Pareto竞争机制和微粒知识库。该知识库用于存储微粒在每次飞行循环后的飞行经验，知识库的更新是通过考虑一个基于地理学的系统，该系统是就每个微粒的目标函数值而言来定义的。这个知识库被微粒用来确定一个指导搜索的领导者。同时非劣解的确定是通过将候选个体与从种群中随机选出的比较集进行比较，因此比较集的参数对算法成功与否有着至关重要的影响。若参数过大，则容易发生早熟收敛的现象，而参数过小，则种群中选出的非劣解的数量可能过少。

PSO模拟的是鸟群的捕食行为。设想场景：一群鸟在随机搜索食物。在这个区域里只有一块食物。所有鸟都不知道食物在哪里。但是他们知道当前的位置距离食物还有多远。那么找到食物的最优策略是什么？最简单有效的就是搜寻目前离食物最近的鸟的周围区域。

鸟群在整个搜寻过程中，通过相互传递各自的信息，让其他鸟知道自己当前的位置，通过这样的协作来判断自己找到的是不是最优解，同时也将最优解的信息传递给整个鸟群，最终整个鸟群都能聚集在食物源的周围，即找到了最优解。

PSO中，每个优化问题的解都是搜索空间的一只鸟，我们称之为“粒子”。所有的粒子都有一个被优化的函数决定的适应值，每个粒子还有一个速度决定他们飞翔的方向和距离。然后粒子们就追随当前的最优粒子在解空间中搜索。

PSO初始化为一群随机粒子（随机解）。然后通过迭代找到最优解，在每一次迭代中，粒子通过跟踪两个“极值”来更新自己。第一个就是粒子本身所找到的最优解，这个解叫做个体极值。另一个极值是整个种群目前找到的最优解，这个极值是全局机制。另外也可以不用整个种群而只是用其中一部分作为粒子的邻居，那么在所有邻居中的极值就是局部极值。
PSO初始化为一群随机粒子（随机解）。然后通过迭代找到最优解。在每一次迭代中，粒子通过跟踪两个“极值（pbest和gbest）”来更新自己。在找到这两个最优值后，粒子通过下面的公式来更新自己的速度和位置。

本仿真算例：

m：加工单元4个

n：工件10个

k：工序数

D: AGV 2个

表格为各工件工序的在加工单元上的加工时间（单位: min），0表示该机床无法加工此道工序

加工单元工件号	工序号	m1	m2	m3		m4
n1	k11	3	0	0		3
	k12	0	4	4		4
	k13	3	3	0		3
n2	k21	2	2	0		2
n2	k22	0	5	0		5
n3	k31	3	3	0		0
	k32	6	0	6		6
	k33	4	4	4		0
n4	k41	5	0	0		5
	k42	0	0	3		3
	k43	2	2	2		0
	k44	4	4	4		0
n5	k51	3	0	0		3
n5	k52	0	2	2		2
n6	k61	2	2	0		2
	k62	4	0	0		4
	k63	3	3	3		3
n7	k71	3	0	3		3
n7	k72	2	2	2		2
n8	k81	1	1	1		1
	k82	3	0	3		0
	k83	0	3	0		3
n9	k91	2	2	2		0
	k92	3	0	3		0
	k93	5	5	5		0
n10	k101	3	3		3	0
	k102	0	0		2	2
	k103	4	0		4	0
	k104	5	5		0	5

AGV搬运时间表（单位:min）

	待加工区	m1	m2	m3	m4	加工完成区
待加工区	0	1	0.5	0.5	1	0
m1	1	0	0.5	1	1.5	1
m2	0.5	0.5	0	0.5	1	0.5
m3	0.5	1	0.5	0	0.5	0.5
m4	1	1.5	1	0.5	0	1
加工完成区	0	1	0.5	0.5	1	0

2.仿真效果预览

matlab2022a仿真结果如下：

3.MATLAB核心程序

dim   = nvar; %变量数量         
Lmin =  lb;
Lmax =  ub;
c1    = 2.5;         %学习因子1
c2    = 2.5;         %学习因子2
w     = 1;         %惯性权重
Iters = 200;         %最大迭代次数
Num   = 400;           

%初始化种群的个体(可以在这里限定位置和速度的范围)
for i = 1:nvar
    for j = 1:Num
        x(j,i) = Lmin(i) + (Lmax(i)-Lmin(i))*rand(1,1);  %随机初始化位置
        v(j,i) = Lmin(i) + (Lmax(i)-Lmin(i))*rand(1,1);  %随机初始化速度
    end
end
y     = [];

%先计算各个粒子的适应度，并初始化Pi和Pg
for i=1:Num
    i
    [p(i),z1,z2,z3] = func_F(x(i,:));
    y(i,:)          = x(i,:);
end
%全局最优
pg = x(1,:);             

for i=2:Num
    i
    [pa(i),z1,z2,z3] = func_F(x(i,:));
    pb(i)            = func_F(pg);
    if pa(i) < pb(i)
       pg=x(i,:);
    end
end

for t=1:Iters
    t
    for i=1:Num
        v(i,:)           = w*v(i,:)+c1*rand*(y(i,:)-x(i,:))+c2*rand*(pg-x(i,:));
        x(i,:)           = x(i,:)+v(i,:);
        [pa(i),z1,z2,z3] = func_F(x(i,:));
        if pa(i)<p(i)
           p(i)  = pa(i);
           y(i,:)= x(i,:);
        end

        [pb(i),z1,z2,z3] = func_F(pg);
        if p(i)<pb(i)
           pg=y(i,:);
        end
        
 
    [Pbest(i),z1_(i),z2_(i),z3_(i)] = func_F(pg);    
    end
    Fs(t)=mean(Pbest);
    Y1       = z1_;
    Y2       = z2_;
    Y3       = z3_;
    tmps     = pg;
end

figure;
plot(Fs);
xlabel('PSO迭代次数');
ylabel('fitness');
grid on

pg
 

II=15;
Total_time = Y1(II);
Total_WT   = Y2(II);
max_WT     = Y3(II);
X          = tmps;
% Z=fitness(X)
%显示甘特图
[V,I] = sort(X(1:end));
 
seq1 = sort(I(1:N1));
seq2 = sort(I(N1+1:N1+N2));
seq3 = sort(I(N1+N2+1:N1+N2+N3));
seq4 = sort(I(N1+N2+N3+1:N1+N2+N3+N4));%选择了任务3的计算机器
seq5 = sort(I(N1+N2+N3+N4+1:N1+N2+N3+N4+N5));%选择了任务3的计算机器

x1   = mach(seq1);
x2   = mach(seq2);
x3   = mach(seq3);
x4   = mach(seq4);
x5   = mach(seq5);


for i = 1:N1
    c1(1,i) = time1(x1(i)); 
end
for i = 1:N2
    c2(1,i) = time2(x2(i)); 
end
for i = 1:N3
    c3(1,i) = time3(x3(i)); 
end
for i = 1:N4
    c4(1,i) = time4(x4(i)); 
end
for i = 1:N5
    c5(1,i) = time5(x5(i)); 
end

xx  = [x1,x2,x3,x4,x5]
seq = [seq1,seq2,seq3,seq4,seq5] 
cc  = [c1,c2,c3,c4,c5] 


 
%[x y w h]
LEN = N1+N2+N3+N4+N5;
rec1(:,LEN)=[0,0,0,0]';
rec2(:,LEN)=[0,0,0,0]';
rec3(:,LEN)=[0,0,0,0]';
rec4(:,LEN)=[0,0,0,0]';
rec5(:,LEN)=[0,0,0,0]';


cnt1=0;
cnt2=0;
cnt3=0;
cnt4=0;
cnt5=0;
color1=[];
color2=[];
color3=[];
color4=[];
color5=[];
Timesave1=[];
Timesave2=[];
Timesave3=[];
Timesave4=[];
Timesave5=[];
P1=[];
P2=[];
P3=[];
P4=[];
P5=[];
for i = 1:LEN
    idx = find(seq == i);
    time= cc(idx);
    if xx(idx)==1;%machine1
        cnt1=cnt1+1;
        cnt2=cnt2;
        cnt3=cnt3;
        cnt4=cnt4;
        cnt5=cnt5;
        if idx<=N1
            color1=[color1,'c'];
            color2=[color2];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P1(:,cnt1)=[1;1];
        end
        if idx>N1&idx<=N1+N2
            color1=[color1,'g'];
            color2=[color2];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P1(:,cnt1)=[1;2];
        end
        
        
        if idx>N1+N2&idx<=N1+N2+N3
            color1=[color1,'y'];
            color2=[color2];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P1(:,cnt1)=[1;3];
        end
        
        if idx>N1+N2+N3&idx<=N1+N2+N3+N4
            color1=[color1,'r'];
            color2=[color2];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P1(:,cnt1)=[1;4];
        end
        
        if idx>N1+N2+N3+N4&idx<=N1+N2+N3+N4+N5
            color1=[color1,'b'];
            color2=[color2];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P1(:,cnt1)=[1;5];
        end
        
        if cnt1==1
           rec1(1,cnt1) = 0;
           rec1(2,cnt1) = 1-0.3;
           rec1(3,cnt1) = cc(idx);
           rec1(4,cnt1) = 0.6; 
        else
           rec1(1,cnt1) = rec1(1,cnt1-1)+rec1(3,cnt1-1);
           rec1(2,cnt1) = 1-0.3;
           rec1(3,cnt1) = cc(idx);
           rec1(4,cnt1) = 0.6; 
        end
        Timesave1(cnt1)=cc(idx);
        
    end
    
    
    
    if xx(idx)==2;%machine2
        cnt1=cnt1;
        cnt2=cnt2+1;
        cnt3=cnt3;
        cnt4=cnt4;
        cnt5=cnt5;
        if idx<=N1
            color1=[color1];
            color2=[color2,'c'];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P2(:,cnt2)=[2;1];
        end
        if idx>N1&idx<=N1+N2
            color1=[color1];
            color2=[color2,'g'];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P2(:,cnt2)=[2;2];
        end
        if idx>N1+N2&idx<=N1+N2+N3
            color1=[color1];
            color2=[color2,'y'];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P2(:,cnt2)=[2;3];
        end
        if idx>N1+N2+N3&idx<=N1+N2+N3+N4
            color1=[color1];
            color2=[color2,'r'];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P2(:,cnt2)=[2;4];
        end
        if idx>N1+N2+N3+N4&idx<=N1+N2+N3+N4+N5
            color1=[color1];
            color2=[color2,'b'];
            color3=[color3];
            color4=[color4];
            color5=[color5];
            P2(:,cnt2)=[2;5];
        end
        
        if cnt2==1
           rec2(1,cnt2) = 0;
           rec2(2,cnt2) = 2-0.3;
           rec2(3,cnt2) = cc(idx);
           rec2(4,cnt2) = 0.6; 
        else
           rec2(1,cnt2) = rec2(1,cnt2-1)+rec2(3,cnt2-1);
           rec2(2,cnt2) = 2-0.3;
           rec2(3,cnt2) = cc(idx);
           rec2(4,cnt2) = 0.6; 
        end
        Timesave2(cnt2)=cc(idx);
    end
.........................................................
        Timesave5(cnt5)=cc(idx);
    end  
end
02_103m