基于多元宇宙算法求解电力系统多目标优化问题附Matlab代码

%-------------------------------------------------------------------------- % Economic and Emission dispatch Data for 140 bus power system % This is the main program. % Main paper: % A. Sundaram, “Multiobjective multi-verse optimization algorithm to solve % combined economic, % heat and power emission dispatch problems,” % Applied Soft Computing Journal, % vol. 91, p. 106195, 2020, % doi: 10.1016/j.asoc.2020.106195. % % Data was developed by Dr.Arunachalam Sundaram and is available in % appendix B of the above paper. % Email:mailtoarunachalam@gmail.com %-------------------------------------------------------------------------- clc; clear; close all; tStart=tic; global costdata emissiondata Pd VarMin VarMax nVar Pd=49342; data140 costdata=[data1(:,1:9);data1(:,10:18)]; emissiondata=1e-2.*[... 1 0.0480 -2.2200 60.0000 1.3100 0.0569 2 0.0480 -2.2200 60.0000 1.3100 0.0569 3 0.0762 -2.3600 100.0000 1.3100 0.0569 4 0.0540 -3.1400 120.0000 0.9142 0.0454 5 0.0850 -1.8900 50.0000 0.9936 0.0406 6 0.0854 -3.0800 80.0000 1.3100 0.0569 7 0.0242 -3.0600 100.0000 0.6550 0.02846 8 0.0310 -2.3200 130.0000 0.6550 0.02846 9 0.0335 -2.1100 150.0000 0.6550 0.02846 10 0.4250 -4.3400 280.0000 0.6550 0.02846 11 0.0322 -4.3400 220.0000 0.6550 0.02846 12 0.0338 -4.2800 225.0000 0.6550 0.02846 13 0.0296 -4.1800 300.0000 0.5035 0.02075 14 0.0512 -3.3400 520.0000 0.5035 0.02075 15 0.0496 -3.5500 510.0000 0.5035 0.02075 16 0.0496 -3.5500 510.0000 0.5035 0.02075 17 0.0151 -2.6800 220.0000 0.5035 0.02075 18 0.0151 -2.6600 222.0000 0.5035 0.02075 19 0.0151 -2.6800 220.0000 0.5035 0.02075 20 0.0151 -2.6800 220.0000 0.5035 0.02075 21 0.0145 -2.2200 290.0000 0.5035 0.02075 22 0.0145 -2.2200 285.0000 0.5035 0.02075 23 0.0138 -2.2600 295.0000 0.5035 0.02075 24 0.0138 -2.2600 295.0000 0.5035 0.02075 25 0.0132 -2.4200 310.0000 0.5035 0.02075 26 0.0132 -2.4200 310.0000 0.5035 0.02075 27 1.8420 -1.1100 360.0000 0.9936 0.0406 28 1.8420 -1.1100 360.0000 0.9936 0.0406 29 1.8420 -1.1100 360.0000 0.9936 0.0406 30 0.0850 -1.8900 50.0000 0.9936 0.0406 31 0.0121 -2.0800 80.0000 0.9142 0.0454 32 0.0121 -2.0800 80.0000 0.9142 0.0454 33 0.0121 -2.0800 80.0000 0.9142 0.0454 34 0.0012 -3.4800 65.0000 0.6550 0.02846 35 0.0012 -3.2400 70.0000 0.6550 0.02846 36 0.0012 -3.2400 70.0000 0.6550 0.02846 37 0.0950 -1.9800 100.0000 1.4200 0.0677 38 0.0950 -1.9800 100.0000 1.4200 0.0677 39 0.0950 -1.9800 100.0000 1.4200 0.0677 40 0.0151 -2.6800 22.0000 0.5035 0.02075 ]; emissiondata=[emissiondata;emissiondata;emissiondata;emissiondata(1:20,:)]; emissiondata(:,1)=1:140; %% Max_time=2500; N=20; ArchiveMaxSize=100; % max_iter=Max_time; nVar=140; % Number of Decision Variables VarSize=[1 nVar]; % Size of Decision Variables Matrix VarMin=max(costdata(:,2),(costdata(:,9)-costdata(:,8))); % Lower Bound of Variables VarMax=min(costdata(:,3),(costdata(:,9)+costdata(:,7))); % Upper Bound of Variables fobj=@(x) IEEE140aobj(x); dim=nVar; lb=VarMin'; ub=VarMax'; obj_no=2; Best_universe=zeros(1,dim); Best_universe_Inflation_rate=inf*ones(1,obj_no); Archive_X=zeros(ArchiveMaxSize,dim); Archive_F=ones(ArchiveMaxSize,obj_no)*inf; Archive_member_no=0; WEP_Max=1; WEP_Min=0.2; for i=1:N Universes(i,:)=lcheck140; end Time=1; while Time<Max_time+1 WEP=WEP_Min+Time*((WEP_Max-WEP_Min)/Max_time); TDR=1-((Time)^(1/6)/(Max_time)^(1/6)); for i=1:size(Universes,1) %Boundary checking (to bring back the universes inside search % space if they go beyoud the boundaries Flag4ub=Universes(i,:)>ub; Flag4lb=Universes(i,:)<lb; Universes(i,:)=(Universes(i,:).*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb; Universes(i,:)=lbcoff140bus(Universes(i,:)); %Calculate the inflation rate (fitness) of universes Inflation_rates(i,:)=fobj(Universes(i,:)); %Elitism if dominates(Inflation_rates(i,:),Best_universe_Inflation_rate) Best_universe_Inflation_rate=Inflation_rates(i,:); Best_universe=Universes(i,:); end end [sorted_Inflation_rates,sorted_indexes]=sort(Inflation_rates); for newindex=1:N Sorted_universes(newindex,:)=Universes(sorted_indexes(newindex),:); end %Normaized inflation rates (NI in Eq. (3.1) in the original MVO paper) normalized_sorted_Inflation_rates=normr(sorted_Inflation_rates); Universes(1,:)= Sorted_universes(1,:); % Universes(1,:)=lchecktf1(Universes(1,:)); [Archive_X, Archive_F, Archive_member_no]=UpdateArchive(Archive_X, Archive_F, Universes, Inflation_rates, Archive_member_no); if Archive_member_no>ArchiveMaxSize Archive_mem_ranks=RankingProcess(Archive_F, ArchiveMaxSize, obj_no); [Archive_X, Archive_F, Archive_mem_ranks, Archive_member_no]=HandleFullArchive(Archive_X, Archive_F, Archive_member_no, Archive_mem_ranks, ArchiveMaxSize); else Archive_mem_ranks=RankingProcess(Archive_F, ArchiveMaxSize, obj_no); end Archive_mem_ranks=RankingProcess(Archive_F, ArchiveMaxSize, obj_no); % to improve coverage index=RouletteWheelSelection(1./Archive_mem_ranks); if index==-1 index=1; end Best_universe_Inflation_rate=Archive_F(index,:); Best_universe=Archive_X(index,:); %Update the Position of universes for i=2:size(Universes,1)%Starting from 2 since the firt one is the elite Back_hole_index=i; for j=1:size(Universes,2) r1=rand(); if r1<normalized_sorted_Inflation_rates(i) White_hole_index=RouletteWheelSelection(-sorted_Inflation_rates);% for maximization problem -sorted_Inflation_rates should be written as sorted_Inflation_rates if White_hole_index==-1 White_hole_index=1; end %Eq. (3.1) in the paper Universes(Back_hole_index,j)=Sorted_universes(White_hole_index,j); % Universes(Back_hole_index,j)=lchecktf1(Universes(Back_hole_index,j)); end if (size(lb',1)==1) %Eq. (3.2) in the original MVO paper if the boundaries are all the same r2=rand(); if r2<WEP r3=rand(); if r3<0.5 Universes(i,j)=Best_universe(1,j)+TDR*((ub-lb)*rand+lb); end if r3>0.5 Universes(i,j)=Best_universe(1,j)-TDR*((ub-lb)*rand+lb); end end end if (size(lb',1)~=1) %Eq. (3.2) in the original MVO paper if the upper and lower bounds are %different for each variables r2=rand(); if r2<WEP r3=rand(); if r3<0.5 Universes(i,j)=Best_universe(1,j)+TDR*((ub(j)-lb(j))*rand+lb(j)); end if r3>0.5 Universes(i,j)=Best_universe(1,j)-TDR*((ub(j)-lb(j))*rand+lb(j)