粒子群优化算法

%% 粒子群算法Particle Swarm Optimization
clc;
clear;
close all;
%% 背景
% 使用单个自变量进行寻优
%% 描述
%        种群数量：  粒子群算法的最大特点就是速度快，因此初始种群取50-1000都是可以的，虽然初始种群越大收敛性会更好，不过太大了也会影响速度；
%        迭代次数：  一般取100~4000，太少解不稳定，太多浪费时间。对于复杂问题，进化代数可以相应地提高；
%        惯性权重：  该参数反映了个体历史成绩对现在的影响，一般取0.5~1；
%        学习因子：  一般取0~4，此处要根据自变量的取值范围来定，并且学习因子分为个体和群体两种；
%        空间维数：  粒子搜索的空间维数即为自变量的个数。
%        位置限制：  限制粒子搜索的空间，即自变量的取值范围，对于无约束问题此处可以省略。
%        速度限制：  如果粒子飞行速度过快，很可能直接飞过最优解位置，但是如果飞行速度过慢，会使得收敛速度变慢，因此设置合理的速度限制就很有必要了。
% 最大值：32.1462
% 变量取值：18.3014
%% 初始化参数值
popsize_x=500;
popsize_y=1;            %列的个数即变量的个数
li_up=20;               %位置限制
li_do=0;
li_velocity=[-1,1];     %速度限制
w = 0.1;                        % 惯性权重,较大时可以全局寻优0
v=zeros(popsize_x,popsize_y);
iT_max=50;             %最大迭代次数
pop = init_pop(popsize_x,popsize_y,li_up,li_do);%初始化位置
c1=0.1;
c2=0.1;
%% 寻找新的粒子位置
for i=1:iT_max
    
    fitvalue=fitness(pop);                         %计算初始的适应度值；
    [individual_best,fitvalue1]=best(pop,fitvalue);           %寻找最佳个体； 
    history(:,i)=[individual_best,fitvalue1];
    history_best=max(history(1,:));
%% 速度更新和位置更新为核心。
    v=v * w/(10*i) + c1 * rand * (individual_best - pop) + c2 * rand * (repmat(history_best,popsize_x,popsize_y) - pop);% 速度更新,此公式关系整体收敛的效果
    pop=pop+v;                              % 位置更新
    
    [pop_x,pop_y]=size(pop);               %限制粒子的坐标范围
    for j=1:pop_y
        if any(pop(:,j)>li_up)==1
            pop(:,j)=li_up;
        elseif any(pop(:,j)<li_do)==1
            pop(:,j)=li_do;
        end
        pop(:,j)=pop(:,j);
    end
    diedai(i)=i;
end
%% 绘图
f=@(x)x .* sin(x) .* cos(2 * x) - 2 * x .* sin(3 * x);
subplot(1,2,1);
ezplot(f,[1,0.1,20]);
hold on;
plot(history_best,fitvalue1,'o');
title('结果展示图');
xlabel('变量');
ylabel('值');
subplot(1,2,2);
plot(diedai,history(1,:))
title('迭代收敛图');
xlabel('迭代次数');
ylabel('变量值');
%% 在迭代结果不稳定的情况时，且结果不准确时，一般都是由于粒子移动速度太快导致错过最优解

function pop = init_pop(popsize_x,popsize_y,li_up,li_do)
%UNTITLED6 此处显示有关此函数的摘要
%   此处显示详细说明
    pop=round(rand(popsize_x,popsize_y).*(li_up-li_do));
%     for i=1:popsize_x
%         for j=1:popsize_y
%             if pop(i,j)>li_up
%                 pop(i,j)=li_up;
%             elseif pop(i,j)<li_do
%                 pop(i,j)=li_do;
%             end
%             pop(i,j)=pop(i,j);
%         end
%     end
             
end

function [best_individual,fitness]=best(pop,fitvalue)
%% 选择最大值
    [ind_x,ind_y]=sort(fitvalue);
    best_individual=pop(ind_y(end),:);      %选择最大值
    fitness=ind_x(end);
%     best_individual=pop(ind_y(1),:);      %选择最小值
%     fitness=ind_x(1);
end

function fitvalue=fitness(x)
%% 计算粒子的适应度值
    fitvalue=x .* sin(x) .* cos(2 * x) - 2 * x .* sin(3 * x);
end