1 主要内容

该程序复现《A Modified Binary PSO to solve the Thermal Unit Commitment Problem》第五章内容，主要做的是一个考虑需求响应的机组组合问题，首先构建了机组组合问题的基本模型，在此基础上，进一步考虑负荷侧管理，也就是需求响应，在调控过程中通过补偿引导负荷侧积极进行需求响应，在模型的求解上，采用了一种基于改进二进制粒子群算法的求解方法，相较于传统的粒子群算法，更加创新，而且求解的效果更好，代码出图效果非常好。该程序函数比较多，主函数为Swarm_generator，运行结果已经保存在Graphs文件夹内部，可以通过运行Graphs.m直接得到出图结果。程序采用matlab编程，注释为英文，适合具有编程经验的同学下载学习！

2 部分代码

N=size(I,1); %%%Number of TGU's
T=24;        %%%Study period duration

%%%%fmincon.m solver options for F_LIM_ED
OPTS = optimoptions('fmincon','Algorithm','sqp','display','off','ConstraintTolerance',1e-2,'OptimalityTolerance',1.0000e-02,'MaxIterations',100, 'StepTolerance',1e-2,'MaxFunctionEvaluations',100);

%%%%%%%%%%%%%%%%%%%%%%%%%%UNIT DATA%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%Generator limits
PGI_MAX=I(:,1); %%%Upper generation limit for each generator
PGI_MIN=I(:,2); %%%Lower generation limit for each generator

%%%I/O curve for each generator is modeled by a smooth quadratic function 
%%%i.e C(Pgi)=ai*Pgi^2+bi*Pgi+ci 
ai=I(:,3);
bi=I(:,4);
ci=I(:,5);

[I_C_SORT_EXS,IDX_EXS,I_C_SORT_INS,IDX_INS] = AFLC(ai,bi,ci,PGI_MAX); %%%Unit priority list assuming same fuel cost for every unit

%%MINIMUM TIME DATA
MUT=I(:,6);     %%%Minimum up-time for each generator
MDT=I(:,7);     %%%Minimum down-time for each generator

%%START-UP AND SHUT-DOWN COSTS
SU_H=I(:,8);
SU_C=I(:,9);
CSH=I(:,10);
init_status=I(:,11);

%%%%OTHER CONSTRAINTS
M_R=I(:,13);    %%MR units 
U_N=I(:,14);    %%Unavailable units

%%%% MR & UN hours & units %%Hr begin  Hr end
MR=[  10          15  ];
UN=[                  ];

%%%% Prohibitive operative zones (POZ)%%%Unit    %%%POZ   
POZ=[  [1  150 165; 1  448,453][2  90  110; 2  240 250][8  20   30; 8  40  45][10 12   17; 10 35  45]];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%MAIN LOOP%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%Initiaize solution (NxT) schedule solution as an empty array
Init_SOL=zeros(N,T);

complete=ones(1,T); %%%Nth unit operation schedule is ON for all t
sum_Pgi_max=sum(PGI_MAX);

%%%%%Define swarm size (Either particle or bee)
SWARM_SZ=20;
PART_BEE=0;
YES_CNT=0;
NO_CNT=0;

population=zeros(SWARM_SZ,N*T);
total_COST=zeros(SWARM_SZ,1);

while PART_BEE0); %%%Find the exception indices for must run units
UN_idx=find(U_N>0); %%%Find the exception indices for unavailable units 

%%%%%%%%%%%%Determine initial status of units %%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

Tdiff_UP=zeros(N,3);%%%Exception indices for coupling constraints Tdiff_UP
Tdiff_DW=zeros(N,3);%%%Exception indices for coupling constraints Tdiff_DW

for k=1:NTdiff_UP(k,1)=k;  Tdiff_DW(k,1)=k;  if ismember(k,MR_idx) || ismember(k,UN_idx) %%%If must run or unavailable unit, continue schedule repair without modifying INIT_SOLcontinue
elseP=1+(PGI_MAX(k)/sum_Pgi_max);if init_status(k)>0 %%Unit Online for ti=init_status hoursTdiff_U=init_status(k)-MUT(k);if Tdiff_U<0Tdiff_UP(k,2)=abs(Tdiff_U);Tdiff_UP(k,3)=1;Init_SOL(k,1:Tdiff_UP(k,2))=1;                else Tdiff_UP(k,2)=0;Tdiff_UP(k,3)=0;end
for j=Tdiff_UP(k,2)+1:TX=-1+(P-(-1))*rand;if X>=0Init_SOL(k,j)=1;elseInit_SOL(k,j)=0;end    end
else             %%Unit is offline for ti=init_status hoursTdiff_D=MDT(k)+init_status(k);if Tdiff_D>0Tdiff_DW(k,2)=abs(Tdiff_D);Init_SOL(k,1:Tdiff_D)=0;Tdiff_DW(k,3)=1;  else Tdiff_DW(k,2)=0;Tdiff_DW(k,3)=0;end
for j=Tdiff_DW(k,2)+1:TX=-1+(P-(-1))*rand;if X>=0Init_SOL(k,j)=1;elseInit_SOL(k,j)=0;end    endendend
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%Check & Repair Coupling constraints %%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%Decompose commitment schedule into T single-hour unit combinations 

Exs_On_cap=zeros(3,T); %%Excess online capacity
Ins_On_cap=zeros(3,T); %%Insufficient online capacity

f_unit_ins_all=zeros(N,T); %%%This unit is not part of the exceptions MR,UN,Initialization

free_unit_ins=zeros(N,T);  %%%These units are offline and available to be turned off 

%%Check Excess online capacity first for full commitment
check_full_exs=ones(1,N)*PGI_MIN;     %%%Check if system has excess online capacity for a complete commitment
check_1=find(check_full_exs>P_D);
if isempty(check_1)clear excp_exs f_unit_exs_all free_unit_exs Exs_On_capelse %%%Excess online capacity schedule repairfor t=1:Tu_t=Init_SOL(:,t); %%Decompose into T single hour unit combinations
Exs_On_cap(1,t)=u_t'*PGI_MIN;Exs_On_cap(2,t)=P_D(t);
if Exs_On_cap(1,t)>P_D(t)Exs_On_cap(3,t)=1;format_exs='Excess capacity at hour %u. Repairing schedule...\n';fprintf(format_ins,t);  endend end     

status_e=0;
status_f=0;
status_g=0;

status=[status_e,status_f,status_g];
trials=0;
max_trials=15;
while ~all(status) trials=trials+1;if trials>max_trialsbreakendif trials==1%%%Insufficient online capacity initial schedule repairfor t=1:Tu_t=Init_SOL(:,t); %%Decompose into T single hour unit combinations
Ins_On_cap(1,t)=u_t'*PGI_MAX;Ins_On_cap(2,t)=P_D(t)+SRREQ(t);
if Ins_On_cap(1,t)