水电站中长期优化调度主要是在满足电力系统出力要求以及下游综合用水要求等的前提下，合理地安排水库各个月份的发电流量，使得调度期内的总发电量最大。即在给定预报入库流量过程线、下游综合用水流量过程线、调度期起始水位和终止水位等条件下，通过总水量的合理运用使调度期内水电站的总发电量最大。本源码由GreenSim团队原创，转载请注明，有意购买源码或代写相关程序，请与GreenSim团队联系(主页http://blog.sina.com.cn/greensim)。

%% 参数设置

A=8.5;%出力系数，常数

Tt=730*ones(12,1);%第t个时段的小时数

%注意：一年按363天*24小时算，均分为12个月

HtLB=55*ones(12,1);%第t时段水位约束的下界，单位：米

HtUB=[65;65;65;61;61;61;61;65;65;65;65;65];%第t时段水位约束的上界，单位：米

VtLB=zeros(12,1);

VtUB=zeros(12,1);

for i=1:12

VtLB(i)=Ht2Vt(HtLB(i));

VtUB(i)=Ht2Vt(HtUB(i));

end

%注意：蓄水量Vt和水位Ht之间有一一对应的关系，单位：立方米

NtLB=260000*ones(12,1);%出力约束的下界，单位：千瓦

NtUB=1400000*ones(12,1);%出力约束的上届，单位：千瓦

%注意：Nt=A*Qt*Ht

QtLB=308*ones(12,1);%泄流量下界，单位：立方米/秒

QtUB=29200*ones(12,1);%泄流量上界，单位：立方米/秒

qt=[373;859;1568;2100;3210;5049;1596;1160;925;781;572;1010];%入库流量，单位：立方米/秒

%注意：以上三个量，时间单位相乘时，小时乘以3600转化成秒

% GreenSim团队原创作品，转载请注明

% Email:greensim@163.com

% GreenSim团队主页：http://blog.sina.com.cn/greensim

%

[color=red]欢迎访问GreenSim——算法仿真团队→[url=http://blog.sina.com.cn/greensim]http://blog.sina.com.cn/greensim[/url][/color]

%% 调用粒子群算法

K=60;

N=80;

w=0.5;

c1=0.3;

c2=0.2;

[BESTX,BESTY,ALLX,ALLY]=PSO(K,N,w,c1,c2,VtLB,VtUB,QtLB,QtUB,NtLB,NtUB,qt,A,Tt);

%%

X=BESTX{K};

[Vt,Qt,St]=DeCode(X);

disp('最佳蓄水量');

disp(Vt);

disp('最佳平均发电流量');

disp(Qt);

disp('最佳平均弃水流量');

disp(St);

disp('最大总发电量');

disp(-BESTY(K));

function

[BESTX,BESTY,ALLX,ALLY]=PSO(K,N,w,c1,c2,VtLB,VtUB,QtLB,QtUB,NtLB,NtUB,qt,A,Tt)

%% 此函数实现粒子群优化算法，用于水电站调度优化

% GreenSim团队原创作品，转载请注明

% Email:greensim@163.com

% GreenSim团队主页：http://blog.sina.com.cn/greensim

%

[color=red]欢迎访问GreenSim——算法仿真团队→[url=http://blog.sina.com.cn/greensim]http://blog.sina.com.cn/greensim[/url][/color]

%% 输入参数列表

%

K 迭代次数

%

N 粒子个数

%

w 惯性因子

%

c1 加速因子，针对历史最优状态

%

c2 加速因子，针对全局最优状态

%

VtLB 蓄水量约束下界，立方米，N*1向量

%

VtUB 蓄水量约束上界，立方米，N*1向量

%

QtLB 平均发电流量约束下界，立方米/秒，N*1向量

%

QtUB 平均发电流量约束上界，立方米/秒，N*1向量

%

NtLB 出力约束下界，千瓦，N*1向量

%

NtUB 出力约束上界，千瓦，N*1向量

%

qt 入库流量，立方米/秒，N*1向量

%

A 出力系数

%

Tt 一个周期的小时数

%% 输出参数列表

% BESTX K×1细胞结构，每一个元素是M×1向量，记录每一代的最优个体

% BESTY K×1矩阵，记录每一代的最优个体的评价函数值

%

ALLX K×1细胞结构，每一个元素是M×N矩阵，记录全部个体

%

ALLY K×N矩阵，记录全部个体的评价函数值

%% 第一步：粒子状态初始化

M=3*length(VtLB);

farm=zeros(N,M);

for i=N

[Vt,Qt,St]=Initialize2(VtLB,VtUB,QtLB,QtUB,NtLB,NtUB,qt,A,Tt);

X=EnCode(Vt,Qt,St);

farm(i,:)=X';

end

%粒子历史最优状态初始化

Pfarm=farm;

%粒子群全局最优状态

for i=1:N

X=Pfarm(i,:);

[Vt,Qt,St]=DeCode(X');

SE=ObjFun(A,Qt,St,Vt,Tt);

Fitness(i)=SE;

end

MinFit=min(Fitness);

POS=find(Fitness==MinFit);

Gfarm=Pfarm(POS(1),:);

%输出变量初始化

ALLX=cell(K,1);%细胞结构，每一个元素是M×N矩阵，记录每一代的个体

ALLY=zeros(K,N);%K×N矩阵，记录每一代评价函数值

BESTX=cell(K,1);%细胞结构，每一个元素是M×1向量，记录每一代的最优个体

BESTY=zeros(K,1);%K×1矩阵，记录每一代的最优个体的评价函数值

k=1;%迭代计数器初始化

%% 第二步：迭代过程

while k<=K

%% 粒子状态更新

newfarm=farm;

for

i=1:N

newfarm(i,:)=w*farm(i,:)+c1*rand*(Pfarm(i,:)-farm(i,:))+c2*rand*(Gfarm-farm(i,:));

X=newfarm(i,:);

[Vt,Qt,St]=DeCode(X');

[Flag,Vt,Qt,St]=Correct(Vt,Qt,St,VtLB,VtUB,QtLB,QtUB,NtLB,NtUB,qt,A,Tt);

X=EnCode(Vt,Qt,St);

newfarm(i,:)=X';

end

%% 历史最优状态和全局最优状态更新

NEWFIT=zeros(1,N);

for

i=1:(N)

aa=newfarm(i,:);

[Vt,Qt,St]=DeCode(aa');

FitA=ObjFun(A,Qt,St,Vt,Tt);

bb=Pfarm(i,:);

[Vt,Qt,St]=DeCode(bb');

FitB=ObjFun(A,Qt,St,Vt,Tt);

NEWFIT(i)=FitA;

if FitA

Pfarm(i,:)=A;

end

if FitA

Gfarm=A;

end

end

%% 记录

ALLX{k}=newfarm;

ALLY(k,:)=NEWFIT;

BESTX{k}=Gfarm;

[Vt,Qt,St]=DeCode(Gfarm');

BESTY(k)=ObjFun(A,Qt,St,Vt,Tt);

disp(k);

k=k+1;

end