电力系统 matlab,利用matlab解决电力系统规划问题

问题

设计思路

先用传统的基荷、腰荷、峰荷策略解决第一问，然后借助第一问的结果解决第二问。最后后将第一问的策略用代码实现，解决第三问。

一.

1. 划分负荷

明显，无论在任何负载率下，A、B、C发电的经济性是依次变差的，再加上C机组的启停都有额外的费用。所以它们出力的优先级也应当是依次递减。

将日负荷曲线分为基荷、腰荷和峰荷。在10%的旋转备用率条件下，只开A、B机组能够达到的开机容量最多能到达10800MW，无法满足所有负荷需求，所以C机组必须作为备用，在峰荷时出力。

使A、B机组满足基荷、腰荷的负荷要求，A、B、C机组满足峰荷的负荷要求。

日负荷直方图

观察日负荷直方图，有两种划分负荷的方案。

方案甲：5000MW以下为基荷，5000~8000MW为腰荷，8000MW以上为峰荷。

方案乙：5000MW以下为基荷，5000~9000MW为腰荷，9000MW以上为峰荷。

方案甲：

约束条件：

式中A、B分别表示A、B电站机组当日开机的数量。

A、B的开机容量必须大于腰荷，然后提供10%的旋转备用率。

$1000A+600B\geq8800$

A、B的最小技术出力必须小于基荷，自动满足10%的旋转备用率。

$1000A\times50 \% + 600 B\times60\%\le5500$

A、B的数量满足条件

$0 \le A \le 6 \\ 0\le B\le 8$

峰荷时，系统的开机容量需满足旋转备用率的要求，假设负荷率为100%时，C机组全开。

$1000A+600B+350\times4\geq10000\times\left(1+10\%\right)$

解不等式组，无解，下面寻找最接近的整数解。

因为A无论在何种条件下都是最经济的选择，所以先假设：

得：

