24小时热门版块排行榜

返回列表

离子木木

金虫 (正式写手)

应助: 0 (幼儿园)
金币: 637.5
散金: 20
帖子: 331
在线: 24小时
虫号: 4551087
注册: 2016-03-29
性别: GG
专业: 制造系统与自动化

[求助] 以信噪比为适应度函数的自适应粒子群优化算法Matlab Code 修改已有1人参与

本人最近写了一篇以信噪比为适应度函数的自适应粒子群优化算法Matlab Code，但是运行总是出错误，恳求哪位Matlab大神帮忙看看修改修改，谢谢。
信噪比SNR适应度函数：
function SNR=fitness(x)
N=20000;
fs=20;
h=1/fs;
T=(0:N-1)*h;
AL=0.1;
fL=0.01;
u=AL*sin(2*pi*fL.*T);
D=0.5; %噪声强度
ff=fs*(0:N/100-1)/N;
sz=sqrt(2*D).*randn(1,N); %产生的白噪声序列
s=u+sz;       %信号
xx(1,:)=sr(x(1),x(2),h,s);
Yf=fft(xx(1,:));
P=Yf.*conj(Yf)/N; %功率谱
G=sum(P);
H=P(11);
GG=(G-H)/10;
SNR=10*log(H/GG)/log(10);  %以SNR作为适应度函数

非线性微分方程离散并求解的数值算法：
function x=sr(a,b,h,x1)
x=zeros(1,length(x1));
for i=1:length(x1)-1
k1=h*(x(1)*x(i)-x(2)*x(i).^3+x1(i));
k2=h*(x(1)*(x(i)+k1/2)-x(2)*(x(i)+k1/2).^3+x1(i));
k3=h*(x(1)*(x(i)+k2/2)-x(2)*(x(i)+(sqrt(2)-1)*k1/2+(2-sqrt(2))*k2/2).^3+x1(i+1));
k4=h*(x(1)*(x(i)+k3)-x(2)*(x(i)-sqrt(2)*k2/2+(2+sqrt(2))*k3/2).^3+x1(i+1));
x(i+1)=x(i)+(1/6)*(k1+(2-sqrt(2))*k2+(2+sqrt(2))*k3+k4);
end

主函数：
clear all;
clc;
format long;
%------给定初始化条件-------------------
c1=2;          %学习因子1
c2=2;          %学习因子2
wmax=0.9;             %惯性权重
wmin=0.6;             %惯性权重
MaxDT=100;          %最大迭代次数
D=2;                %搜索空间维数（未知数个数）
N=40;                %初始化群体个体数目
eps=10^(-6);          %设置精度(在已知最小值时候用)
%------初始化种群的个体(可以在这里限定位置和速度的范围)------------
for i=1:N
for j=1:D
      x(i,j)=rand(1);  %随机初始化位置
      v(i,j)=rand(1);  %随机初始化速度
end
end
for i=1:N
p(i)=fitness(x(i,:));
y(i,:)=x(i,:);
end
pg=x(1,:);          %Pg为全局最优
for i=2:N
if fitness(x(i,:))<fitness(pg)
      pg=x(i,:);
end
end
for t=1:MaxDT
for j=1:N
      fv(j)=fitness(x(j,:));
end
   fvag=sum(fv)/N;
   fmin=min(fv);
for i=1:N
      if fv(i)<=fvag
         w=wmin-(fv(i)-fmin)*(wmax-wmin)/(fvag-fmin);
      else
         w=wmax;
      end
      v(i,:)=w*v(i,:)+c1*rand*(y(i,:)-x(i,:))+c2*rand*(pg-x(i,:));  %实现速度的更新
      x(i,:)=x(i,:)+v(i,:);  %实现位置的更新
      if x(i,1)>1||x(i,1)<0
         x(i,1)=rand(1);
      end
      if x(i,2)>1||x(i,2)<0
         x(i,2)=rand(1);
      end
      if fitness(x(i,:))<p(i) %判断当此时的位置是否为最优的情况，当不满足时继续更新
         p(i)=fitness(x(i,:));
         y(i,:)=x(i,:);
      end
      if p(i)<fitness(pg)
         pg=y(i,:);
      end
end
Pbest(t)=fitness(pg);
end
%------最后从所得到的结果中取出最优的解
disp('*************************************************************')
disp('函数的全局最优位置为：')
Solution=pg'
disp('最后得到的优化极值为：')
Result=fitness(pg)
disp('*************************************************************')
figure(1)
r=[1:1:100];
plot(r,Pbest,'r--','linewidth',1);
xlabel('迭代次数')
ylabel('适应度值')
title('自适应权重PSO算法')
hold on
Solution=pg'
Result=fitness(pg)
figure(2)
plot(x,'DisplayName','x')
axis([0,40,-5,5])
%------算法结束---DreamSun GL & HF,适应度函数源程序（fitness.m）------------

回复此楼