基于matlab的多目标粒子群算法（开源）

代码仅供参考，代码如下：

1.基于最小角度的gbest引导；

2.基于概率和支配关系的pbest引导；

3.基于支配关系和粒子密度的最优集选择；

4.增加克隆变异，对处于目标函数最值得粒子进行克隆与变异，当找到比最值点还要小或者大的点时，替换原来的点；（可以注释这部分）

（代码仅供参考，能运行，如有错误，欢迎留言）

xlimit行数要与适应度函数变量数目相等！！！！！！！！！！！！

结果说明：[变量1,变量2...变量n,适应度函数1,适应度函数2...,适应度函数m];

-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------2022/04/13更新----------------------------------------------------------

clc;
close all;
format long

%%%%%%%%%%%%%%%%%%%多目标粒子群算法%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%多目标多变量%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%基于最小角度的全局引导，基于概率的个体引导%%%%%%%%%%%%%%%
%%%%%基于支配关系和粒子密度的A集更新，基于克隆变异的多样性更新%%%%%%%

%% 初始设置
%1. --决策空间的设置-- %
xlimit = [-100,100];
% xlimit = [-5,5;-5,5];   %% 确定位置变量范围，维数与变量个数相等，左小右边大变量名称从左往右为
% xlimit = [0,pi/2;0,pi/2];
% xlimit = [1,1.6;0.3,0.8];
% xlimit = [75,95;4,6;41.5,49];
% xlimit = [1,1.6;0.3,0.8;0.2,0.7;-2.3,2.6];
% xlimit = [1,1.6;0.3,0.8;0.2,0.7;-4,2.6;-5.1,10.25];
c2 = 2;                            %%社会学习系数
c1 = 0.75;                         %%自我学习系数
w = 0.2;                           %%惯性系数
pointnum = 4;                      %几位小数运算
maxgen = 80;                       %迭代步数     100-1000
popsize = 600;                     %种群数    
vrange = 0.3;                      %%速度范围                     
achievemaxnum = 200;               %achieve集最大个数
copymode = 1;                      %是否变异
copynum = 500;                     %克隆个数
colrag = [0.01,0.35];               %变异范围
popdimension = size(xlimit,1);     %变量维数目
vlimit = vrange*[-xlimit(:,2)+xlimit(:,1),xlimit(:,2)-xlimit(:,1)];        %0.1-0.5

%% 2. --种群位置与速度初始化-- %
pop = repmat(xlimit(:,1)',popsize,1)+repmat(diff(xlimit'),popsize,1).*rand(popsize,popdimension);  %一行为一个点
v = repmat(vlimit(:,1)',popsize,1)+repmat(diff(vlimit'),popsize,1).*rand(popsize,popdimension);

%% 3. --计算初始种群的适应度值-- %
fitness = round(fcn(pop)*10^pointnum)/(10^pointnum);       %行向量
fcnnum = size(fitness,2);  %目标函数个数

%% 4. --初始pbest赋值-- %
pbest = pop;
pbestfitness = fitness;
tic;

%% 5. --外部归档集archive初始化-- %
archive = [];
archivefitness = [];

%% 6. --A集更新-- %
for i=1:1:popsize   %每次只考虑将一个粒子加到外部归档集
    %{
    输入-pop、popsize、fitness、archive、archivefitness、achievemaxnum
    输出 archive archivefitness
    %}
    if isempty(archive)       %第一个粒子直接存入archive
        archive(1,:) = pop(i,:);
        archivefitness(1,:) = fitness(i,:);
    else
        bijiaojvzhen = repmat(fitness(i,:),size(archivefitness,1),1);
        bigcop = bijiaojvzhen >= archivefitness;
        equcop = bijiaojvzhen == archivefitness;
        smallcop = bijiaojvzhen <= archivefitness;
        maybedeleterow = all(bigcop > 0,2);   %全大于等于0的行为1
        notdeleterow = all(equcop > 0,2);     %全等于0的行为1
        dominaterow = all(smallcop>0,2);         %被A集支配的行为1
        todeleterow = logical(maybedeleterow - notdeleterow);   %要删除的行为1
        if sum(dominaterow,'all')>0                                        %pop被A集支配
            archive(todeleterow,:) = [];                                    %删除被pop支配A集粒子
            archivefitness(todeleterow,:) = [];                             %删除被pop支配A集粒子
        else                                                                %pop不能被A集支配
            archive(todeleterow,:) = [];                                    %删除被pop支配A集粒子
            archivefitness(todeleterow,:) = [];                             %删除被pop支配A集粒子
            extremum = [min(archivefitness,[],1);...                        %算出适应度端点
                max(archivefitness,[],1)];
            if size(archive,1) < 2                                          %如果空了就直接加进去
                archive(end+1,:) = pop(i,:);
                archivefitness(end+1,:) = fitness(i,:);
            elseif any((fitness(i,:)>extremum(2,:))>0,'all') ||...
                    any((fitness(i,:)<extremum(1,:))>0,'all')               %如果新的点不在旧的区间直接加进去
                archive(end+1,:) = pop(i,:);
                archivefitness(end+1,:) = fitness(i,:);
            else                                                            %如果新的点在旧的区间网格更新
                archivefitnessdist = diff(extremum,1,1);                    %适应度距离 (1行向量)
                elementdist = archivefitnessdist/achievemaxnum;             %单元距离 （1行向量）
                updatacol = (elementdist > 0);                              %对区间大于0的目标
                popmeshnum = zeros(1,fcnnum);                               %粒子网格号码（1行向量）
                archivemeshnum = zeros(size(archive,1),fcnnum);             %A集网格号码（多行向量）
                popmeshnum(:,updatacol) = ceil((fitness(i,updatacol)-extremum(1,updatacol))...
                    ./elementdist(:,updatacol));                            %更新粒子网格号码（1行向量）
                archivemeshnum(:,updatacol) = ceil((archivefitness(:,updatacol)-extremum(1,updatacol))...
                    ./elementdist(:,updatacol));                            %更新A集网格号码（多行向量）
                archivemeshnumuiq = unique(archivemeshnum,'rows');          %A集唯一网格号码（多行向量）
                if ~ismember(popmeshnum,archivemeshnumuiq)                   %如果号码没有粒子，加入
                    archive(end+1,:) = pop(i,:);
                    archivefitness(end+1,:) = fitness(i,:);
                else                                                         %如果号码有粒子了
                    reshapearchivemeshnumuiq = reshape(archivemeshnumuiq.',1,...
                        size(archivemeshnumuiq,2),size(archivemeshnumuiq,1));    %A集唯一网格号码转为三维矩阵（1行）
                    reparchivemeshnum = repmat(archivemeshnum,1,1,...
                        size(archivemeshnumuiq,1));                              %A集网格号码转为三维矩阵&#