$\left\{\begin {matrix} A=6\\ B=6 \end{matrix}\right.$

此时，不满足约束条件2：

$1000A\times50\%+600B\times60\%=5160\approx5000$

实际系统中，发电量可以略高于用电量。

方案乙：

约束条件：

式中A、B分别表示A、B电站机组当日开机的数量。

A、B的开机容量必须大于腰荷，然后提供10%的旋转备用率。

$1000A+600B\geq9900$

A、B的最小技术出力必须小于基荷，自动满足10%的旋转备用率。

$1000A\times50\%+600B\times60\%\le5000$

A、B的数量满足条件。

$0\le A\le6\\ 0\le B\le8$

峰荷时，系统的开机容量需满足旋转备用率的要求，假设负荷率为100%时，C机组全开。

$1000A+600B+350×4≥10000×(1+10\%)$

解不等式组(1，2，3)，无解。

所以选择方案甲：

$\left\{\begin {matrix} A=6\\ B=6 \end{matrix}\right.$

2. 计算A、B运行成本

下面计算基荷和腰荷的各个工况下A、B每台发电机的出力，已解决在一定条件下成本最小的规划问题。

约束条件：

$P=1000X_{a1}+1000X_{a2}+1000X_{a3}+1000X_{a4}+1000X_{a5}+1000X_{a6}+ 600X_{b1}+600X_{b2}+600X_{b3}+600X_{b4}+600X_{b5} \\ 500\le\ X_{a1},X_{a2},X_{a3},X_{a4},X_{a5},X_{a6}\le1000 \\ {360\le X}_{b1},X_{b2},X_{b3},X_{b4},X_{b5}\le600 \\ minF=\sum_{i=1}^{6}{X_{ai}\left(300-10\frac{X_{ai}}{1000}\right)\times\frac{600}{1000}+\sum_{i=1}^{5}{X_{bi}\left(322-12\frac{X_{bi}}{600}\right)\times\frac{600}{1000}}}$

可以发现，此规划问题的目标函数是二次函数，此问题为二次规划问题，无法用线性规划的单纯形法。所以需采取其他方法。

B机组出力的调节范围为：

${minP}_b=600\times5\times60\%=1800\mathrm{MW} \\ {maxP}_b=600\times5\times100\%=3000\mathrm{MW}$

A机组出力的调节范围为：

$minP_a=1000\times6\times50\%=3000\mathrm{MW}\\ maxP_a=1000\times6\times100\%=6000\mathrm{MW}$

由于B机组的单位成本始终高于A机组，所以尽量减小B机组的出力即可减小整体发电成本。于是，可以将基荷和腰荷时期不同的负荷率对应到A、B不同的出力水平，如下：

负荷/MW

5500

6500

8000

A出力功率/MW

3700

4700

6000

B出力功率/MW

1800

2000

使B的出力始终最小即可。

确定了A、B机组的出力之和，下一步即需确定A、B机组的每台出力，即：在总出力确定的情况下，如何为每台机组分配出力，使得总成本最小。

对于A机组：

如某一小时的A机组运行成本

$F_a=\sum_{i=1}^{6}{X_{ai}\left(300-10\frac{X_{ai}}{1000}\right)\times\frac{600}{1000}} =\sum_{i=1}^{6}\left(-\frac{6}{1000}{X_{ai}}^2+180X_{ai}\right)$

是

$X_{ai}$ 的开口向下的二次函数，所以

是

$X_{ai}$ 的上凸函数。

由琴生不等式及其推论可知，应使机组的出力尽量集中于最大出力和最小技术出力处，仅有一台机组处于最大出力和最小技术出力之间时，成本最小。

琴生不等式：

当

$a_1+a_2+a_3+\ldots+a_n=1，且f\left(X\right)$ 为凸函数时：

$a_1f\left(X_1\right)+a_2f\left(X_2\right)+a_3f\left(X_3\right)+\ldots+a_nf\left(X_n\right)\le \ f\left(X_1+X_2+X_3+\ldots+X_n\right)$

B机组同理

所以，可以得出基荷和腰荷时所有发电机的出力水平表如下：

负荷

5500

6500

8000

$X_{a1}$

840

1000

$X_{a2}$

500

1000

$X_{a3}$

500

840

1000

$X_{a4}$

500

1000

$X_{a5}$

500

1000

$X_{a6}$

500

840

$X_{b1}$

360

$X_{b2}$

360

$X_{b3}$

360

$X_{b4}$

360

$X_{b5}$

360

$X_{b6}$

360

总成本

997447.2

1168447

1424947

3. 计算C启停及运行成本

下面考虑峰荷时即负荷率为90%和100%的的情况：

负荷率90%时

为了保证10%的旋转备用率，必须使总开机容量到达

$9900\mathrm{MW}$

$1000A+600B+350C=9900\mathrm{MW}$

所以最少开启2台C机组

$350\mathrm{MW}$ ，总开机容量达到

$10600\mathrm{MW}$ 。

C机组出力调整为最小：

$20\%\times350=70\mathrm{MW}$

则C机组总共出力

$210\mathrm{MW}$ ，A、B需出力：

$9000-2\times70=8860\mathrm{MW}$

负荷率100%时

为了保证10%的旋转备用率，必须使总开机容量到达

$11000\mathrm{MW}$ ，

$1000A+600B+350C=11000\mathrm{MW}$

所以最少开启4台C机组

$350\mathrm{MW}$ ，总开机容量达到

$11000\mathrm{MW}$ 。

目标出力

5500

6500

8000

9000

10000

$X_{a1}$

840

1000

$X_{a2}$

500

1000

$X_{a3}$

500

840

1000

$X_{a4}$

500

1000

$X_{a5}$

500

1000

$X_{a6}$

500

1000

$X_{b1}$

360

600

$X_{b2}$

360

600

$X_{b3}$

360

600

$X_{b4}$

360

410

600

$X_{b5}$

360

600

$X_{b6}$

360

600

$X_{c1}$

190

$X_{c2}$

$X_{c3}$

$X_{c4}$

总出力

5500

6500

8000

9000

10000

开机容量

9360

9950

11000

总成本

997447.2

1168447

1424947

1625101

1900940

4. 结论

得出电站日运行规划

时间(h)

目标出力(MW)

各机组出力(MW)

实际出力(MW)

开机容量(MW)

单位小时出力成本(元)

开停机成本(元)

$X_{a1}$

$X_{a2}$

$X_{a3}$

$X_{a4}$

$X_{a5}$

$X_{a6}$

$X_{b1}$

$X_{b2}$

$X_{b3}$

$X_{b4}$

$X_{b5}$

$X_{b6}$

$X_{c1}$

$X_{c2}$

$X_{c3}$

$X_{c4}$

6500

1000

840

500

360

6500

9360

1168447

5000

840

500

360

5500

9360

997447.2

5000

840

500

360

5500

9360

997447.2

5000

840

500

360

5500

9360

997447.2

5000

840

500

360

5500

9360

997447.2

6500

1000

840

500

360

6500

9360

1168447

6500

1000

840

500

360

6500

9360

1168447

8000

1000

840

360

8000

9360

1424947

8000

1000

840

360

8000

9360

1424947

9000

1000

600

410

360

9000

9950

1625101

36000

9000

1000

600

410

360

9000

9950

1625101

8000

1000

840

360

8000

9360

1424947

36000

8000

1000

840

360

8000

9360

1424947

8000

1000

840

360

8000

9360

1424947

8000

1000

840

360

8000

9360

1424947

8000

1000

840

360

8000

9360

1424947

9000

1000

600

410

360

9000

9950

1625101

36000

9000

1000

600

410

360

9000

9950

1625101

10000

1000

600

190

10000

11000

1900940

108000

10000

1000

600

190

10000

11000

1900940

9000

1000

600

410

360

9000

9950

1625101

108000

8000

1000

840

360

8000

9360

1424947

3600

8000

1000

840

360

8000

9360

1424947

6500

1000

840

500

360

6500

9360

1168447

成本(元)：

33415488.17

396000

总成本(元)：

33811488.17

二.

1.建立停运容量概率模型

对于确切状态概率，有递推公式

为新增一台机组(容量

$C\ \mathrm{MW}$ ，

)后，系统停运容量为

的确切概率，

对于第一台机组，

当

时，

累积状态概率公式

取步长为

$50\mathrm{MW}$ 编程迭代计算停运概率。

2.matlab程序代码

\main.m

%main.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

clear

clc

cmat=[ones(6,1)*1000;ones(8,1)*600;ones(4,1)*350]';

rmat=[ones(6,1)*0.1;ones(8,1)*0.08;ones(4,1)*0.12]';

p=[];

p(1)=1;

for i=1:18

p2=[];%p2代表p(x)

c=cmat(i);

r=rmat(i);

p1=p;%p1代表p'(x)

p1(x2i(sum(cmat(1:i))))=0;

for x=0:50:sum(cmat(1:i))%x从0到装机容量

if x

pxc=0;

else

pxc=p1(x2i(x-c));

end

px=p1(x2i(x));

p2=[p2,(1-r)*px+r*pxc];

end

p=p2;

end

P=p*(tril(ones(245)));

function i= x2i(x)

%此函数用于将x值转换为i值

i=x/50+1;

end

三.

筛选出所有的可行的装机组合，然后将其对应的待机组合求出。

将“一”中的策略用代码实现，输入待机组合，输出机组出力表和日运行成本，

同时根据装机组合求其装机成本。

比较各种方案的总成本。

程序代码

%main.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

clc

clear

global PowerPlantCombine;%存放A、B、C电站的个数，(1x3)矩阵

LoadRate=[65,55,55,55,55,65,65,80,80,90,90,80,80,80,80,80,90,90,100,100,90,80,80,65]/100;

MaxLoad=12500;

Load=LoadRate*MaxLoad;

%PowerPlantNumberMat=[6;6;4];

PowerPlantCombineMat=candidatePowerPlantCombine();

dayCostmat=[];

for i = 1:length(PowerPlantCombineMat)%遍历所有可能的电厂组合方式

PowerPlantCombine=PowerPlantCombineMat(:,i);

dayCost=0;

Outputmat=[];

for j=1:24%遍历24小时，根据负载求每个电厂的出力组合

load=Load(j);

Output=output(load);%出力组合

Cost=cost(Output);%出力费用

if j==1

dayCost=dayCost+Cost;

elseif j==24

dayCost=dayCost+Cost+CSwitchCost(Outputmat(1,:),Output);%考虑24点开1点关的C机组的开停机费用

else

dayCost=dayCost+Cost+CSwitchCost(Outputmat(end,:),Output);%考虑C机组的开停机费用

end

Outputmat=[Outputmat;Output];%24小时的出力表

end

dayCostmat=[dayCostmat,dayCost];

end

powerPlantBuildCombineMat=convert(PowerPlantCombineMat);%根据待机容量求相应的装机容量。(待机容量即为当日开机了的所有机组的总容量)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

buildCostMat=([44067.24,29745.387,12333.3]*10000)*(powerPlantBuildCombineMat-PowerPlantCombineMat);%计算装机成本

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

costmat=dayCostmat*365+buildCostMat;%计算总成本

[min,index]=min(costmat);%找到所有费用的最小值

%还需考虑如果最小值不满足C机组的启停约束，需退而求其次

aim=[costmat;powerPlantBuildCombineMat;PowerPlantCombineMat];

aim(:,index)

%cost.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

function cost = cost(output)

%此函数输入每小时的output(出力)矩阵，输出当小时的发电成本

global PowerPlantCombine;

Anumber=PowerPlantCombine(1);

Bnumber=PowerPlantCombine(2);

Cnumber=PowerPlantCombine(3);

Aoutput=output(1:Anumber);

Boutput=output(Anumber+1:Anumber+Bnumber);

Coutput=output(Anumber+Bnumber+1:end);

Acost=sum((300-10.*Aoutput./1000).*600.*Aoutput./1000);

Bcost=sum((322-12.*Boutput./1000).*600.*Boutput./1000);

Ccost=sum((266-16.*Coutput./1000).*1800.*Coutput./1000);

cost=Acost+Bcost+Ccost;

end

%convert.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

function powerPlantbuildCombineMat = convert(powerPlantonCombineMat)

%通过此函数输入待机组合输出对应的装机组合

powerPlantbuildCombineMat=[];

for i=1:length(powerPlantonCombineMat)

powerPlantOnCombine=powerPlantonCombineMat(:,i);

Anumber=powerPlantOnCombine(1);

Bnumber=powerPlantOnCombine(2);

Cnumber=powerPlantOnCombine(3);

if Anumber<6

Anumber=6;

end

if Bnumber<8

Bnumber=8;

end

if Cnumber<4

Cnumber=4;

end

powerPlantOnLoad=Anumber*1000+Bnumber*600+Cnumber*350;

if powerPlantOnLoad<15000

Cnumber=Cnumber+ceil((15000-powerPlantOnLoad)/350);

end

powerPlantbuildCombine=[Anumber;Bnumber;Cnumber];

powerPlantbuildCombineMat=[powerPlantbuildCombineMat,powerPlantbuildCombine];

end

%condidatePowerPlant.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

function A = candidatePowerPlantCombine()

%此函数输出可行的候选待机发电机组合

% 此处显示详细说明

T=[];

J=[];

K=[];

for j=1:8

for k=1:10

for t=1:7

if (1000*j+600*k>=11000)&(500*j+360*k<=6875)&(1000*j+600*k>=13750-350*t)

J=[J,j];

K=[K,k];

T=[T,t];

end

A=[J;K;T];

end

%CSwitchCost.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

function switchCost = CSwitchCost(outputFormer,outputNow)

% 输入两个相邻两小时的机组出力状态，通过比较C机组启停状态，计算C机组起停的费用

% 此处显示详细说明

global PowerPlantCombine;

Anumber=PowerPlantCombine(1);

Bnumber=PowerPlantCombine(2);

Cnumber=PowerPlantCombine(3);

CoutputFormer=outputFormer(Anumber+Bnumber+1:end);

CoutputNow=outputNow(Anumber+Bnumber+1:end);

CBoolFormer=logical(CoutputFormer);

CBoolNow=logical(CoutputNow);

Switch=(CBoolFormer-CBoolNow).^2;%是否开/关，若为“1”则开或关

switchNumber=sum(Switch);%总开/关次数

switchCost=switchNumber*10*1800;%费用

end

%output.m

%华中科技大学电气与电子工程学院

%强电磁与新技术国家重点实验室

%季天泽

%M201877109

function output = output(load)

%输入负荷，输出所有电站的出力水平矩阵，也就是“一”中的部分

%将出力大小分为四个范围，以0，M1，M2，M3，M4为界，分别考虑出力情况

%M1~M2为基荷

%M2~M3为腰荷

%M3~M4为峰荷

%其他范围无解

global PowerPlantCombine

Anumber=PowerPlantCombine(1);%A机组的数量

Bnumber=PowerPlantCombine(2);%B机组的数量

Cnumber=PowerPlantCombine(3);%C机组的数量

MaxA=Anumber*1000;%A、B、C所有机组的出力范围

MaxB=Bnumber*600;

MaxC=Cnumber*350;

MinA=Anumber*500;

MinB=Bnumber*360;

MinC=Cnumber*70;

M1=MinA+MinB;%0~M1无解

M2=MaxA+MinB;%M1~M2,C停机，B最小出力

M3=MaxA+MaxB;%M2~M3，C停机，A最大出力

M4=MaxA+MaxB+MaxC;%M3~M4，C开机，A、B最大出力(B可能有一台机组需减小出力)

%M4以上，无解

%通过总负荷来得到A、Bnumber、C机组分别的总负荷，对于相同机组，使更多的机组处于出力区间的端点处，一台机组处于出力区间之内(琴生不等式)

if (M1<=load)&&(load

Aload=load-MinB;

Bload=MinB;

AmaxNumber=floor((Aload-500*Anumber)/500);

AchangableLoad=(Aload-500*Anumber)-500*AmaxNumber+500;

outputA=[ones(1,AmaxNumber)*1000,[AchangableLoad],ones(1,Anumber-AmaxNumber-1)*500];

outputB=ones(1,Bnumber)*360;

outputC=ones(1,Cnumber)*0;

elseif (M2<=load)&&(load

Aload=MaxA;

Bload=load-MaxA;

BmaxNumber=floor((Bload-360*Bnumber)/240);

BchangableLoad=(Bload-360*Bnumber)-240*BmaxNumber+360;

outputB=[ones(1,BmaxNumber)*600,BchangableLoad,ones(1,Bnumber-BmaxNumber-1)*360];

outputA=ones(1,Anumber)*1000;

outputC=ones(1,Cnumber)*0;

elseif (M3<=load)&&(load

Aload=MaxA;

Bload=MaxB;

Cload=load-MaxA-MaxB;

if Cload<70%如果A、B都为最大出力将导致C的目标出力小于其可能的最小出力，那么需将C设置为最小出力(一台机为70，其余为0)，重新设置B的出力。

Cload=70;

Bload=load-Aload-Cload;

BmaxNumber=floor((Bload-360*Bnumber)/240);

BchangableLoad=(Bload-360*Bnumber)-240*BmaxNumber+360;

outputB=[ones(1,BmaxNumber)*600,BchangableLoad,ones(1,Bnumber-BmaxNumber-1)*360];

outputA=ones(1,Anumber)*1000;

outputC=ones(1,Cnumber)*0;

outputC(1)=70;

% outputB=ones(1,Bnumber);

% outputB(end)=Bload-(Bnumber-1)*600;

else

ConNumber=ceil(Cload/350);%确认C开机的数量

CmaxNumber=ceil((Cload-70*ConNumber)/280-1);

CchangableLoad=(Cload-70*ConNumber)-280*CmaxNumber+70;

outputC=[ones(1,CmaxNumber)*350,CchangableLoad,ones(1,ConNumber-CmaxNumber-1)*70,ones(1,Cnumber-ConNumber)*0];

outputB=ones(1,Bnumber)*600;

outputA=ones(1,Anumber)*1000;

end

else

'此发电机组合错误'

PowerPlantCombine

end

output=[outputA,outputB,outputC];

end

你可以随意更改我的代码的任何部分，如果我的思路和代码对你有所帮助，请点赞或评论，这对我非常重要，谢谢

电力系统 matlab,利用matlab解决电力系统规划问题相关推荐

环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序
环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序,程序通用性强,通过修改参数可以得到任意节点和网络的环形网络牛拉法潮流计算. YID:856064261479453 ...
环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序，程序通用性强
环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序,程序通用性强,通过修改参数可以得到任意节点和网络的环形网络牛拉法潮流计算. 现有:6960642614794538 ...
环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序，程序通用性强，通过修改参数
环形网络潮流计算matlab 利用matlab编程计算任意环形网络牛拉法潮流计算程序,程序通用性强,通过修改参数可以得到任意节点和网络的环形网络牛拉法潮流计算. 现有:6960642614794538 ...
热传递 matlab,利用matlab程序解决热传导问题
利用matlab程序解决热传导问题 1哈佛大学能源与环境学院课程作业报告作业名称:传热学大作业--利用 matlab 程序解决热传导问题院系:能源与环境学院专业:建筑环境与设备工程学号:5201314 ...
神经网络人口预测matlab,利用MATLAB工具模拟实现人口预测实例分析.doc
摘要:人口预测对于国民经济计划的制定和社会战略目标的决策具有重要参考价值．本文主要利用MATLAB软件对中国1974-2015年的人口总数数据进行统计分析,建立BP神经网络预测模型.把它与一元非线性回 ...
二维静电场有限元matlab,利用MATLAB计算电磁场有关分布
<利用MATLAB计算电磁场有关分布>由会员分享,可在线阅读,更多相关<利用MATLAB计算电磁场有关分布(17页珍藏版)>请在人人文库网上搜索. 1.电磁场实验报告实验一模 ...
物联网控制原理与技术--基于Matlab/利用MATLAB进行频域分析（伯德图）的应用（超详细/设计/实验/作业/练习）
目录课程名:物联网控制原理与技术内容/作用:设计/实验/作业/练习学习:利用MATLAB进行频域分析(伯德图) 一.前言二.环境与设备三.原理四.内容五.结果与分析课程名:物联网控制原 ...
投入产出实例matlab,利用MATLAB进行基础的投入产出表分析实验报告.doc
利用MATLAB进行基础的投入产出表分析实验报告学生实验报告学院: XX学院课程名称:宏观经济数量分析方法与模型专业班级: 姓名: 学号: 学生实验报告学生姓名学号同组人实验 ...
单闭环调速仿真matlab,利用Matlab仿真平台设计单闭环直流调速系统
内容简介: 毕业论文利用Matlab仿真平台设计单闭环直流调速系统共21页,3989字. 目录一.摘要 --------------------------2 二.总体方案设计 ------- ...
多普勒谱matlab,利用MATLAB仿真多普勒效应.doc
利用MATLAB仿真多普勒效应利用MATLAB仿真多普勒效应某某某摘要:分析多普勒效应特性,建立数学模型,利用MATLAB软件对其进行仿真试验,进行定量分析,根据仿真试验结果绘制出听者接收到的 ...

电力系统 matlab,利用matlab解决电力系统规划问题

电力系统 matlab,利用matlab解决电力系统规划问题相关推荐

最新文章

热门文章