Getting Infeasibility row when the constraints function is multiplied by a constant parameter.

LIESEA LEUNG

unread,

Mar 28, 2019, 10:11:38 AM3/28/19

to YALMIP

Dear Dr.Löfberg

I can get the answer after adding this constriants.

for x3=1:6
    constraints=[constraints,10<=(sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12))))<=Z3(x3,1)];
end

but when I multiply a parameter in this constraint function. It informs that 'Infeasibility row 'c697': 0 <= -10.'

para=T*0.92*P_c_peb
for x3=1:6
    constraints=[constraints,10<=para*(sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12))))<=Z3(x3,1)];
end

Thank you very much if you could reply.

Johan Löfberg

unread,

Mar 28, 2019, 1:32:01 PM3/28/19

to YALMIP

sounds like the rhs evaluates to 0, i.e T*0.92*P_c_peb=0

Message has been deleted

LIESEA LEUNG

unread,

Mar 28, 2019, 8:54:15 PM3/28/19

to YALMIP

Thx, i define para=30*0.92*117, it could get the answer. But when I want to simplify my code, I rewrite the code

Sn=324;
SOC_use=normrnd(70,7,[6,8]);
Z1=(sum(SOC_use,2)/Sn-0.3+0.2)*Sn;
Z3=(sum(SOC_use(:,1:7),2)/Sn+1-0.3)*Sn;
para=30*0.92*20;
Z2=sdpvar(6,1);
for x3=1:6
    for y3=1:6
        Z2(x3,1)=para*sum(C(x3,P2(y3*2-1):P2(y3*2)))+Z2(x3,1);
    end
end
constraints=[constraints, 10<=Z2<=Z3];

to replace the code

for x3=1:6
    constraints7=[constraints7,10<=30*0.92*20*(sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12))))<=Z3(x3,1)];
end

These two different code get different Optimal solution.

Johan Löfberg

unread,

Mar 29, 2019, 2:52:32 AM3/29/19

to YALMIP

way too messy for anyone but you to see where your bug is (unless you've run this twice, then of course you get different solutions as you randomize)

LIESEA LEUNG

unread,

Mar 29, 2019, 3:50:19 AM3/29/19

to YALMIP

My idea is that the code

Z2=sdpvar(6,1)
for x3=1:6
    for y3=1:6
        Z2(x3,1)=para*sum(C(x3,P2(y3*2-1):P2(y3*2)))+Z2(x3,1);
    end
end
constraints=[constraints, 10<=Z2<=Z3];

is equivalent to the second one.

for x3=1:6
    constraints7=[constraints7,10<=para*(sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):  P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12))))<=Z3(x3,1)];
end

So the solution in File 1 and File 2 (different in #CONSTRAINTS6) should be smiliar. In fact, only the File 1 can get the result I want.

test_n6_feasible.m

test_n6_no_feasible.m

Johan Löfberg

unread,

Mar 29, 2019, 4:02:10 AM3/29/19

to YALMIP

You're simply not defining the same thing, i.e., your vectorized code i flawed

Add

themiddle = [];
theright = [];
for x3=1:6
    themiddle = [themiddle;para*(sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12))))];
    theright = [theright;Z3(x3,1)];    
end

and then compare those with your vectorization, and the expression in the middle is not the same (the first coefficient differens for instance)

K>> sdisplay(Z2)

ans =

  6×1 cell array

    {'1656*C(1)+552*C(67)+552*C(73)+552*C(79)+552*C(85)+552*C(91)+552*C(97)+internal(6080)' }
    {'1656*C(2)+552*C(68)+552*C(74)+552*C(80)+552*C(86)+552*C(92)+552*C(98)+internal(6081)' }
    {'1656*C(3)+552*C(69)+552*C(75)+552*C(81)+552*C(87)+552*C(93)+552*C(99)+internal(6082)' }
    {'1656*C(4)+552*C(70)+552*C(76)+552*C(82)+552*C(88)+552*C(94)+552*C(100)+internal(6083)'}
    {'1656*C(5)+552*C(71)+552*C(77)+552*C(83)+552*C(89)+552*C(95)+552*C(101)+internal(6084)'}
    {'1656*C(6)+552*C(72)+552*C(78)+552*C(84)+552*C(90)+552*C(96)+552*C(102)+internal(6085)'}

K>> sdisplay(themiddle)

