Re: Some issues with 'implies'

[Johan Löfberg]

In your "if x==1 then y == 0" it is not clear if you intend "if all(x==1) then y==0" or if it is a vectorized notation for "if x(i)==1 then y(i)==0"

...and then I don't see (or understand?) any actual question either

[ran tao]

thanks for your reply

i am sorry i did not make it clear.(may due to my poor English)

i want to generator the four binary variables as shown in the picture.

i think it may be a vectorized notation for "if x(i)==1 then y(i)==0"

the description of constraints is:

1, if u_g_dn is 1 to 0 from time t-1 to t, the u_g_cd equal to 1 at time t. if not, the u_g_cd equal to 0

2, if u_g_cd is 1 to 0 from time t-1 to t, the u_g_ht or u_g_dn equal to 1 at time t. if not, they equal to 0

3, if u_g_ht is 1 to 0 from time t-1 to t, the u_g_cd or u_g_ss equal to 1 at time t. if not, they equal to 0

4, if u_g_ss is 1 to 0 from time t-1 to t, the u_g_ht  equal to 1 at time t. if not, they equal to 0

and, that is my code

clc,clear;

%% 预定义参数

Horizon = 96;

Load = 100 + 60*sin((1:Horizon)*2*pi/96);

ngen = 5;

% 各状态成本

c_g_dn = [50;20;10;45;60];

c_g_cd = [0;0;0;0;0];

c_g_ht = c_g_dn + 20;

c_g_ss = c_g_dn + 30;

% 各状态最小持续时间

minup_cd = [4;3;5;3;5];

minup_ht = [1;1;2;1;2];

minup_ss = [6;6;8;6;8];

%% 决策变量

% 就绪状态（停机状态）

u_g_dn = binvar(ngen, Horizon);

% 冷状态

u_g_cd = binvar(ngen, Horizon);

% 热状态

u_g_ht = binvar(ngen, Horizon);

% 旋转状态（启动状态）

u_g_ss = binvar(ngen, Horizon);

% 机组出力

x_P_g = sdpvar(ngen, Horizon);

%% 约束条件

cons = [];

% cons = [cons, u_g_ss(4,20)==1];

cons = [cons, sum(x_P_g, 1) == Load];

cons = [cons, u_g_ss .* 10 <= x_P_g <= u_g_ss .* 32];

for t = 2:Horizon

    for k = 1:ngen

        % 升级

        cons = [cons, implies(u_g_ht(k,t-1)-u_g_ht(k,t)==1, u_g_ss(k,t) + u_g_ss(k,t)==1)];

        cons = [cons, implies(u_g_cd(k,t-1)-u_g_cd(k,t)==1, u_g_ht(k,t) + u_g_dn(k,t)==1)];

        cons = [cons, implies(u_g_dn(k,t-1)-u_g_dn(k,t)==1, u_g_cd(k,t)==1)];

        cons = [cons, implies(u_g_ss(k,t-1)-u_g_ss(k,t)==1, u_g_ht(k,t)==1)];

%         % 降级

%         cons = [cons, implies(u_g_ss(:,t-1)-u_g_ss(:,t)==1, u_g_ht(:,t)==1)];

%         cons = [cons, implies(u_g_ht(:,t-1)-u_g_ht(:,t)==1, u_g_cd(:,t)==1)];

%         cons = [cons, implies(u_g_cd(:,t-1)-u_g_cd(:,t)==1, u_g_dn(:,t)==1)];

%          % 升级

%         cons = [cons, implies(u_g_ss(:,t)==1, u_g_ht(:,t-1)-u_g_ht(:,t)==1)];

%         cons = [cons, implies(u_g_ht(:,t)==1, u_g_cd(:,t-1)-u_g_cd(:,t)==1)];

%         cons = [cons, implies(u_g_cd(:,t)==1, u_g_dn(:,t-1)-u_g_dn(:,t)==1)];

%         % 降级

%         cons = [cons, implies(u_g_ht(:,t)==1, u_g_ss(:,t-1)-u_g_ss(:,t)==1)];

%         cons = [cons, implies(u_g_cd(:,t)==1, u_g_ht(:,t-1)-u_g_ht(:,t)==1)];

%         cons = [cons, implies(u_g_dn(:,t)==1, u_g_cd(:,t-1)-u_g_cd(:,t)==1)];

    end

end

cons = [cons, u_g_dn + u_g_cd + u_g_ht + u_g_ss <= 1];

%冷状态最小持续时间

for k = 2:Horizon

     for unit = 1:ngen

      indicator = u_g_cd(unit,k) - u_g_cd(unit,k-1);

      range = k:min(Horizon,k+minup_cd(unit)-1);

      cons = [cons, u_g_cd(unit,range) >= indicator];

     end

end

%热状态最小持续时间

for k = 2:Horizon

     for unit = 1:ngen

      indicator = u_g_ht(unit,k) - u_g_ht(unit,k-1);

      range = k:min(Horizon,k+minup_ht(unit)-1);

      cons = [cons, u_g_ht(unit,range) >= indicator];

     end

end

%旋转状态最小持续时间

for k = 2:Horizon

     for unit = 1:ngen

      indicator = u_g_ss(unit,k) - u_g_ss(unit,k-1);

      range = k:min(Horizon,k+minup_ss(unit)-1);

      cons = [cons, u_g_ss(unit,range) >= indicator];

     end

end

%% 目标函数

F = u_g_dn .* repmat(c_g_dn,1, Horizon) + u_g_cd .* repmat(c_g_cd, 1, Horizon) + u_g_ht .* repmat(c_g_ht, 1, Horizon) + x_P_g .* repmat(c_g_ss, 1, Horizon);

F = sum(sum(F));

%%优化

optimize(cons, F, sdpsettings('solver', 'gurobi'));

u_g_dn = value(u_g_dn);

u_g_cd = value(u_g_cd);

u_g_ht = value(u_g_ht);

u_g_ss = value(u_g_ss);

x_P_g = value(x_P_g);

%% 画图

for ii = 1:ngen

    subplot(ngen,1,ii)

    stairs(u_g_dn(ii,:),'k','LineWidth',2);

    hold on

    stairs(u_g_cd(ii,:),'b','LineWidth',2);

    hold on

    stairs(u_g_ht(ii,:),'y','LineWidth',2);

    hold on

    stairs(u_g_ss(ii,:),'r','LineWidth',2);

    hold on

    ylim([-0.5,1.5]);

    xlim([0,97]);

%     xticks(1:97);

%     xticklabels([1:3:97]);

end

[Johan Löfberg]

implications involving equalities are nasty to model automatically

since this is such an easy model if x-y == 1 then z == 1, you should do it manually, z >= 1 + (x - y - 1)