电力系统中机组组合强大的Yalmip+Cplex（入门学习+机组组合问题Matlab实现）

目录

​​1、概述​​

​​（1）电力调度优化理论及其应用​​

​​（2）电力调度的机组组合UC ​​

​​（3）组合模型及算法 ​​

​​2、人性化的Yalmip+Cplex ​​

​​（1）热身运动——线性规划​​

​​3、Yalmip+Cplex 解决机组组合（YYDS）​​

​​3.1 算例​​

​​3.2 目标函数​​

​​ 3.3 约束条件​​

​​3.4 模型约束条件处理​​

​​3.5 Matlab代码实现​​

1、概述

（1）电力调度优化理论及其应用

（2）电力调度的机组组合UC 

（3）组合模型及算法 

2、人性化的Yalmip+Cplex 

（1）热身运动——线性规划

         尽管 python 比较火，cplex 对 python 的支持目前还不是太全，相关的学习资料比较少，ibm 自己 出的资料对 python 包的介绍也很简略，例子及相关类方法的介绍也不详细，这一点远没有对 java 或 c++ 支持地好。cplex 在 python 中没有重载加减乘除符号吗？目前给的例子都是输入系数矩阵这种形式。使用起来非常不方便的。

         所以推出了Yalmip+Cplex（matlab），matlab中这个库的下载方法，后面的文章会讨论matlab中Yalmip+Cplex这个强大的库在我们电力系统中的强大作用（噢头曼！！!)，虽然我很喜欢用Python，但是python是这几年才火起来，所以这个Cplex库还不是很完善，所以对于这个库强烈推荐matlab。下面先来点简单的线性规划，让大家热身。

                                  

%（1）设定决策变量X（1）、X（2）
%（2）sdpvar：实数变量；binvar：0—1变量；intvar：整型变量
%（3）Yalmip默认是对称的，要求非对称用full
 
x=sdpvar(2,1,'full');   
z=-2*x(1)+4*x(2);         %目标函数，默认最小
st=[];                    %设定一系列约束
st=[st,-3*x(1)+x(2)<=6];
st=[st,x(1)+2*x(2)>=4];
st=[st,x(1)+3*x(2)==4];
st=[st,x(2)>=-3];
ops=sdpsettings('solver','cplex');  %设定求解器
r=optimize(st,z);                   %求解，如果最大值用-z
x=value(x);                         %查看求解结果x的值
z=value(z);                         %查看目标函数的最优解

%结果
Optimal solution found.
 
>> x
 
x =
 
    13
    -3
 
>> z
 
z =
 
   -38

 ​

3、Yalmip+Cplex解决机组组合（YYDS）

3.1 算例

基于已知的系统数据，求解计划时间内的6机30节点41支路的功率情况与机组的开停机情况，使得系统总成本达到最小。该问题的决策变量由两类，第一类是各时段机组的出力，为连续变量。第二类是各时段机组的启停状态,为整数变量，0表示关停,1表示启动。
本问题属于混合整数规划(MIP）问题，即要在决策变量的可行解空间里找到一组最优解，使得目标函数尽可能取得极值。对于混合整数规划，CPLEX提供了快速的MIP求解方法。

3.2 目标函数

 3.3 约束条件

 

3.4 模型约束条件处理

模型中有许多约束条件，这些约束条件不能用程序直接表达，只能用变换来表达。

3.5 Matlab代码实现

本文只给出部分代码​

%=====原二次函数式的最优成本目标==========
% for i=1:gennum
%     for t=1:T
%     totalcost=totalcost+para(i,1)*p(i,t).^2+para(i,2)*p(i,t)+para(i,3)*u(i,t);  %煤耗成本
%     totalcost=totalcost+costH(i,t);                                %启动成本
%     totalcost=totalcost+costJ(i,t);                                %关停成本
%     end
% end
%% 负荷平衡约束
for t=1:T
    st=st+[sum(p(:,t))==PL(1,t)];%负荷平衡约束;
end
%% 机组出力上下限约束
for t=1:T
    for i=1:gennum
        st=st+[u(i,t)*limit(i,2)<=p(i,t)<=u(i,t)*limit(i,1)];%机组出力上下限约束
    end
end
%% 机组爬坡约束
%按下式进行推导编程
% %启动最大升速率
% Su=(Pmax+Pmin)/2;
% %停机最大降速率
% Sd=(Pmax+Pmin)/2;
%Ru=Rud;Rd=Rud;
% %上爬坡约束
% for t=2:T
% st=st+[p(:,t)-p(:,t-1)<=u(:,t-1).*(Ru-Su)+Su];
% end
% %下爬坡约束
% for t=2:T
%st=st+[p(:,t-1)-p(:,t)<=u(:,t).*(Rd-Sd)+Sd];
% end
%展开表达式：
for t=2:T
    for i=1:gennum
        % st=st+[-Rud(i,1)*u(i,t)+(u(i,t)-u(i,t-1))*limit(i,2)-limit(i,1)*(1-u(i,t))<=p(i,t)-p(i,t-1)];
        % st=st+[p(i,t)-p(i,t-1)<=Rud(i,1)*u(i,t-1)+(u(i,t)-u(i,t-1))*limit(i,2)+limit(i,1)*(1-u(i,t))];
        %由于原式可能关机以后就无法再开动了，改用下式
        st=st+[p(i,t-1)-p(i,t)<=Rud(i,1)*u(i,t)+(1-u(i,t))*(limit(i,2)+limit(i,1))/2];%下坡
        st=st+[p(i,t)-p(i,t-1)<=Rud(i,1)*u(i,t-1)+(1-u(i,t-1))*(limit(i,2)+limit(i,1))/2];%上坡
    end
end
%% 热备用约束
hp=0.05;%热备用系数
for t=1:T
    st=st+[sum(u(:,t).*limit(:,1)-p(:,t))>=hp*PL(1,t)];
end
%% 启停时间约束

Tried aggregator 2 times.
MIP Presolve eliminated 3840 rows and 23628 columns.
MIP Presolve modified 5356 coefficients.
Aggregator did 10 substitutions.
Reduced MIP has 1142 rows, 1130 columns, and 3962 nonzeros.
Reduced MIP has 144 binaries, 0 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.05 sec. (24.83 ticks)
Probing time = 0.00 sec. (0.39 ticks)
Tried aggregator 1 time.
Reduced MIP has 1142 rows, 1130 columns, and 3962 nonzeros.
Reduced MIP has 144 binaries, 230 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.02 sec. (2.70 ticks)
Probing time = 0.00 sec. (0.37 ticks)
MIP emphasis: balance optimality and feasibility.
MIP search method: dynamic search.
Parallel mode: deterministic, using up to 16 threads.
 
Node log . . .
Best integer =   1.527930e+07  Node =       0  Best node =   5.968625e+04
Best integer =   1.080804e+07  Node =       0  Best node =   5.968625e+04
Best integer =   6.975152e+06  Node =       0  Best node =   4.781375e+06
Best integer =   6.024695e+06  Node =       0  Best node =   6.024695e+06
GUB cover cuts applied:  10
Cover cuts applied:  13
Implied bound cuts applied:  32
Flow cuts applied:  23
Mixed integer rounding cuts applied:  54
Zero-half cuts applied:  7
Gomory fractional cuts applied:  4
 
ans =
 
   6.0247e+06