ans =

  6×1 cell array

    {'552*C(1)+552*C(7)+552*C(13)+552*C(19)+552*C(25)+552*C(31)+552*C(37)+552*C(43)+552*C(49)+552*C(55)+552*C(61)+552*C(67)+552*C(85)+552*C(91)+552*C(97)+552*C(103)+552*C(127)+552*C(133)+552*C(139)+552*C(163)+552*C(169)+552*C(175)+552*C(193)+552*C(199)+552*C(205)+552*C(211)+552*C(235)+552*C(241)+552*C(247)+552*C(253)+552*C(259)+552*C(265)+552*C(271)+552*C(277)+552*C(283)'            }
    {'552*C(2)+552*C(8)+552*C(14)+552*C(20)+552*C(26)+552*C(32)+552*C(38)+552*C(44)+552*C(50)+552*C(56)+552*C(62)+552*C(68)+552*C(74)+552*C(98)+552*C(104)+552*C(110)+552*C(134)+552*C(140)+552*C(146)+552*C(164)+552*C(170)+552*C(176)+552*C(182)+552*C(206)+552*C(212)+552*C(218)+552*C(242)+552*C(248)+552*C(254)+552*C(260)+552*C(266)+552*C(272)+552*C(278)+552*C(284)'                      }
    {'552*C(3)+552*C(9)+552*C(15)+552*C(21)+552*C(27)+552*C(33)+552*C(39)+552*C(45)+552*C(51)+552*C(57)+552*C(63)+552*C(69)+552*C(75)+552*C(81)+552*C(105)+552*C(111)+552*C(117)+552*C(141)+552*C(147)+552*C(153)+552*C(177)+552*C(183)+552*C(189)+552*C(207)+552*C(213)+552*C(219)+552*C(225)+552*C(249)+552*C(255)+552*C(261)+552*C(267)+552*C(273)+552*C(279)+552*C(285)'                      }
    {'552*C(4)+552*C(10)+552*C(16)+552*C(22)+552*C(28)+552*C(34)+552*C(40)+552*C(46)+552*C(52)+552*C(58)+552*C(64)+552*C(70)+552*C(76)+552*C(82)+552*C(88)+552*C(106)+552*C(112)+552*C(118)+552*C(124)+552*C(148)+552*C(154)+552*C(160)+552*C(178)+552*C(184)+552*C(190)+552*C(196)+552*C(214)+552*C(220)+552*C(226)+552*C(232)+552*C(256)+552*C(262)+552*C(268)+552*C(274)+552*C(280)+552*C(286)'}
    {'552*C(5)+552*C(11)+552*C(17)+552*C(23)+552*C(29)+552*C(35)+552*C(41)+552*C(47)+552*C(53)+552*C(59)+552*C(65)+552*C(71)+552*C(77)+552*C(83)+552*C(89)+552*C(95)+552*C(119)+552*C(125)+552*C(131)+552*C(149)+552*C(155)+552*C(161)+552*C(167)+552*C(191)+552*C(197)+552*C(203)+552*C(221)+552*C(227)+552*C(233)+552*C(239)+552*C(263)+552*C(269)+552*C(275)+552*C(281)+552*C(287)'            }
    {'552*C(6)+552*C(12)+552*C(18)+552*C(24)+552*C(30)+552*C(36)+552*C(42)+552*C(48)+552*C(54)+552*C(60)+552*C(66)+552*C(72)+552*C(78)+552*C(84)+552*C(90)+552*C(96)+552*C(102)+552*C(120)+552*C(126)+552*C(132)+552*C(138)+552*C(162)+552*C(168)+552*C(174)+552*C(198)+552*C(204)+552*C(210)+552*C(234)+552*C(240)+552*C(246)+552*C(270)+552*C(276)+552*C(282)+552*C(288)'                       }

LIESEA LEUNG

unread,

Mar 29, 2019, 4:49:25 AM3/29/19

to YALMIP

Thank you, I figured it out. But if I want this simplify this code,

sum(C(x3,P2(x3,1):P2(x3,2)))+sum(C(x3,P2(x3,3):P2(x3,4)))+sum(C(x3,P2(x3,5):P2(x3,6)))+sum(C(x3,P2(x3,7):P2(x3,8)))+sum(C(x3,P2(x3,9):P2(x3,10)))+sum(C(x3,P2(x3,11):P2(x3,12)))

How should I write.

Johan Löfberg

unread,

Mar 29, 2019, 4:52:19 AM3/29/19

to YALMIP

Using MATLAB indexing acrobatics as with any MATLAB code. How: not obvious

LIESEA LEUNG

unread,

Mar 29, 2019, 4:57:48 AM3/29/19

to YALMIP

Ok, again, thank you very much for your reply.

LIESEA LEUNG

