📄 soa_wgyh_zyk_dai_57.m
字号:
% FDA_RUN
%
% Filter Design Algorithm that uses a seeker optimization algorithm to design and
% optimizes an IIR transfer function that fits a desired magnitude
% response. Adjustable inputs include the desired magnitude response
% (H_d), frequency specification points (W), frequnecy weighting values
% (Q), transfer fucntion order (alpha), population size (N), probability
% of crossover (p_c), maximimum number of generation executions
% (gen_max), minimum fitness stopping criteria (fit_min), and number of
% optimization attempts (attempt_max).
%
% SELECT_RANK, CROSS_NORMAL, and REPLACE_TRUNCATE functions are required by
% FDA_RUN. Function FDA_ANALYZE is needed to analyze and plot the results.
%
% Author: Chaohua Dai, 05-19-07
% Southwest Jiaotong University, Chengdu, 610031, China
% Housekeeping
clear
close
rand('state',sum(100*clock)) % create a random seed (VERY IMPORTANT!!!!)
tic
home
Size = 24;
% max number of generations per problem
gen_max =200;
% min fitness level
fit_min = 0;
thresholdFun=1e-3;
% number of design attempts
trnum = 50;
% number of vectors per element
Variables=31;
%[baseMVA, bus, gen, branch] = loadcase(casename);
casename='case57a';
opt=mpoption;
opt=mpoption(opt,'PF_ALG',2);
Variablesmax=[40.9 40.9 40.9 40.9 40.9 40.9 40.9 40.9 40.9 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.124 1.10 1.10 1.10 1.10 1.10 1.10 1.10];
Variablesmin=[-20 -20 -20 -20 -20 -20 -20 -20 -20 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9];
% Variablesmax=[40 40 40 40 40 40 40 40 40 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10 1.10];
% Variablesmin=[-20 -20 -20 -20 -20 -20 -20 -20 -20 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9];
% Umax=2;
% Umin=-2;
%
% Variablesmax=Umax.*ones(1,Variables);
% Variablesmin=Umin.*ones(1,Variables);
% data1=load('2006data.txt');
% ttt=data1(:,1);
% aaa=data1(:,3);
% hhh=data1(:,4);
%
% data=aaa;
% tau = 3;
% m=6;
% tsmax=max(data);
% tsmin=min(data);
% data=(data-tsmin)/(tsmax-tsmin);
%
% data11=data(9600:9744);
% N=length(data11);
% [n1 x1]=reconstitution(data11,N,m,tau);%reconstitute state space
%
% % fcmdata=[n11 x11'];
% % [center, U, obj_fcn] = fcm(fcmdata, 100,[NaN,150,NaN,NaN]);
% % n1=center(:,1:m);
% % x1=center(:,m+1)';
%
%
% data11=data(9729:9792);%(1041:1104);
% N=length(data11);
% [n2 x2]=reconstitution(data11,N,m,tau);%reconstitute state space
% %-------------------my program
%
% Variablesmax=[10 0.02 100000000 0.1];%[20.133 0.001766 8.6698e+005][40.472 0.01 3.0316e+005][10 0.011705 5.8847e+005]%[2.0926 0.0089522 2729.7 0.010215 1][1.6323 0.0005266 2929.4 0.010799 1.3103][0.5 0.0001 3869.1 0.0051837 3.0148][1.75 2901.5 100 0.025];%[7 10000 100];
% Variablesmin=[1 0.0001 10000000 0.01];%[1.75 2901.5 0.01 0.01];%[0.01 10 1 ];1.8939 0.001[2.7579 0.001097 3748.2 0.022802]
delta_xmax=0.3*(Variablesmax-Variablesmin);
delta_xmin=-delta_xmax;
w_max=0.9;w_min=0.4;
mu_max=0.9;mu_min=0.0111;
% number of elements in population
%寻找邻域
numregion=3;
for region=1:1:numregion
startreg(region)=floor((region-1)*Size/numregion+1);
endreg(region)=floor(region*Size/numregion);
Sizereg(region)=endreg(region)-startreg(region)+1;
end
fFun=zeros(trnum,gen_max);
average_k=zeros(1,trnum); %
average_time=zeros(1,trnum);
stop_num=gen_max*ones(1,trnum);
num_stop=gen_max;
failnum=0;
fBestFun=ones(1,trnum); %
F=10000.*ones(1,Size);%当代函数值
F_t_1=F; %上一代函数值
F_t_2=F; %上上代函数值
fBestS=ones(trnum,Variables);% 记录trnum次实验的最佳个体
lBestS=ones(Size,Variables); %局部最佳个体
gBestS=ones(1,Variables);
gBestSreg=ones(numregion,Variables); %全局最佳个体
% Sign_x=ones(Size,Variables);
Sign_x=ones(1,Variables);
stall = cputime;
for tr=1:trnum
% Display the attempt in progress
disp(['Attempt = ',num2str(tr)]);
num_stop=gen_max;
% n1=n11;
% x1=x11;
% fcmnum = 100;%50+tr*25;
% %fcm
% fcmdata=[n11 x11'];
% [center, U, obj_fcn] = fcm(fcmdata, fcmnum ,[NaN,200,NaN,NaN]);
% n1=center(:,1:m);
% x1=center(:,m+1)';
% create the initial population of filters
% with a uniform PDF in the unit circle
% build a structure for handling the population of elements,
% corresponding fitness, the selection subset of elements, and
% generated offspring elements
E = ones(Size,1)*Variablesmin+(ones(Size,1)*(Variablesmax-Variablesmin)).*rand(Size,Variables);
% for i=1:Size
% for j=1:3
% E(i,j) = floor(E(i,j));
% end
% for j=4:7
% E(i,j) = floor(E(i,j).*40)/40;
% end
% end
E_t_2=E; %上代个体
E_t_1=E; %上上代个体
lBestFun=inf*ones(1,Size);
gBestFun=inf;
gBestFunreg(1:1:numregion)=inf;
flag=0;
time1=cputime;
% Start the FDA
for k = 1:gen_max
% Evaluate fitness of current generation
for i=1:Size
for j=1:9
E(i,j) = floor(E(i,j));
end
for j=10:24
E(i,j) = floor(E(i,j).*40)/40;
end
end
for s=1:Size
% calculate fitness
F_t_2(s)=F_t_1(s);
F_t_1(s)=F(s);
F(s)=runpf_test57(casename,E(s,:),opt);
GS=[tr,k,s]
end
par_xx(tr,k,:,:) = E;
[BestFun,indexBestF]=min(F); %
[WorstFun,indexWorstF]=max(F); %
%weight=(G-k)/gen_max; %inertia weight
weight=w_max-k*(w_max-w_min)/gen_max;
bFun(k)=BestFun;% 不同于gBestFun,当算法不收敛时，两者不相等
if BestFun<=gBestFun,% 最小问题，所以gBestFun大时需更新
gBestFun=BestFun; %全局最佳函数值
gBestS=E(indexBestF(1),:); %全局最佳个体
end
gBestS_xx(tr,k,:)=gBestS;
gBestFun_xx(tr,k)=gBestFun;
for region=1:1:numregion,
[BestFunreg,IndexBestreg]=min(F(startreg(region):endreg(region)));
IndexBestreg=IndexBestreg+startreg(region)-1;
if gBestFunreg(region)>BestFunreg,% 最小问题，所以gBestFun大时需更新
gBestFunreg(region)=BestFunreg;%全局最佳函数值
gBestSreg(region,:)=E(IndexBestreg(1),:); %全局最佳个体
end
[WorstFunreg,IndexWorstreg]=max(F(startreg(region):endreg(region)));
IndexWorstreg=IndexWorstreg+startreg(region)-1;
% if gWorstFunreg(region) < WorstFunreg,% 最小问题，所以gBestFun大时需更新
% gWorstFunreg(region)=WorstFunreg;%全局最佳函数值
gWorstSreg(region,:)=E(IndexWorstreg(1),:); %全局最佳个体
% end
end
for s=1:1:Size,
if lBestFun(s)>F(s),% 最小问题，所以lBestFun大时需更新
lBestFun(s)=F(s);
lBestS(s,:)=E(s,:); %局部最佳个体
end
end
if gBestFun<=thresholdFun & flag==0 %检查是否满足 求最小值问题
flag=1;
average_k(tr)=k;
elseif gBestFun>thresholdFun & k==gen_max % 检查是否失败
average_k(tr)=gen_max;
failnum=failnum+1;
end
% check for fit_min stopping criteria
if(gBestFun <= fit_min)
disp('Minimum fitness achieved.');
num_stop=k;
break ;
end
for region=1:1:numregion,
[OderF,IndexF]=sort(F(startreg(region):endreg(region)),'descend');
IndexF=IndexF+startreg(region)-1;
%% 邻域间交流
si=0;
for region1=1:1:numregion,
if region1~=region,
si=si+1;
E(IndexF(si),:)=gBestSreg(region1,:);
F(IndexF(si))=gBestFunreg(region1);
end
end
% En=weight*abs(E(IndexF(end),:)-E(IndexF(1),:)); %weight*
%
% aver_E=sum(E(startreg(region):endreg(region),:),1)/(startreg(region)-endreg(region)+1);
En=weight*abs(gBestSreg(region,:)-gWorstSreg(region,:)); %
% En=weight*abs(gBestSreg(region,:)-aver_E); %
%
% radius=(Variablesmax-Variablesmin)/3;
% En=weight*radius.*rand(1,Variables);
% En=weight*radius;
for s=startreg(region):1:endreg(region), % kk+old_num-1,% %开始对种群逐个个体进行进化操作
IndexF_local=find(IndexF==s);%+startreg(region)-1;
mu=(mu_max-(mu_max-mu_min)*(Sizereg(region)-IndexF_local)/(Sizereg(region)-1));%*rand(1,Variables);
% IndexF_local=find(IndexF==s)+startreg(region)-1;
% mu=(mu_max-(mu_max-mu_min)*(Size-IndexF_local)/(Size-1));
mu=mu+(1.0-mu).*rand(1,Variables);
x_ldirect=lBestS(s,:)-E(s,:); %局部最优个体的方向
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -