0
0

Matlab|含风电-光伏-光热电站电力系统N-k安全优化调度模型

**1 **主要内容

该程序参考《光热电站促进风电消纳的电力系统优化调度》光热电站模型，主要做的是考虑N-k安全约束的含义风电-光伏-光热电站的电力系统优化调度模型，从而体现光热电站在调度灵活性以及经济性方面的优势。同时代码还考虑了光热电站对风光消纳的作用，对比了含义光热电站和不含光热电站下的弃风弃光问题，同时还对比了考虑N-k约束下的调度策略区别。以14节点和118节点算例为例，对模型进行了系统性的测试，复现效果良好，是学习N-k约束以及光热电站调度的必备程序！程序采用matlab+cplex（mosek/gurobi）进行求解，可以选择已经安装的求解器进行求解。

  • 程序算例

程序对于118节点系统采用了四个算例进行对比，14节点系统有3种算例对比，并增加了弃风量的对比程序。

  • 程序模型

  • 程序亮点

  1. 采用光热电站模型，也是最近研究比较热的一个方向。
  2. 采用转移分布因子矩阵处理潮流问题，这也是很多文献中都采用的方法。​

**2 **部分程序

  1. clc; clear; close all; % 关闭所有已打开的绘图窗口
  2. %% 参数设定
  3. NT = 24; % 时间范围
  4. CoeffReseve_load = 0.03; 
  5. CoeffReserve_VRE = 0.05; 
  6. yita_TES = 0.98;  
  7. yita_PB = 0.415;  
  8. % 文章里Table 2的数据
  9. Capacity_TES_CSP = 0; 
  10. initial_TES_t0 = 0.5;  
  11. initial_TES_t1 = 0.78;
  12. TES_initial = 0.5;         
  13. beta_Load = 3*10e3;  
  14.     mpc = case14_1; % 载入数据 matpower 数据格式
  15. %% 有功负荷 24h所有节点总的
  16. %    mpc.load = [
  17. %        2842.42  3020.2  3296.96  3444.44  3607.07  3891.91  4070.7  4295.95  4476.76  4661.61  4859.59  5077.77  ...
  18. %        4717.17  4519.19  4301.01  3995.95  3703.03  3806.06  4037.37  4063.63  3721.21  3245.45  3097.97  2827.27
  19. %    ]/6.3; 
  20. ​
  21.    mpc.load = [
  22.        683.42  792.2  896.96  1044.44  1087.07  1121.91  1200.7  1235.95  1326.76  1461.61  1489.59  1577.77  ...
  23.        1417.17  1219.19  1101.01  1075.95  903.03  1186.06  1237.37  1463.63  1221.21  1005.45  827.97  807.27
  24.     ]/2; 
  25. ​
  26. ​
  27.     mpc.P_RE = [0.00   0.00   0.00   0.00   0.00   0.00   15.76   43.17   82.35   109.44   122.55   146.10   ...% PV
  28.                 126.66   86.05   60.05   52.82   25.78   4.28   0.00   0.00   0.00   0.00   0.00   0.00  
  29.                 100.26   133.95   147.28   134.11   170.52   159.44   138.55   72.83   58.83   73.37   79.90   80.54 ...  % Wind
  30.                 91.96   101.68   121.49   122.93   133.11   162.44   130.95   133.25   151.26   139.33   120.60   90.33
  31.                 ]*1; % 可再生能源 24小时数据（实际发电量）
  32. %% 电网相关名称
  33.     baseMVA = mpc.baseMVA;
  34.     bus = mpc.bus;
  35.     gen = mpc.gen;
  36.     branch = mpc.branch;
  37.     gencost = mpc.gencost;
  38.     RE = mpc.RE;
  39.     CSP = mpc.CSP;
  40.     P_RE = mpc.P_RE;
  41. ​
  42. N = length(bus(:,1));      % 网络中所有节点数
  43. N_Br = length(branch(:,1));% 线路数
  44. N_Gen = length(gen(:,1));  % 火电发电机组数
  45. N_RE = length(RE(:,1));    % 可再生能源节点机组数
  46. N_CSP = length(CSP(:,1));  % CSP发电站数
  47. ​
  48. % 常规机组相关数据提取, 取数据矩阵中的列向量 和功率有功的项，均需标幺值化，以便运算和求解
  49. P_Gen_max = gen(:,9)/baseMVA; 
  50. P_Gen_min = gen(:,10)/baseMVA; 
  51. type_Gen = gen(:,22); 
  52. P_Gen_up = gen(:,23) /baseMVA;  
  53. P_Gen_down = gen(:,24) /baseMVA;
  54. T_Gen_min_on = gen(:,25); 
  55. T_Gen_min_off = gen(:,26); 
  56. c_ST_g = gen(:,28);
  57. c_G_g = gen(:,30); 
  58. ​
  59. % CSP机组相关数据提取
  60. P_CSP_max = CSP(:,9)/baseMVA; 
  61. P_CSP_min = CSP(:,10)/baseMVA; 
  62. P_CSP_up = CSP(:,23)/baseMVA;   
  63. P_CSP_down = CSP(:,24)/baseMVA; 
  64. T_CSP_min_on = CSP(:,25); 
  65. T_CSP_min_off = CSP(:,26);
  66. c_CSP_g = CSP(:,30);       
  67. ​
  68. PtCSP_fore = [ % 可用的太阳能热功率向量 
  69.     0.00   0.00   0.00   0.00   0.00   0.00   190.57   390.57   790.57 990.57   1390.57   1891.03 ...
  70.     2111.64   2200.92   2202.36   2118.26   1895.37   1408.35   0.00   0.00   0.00   0.00   0.00   0.00 ]/20;
  71. PtCSP_fore = PtCSP_fore/baseMVA; 
  72. P_RE = P_RE/baseMVA; % 可再生能源PV WT机组出力
  73. ​
  74. beta_Load = beta_Load*baseMVA^2; % $/MWh -> $/p.u.
  75. ​
  76. M_bus_G = zeros(N,N_Gen); % 发电机机组-索引矩阵
  77. for row = 1:N
  78.     if abs(find(mpc.gen(:,1) == row)) > 0  % 发电机节点号 与 行号对应
  79.         M_bus_G(row,find(mpc.gen(:,1) == row)) = 1; % M_bus_G相应处置1
  80.     end
  81. end
  82. ​
  83. M_bus_RE = zeros(N,N_RE); % 可再生能源机组-索引矩阵
  84. for row = 1:N
  85.     if abs(find(mpc.RE(:,1) == row))>0
  86.         M_bus_RE(row,find(mpc.RE(:,1) == row)) = 1;
  87.     end
  88. end
  89. ​
  90. M_bus_CSP = zeros(N,N_CSP); % CSP机组-索引矩阵
  91. for row = 1:N
  92.     if abs(find(mpc.CSP(:,1) == row))>0
  93.         M_bus_CSP(row,find(mpc.CSP(:,1) == row)) = 1;
  94.     end
  95. end
  96. GSDF = makePTDF(mpc); % 发电转移分布因子矩阵，表征节点注入功率在全网络的分布
  97. ​
  98. %% 负荷矩阵数据，按照 算例数据mpc.bus(:,3) 中各节点负荷的比例分配
  99.     PD = bus(:,3)/baseMVA; 
  100.     P_factor = PD/sum(PD);
  101.     P_sum = mpc.load/baseMVA; 
  102.     PD = P_factor*P_sum;      
  103. ​
  104. %% 决策变量命名
  105.     PG_G = sdpvar(N_Gen,NT,'full');  
  106.     PG_RE = sdpvar(N_RE,NT,'full');   % （风光并网量）
  107.     PG_CSP = sdpvar(N_CSP,NT,'full'); 
  108.     PC_Load = sdpvar(N,NT,'full');   
  109.     
  110.     onoff_gen = binvar(N_Gen,NT,'full');
  111.     onoff_CSP = binvar(N_CSP,NT,'full'); 
  112.     
  113.     Branch = sdpvar(N_Br,NT,'full');   
  114.     Cost_StartUp  = sdpvar(N_Gen,NT-1,'full');
  115.     
  116.     Pt_TES_charge = sdpvar(N_CSP,NT,'full');  
  117.     Pt_TES_discharge= sdpvar(N_CSP,NT,'full');
  118.     Et_TES = sdpvar(N_CSP,NT,'full');         
  119.     
  120. %% 约束条件列写   
  121.     Cons = [];
  122.     for t = 1:NT
  123.         if t >= 2 % type(1-水电, 2-火电机组)
  124.             for i = 1:N_Gen % 火电机组-最小启/停时间约束 式(8-9)
  125.                 if (type_Gen(i,1)==2) || (type_Gen(i,1)==5) 
  126.                     for tao = t + 1:min(t+T_Gen_min_on(i,1)-1,NT)   
  127.                         Cons = [Cons, onoff_gen(i,t)-onoff_gen(i,t-1) <= onoff_gen(i,tao)];
  128.                     end
  129.                     for tao = t + 1:min(t+T_Gen_min_off(i,1)-1,NT) 
  130.                         Cons = [Cons, onoff_gen(i,t-1)-onoff_gen(i,t) <= 1-onoff_gen(i,tao)];
  131.                     end
  132.                 end
  133.             end
  134.             for i = 1:N_CSP  
  135.                 for tao = t+1:min(t+T_CSP_min_on(i,1)-1,NT)
  136.                     Cons = [Cons, onoff_CSP(i,t)-onoff_CSP(i,t-1) <= onoff_CSP(i,tao)]; % CSP机组最小启/停时间约束
  137.                 end
  138.                 for tao = t+1:min(t+T_CSP_min_off(i,1)-1,NT)
  139.                     Cons = [Cons, onoff_CSP(i,t-1)-onoff_CSP(i,t) <= 1-onoff_CSP(i,tao)];
  140.                 end
  141.             end
  142.         end 
  143.         if t >= 2 % 火电机组 爬坡约束 式(6-7)
  144.             Cons = [Cons,  PG_G(:,t) - PG_G(:,t-1) <= ...
  145.                      onoff_gen(:,t-1).* P_Gen_up*60 + ... 
  146.                     (onoff_gen(:,t)-onoff_gen(:,t-1)) .* P_Gen_min + ... 
  147.                     (1-onoff_gen(:,t)) .* P_Gen_max];  
  148.             Cons = [Cons, -PG_G(:,t) + PG_G(:,t-1) <= ...
  149.                     onoff_gen(:,t) .* P_Gen_down*60 + ...
  150.                    (onoff_gen(:,t-1)-onoff_gen(:,t)) .* P_Gen_min + ...  
  151.                    (1-onoff_gen(:,t-1)) .* P_Gen_max];
  152.                
  153.             % CSP 机组 爬坡约束 式(6-7)
  154.             Cons = [Cons,  PG_CSP(:,t) - PG_CSP(:,t-1) <= ...
  155.                      onoff_CSP(:,t-1).* P_CSP_up*60 + ... %  
  156.                      (onoff_CSP(:,t)-onoff_CSP(:,t-1)) .* P_CSP_min + ...
  157.                      (1-onoff_CSP(:,t)) .* P_CSP_max]; 
  158.             Cons = [Cons, -PG_CSP(:,t) + PG_CSP(:,t-1) <= onoff_CSP(:,t) .* P_CSP_down*60 + ...  
  159.                     (onoff_CSP(:,t-1)-onoff_CSP(:,t)) .* P_CSP_min + ...  
  160.                     (1-onoff_CSP(:,t-1)) .* P_CSP_max];
  161.         end
  162.     end
  163.     % 机组出力的上下边界约束-式(3) % t(1-水电,2-火电, 5-燃气发电机组 6-CSP)
  164.   Ind_2_5 = union(find(type_Gen(:,1) == 2),find(type_Gen(:,1) == 5)); 
  165.    Cons = [Cons, onoff_gen(Ind_2_5,:) .* (P_Gen_min(Ind_2_5,1) * ones(1,NT)) ...    
  166.            <= PG_G(Ind_2_5,:) <= ...
  167.           onoff_gen(Ind_2_5,:) .* (P_Gen_max(Ind_2_5,1) * ones(1,NT))];  
  168.         
  169.    
  170.         
  171.     Cons = [Cons, onoff_CSP.*(P_CSP_min*ones(1,NT)) <= PG_CSP <= onoff_CSP.*(P_CSP_max*ones(1,NT))]; % CSP机组出力-边界约束
  172. %     Cons = [Cons, onoff_CSP == ones(1,24)]; % CSP机组 
  173.   
  174.     Cons = [Cons, sum(PG_G,1) + sum(PG_RE,1) + sum(PG_CSP,1) == sum(PD - PC_Load,1)]; % 式(2)
  175.    
  176.     Cons = [Cons, Branch == GSDF*(M_bus_G*PG_G + M_bus_RE*PG_RE + M_bus_CSP*PG_CSP - (PD-PC_Load))]; % 
  177. %     Cons = [Cons, -branch(:,6)*ones(1,NT) <= GSDF*(M_bus_G*PG_G+M_bus_RE*PG_RE+M_bus_CSP*PG_CSP-(PD- PC_Load)) <= branch(:,6)*ones(1,NT)]; % 
  178.     Cons = [Cons, -999*ones(N_Br,NT) <= GSDF*(M_bus_G*PG_G+M_bus_RE*PG_RE+M_bus_CSP*PG_CSP-(PD-PC_Load)) <= 999*ones(N_Br,NT)]; % 118系统有186条线路
  179.    
  180.     Cons = [Cons, 0 <= PG_RE <= P_RE]; % 可再生出力
  181.  
  182.     Cons = [Cons, [60;50;100;80;40]/baseMVA * ones(1,24) <= PG_G ];
  183.   
  184.     Cons = [Cons, 0 <= PC_Load <= PD]; % 式(22)    
  185.   
  186.     Cons = [Cons, sum(onoff_gen .* (P_Gen_max*ones(1,NT)) - PG_G,1) + ...
  187.             sum(onoff_CSP .* (P_CSP_max*ones(1,NT)) - PG_CSP,1) >= ...
  188.             sum(CoeffReseve_load*PD,1) + sum(CoeffReserve_VRE*PG_RE,1) ];
  189.    
  190.     Cons = [Cons, Cost_StartUp >= (onoff_gen(:,2:NT) - onoff_gen(:,1:NT-1)) .* (c_ST_g*ones(1,NT-1))]; % 传统机组启动成本
  191.     Cons = [Cons, Cost_StartUp >= 0];
  192.     
  193. %%%%%% CSP电站运转内部约束 %%%%%%
  194.     E_TES_max = Capacity_TES_CSP * P_CSP_max; 
  195.    
  196.     Cons = [Cons, PG_CSP/yita_PB + Pt_TES_charge - Pt_TES_discharge <= PtCSP_fore]; % CSP输出电功率与TES充/放热功率,预测光热功率关系
  197.   
  198.     Cons = [Cons, Et_TES(:,2:NT) == Et_TES(:,1:NT-1) + Pt_TES_charge(:,1:NT-1)*yita_TES - Pt_TES_discharge(:,1:NT-1)/yita_TES];
  199.   Cons = [Cons, Et_TES(:,1) == TES_initial * E_TES_max]; 
  200.     Cons = [Cons, Et_TES(:,1) == Et_TES(:,NT)];          
  201.    
  202.     Cons = [Cons, 0 <= Pt_TES_charge    <= Capacity_TES_CSP*ones(N_CSP,NT)]; 
  203.     Cons = [Cons, 0 <= Pt_TES_discharge <= Capacity_TES_CSP*ones(N_CSP,NT)];
  204.   
  205.     Cons = [Cons, 0 <= Et_TES <= E_TES_max * ones(1,NT)];
  206. ​
  207. %% 目标函数 
  208.     obj = sum(c_G_g'*PG_G) + sum(c_CSP_g'*PG_CSP) + sum(sum(Cost_StartUp) + beta_Load*sum(sum(PC_Load)) ); 
  209.     % 机组的边际发电成本 + 启动成本 + 负荷削减成本
  210.     
  211.     % 运行调度 
  212.     ops = sdpsettings('solver','cplex'); %  gurobi
  213.     ans = optimize(Cons,obj,ops)
  214.     
  215. %% 求解成功后取值
  216.   PG_G = value(PG_G)  ; 
  217.     PG_RE = value(PG_RE) ;  
  218.     PG_CSP = value(PG_CSP) ; 
  219.     PC_Load = value(PC_Load) ;   
  220.     onoff_gen = value(onoff_gen) ; 
  221.     onoff_CSP = value(onoff_CSP) ; 
  222.     Branch = value(Branch) ;   
  223.     Cost_StartUp  = value(Cost_StartUp);
  224.     obj = value(obj); % 总成本
  225.     Pt_TES_charge = value(Pt_TES_charge);   
  226.     Pt_TES_discharge = value(Pt_TES_discharge); 
  227.     Et_TES = value(Et_TES);                 
  228.     
  229. disp(['IEEE14 不考虑N-k的和无CSP的经济调度情况,运行成本为 ', num2str(obj)])
  230. %% 绘图 
  231. % 已知的相关输入数据
  232.     figure
  233.     subplot(3,1,1)
  234.     plot(PtCSP_fore * baseMVA,'m-o');
  235.   title('CSP预测功率值')
  236.   xlabel('时间(h)');
  237.     ylabel('功率(MW)');
  238.     
  239.     subplot(3,1,2)
  240.     plot(P_RE(1,:) * baseMVA,'m-o'); hold on
  241.     plot(P_RE(2,:) * baseMVA,'b-s');
  242.   title('可再生能源预测出力值')
  243.   xlabel('时间(h)');
  244.     ylabel('功率(MW)');
  245.     legend('光伏','风电')
  246.     
  247.     subplot(3,1,3)
  248.     plot(sum(PD) * baseMVA,'r-v');
  249.   title('24h负荷值')
  250.   xlabel('时间(h)');
  251.     ylabel('功率(MW)');
  252.     
  253.     
  254. ​
  255.     
  256. %    subplot(2,1,2)
  257. %  bar(baseMVA*PG_RE',0.75,'stack'); hold on; % 各PV、Wind机组出力
  258. %    legend('PV','Wind')
  259. %    title('电网中可再生能源机组出力')
  260. %  xlabel('时间(h)');
  261. %    ylabel('功率(MW)');
  262.     
  263. %    figure
  264. %    surf(baseMVA*PC_Load);
  265. %    title('负荷削减量')
  266. %  xlabel('时间(h)');
  267. %    ylabel('功率(MW)');
  268. ​
  269. ​

**3 **部分结果

4 下载链接

标签： matlab 光热电站 优化调度
本文转载自: https://blog.csdn.net/superone89/article/details/135855432
版权归原作者 科研工作站 所有， 如有侵权，请联系我们删除。
登录后发布评论

“Matlab|含风电-光伏-光热电站电力系统N-k安全优化调度模型”的评论:

还没有评论
关于作者
...
overfit同步小助手
文章同步
相关阅读

网络安全法-网络运行安全

使用selenium/drissionpage时如何阻止chrome自动跳转http到https

docker desktop 里部署的Open WebUI 管理员密码忘记了的处理方法

在ubuntu20.04中搭建onsite比赛运行环境

利用开源的低代码表单设计器FcDesigner高效管理和渲染复杂表单结构

Kafka学习笔记

【前端】浏览器输入url到页面呈现发生了什么？