一、粒子群演算法的概述

    粒子群演算法(PSO)屬於群智慧演算法的一種，是通過模擬鳥群捕食行為設計的。假設區域裡就只有一塊食物（即通常優化問題中所講的最優解），鳥群的任務是找到這個食物源。鳥群在整個搜尋的過程中，通過相互傳遞各自的資訊，讓其他的鳥知道自己的位置，通過這樣的協作，來判斷自己找到的是不是最優解，同時也將最優解的資訊傳遞給整個鳥群，最終，整個鳥群都能聚集在食物源周圍，即我們所說的找到了最優解，即問題收斂。

二、粒子群演算法的流程

    粒子群演算法通過設計一種無質量的粒子來模擬鳥群中的鳥，粒子僅具有兩個屬性：速度和位置，速度代表移動的快慢，位置代表移動的方向。每個粒子在搜尋空間中單獨的搜尋最優解，並將其記為當前個體極值

（PSO流程）

下面我們具體解釋下流程圖裡面的每一個步驟：

1、初始化

   首先，我們需要設定最大的速度區間，防止超出最大的區間。位置資訊即為整個搜尋空間，我們在速度區間和搜尋空間上隨機初始化速度和位置。設定群體規模

2、個體極值與全域性最優解

   個體極值為每個粒子找到的歷史上最優的位置資訊，並從這些個體歷史最優解中找到一個全域性最優解，並與歷史最優解比較，選出最佳的作為當前的歷史最優解。

3、更新速度和位置的公式

   更新公式為：

其中，稱為慣性因子，和稱為加速常數，一般取。表示區間上的隨機數。表示第個變數的個體極值的第維。表示全域性最優解的第維。

4、終止條件

有兩種終止條件可以選擇，一是最大代數：；二是相鄰兩代之間的偏差在一個指定的範圍內即停止。我們在實驗中選擇第一種。

三、實驗

    我們選擇的測試函式是：Griewank。其基本形式如下：

影象為：

（Griewank函式影象）

在實驗中我們選擇的維數是20；MATLAB程式程式碼如下：

主程式：

1. c1=2;%學習因子

2. c2=2;%學習因子

3. Dimension=20;

4. Size=30;

5. Tmax=500;

6. Velocity_max=1200;%粒子最大速度

7. F_n=2;%測試函式名

8. Fun_Ub=600;%函式上下界

9. Fun_Lb=-600;

10. Position=zeros(Dimension,Size);%粒子位置

11. Velocity=zeros(Dimension,Size);%粒子速度

12. Vmax(1:Dimension)=Velocity_max;%粒子速度上下界

13. Vmin(1:Dimension)=-Velocity_max;

14. Xmax(1:Dimension)=Fun_Ub;%粒子位置上下界，即函式自變數的上下界

15. Xmin(1:Dimension)=Fun_Lb;

16. [Position,Velocity]=Initial_position_velocity(Dimension,Size,Xmax,Xmin,Vmax,Vmin);

17. Pbest_position=Position;%粒子的歷史最優位置，初始值為粒子的起始位置，儲存每個粒子的歷史最優位置

18. Gbest_position=zeros(Dimension,1);%全域性最優的那個粒子所在位置，初始值認為是第1個粒子

19. for j=1:Size

20. Pos=Position(:,j);%取第j列，即第j個粒子的位置

21. fz(j)=Fitness_Function(Pos,F_n,Dimension);%計算第j個粒子的適應值

22. end

23. [Gbest_Fitness,I]=min(fz);%求出所有適應值中最小的那個適應值，並獲得該粒子的位置

24. Gbest_position=Position(:,I);%取最小適應值的那個粒子的位置，即I列

25. for itrtn=1:Tmax

26. time(itrtn)=itrtn;

27. Weight=1;

28. r1=rand(1);

29. r2=rand(1);

30. for i=1:Size

31. Velocity(:,i)=Weight*Velocity(:,i)+c1*r1*(Pbest_position(:,i)-Position(:,i))+c2*r2*(Gbest_position-Position(:,i));

32. end

33. %限制速度邊界

34. for i=1:Size

35. for row=1:Dimension

36. if Velocity(row,i)>Vmax(row)

37. Veloctity(row,i)=Vmax(row);

38. elseif Velocity(row,i)<Vmin(row)

39. Veloctity(row,i)=Vmin(row);

40. else

41. end

42. end

43. end

44. Position=Position+Velocity;

45. %限制位置邊界

46. for i=1:Size

47. for row=1:Dimension

48. if Position(row,i)>Xmax(row)

49. Position(row,i)=Xmax(row);

50. elseif Position(row,i)<Xmin(row)

51. Position(row,i)=Xmin(row);

52. else

53. end

54. end

55. end

56. for j=1:Size

57. P_position=Position(:,j)';%取一個粒子的位置

58. fitness_p(j)=Fitness_Function(P_position,F_n,Dimension);

59. if fitness_p(j)< fz(j) %粒子的適應值比運動之前的適應值要好，更新原來的適應值

60. Pbest_position(:,j)=Position(:,j);

61. fz(j)=fitness_p(j);

62. end

63. if fitness_p(j)<Gbest_Fitness

64. Gbest_Fitness=fitness_p(j);

65. end

66. end

67. [Gbest_Fitness_new,I]=min(fz);%更新後的所有粒子的適應值，取最小的那個，並求出其編號

68. Best_fitness(itrtn)=Gbest_Fitness_new; %記錄每一代的最好適應值

69. Gbest_position=Pbest_position(:,I);%最好適應值對應的個體所在位置

70. end

71. plot(time,Best_fitness);

72. xlabel('迭代的次數');ylabel('適應度值P_g');

初始化：

1. function [Position,Velocity] = Initial_position_velocity(Dimension,Size,Xmax,Xmin,Vmax,Vmin)

2. for i=1:Dimension

3. Position(i,:)=Xmin(i)+(Xmax(i)-Xmin(i))*rand(1,Size); % 產生合理範圍內的隨機位置，rand(1,Size)用於產生一行Size個隨機數

4. Velocity(i,:)=Vmin(i)+(Vmax(i)-Vmin(i))*rand(1,Size);

5. end

6. end

適應值計算：

1. function Fitness=Fitness_Function(Pos,F_n,Dimension)

2. switch F_n

3. case 1

4. Func_Sphere=Pos(:)'*Pos(:);

5. Fitness=Func_Sphere;

6. case 2

7. res1=Pos(:)'*Pos(:)/4000;

8. res2=1;

9. for row=1:Dimension

10. res2=res2*cos(Pos(row)/sqrt(row));

11. end

12. Func_Griewank=res1-res2+1;

13. Fitness=Func_Griewank;

14. end

最終的收斂曲線：

（收斂曲線）