MATLAB: Solving the Unit Commitment Problem using optimproblem

economic load dispatchintlinprogMATLABmilpoptimizationoptimproblemunit commitment

I am trying to solve the Unit Commitment problem using the optimproblem function. I will explain my problem with an example. Lets say I have three generators which will be operated to meet various loads for a 24 hours period. Now the cost curves for these units are quadratic, so piecewise linearization is required if I am going to use the intlinprog function. I have successfully linearized the system and the problem now involves minimizing:
     Total C(P) = C_1(P_1min)T_1 + s_11*P_11*T_1 + s_12*P_12*T_1 + … + s_1n*P_1n*T_1
            + C_2(P_2min)T_2 + s_21*P_21*T_2 + s_22*P_22*T_2 + … + s_2n*P_2n*T_2
            + C_3(P_3min)T_3 + s_31*P_31*T_3 + s_32*P_32*T_3 + … + s_3n*P_3n*T_3
Where Sik is the slope for each segment of the linearized curve, Pik are the power increments (which I need to determine) for each segment. T_i is a binary variable which keeps track of the state of unit and Ci(Pimin) is the cost at Pimin. Now based on my research, the problem should now be formulated as minimize:
    Total C(P) = C_1(P_1min)T_1 + s_11*Z_11 + s_12*Z_12 + … + s_1n*Z_1n
            + C_2(P_2min)T_2 + s_21*Z_21 + s_22*Z_22 + … + s_2n*Z_2n
            + C_3(P_3min)T_3 + s_31*Z_31 + s_32*Z_32 + … + s_3n*Z_3n
Where Ti = 0 implies that Zik = 0 and Ti = 1 implies that Zik = Pik.
Now I have two optimization variable:
  1. Pik
  2. Ti
With the optimization problem being made up of two part:
  1. Sum of ALL Sik*Pik values
  2. Sum of Ci(Pimin)*Ti values.
The following code illustrates this (with the addition of another optimization variable which isn't so important):
 %Amount of power generated in an hour by a plant
power = optimvar('power',nHours,plant,'LowerBound',0,'UpperBound',maxGenConst);
%Indicator if plant is operating during an hour 
isOn = optimvar('isOn',nHours,actual_plant,'Type','integer','LowerBound',0,'UpperBound',1);
%Indicator if plant is starting up during an hour
startup = optimvar('startup',nHours,actual_plant,'Type','integer','LowerBound',0,'UpperBound',1);
    %Costs
powerCost = sum(power*f,1);
isOnCost = sum(isOn*OperatingCost',1);
startupCost = sum(startup*StartupCost',1);
 % set objective
powerprob.Objective = powerCost + isOnCost + startupCost;
   % power generation over all plans meets hourly demand
powerprob.Constraints.isDemandMet = sum(power,2) >= Load;
    % only gen power when plant is on
% if isOn=0 power must be zero, if isOn=1 power must be less than maxGenConst
powerprob.Constraints.powerOnlyWhenOn = power <= maxGenConst.*isOn;   
    % if on, meet MinGenerationLevel
% if isOn=0 power >= 0, if isOn=1 power must be more than minGenConst
powerprob.Constraints.meetMinGenLevel = power >= maxGenConst.*isOn;
My problem lies in this calculation: maxGenConst.*isOn. These have different sizes as maxGenConst is a nHours * (num_of_gen*num_of_segments) whereas isOn is a nHour * num_of_gen. The reason this occurs is because after linearization, the cost equation becomes (num_of_segments + 1) equations. I would appreciate some help in formulating the problem so I can use the optimproblem function. I have attached my full code and data.

Best Answer

I was able to solve my problem. I made maxGenConst the same size as isOn. Then I did the dot multiplication with isOn and was able to index the optimization expression and reshape into the proper dimensions. I have attached my code with data for a test system for anyone who will need it. This code considers power generation constraints, minimum up and down times, start up costs and the basic load demand constraints. I have included the data for a 6 unit system with the necessary data and the load demand over a 24hr period. This code does not consider losses in the system.
Related Question