1.算法仿真效果
matlab2022a仿真结果如下：
2.算法涉及理论知识概要
一般来说，海鸥是群居性的，它们用自己的智慧来寻找并攻击猎物。海鸥最重要的是它们的迁徙和攻击行为。迁徙被定义为海鸥从一个地方到另一个地方的季节性迁移，以寻找最丰富食物来源，以提供足够的能量。该行为可描述为：
在迁徙过程中，它们成群结队地迁徙。海鸥的初始位置不同，以避免相互碰撞；
在一个群体中，海鸥可以朝着最适合生存的海鸥的方向前进；
其他海鸥可以根据最适合的海鸥更新它们的初始位置。
海鸥在迁徙过程中，会经常攻击其他候鸟，在攻击过程中，它们可以做出螺旋形的自然运动，如图1所示。
海鸥的迁徙和攻击行为。海鸥根据季节更替进行迁徙，迁徙飞行时海鸥会避免相互碰撞；海鸥会攻击猎物，攻击时呈螺旋形的运动形态；在一个群体中，海鸥朝着最佳位置的方向前进。
1.初始化参数；
2.初始化种群位置；
3.计算适应度值并保留全局最优位置；
4.迁徙，全局搜索：
抽象海鸥的迁徙行为主要有三步，第一要满足避免碰撞条件（这样可以确保种群多样性）；第二要计算最佳位置的方向；第三要根据该方向移动到新的位置。
SOA的迭代伪代码如下图6所示：
3.MATLAB核心程序
%%参数设置
Umax=0.9500;%最大隶属度值
Umin=0.0111;%最小隶属度值
Wmax=0.9;%权重最大值
Wmin=0.1;%权重最小值
Dim = 3; % 维数
SwarmSize =30; % 粒子群规模
MaxIter = 100; % 最大迭代次数
MinFit = 10; % 最小适应值
Ub = [100 100 100];
Lb = [0 0 0];
%%种群初始化
Range = ones(SwarmSize,1)*(Ub-Lb);
Swarm = rand(SwarmSize,Dim).*Range + ones(SwarmSize,1)*Lb; % 初始化粒子群
fSwarm = zeros(SwarmSize,1);
for i=1:SwarmSize
fSwarm(i,:) = PID_SOA(Swarm(i,:)); % 粒子群的适应值
end
%%个体极值和群体极值
[bestf bestindex]=min(fSwarm);
zbest=Swarm(bestindex,:); % 全局最佳
gbest=Swarm; % 个体最佳
fgbest=fSwarm; % 个体最佳适应值
fzbest=bestf; % 全局最佳适应值
%迭代寻优
Di=0*rand(SwarmSize,Dim);
Buchang=0*rand(SwarmSize,Dim);
C=0*rand(SwarmSize,Dim);
Diego=0*rand(SwarmSize,Dim);
Dialt=0*rand(SwarmSize,Dim);
Dipro=0*rand(SwarmSize,Dim);
iter = 0;
y_fitness = zeros(1,MaxIter); % 预先产生4个空矩阵
K_p = zeros(1,MaxIter);
K_i = zeros(1,MaxIter);
K_d = zeros(1,MaxIter);
while( (iter < MaxIter) && (fzbest > MinFit) )
for i=1:SwarmSize
W=Wmax-iter*(Wmax-Wmin)/MaxIter;
Diego(i,:)=sign(gbest(i,:) -Swarm(i,:));%确定利己方向
Dialt(i,:)=sign(zbest -Swarm(i,:));%确定利他方向
if PID_SOA(gbest(i,:))>=PID_SOA(Swarm(i,:))%确定预动方向
Dipro(i,:)=-Di(i,:);
else
Dipro(i,:)=Di(i,:);
end
Di(i,:)=sign(W* Dipro(i,:)+0.5*Diego(i,:)+0.5*Dialt(i,:));%确定经验梯度方向
[Orderfgbest,Indexfgbest]=sort(fgbest,'descend');
u=Umax-(SwarmSize-Indexfgbest(i))*(Umax-Umin)/(SwarmSize-1);
U=u+(1-u)*rand;
H=(MaxIter-iter)/MaxIter;%迭代过程中权重的变化
C(i,:)=H*abs(zbest-10*rand(1,3));%确定高斯函数的参数
T=sqrt(-log(U));
Buchang(i,:)=C(i,:)*T;%确定搜索不长的大
Buchang(i,find(Buchang(i,:)>3*max(C(i,:))))=3*max(C(i,:));
%更新位置
Swarm(i,:)=Swarm(i,:)+Di(i,:).*Buchang(i,:);
Swarm(i,find(Swarm(i,:)>100))=100;
Swarm(i,find(Swarm(i,:)<0))=0;
% 适应值
fSwarm(i,:) =PID_SOA(Swarm(i,:));
% 个体最优更新
if fSwarm(i) < fgbest(i)
gbest(i,:) = Swarm(i,:);
fgbest(i) = fSwarm(i);
end
% 群体最优更新
if fSwarm(i) < fzbest
zbest = Swarm(i,:);
fzbest = fSwarm(i);
end
end
iter = iter+1; % 迭代次数更新
y_fitness(1,iter) = fzbest; % 为绘图做准备
K_p(1,iter) = zbest(1);
K_i(1,iter) = zbest(2);
K_d(1,iter) = zbest(3);
end
%%绘图
figure(1) % 绘制性能指标ITAE的变化曲线
plot(y_fitness,'LineWidth',2)
title('最优个体适应值');
xlabel('迭代次数');ylabel('适应值');
grid on
figure(2) % 绘制PID控制器参数变化曲线
plot(K_p,'LineWidth',2)
hold on
plot(K_i,'k','LineWidth',2)
hold on
plot(K_d,'--r','LineWidth',2)
title('Kp、Ki、Kd 优化曲线');
xlabel('迭代次数');
ylabel('参数值');
legend('Kp','Ki','Kd');
grid on