%% I. 清空环境
clc
clear

%% II. 绘制目标函数曲线
figure
[x,y] = meshgrid(1:1:100,1:1:100);%x，y 定义域
global z; %全局变量 存放z

for i=1:100
for j = 1:100
z(i,j) = x(i,j)+y(i,j);%z=x+y
end
end

mesh(x,y,z)
hold on

%% III. 参数初始化
c1 = 1.49445;
c2 = 1.49445;

maxgen = 100; % 进化次数
sizepop = 100; %种群规模
%设置最大最小速度，最大最小种群位置
Vmax = 1;
Vmin = -1;
pop1max = 1;
pop1min = -1;
pop2max = 1;
pop2min = -1;

%% IV. 产生初始粒子和速度

for i = 1:sizepop
% 产生一个种群
pop(i,:) = [50rand+51,50rand+51]; %初始种群
V(i,:) = [5rand+51,5rand+51]; %初始化速度
% 计算适应度
fitness(i) = fun(pop(i,:)); %染色体的适应度
end

%% V. 个体极值和群体极值
[bestfitness bestindex] = max(fitness);
zbest = pop(bestindex,:); %全局最佳
gbest = pop; %个体最佳
fitnessgbest = fitness; %个体最佳适应度值
fitnesszbest = bestfitness; %全局最佳适应度值

%% VI. 迭代寻优
for i = 1:maxgen
for j = 1:sizepop
% 速度更新
V(j,:) = V(j,:) + c1rand(gbest(j,:) - pop(j,:)) + c2rand(zbest - pop(j,:));
V(j,find(V(j,:)>Vmax)) = Vmax;
V(j,find(V(j,:)<Vmin)) = Vmin;
% 种群更新
pop(j,:) = pop(j,:) + V(j,:);
if (pop(j,1)>pop1max)||(pop(j,2)>pop2max)
pop(j,:) = [pop1max,pop2max];
elseif (pop(j,1)<pop1min)||(pop(j,2)<pop2min)
pop(j,:) = [pop1min,pop2min];
end
% 适应度值更新
fitness(j) = fun(pop(j,:));
end
for j = 1:sizepop
% 个体最优更新
if fitness(j) > fitnessgbest(j)
gbest(j,:) = pop(j,:);
fitnessgbest(j) = fitness(j);
end
% 群体最优更新
if fitness(j) > fitnesszbest
zbest = pop(j,:);
fitnesszbest = fitness(j);
end
end
s(i) = fitnesszbest;
end
%% VII.输出结果
%输出最优适应度、最优第IV层厚度、最优第II层厚度，并在图中标出
[fitnesszbest, zbest]
plot3(zbest(2), zbest(1),fitnesszbest,‘bo’,‘linewidth’,1.5)
axis tight
shading interp
colormap(jet)
figure
plot(s)
title(‘最优个体适应度’,‘fontsize’,12);
xlabel(‘进化代数’,‘fontsize’,12);ylabel(‘适应度’,‘fontsize’,12);

function g = fun(n)
%函数用于计算粒子适应度值
%x input 输入粒子
%y output 粒子适应度值
g=n1+n2;
end