多目标粒子群算法matlab

最新推荐文章于 2024-09-05 07:30:46 发布

电磁MATLAB

最新推荐文章于 2024-09-05 07:30:46 发布

阅读量2.5k

点赞数 2

文章标签：算法 matlab 机器学习

本文链接：https://blog.csdn.net/weixin_47365903/article/details/130259727

版权

多目标粒子群算法是一种结合PSO和多目标优化的算法，用于解决复杂优化问题，特别是多目标问题。该算法通过加权平均处理多个目标，避免局部最优解。在Matlab中，算法通过初始化粒子群、计算适应度值、更新速度和位置等步骤进行迭代优化。程序示例展示了如何在Matlab中实现这一过程，包括目标函数计算、粒子更新和种群管理等关键步骤。

摘要由CSDN通过智能技术生成

多目标粒子群算法matlab
1 算法介绍
多目标粒子群算法是一种常用于多目标优化问题的启发式算法。该算法结合了粒子群算法（PSO）和多目标优化技术，能够在不同的目标函数中寻找最佳解，解决了单目标优化算法无法解决的问题，同时也可以避免在多目标问题中出现的局部最优解的问题。
多目标粒子群算法在处理多目标问题时，需要将多个目标融合，通过加权平均方式提高粒子的全局最优性。具体来说，该算法在每代中更新所有粒子的位置和速度，将每个粒子看作是在目标函数空间中搜索的一个点。扩展标准粒子群算法，每个粒子都有多个的位置，对应了多目标函数的结果。每个个体根据当前位置和速度计算适应度值来调整粒子位置，并更新粒子速度。此外，在更新过程中，通过粒子之间的交流来分享信息，以帮助寻找更好的解。
多目标粒子群算法的主要流程可以总结为以下几个步骤：

1.初始化粒子群——随机初始化一定数量的粒子位置和速度，可以采用随机分布函数。同时，设置最大迭代次数和收敛容差。
2.计算适应度值——对于每个粒子位置，按照各个目标函数计算相应的适应度值，并将其存储。
3.更新粒子速度和位置——根据每个粒子的当前位置和速度，通过标准粒子群算法的迭代公式计算更新后的速度和位置。
4.更新最优解——根据每个粒子的当前位置和目标函数值，更新全局最优解。
5.结束迭代——如果达到最大迭代次数或收敛容差，则停止算法运行，输出最优解。
在这里插入图片描述

多目标粒子群算法具有的优点包括全局搜索能力强、易于实现、收敛速度较快、适用于多种类型的问题等。但是，该算法也存在一些问题，如初始参数的敏感性、容易陷入局部最优解等。为了克服这些问题，通常需要结合其他优化算法，如模拟退火和遗传算法等。
总之，多目标粒子群算法是一种强大的优化算法，不仅可以用于解决多目标问题，而且具有可扩展性和灵活性。

2 算例
在这里插入图片描述

3 程序运行结果
在这里插入图片描述

4 matlab程序
1）主函数

clc;
clear;
close all;

%% 问题定义
CostFunction=@(x)fun_obj(x);
nVar=50;%变量个数
VarMin=0;%变量下限
VarMax=1;%变量上限
VarSize=[1 nVar];
VelMax=(VarMax-VarMin)/10;

%% MOPSO 设置
nPop=300; % Population Size
nRep=60; % Repository Size
MaxIt=200; % Maximum Number of Iterations
phi1=2.05;%学习因子
phi2=2.05;
phi=phi1+phi2;
chi=2/(phi-2+sqrt(phi^2-4phi));
w=chi; % Inertia Weight
wdamp=1; % Inertia Weight Damping Ratio
c1=chiphi1; % Personal Learning Coefficient
c2=chi*phi2; % Global Learning Coefficient
alpha=0.1; % Grid Inflation Parameter
nGrid=10; % Number of Grids per each Dimension
beta=4; % Leader Selection Pressure Parameter
gamma=2; % Extra (to be deleted) Repository Member Selection Pressure

