💥💥💥💞💞💞欢迎来到本博客❤️❤️❤️💥💥💥
📝目前更新：🌟🌟🌟电力系统相关知识，期刊论文，算法，机器学习和人工智能学习。
🚀支持：🎁🎁🎁如果觉得博主的文章还不错或者您用得到的话，可以关注一下博主，如果三连收藏支持就更好啦！这就是给予我最大的支持！
👨🎓博主课外兴趣：中西方哲学，送予读者：
👨💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学，什么是电的时候，不要觉得这些问题搞笑，哲学就是追究终极问题，寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能让人胸中升起一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它居然给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“真理”上的尘埃吧。
或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎
📋📋📋本文目录如下：⛳️⛳️⛳️
目录
0 概述
多元宇宙算法求解电力系统多目标优化算法有很好的效果，代码换成自己的目标函数，加上约束和惩罚项等。本文用多元宇宙算法求解电力系统多目标优化问题——电力系统环境经济调度问题。
提出了一种求解电力系统环境经济调度的新方法,该方法利用宇宙空间在随机创建过程中高膨胀率的物体随虫洞在空间移动物体的规律,通过对白洞和黑洞间随机传送物体来实现最优搜索. 算法具有运算速度快,收敛性强,适用于高维计算等特点.以总燃料费用最低和总污染排放最少为多目标建立环境经济调度模型,最后,通过发电厂传统10机组和40机组算例进行仿真.结果表明:本文所提算法具有经济性和有效性.
1 环境经济调度数学模型
2 多元宇宙算法
3 讲解
3.1 主函数：
-
clc;
-
clear;
-
close
all;
-
tStart
=tic;
-
%
global costdata emissiondata B B0 B00 Pd VarMin VarMax nVar
-
global data B B0 B00 Pd VarMin VarMax nVar
-
Pd
=
2000;
-
data
=xlsread(
'IEEE10.xls');
-
-
B1
=xlsread(
'B10.xls');
-
B
=B1(
1:
10,
1:
10);
-
B0
=[
0
0
0
0
0
0
0
0
0
0];
-
B00
=
0;
-
% B
=B1(
1:
3,
1:
3);
-
% B0
=[
0
0
0
0
0
0
0
0
0
0];
-
% B00
=
0;
-
-
-
%
%
-
-
Max_time
=
600;
%迭代次数
-
N
=
100;
-
ArchiveMaxSize
=
100;
-
% max_iter
=Max_time;
-
nVar
=
10;
% 机组个数
-
-
VarSize
=[
1 nVar];
% 决策变量矩阵的大小
-
VarMin
=data(:,
2);
%机组出力下限
-
VarMax
= data(:,
3);
% 机组出力上限
-
-
fobj
=@(x) IEEE3aobj(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,:)
=lcheck3;
-
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)
-
-
%边界检查（如果宇宙超出边界，则将它们带回搜索空间内）
-
Flag4ub
=Universes(i,:)
>ub;
-
Flag4lb
=Universes(i,:)
<lb;
-
Universes(i,:)
=(Universes(i,:).
*(
~(Flag4ub
+Flag4lb)))
+ub.
*Flag4ub
+lb.
*Flag4lb;
-
Universes(i,:)
=lbcoff3bus(Universes(i,:));
-
%计算宇宙的通货膨胀率(适合度)
-
Inflation_rates(i,:)
=fobj(Universes(i,:));
-
-
%精英主义
-
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
-
-
%原始MVO论文中的标准化通货膨胀率
-
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);
-
% 提高复盖率
-
index
=RouletteWheelSelection(
1.
/Archive_mem_ranks);
-
if index
=
=
-1
-
index
=
1;
-
end
-
Best_universe_Inflation_rate
=Archive_F(index,:);
-
Best_universe
=Archive_X(index,:);
-
-
-
-
%更新宇宙的位置
-
for i
=
2:size(Universes,
1)
%从
2开始，因为第
1位是精英
-
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);
% 对于最大化问题，排序的通货膨胀率应该写成排序的通货膨胀率
-
if White_hole_index
=
=
-1
-
White_hole_index
=
1;
-
end
-
%Eq. (
3.1)
-
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)
-
%如果边界都是一样的，那么原MVO论文中的公式(3.2)就会出现
-
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)
-
%公式(
3.2 )在原始MVO论文中，如果对每个变量的上下界不同
-
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));
-
end
-
end
-
end
-
-
end
-
Universes(i,:)
=lbcoff3bus(Universes(i,:));
-
end
-
display([
'At the iteration ', num2str(
Time),
' there are ', num2str(Archive_member_no),
' non-dominated solutions in the archive']);
-
Time
=
Time
+
1;
-
%
-
end
-
plot(Archive_F(:,
1),Archive_F(:,
2),
'Ro',
'LineWidth',
2,...
-
'MarkerEdgeColor',
'r',...
-
'MarkerFaceColor',
'r',...
-
'MarkerSize',
2);
-
xlabel(
'污染排放量')
-
ylabel(
'煤耗量')
-
title(
'Pareto最前沿')
-
% Universes
-
Archive_F(:,
1)
-
Archive_F(:,
2)
-
Best_universe
clc;
clear;
close all;
tStart=tic;
% global costdata emissiondata B B0 B00 Pd VarMin VarMax nVar
global data B B0 B00 Pd VarMin VarMax nVar
Pd=2000;
data=xlsread('IEEE10.xls');
B1=xlsread('B10.xls');
B=B1(1:10,1:10);
B0=[0 0 0 0 0 0 0 0 0 0];
B00=0;
% B=B1(1:3,1:3);
% B0=[0 0 0 0 0 0 0 0 0 0];
% B00=0;
%%
Max_time=600; %迭代次数
N=100;
ArchiveMaxSize=100;
% max_iter=Max_time;
nVar=10; % 机组个数
VarSize=[1 nVar]; % 决策变量矩阵的大小
VarMin=data(:,2); %机组出力下限
VarMax= data(:,3); % 机组出力上限
fobj=@(x) IEEE3aobj(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,:)=lcheck3;
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)
%边界检查（如果宇宙超出边界，则将它们带回搜索空间内）
Flag4ub=Universes(i,:)>ub;
Flag4lb=Universes(i,:)<lb;
Universes(i,:)=(Universes(i,:).*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb;
Universes(i,:)=lbcoff3bus(Universes(i,:));
%计算宇宙的通货膨胀率(适合度)
Inflation_rates(i,:)=fobj(Universes(i,:));
%精英主义
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
%原始MVO论文中的标准化通货膨胀率
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);
% 提高复盖率
index=RouletteWheelSelection(1./Archive_mem_ranks);
if index==-1
index=1;
end
Best_universe_Inflation_rate=Archive_F(index,:);
Best_universe=Archive_X(index,:);
%更新宇宙的位置
for i=2:size(Universes,1)%从2开始，因为第1位是精英
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);% 对于最大化问题，排序的通货膨胀率应该写成排序的通货膨胀率
if White_hole_index==-1
White_hole_index=1;
end
%Eq. (3.1)
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)
%如果边界都是一样的，那么原MVO论文中的公式(3.2)就会出现
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)
%公式( 3.2 )在原始MVO论文中，如果对每个变量的上下界不同
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));
end
end
end
end
Universes(i,:)=lbcoff3bus(Universes(i,:));
end
display(['At the iteration ', num2str(Time), ' there are ', num2str(Archive_member_no), ' non-dominated solutions in the archive']);
Time=Time+1;
%
end
plot(Archive_F(:,1),Archive_F(:,2),'Ro','LineWidth',2,...
'MarkerEdgeColor','r',...
'MarkerFaceColor','r',...
'MarkerSize',2);
xlabel('污染排放量')
ylabel('煤耗量')
title('Pareto最前沿')
% Universes
Archive_F(:,1)
Archive_F(:,2)
Best_universe
4 运行结果
10机组运行结果：
本文提出了一种求解电力系统环境经济调度的新方法,计及阀点效应和污染排放因素,建立多目标规划模型,利用PPF定价原则权衡多重因素. 多元宇宙算法在求解EED问题时具有计算精度高,收敛速度快等特点,在求解高维度问题表现更佳,适用于其他工程问题研究.
4 IEEE140节点仿真算例及Matlab代码【创新】
🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。
5 参考文献
[1]刘世宇,王孜航,杨德友.多元宇宙算法及其在电力系统环境经济调度的应用[J].东北电力大学学报,2018,38(04):19-26.
6 完整Matlab代码实现
转载：https://blog.csdn.net/m0_73907476/article/details/127964378