【优化算法】多目标粒子群优化算法（MOPSO）【含Matlab源码 033期】

一、获取代码方式

获取代码方式1：

完整代码已上传我的资源：【优化算法】多目标粒子群优化算法（MOPSO）【含Matlab源码 033期】

获取代码方式2：

通过订阅紫极神光博客付费专栏，凭支付凭证，私信博主，可获得此代码。

备注：

订阅紫极神光博客付费专栏，可免费获得1份代码（有效期为订阅日起，三天内有效）；

二、粒子群算法简介

1 粒子群算法的概念

粒子群优化算法(PSO：Particle swarm optimization) 是一种进化计算技术（evolutionary computation）。源于对鸟群捕食的行为研究。粒子群优化算法的基本思想：是通过群体中个体之间的协作和信息共享来寻找最优解．

PSO的优势：在于简单容易实现并且没有许多参数的调节。目前已被广泛应用于函数优化、神经网络训练、模糊系统控制以及其他遗传算法的应用领域。

2 粒子群算法分析

2.1基本思想

粒子群算法通过设计一种无质量的粒子来模拟鸟群中的鸟，粒子仅具有两个属性：速度和位置，速度代表移动的快慢，位置代表移动的方向。每个粒子在搜索空间中单独的搜寻最优解，并将其记为当前个体极值，并将个体极值与整个粒子群里的其他粒子共享，找到最优的那个个体极值作为整个粒子群的当前全局最优解，粒子群中的所有粒子根据自己找到的当前个体极值和整个粒子群共享的当前全局最优解来调整自己的速度和位置。下面的动图很形象地展示了PSO算法的过程：

2 更新规则

PSO初始化为一群随机粒子(随机解)。然后通过迭代找到最优解。在每一次的迭代中，粒子通过跟踪两个“极值”(pbest，gbest)来更新自己。在找到这两个最优值后，粒子通过下面的公式来更新自己的速度和位置。

公式(1)的第一部分称为【记忆项】，表示上次速度大小和方向的影响；公式(1)的第二部分称为【自身认知项】，是从当前点指向粒子自身最好点的一个矢量，表示粒子的动作来源于自己经验的部分；公式(1)的第三部分称为【群体认知项】，是一个从当前点指向种群最好点的矢量，反映了粒子间的协同合作和知识共享。粒子就是通过自己的经验和同伴中最好的经验来决定下一步的运动。以上面两个公式为基础，形成了PSO的标准形式。

公式(2)和 公式(3)被视为标准PSO算法。

3 PSO算法的流程和伪代码

三、部分源代码

%%%mopso;clc;clear;close all;%% Problem DefinitionCostFunction=@(x) ZDT(x);% Cost FunctionnVar=5; % Number of Decision VariablesVarSize=[1 nVar]; % Size of Decision Variables MatrixVarMin=0;% Lower Bound of VariablesVarMax=1;% Upper Bound of Variables%% MOPSO ParametersMaxIt=200; % Maximum Number of IterationsnPop=200; % Population SizenRep=100; % Repository Size库的大小w=0.5; % Inertia Weightwdamp=0.99; % Intertia Weight Damping Rate惯性重量阻尼率c1=1;% Personal Learning Coefficientc2=2;% Global Learning CoefficientnGrid=7; % Number of Grids per Dimension每个维度的网格数alpha=0.1;% Inflation Rate通货膨胀率beta=2; % Leader Selection Pressure领导人选择压力gamma=2; % Deletion Selection Pressure 删除选择压力mu=0.1; % Mutation Rate 突变率%% Initialization 初始化empty_particle.Position=[];empty_particle.Velocity=[];empty_particle.Cost=[];empty_particle.Best.Position=[];empty_particle.Best.Cost=[];empty_particle.IsDominated=[];empty_particle.GridIndex=[];empty_particle.GridSubIndex=[];%网格的子索引pop=repmat(empty_particle,nPop,1);for i=1:nPoppop(i).Position=unifrnd(VarMin,VarMax,VarSize);pop(i).Velocity=zeros(VarSize);pop(i).Cost=CostFunction(pop(i).Position);% Update Personal Best 更新个人最好pop(i).Best.Position=pop(i).Position;pop(i).Best.Cost=pop(i).Cost;end% Determine Domination确定统治pop=DetermineDomination(pop);rep=pop(~[pop.IsDominated]);Grid=CreateGrid(rep,nGrid,alpha);for i=1:numel(rep)rep(i)=FindGridIndex(rep(i),Grid);end%% MOPSO Main LoopMOPSO主循环for it=1:MaxItfor i=1:nPopleader=SelectLeader(rep,beta);pop(i).Velocity = w*pop(i).Velocity ...+c1*rand(VarSize).*(pop(i).Best.Position-pop(i).Position) ...+c2*rand(VarSize).*(leader.Position-pop(i).Position);pop(i).Position = pop(i).Position + pop(i).Velocity;pop(i).Position = max(pop(i).Position, VarMin);pop(i).Position = min(pop(i).Position, VarMax);pop(i).Cost = CostFunction(pop(i).Position);% Apply Mutation 应用突变pm=(1-(it-1)/(MaxIt-1))^(1/mu);if rand<pmNewSol.Position=Mutate(pop(i).Position,pm,VarMin,VarMax);NewSol.Cost=CostFunction(NewSol.Position);if Dominates(NewSol,pop(i))pop(i).Position=NewSol.Position;pop(i).Cost=NewSol.Cost;elseif Dominates(pop(i),NewSol)% Do Nothingelseif rand<0.5pop(i).Position=NewSol.Position;pop(i).Cost=NewSol.Cost;endendendif Dominates(pop(i),pop(i).Best)pop(i).Best.Position=pop(i).Position;pop(i).Best.Cost=pop(i).Cost;elseif Dominates(pop(i).Best,pop(i))% Do Nothingelseif rand<0.5pop(i).Best.Position=pop(i).Position;pop(i).Best.Cost=pop(i).Cost;endendend% Add Non-Dominated Particles to REPOSITORY 向存储库添加非支配粒子rep=[reppop(~[pop.IsDominated])]; %#ok% Determine Domination of New Resository Members 决定新成员的支配地位rep=DetermineDomination(rep);% Keep only Non-Dminated Memebrs in the Repository 在存储库中只保留非支配成员rep=rep(~[rep.IsDominated]);% Update Grid 更新网格Grid=CreateGrid(rep,nGrid,alpha);% Update Grid Indices 更新网格索引for i=1:numel(rep)rep(i)=FindGridIndex(rep(i),Grid);end% Check if Repository is Full 检查存储库是否已满if numel(rep)>nRepExtra=numel(rep)-nRep;for e=1:Extrarep=DeleteOneRepMemebr(rep,gamma);endend% Plot Costs 绘制成本图figure(1);PlotCosts(pop,rep);pause(0.01);% Show Iteration Informationdisp(['Iteration ' num2str(it) ': Number of Rep Members = ' num2str(numel(rep))]);% Damping Inertia Weight 显示迭代信息w=w*wdamp; end%% Resluts

四、运行结果

五、matlab版本及参考文献

1 matlab版本

a

2 参考文献

[1] 包子阳,余继周,杨杉.智能优化算法及其MATLAB实例（第2版）[M].电子工业出版社，.

[2]张岩,吴水根.MATLAB优化算法源代码[M].清华大学出版社，.