%% 初始化
particle=CreateEmptyParticle(nPop);
for i=1:nPop
particle(i).Velocity=0;
particle(i).Position=unifrnd(VarMin,VarMax,VarSize); %unifrnd在[VarMin,VarMax]区间模拟VarSize的数组
% 初始化了自变量X的位置
particle(i).Cost=CostFunction(particle(i).Position); % 求解目标函数
particle(i).Best.Position=particle(i).Position;
particle(i).Best.Cost=particle(i).Cost;
end
particle=DetermineDomination(particle);
rep=GetNonDominatedParticles(particle);
rep_costs=GetCosts(rep);
G=CreateHypercubes(rep_costs,nGrid,alpha);
for i=1:numel(rep)
[rep(i).GridIndex rep(i).GridSubIndex]=GetGridIndex(rep(i),G);
end

%% MOPSO 主循环
for it=1:MaxIt
for i=1:nPop
rep_h=SelectLeader(rep,beta);

    particle(i).Velocity=w*particle(i).Velocity ...
                         +c1*rand*(particle(i).Best.Position - particle(i).Position) ...
                         +c2*rand*(rep_h.Position -  particle(i).Position);

    particle(i).Velocity=min(max(particle(i).Velocity,-VelMax),+VelMax);

    particle(i).Position=particle(i).Position + particle(i).Velocity;

    flag=(particle(i).Position<VarMin | particle(i).Position>VarMax);    %非劣检查:flag
    particle(i).Velocity(flag)=-particle(i).Velocity(flag);
    
    particle(i).Position=min(max(particle(i).Position,VarMin),VarMax);

    particle(i).Cost=CostFunction(particle(i).Position);

    if Dominates(particle(i),particle(i).Best)
        particle(i).Best.Position=particle(i).Position;
        particle(i).Best.Cost=particle(i).Cost;
        
    elseif ~Dominates(particle(i).Best,particle(i))
        if rand<0.5
            particle(i).Best.Position=particle(i).Position; 
            particle(i).Best.Cost=particle(i).Cost;
        end
    end

end

particle=DetermineDomination(particle);
nd_particle=GetNonDominatedParticles(particle);

rep=[rep
     nd_particle];

rep=DetermineDomination(rep);
rep=GetNonDominatedParticles(rep);

for i=1:numel(rep)
    [rep(i).GridIndex rep(i).GridSubIndex]=GetGridIndex(rep(i),G);
end

if numel(rep)>nRep
    EXTRA=numel(rep)-nRep;
    rep=DeleteFromRep(rep,EXTRA,gamma);
    
    rep_costs=GetCosts(rep);
    G=CreateHypercubes(rep_costs,nGrid,alpha);
    
end

% disp(['Iteration ’ num2str(it) ': Number of Repository Particles = ’ num2str(numel(rep))]);

w=w*wdamp;
plot(rep_costs(1,:),rep_costs(2,:),'rx');
xlabel('f1');
ylabel('f2');
drawnow

end

%% 结果

2）子函数

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  MATLAB Code for                                                  %
%                                                                   %
%  Multi-Objective Particle Swarm Optimization (MOPSO)              %
%  Version 1.0 - Feb. 2011                                          %
%                                                                   %
%  According to:                                                    %
%  Carlos A. Coello Coello et al.,                                  %
%  "Handling Multiple Objectives with Particle Swarm Optimization," %
%  IEEE Transactions on Evolutionary Computation, Vol. 8, No. 3,    %
%  pp. 256-279, June 2004.                                          %
%                                                                   %
%  Developed Using MATLAB R2009b (Version 7.9)                      %
%                                                                   %
%  Programmed By: S. Mostapha Kalami Heris                          %
%                                                                   %
%         e-Mail: sm.kalami@gmail.com                               %
%                 kalami@ee.kntu.ac.ir                              %
%                                                                   %
%       Homepage: http://www.kalami.ir                              %
%                                                                   %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function rep_h=SelectLeader(rep,beta)
    if nargin<2
        beta=1;
    end

    [occ_cell_index occ_cell_member_count]=GetOccupiedCells(rep);
    
    p=occ_cell_member_count.^(-beta);
    p=p/sum(p);
    
    selected_cell_index=occ_cell_index(RouletteWheelSelection(p));
    
    GridIndices=[rep.GridIndex];
    
    selected_cell_members=find(GridIndices==selected_cell_index);
    
    n=numel(selected_cell_members);
    
    selected_memebr_index=randi([1 n]);
    
    h=selected_cell_members(selected_memebr_index);
    
    rep_h=rep(h);