unread,

Mar 29, 2019, 12:12:45 PM3/29/19

to YALMIP

Dear Löfberg

I have a new problem. I test the same code but have different parameters.

The first is

options=sdpsettings('solver','cplex');
N=6;T=30;K=24*60/T;M=3;a=5/100;y=50;
C=binvar(N,K);
P_peak=sdpvar;
P_c_peb=54;
for time=1:K
    if time<=16;
        load(1,time)=400;
    elseif time>16&&time<=26
        load(1,time)=900;
    elseif time>26&&time<=32
        load(1,time)=1000;
    elseif time>32&&time<=34
        load(1,time)=600;
    elseif time>34&&time<=38
        load(1,time)=1000;
    elseif time>38&&time<=44
        load(1,time)=1350; 
    elseif time>44&&time<=46
        load(1,time)=800;
    elseif time>46
        load(1,time)=400;
    end
end
for k =1:K       %eletricity price
    if k<=14||(k<=K&&k>46)
       price(k,1)=0.3818;
    elseif (k>14&&k<=20)||(k>30&&k<=36)||(k>42&&k<=46)
       price(k,1)=0.8395;
    elseif (k>20&&k<=30)||(k>36&&k<=42)
       price(k,1)=1.3222;
    end
end
constraints=[P_peak>=sum(C*P_c_peb)+load];
constraints=[constraints,sum(C)<=M];
%objective function
EPC=sum(C*P_c_peb*T*price)/60;   %EPC
f1=K*T/365/24/60;f2=a*(1+a)^y/((1+a)^y-1);
ECC=f1*f2*P_peak*14847;
F=EPC+ECC;  

optimize(constraints,F);
value(C)

and above code get the answer I want.

i.e. value(C)= 0 0 0 0 0 0 0

but when I change the N,K

options = sdpsettings('solver','cplex');   %cplex
N=24;T=5;K=24*60/T;a=5/100;y=50;M=10;
C=binvar(N,K);
P_peak=sdpvar;
P_c_peb=100;
for k =1:K       
    if k<=84||(k>276&&k<=288)
       price(k,1)=0.3818;
    elseif (k>84&&k<=120)||(k>180&&k<=216)||(k>252&&k<=276)
       price(k,1)=0.8395;
    elseif (k>120&&k<=180)||(k>216&&k<=252)
       price(k,1)=1.3222;
    end
end
for j=1:K        %load
    if j<=48
        load(1,j)=200;
    elseif j>48&&j<=108
        load(1,j)=400;
    elseif j>108&&j<=144
        load(1,j)=800;
    elseif j>144&&j<=156
        load(1,j)=1000;
    elseif j>156&&j<=192
        load(1,j)=800;
    elseif j>192&&j<=228
        load(1,j)=1000;        
    elseif j>228&&j<=264
        load(1,j)=1400;
    elseif j>264&&j<=288
        load(1,j)=400;
    end
end
constraints=[P_peak>=sum(C*P_c_peb)+load];  
constraints=[constraints,sum(C)<=M]; 

%objective function
EPC=sum(C*P_c_peb*T*price)/60;   %EPC
f1=K*T/365/24/60;f2=a*(1+a)^y/((1+a)^y-1);
ECC=f1*f2*P_peak*14847;
F=EPC+ECC;

optimize(constraints,F);
value(C)

the answer are

value(C)= NaN NaN NaN NaN NaN NaN......

why after changing the parameters, the answer are different.

Thanks a lot for your reply.

Johan Löfberg

unread,

Mar 29, 2019, 12:16:55 PM3/29/19

to YALMIP

well it looks like you are not even looking at the diagnostic output from optimize,so the obvious guess is that it is infeasible and you thus haven't detected that

LIESEA LEUNG

unread,

Mar 29, 2019, 12:31:14 PM3/29/19

to YALMIP

Sorry for my poor knowledge, I just start learning the yalmip, and I don't understand what diagnostic output from optimize means.

And I mean, under these constraints function,

constraints=[P_peak>=sum(C*P_c_peb)+load];  
constraints=[constraints,sum(C)<=M];

my objective solution are value(C)=0 0 0 0 0.......

But only the first one can get the objective solution. And the difference between the two is only on the parameters.

Johan Löfberg

unread,

Mar 29, 2019, 12:48:28 PM3/29/19

to YALMIP

First example here

https://yalmip.github.io/tutorial/basics/

LIESEA LEUNG

unread,

Mar 29, 2019, 12:56:35 PM3/29/19

to YALMIP

Thanks！

Reply all

Reply to author

Forward