end


function i=RouletteWheelSelection(p)

    r=rand;
    c=cumsum(p);
    i=find(r<=c,1,'first');

end


function costs=GetCosts(pop)

    nobj=numel(pop(1).Cost);
    costs=reshape([pop.Cost],nobj,[]);

end


function [Index SubIndex]=GetGridIndex(particle,G)

    c=particle.Cost;
    
    nobj=numel(c);
    ngrid=numel(G(1).Upper);
    
    str=['sub2ind(' mat2str(ones(1,nobj)*ngrid)];

    SubIndex=zeros(1,nobj);
    for j=1:nobj
        
        U=G(j).Upper;
        
        i=find(c(j)<U,1,'first');
        
        SubIndex(j)=i;
        
        str=[str ',' num2str(i)];
    end
    
    str=[str ');'];
    
    Index=eval(str);
    
end

% 目标函数测试
function z=fun_obj(x)
  z=[0 0];
  z(1)=x(1);
  sum=0;
  for i=2:length(x)
      sum = sum+x(i);
  end
  g=1+9*(sum/(length(x)-1));
  z(2)=g*(1-(z(1)/g)^0.5);
end


function nd_pop=GetNonDominatedParticles(pop)

    ND=~[pop.Dominated];
    
    nd_pop=pop(ND);

end

function [occ_cell_index occ_cell_member_count]=GetOccupiedCells(pop)

    GridIndices=[pop.GridIndex];
    
    occ_cell_index=unique(GridIndices);
    
    occ_cell_member_count=zeros(size(occ_cell_index));

    m=numel(occ_cell_index);
    for k=1:m
        occ_cell_member_count(k)=sum(GridIndices==occ_cell_index(k));
    end
    
end

function dom=Dominates(x,y)

    if isstruct(x)
        x=x.Cost;
    end

    if isstruct(y)
        y=y.Cost;
    end
    
    dom=all(x<=y) && any(x<y);

end


function pop=DetermineDomination(pop)

    npop=numel(pop);
    
    for i=1:npop
        pop(i).Dominated=false;
        for j=1:i-1
            if ~pop(j).Dominated
                if Dominates(pop(i),pop(j))
                    pop(j).Dominated=true;
                elseif Dominates(pop(j),pop(i))
                    pop(i).Dominated=true;
                    break;
                end
            end
        end
    end

end

function rep=DeleteFromRep(rep,EXTRA,gamma)

    if nargin<3
        gamma=1;
    end

    for k=1:EXTRA
        [occ_cell_index occ_cell_member_count]=GetOccupiedCells(rep);

        p=occ_cell_member_count.^gamma;
        p=p/sum(p);

        selected_cell_index=occ_cell_index(RouletteWheelSelection(p));

        GridIndices=[rep.GridIndex];

        selected_cell_members=find(GridIndices==selected_cell_index);

        n=numel(selected_cell_members);

        selected_memebr_index=randi([1 n]);

        j=selected_cell_members(selected_memebr_index);
        
        rep=[rep(1:j-1); rep(j+1:end)];
    end
    
end

function G=CreateHypercubes(costs,ngrid,alpha)

    nobj=size(costs,1);
    
    empty_grid.Lower=[];
    empty_grid.Upper=[];
    G=repmat(empty_grid,nobj,1);
    
    for j=1:nobj
        
        min_cj=min(costs(j,:));
        max_cj=max(costs(j,:));
        
        dcj=alpha*(max_cj-min_cj);
        
        min_cj=min_cj-dcj;
        max_cj=max_cj+dcj;
        
        gx=linspace(min_cj,max_cj,ngrid-1);
        
        G(j).Lower=[-inf gx];
        G(j).Upper=[gx inf];
        
    end

end


function particle=CreateEmptyParticle(n)
    
    if nargin<1
        n=1;
    end

    empty_particle.Position=[];
    empty_particle.Velocity=[];
    empty_particle.Cost=[];
    empty_particle.Dominated=false;
    empty_particle.Best.Position=[];
    empty_particle.Best.Cost=[];
    empty_particle.GridIndex=[];
    empty_particle.GridSubIndex=[];
    
    particle=repmat(empty_particle,n,1);
    
end