Error using svd

1

I am using SDP relaxation to solve a non linear problem.

However, every time I am running my MATLAB I am getting the “Error using svd. SVD does not support sparse matrices. Use SVDS to compute a subset of the singular values and vectors of a sparse matrix.”

After reading a few queries on forum, I believe CVX does not allow to change the SVD to SVDS. Is there any way I can avoid this error? Will moving to cvx python version be easier to solve? I am stuck at that part.

Thank you very much in advance.

2

Please show a reproducible, but minimal (simple) program with input data which illustrates what you want to do and your difficulty in doing that.

3 
NumN = 32;

%load bus 2 & 3 demand == source ieee 4
Pd2 = 100;
Pd3 = 300;
Qd2 = 0;
Qd3 = 100;
%Gen bus 2 & 3 demand
Pd1 =  0;
Pd4 =  0;
Qd1 =  0;
Qd4 = 0;
%Generator max and min limits
Pg1max = 1000;
Pg2max =  1000;
Pg1min = 0;
Pg2min = 0;
Qr1max = 1000;
Qr2max = 1000;
Qr1min = 0;
Qr2min = 0;
% allowable voltage limits V^2 max and min
V1max = 1.50;
V2max =  1.50;
V3max = 1.50;
V4max = 1.50;
V1min = 0;
V2min = 0;
V3min = 0;
V4min = 0;
%constant values 
B = [Pd2 Pd3 Qd2 Qd3 Pd1 Pd4 Qd1 Qd4 0 1.05^2 ...
    Pg1max Pg2max Pg1min Pg2min Qr1max Qr2max Qr1min Qr2min ...
    V1max  V2max V3max V4max V1min V2min V3min V4min 1 1 1 1]';
% Y bus imaginery part of transfer admittance 
r = 0.005+j*0.05;
a = 1/r;
Y = [a -a 0 0 ;
    -a 3*a -a -a;
    0 -a 2*a -a;
    0 -a -a 2*a];
b = imag(Y);
g = real(Y);

% the cost ocefficients == source example 2 ECE 585
A0 = zeros(NumN,NumN);
A0(1,1) = 0.004;
A0(1,2) = 5.3/2;
A0(2,1) = 5.3/2;
A0(2,2) = 500;
A0(3,3) = 0.006;
A0(3,4) = 5.5/2;
A0(4,3) = 5.5/2;
A0(4,4) = 400;

A1 =  zeros(NumN,NumN);
A1(17,19) = -0.5*g(2,1);
A1(19,17) = -0.5*g(2,1);
A1(18,19) = 0.5*b(2,1);
A1(19,18) = 0.5*b(2,1);
A1(18,20) = -0.5*g(2,1);
A1(20,18) = -0.5*g(2,1);
A1(17,20) = -0.5*b(2,1);
A1(20,17) = -0.5*b(2,1);
A1(19,19) = -g(2,2);
A1(20,20) = -g(2,2);
A1(19,21) = -0.5*g(2,3);
A1(21,19) = -0.5*g(2,3);
A1(19,22) = 0.5*b(2,3);
A1(22,19) = 0.5*b(2,3);
A1(20,22) = -0.5*g(2,3);
A1(22,20) = -0.5*g(2,3);
A1(20,21) = -0.5*b(2,3);
A1(21,20) = -0.5*b(2,3);
A1(19,23) = -0.5*g(2,4);
A1(23,19) = -0.5*g(2,4);
A1(19,24) = 0.5*b(2,4);
A1(24,19) = 0.5*b(2,4);
A1(20,24) = -0.5*g(2,4);
A1(24,20) = -0.5*g(2,4);
A1(20,23) = -0.5*b(2,4);
A1(23,20) = -0.5*b(2,4);

A2 =  zeros(NumN,NumN);
A2(19,21) = -0.5*g(3,2);
A2(20,21) = 0.5*b(3,2);
A2(20,22) = -0.5*g(3,2);
A2(19,22) = -0.5*b(3,2);
A2(21,19) = -0.5*g(3,2);
A2(21,20) = 0.5*b(3,2);
A2(22,20) = -0.5*g(3,2);
A2(22,19) = -0.5*b(3,2);
A2(21,21) = -g(3,3);
A2(22,22) = - g(3,3);
A2(23,21) = -0.5*g(3,4);
A2(24,21) = 0.5*b(3,4);
A2(24,22) = -0.5*g(3,4);
A2(23,22) = -0.5*b(3,4);
A2(21,23) = -0.5*g(3,4);
A2(21,24) = 0.5*b(3,4);
A2(22,24) = -0.5*g(3,4);
A2(22,23) = -0.5*b(3,4);

A3 =  zeros(NumN,NumN);
A3(20,17) = -0.5*g(2,1);
A3(20,18) = 0.5*b(2,1);
A3(18,19) = 0.5*g(2,1);
A3(19,17) = 0.5*b(2,1);
A3(17,20) = -0.5*g(2,1);
A3(18,20) = 0.5*b(2,1);
A3(19,18) = 0.5*g(2,1);
A3(17,19) = 0.5*b(2,1);
A3(20,20) = b(2,2);
A3(19,19) = b(2,2);
A3(20,21) = -0.5*g(2,3);
A3(20,22) = 0.5*b(2,3);
A3(22,19) = 0.5*g(2,3);
A3(19,21) = 0.5*b(2,3);
A3(21,20) = -0.5*g(2,3);
A3(22,20) = 0.5*b(2,3);
A3(19,22) = 0.5*g(2,3);
A3(21,19) = 0.5*b(2,3);
A3(20,23) = -0.5*g(2,4);
A3(20,24) = 0.5*b(2,4);
A3(24,19) = 0.5*g(2,4);
A3(19,23) = 0.5*b(2,4);
A3(23,20) = -0.5*g(2,4);
A3(24,20) = 0.5*b(2,4);
A3(19,24) = 0.5*g(2,4);
A3(23,19) = 0.5*b(2,4);

A4 =  zeros(NumN,NumN);
A4(19,22) = -0.5*g(3,2);
A4(20,22) = 0.5*b(3,2);
A4(20,21) = 0.5*g(3,2);
A4(19,21) = 0.5*b(3,2);
A4(23,22) = -0.5*g(3,4);
A4(24,22) = 0.5*b(3,4);
A4(24,21) = 0.5*g(3,4);
A4(23,21) = 0.5*b(3,4);
A4(22,19) = -0.5*g(3,2);
A4(22,20) = 0.5*b(3,2);
A4(21,20) = 0.5*g(3,2);
A4(21,19) = 0.5*b(3,2);
A4(22,23) = -0.5*g(3,4);
A4(22,24) = 0.5*b(3,4);
A4(21,24) = 0.5*g(3,4);
A4(21,23) = 0.5*b(3,4);
A4(22,22) = b(3,3);
A4(21,21) = b(3,3);

A5 =  zeros(NumN,NumN);
A5(1,2) = 0.5;
A5(2,1) = 0.5;
A5(17,17) = -g(1,1);
A5(18,18) = -g(1,1);
A5(19,17) = -0.5*g(1,2);
A5(20,17) = 0.5*b(1,2);
A5(20,18) = -0.5*g(1,2);
A5(19,18) = -0.5*b(1,2);
A5(17,19) = -0.5*g(1,2);
A5(17,20) = 0.5*b(1,2);
A5(18,20) = -0.5*g(1,2);
A5(18,19) = -0.5*b(1,2);

A6 =  zeros(NumN,NumN);
A6(3,4) = 0.5;
A6(4,3) = 0.5;
A6(19,24) = -0.5*b(4,2);
A6(20,24) = -0.5*g(4,2);
A6(20,23) = 0.5*b(4,2);
A6(19,23) = -0.5*g(4,2);
A6(21,24) = -0.5*b(4,3);
A6(22,24) = -0.5*g(4,3);
A6(22,23) = 0.5*b(4,3);
A6(21,23) = -0.5*g(4,3);
A6(24,24) = -g(4,4);
A6(23,23) = -g(4,4);
A6(24,19) = -0.5*b(4,2);
A6(24,20) = -0.5*g(4,2);
A6(23,20) = 0.5*b(4,2);
A6(23,19) = -0.5*g(4,2);
A6(24,21) = -0.5*b(4,3);
A6(24,22) = -0.5*g(4,3);
A6(23,22) = 0.5*b(4,3);
A6(23,21) = -0.5*g(4,3);

A7 =  zeros(NumN,NumN);
A7(9,10) = 0.5;
A7(10,9) = 0.5;
A7(18,18) = b(1,1);
A7(17,17) = b(1,1);
A7(19,18) = -0.5*g(1,2);
A7(20,18) = 0.5*b(1,2);
A7(20,17) = 0.5*g(1,2);
A7(19,17) = 0.5*b(1,2);
A7(18,19) = -0.5*g(1,2);
A7(18,20) = 0.5*b(1,2);
A7(17,20) = 0.5*g(1,2);
A7(17,19) = 0.5*b(1,2);

A8 =  zeros(NumN,NumN);
A8(11,12) = 0.5;
A8(12,11) = 0.5;
A8(19,24) = -0.5*g(4,2);
A8(20,24) = 0.5*b(4,2);
A8(20,23) = 0.5*g(4,2);
A8(19,23) = 0.5*b(4,2);
A8(21,24) = -0.5*g(4,3);
A8(22,24) = 0.5*b(4,3);
A8(22,23) = 0.5*g(4,3);
A8(21,23) = 0.5*b(4,3);
A8(24,24) = b(4,4);
A8(23,23) = b(4,4);
A8(24,19) = -0.5*g(4,2);
A8(24,20) = 0.5*b(4,2);
A8(23,20) = 0.5*g(4,2);
A8(23,19) = 0.5*b(4,2);
A8(24,21) = -0.5*g(4,3);
A8(24,22) = 0.5*b(4,3);
A8(23,22) = 0.5*g(4,3);
A8(23,21) = 0.5*b(4,3);

A9 =  zeros(NumN,NumN);
A9(24,24) = 0; %f4^2

A10 =  zeros(NumN,NumN);
A10(23,23) = 1.05^2; %e4^2

A11 =  zeros(NumN,NumN);
A11(1,2) = 0.5;
A11(2,1) = 0.5;
A11(5,5) = 1;

A12 =  zeros(NumN,NumN);
A12(4,3) = 0.5;
A12(3,4) = 0.5;
A12(7,7) = 1;

A13 =  zeros(NumN,NumN);
A13(1,2) = 0.5;
A13(2,1) = 0.5;
A13(6,6) = -1;

A14 =  zeros(NumN,NumN);
A14(3,4) = 0.5;
A14(4,3) = 0.5;
A14(8,8) = -1;

A15 =  zeros(NumN,NumN);
A15(9,10) = 0.5;
A15(10,9) = 0.5;
A15(13,13) = 1;

A16 =  zeros(NumN,NumN);
A16(11,12) = 0.5;
A16(12,11) = 0.5;
A16(15,15) = 1;

A17 =  zeros(NumN,NumN);
A17(9,10) = 0.5;
A17(10,9) = 0.5;
A17(14,14) = -1;

A18 =  zeros(NumN,NumN);
A18(11,12) = 0.5;
A18(12,11) = 0.5;
A18(16,16) = -1;

A19 =  zeros(NumN,NumN);
A19(17,17) = 1;
A19(18,18) = 1;
A19(25,25) = 1;

A20 =  zeros(NumN,NumN);
A20(19,19) = 1;
A20(20,20) = 1;
A20(27,27) = 1;

A21 =  zeros(NumN,NumN);
A21(21,21) = 1;
A21(22,22) = 1;
A21(29,29) = 1;

A22 =  zeros(NumN,NumN);
A22(23,23) = 1;
A22(24,24) = 1;
A22(31,31) = 1;

A23 =  zeros(NumN,NumN);
A23(17,17) = 1;
A23(18,18) = 1;
A23(26,26) = -1;

A24 =  zeros(NumN,NumN);
A24(19,19) = 1;
A24(20,20) = 1;
A(28,28) = -1;

A25 =  zeros(NumN,NumN);
A25(21,21) = 1;
A25(22,22) = 1;
A25(30,30) = -1;

A26 =  zeros(NumN,NumN);
A26(23,23) = 1;
A26(24,24) = 1;
A26(32,32) = -1;

A27 =  zeros(NumN,NumN);
A27(2,2) = 1; %dg1^2

A28 =  zeros(NumN,NumN);
A28(4,4) = 1; %dg2^2

A29 =  zeros(NumN,NumN);
A29(10,10) = 1; %dr1^2

A30 =  zeros(NumN,NumN);
A30(12,12) = 1; %dr2^2

cvx_begin
  cvx_solver sedumi;
  variable X(NumN,NumN) symmetric;
  minimize trace(A0*X);
  subject to        
    A1*X - B(1) == 0;
    A2*X - B(2) == 0;
    A3*X - B(3) == 0;
    A4*X - B(4) == 0;
    A5*X - B(5) == 0;
    A6*X - B(6) == 0;
    A7*X - B(7) == 0;
    A8*X - B(8) == 0;
    A9*X - B(9) == 0;
    A10*X - B(10) == 0;
    A11*X - B(11) == 0;
    A12*X - B(12) == 0;
    A13*X - B(13) == 0;
    A14*X - B(14) == 0;
    A15*X - B(15) == 0;
    A16*X - B(16) == 0;
    A17*X - B(17) == 0;
    A18*X - B(18) == 0;
    A19*X - B(19) == 0;
    A20*X - B(20) == 0;
    A21*X - B(21) == 0;
    A22*X - B(22) == 0;
    A23*X - B(23) == 0;
    A24*X - B(24) == 0;
    A25*X - B(25) == 0;
    A26*X - B(26) == 0;
    A27*X - B(27) == 0;
    A28*X - B(28) == 0;
    A29*X - B(29) == 0;
    A30*X - B(30) == 0;

     X == semidefinite(NumN);  
 cvx_end
Rango = rank(X)
4

Your problem is reported infeasible. Therefore, X is a matrix having NaN entries. It appears that due to the constraints, CVX has determined the value for some of the elements of X as being 0, while other entries are NaN, Hence X is returned by CVX as a sparse matrix, some of whose entries are NaN. Hence the error message from MATLAB in your command rank(X)

5

clc
clear all

Vbase = 1; % set to 1
Sbase = 1;
Zbase = Vbase*Vbase/Sbase;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%NETWORK PARAMETERS____%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% GENERATOR COST DATA
%[c2(/MWh^2) c1(/Mwh) c0($)]
Cost = [0.11 5 0; %Gen 1
0.085 1.2 0; %Gen 2
0 0 0]; %Gen 3

% NETWORK DATA
%[From To R X b]
branch = [ 1 3 0.065 0.620 0.45;
3 2 0.025 0.750 0.700;
1 2 0.042 0.900 0.300];

% POWER DEMAND
% [ P(MW) Q(MVAR)]
demand1 = [110 40; %bus 1
110 40; %bus 2
95 50]; %bus 3
demand = demand1/Sbase;

% NETWORK LIMITS
% min_P max_P min_Q max_Q min_V max_V
MinMax = [ -10000 10000 -10000 10000 0.9 1.1; %gen 1
-10000 10000 -10000 10000 0.9 1.1; %gen 2
0000 0 -10000 10000 0.9 1.1]; %gen 3

% Number of buses : maxmum of the From - To in branch matrix
nbus = max(max(branch(:,1:2))); %N is paper

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%__NETWORK ADMITTANCE MATRIX (Ybus)%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Y_bus=zeros(nbus,nbus);

from_b = branch(:,1);
to_b = branch(:,2);
y = 1./((branch(:,3)+1i.*branch(:,4)));
b_sus = branch(:,5);
% off diagonal element = -yij
for k=1:nbus
Y_bus(from_b(k),to_b(k))= -y(k);
Y_bus(to_b(k),from_b(k))=Y_bus(from_b(k),to_b(k));
end

%diagonal elements = all connections to node i : yij + bi/2
for m=1:nbus
for n=1:nbus
if from_b(n)==m
Y_bus(m,m)=Y_bus(m,m)+y(n)+1ib_sus(n)/2;
elseif to_b(n)==m
Y_bus(m,m)=Y_bus(m,m)+y(n)+1ib_sus(n)/2;
end
end
end

Y_bus = Y_bus*Zbase;
% [0.2190 - 2.3291i -0.0517 + 1.1087i -0.1673 + 1.5954i
% -0.0517 + 1.1087i 0.0961 - 1.9406i -0.0444 + 1.3319i
% -0.1673 + 1.5954i -0.0444 + 1.3319i 0.2117 - 2.3522i]

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%AC OPF FORMUATION IN SDP FORM_%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Eye = eye(nbus);
for k = 1:nbus
e{k} = Eye(:,k); % standard basis vector
end

for k = 1:nbus
Ynode{k} = e{k}*e{k}’*Y_bus; %equation 3(a)
end

Zero_3 = zeros(3);
for k=1:nbus %all are 6x6 matrix

Yone{k} = 0.5* [ real(Ynode{k}+(Ynode{k}.'))  imag((Ynode{k}.')-Ynode{k});  %equation 3(b)
                 imag(Ynode{k}-(Ynode{k}.'))  real(Ynode{k}+(Ynode{k}.'))];

Ythree{k} = -0.5* [ imag(Ynode{k}+(Ynode{k}.'))  real(Ynode{k}-(Ynode{k}.'));   %equation 3(c)
                    real((Ynode{k}.')-Ynode{k})  imag(Ynode{k}+(Ynode{k}.'))];

Mnode{k} = [  e{k}*e{k}'   Zero_3;                       %equation 3(d)
               Zero_3    e{k}*e{k}'];  
           
Nnode{k} = [ Zero_3     Zero_3;                          %equation 4(a)    
             Zero_3     e{k}*e{k}'];

end

% variable X = vector with Vd then Vq

cvx_begin
%cvx_solver sedumi
variable W(6,6) symmetric %equation 3(f)
minimize Cost(1,1)*(trace(Yone{1}W) + demand(1,1))^2 … % equation 1 and 2© and 3(i)
+ Cost(1,2)(trace(Yone{1}W) + demand(1,1)) + Cost(1,3) …
+ Cost(2,1)(trace(Yone{2}W) + demand(2,1))^2 …
+ Cost(2,2)(trace(Yone{2}*W) + demand(2,1)) + Cost(1,3)

subject to
W >= semidefinite(6); %equation 3(g)

 trace(Yone{1}*W) <= (MinMax(1,2) - demand(1,1))/Sbase;   %equation 3(k)
 trace(Yone{2}*W) <= (MinMax(2,2) - demand(2,1))/Sbase;
 trace(Yone{3}*W) <= (MinMax(3,2) - demand(3,1))/Sbase;
 trace(Yone{1}*W) >= (MinMax(1,1) - demand(1,1))/Sbase;
 trace(Yone{2}*W) >= (MinMax(2,1) - demand(2,1))/Sbase;
 trace(Yone{3}*W) >= (MinMax(3,1) - demand(3,1))/Sbase;

 trace(Ythree{1}*W) <= (MinMax(1,4) - demand(1,2))/Sbase; %equation 3(l)
 trace(Ythree{2}*W) <= (MinMax(2,4) - demand(2,2))/Sbase;
 trace(Ythree{3}*W) <= (MinMax(3,4) - demand(3,2))/Sbase;
 trace(Ythree{1}*W) >= (MinMax(1,3) - demand(1,2))/Sbase;
 trace(Ythree{2}*W) >= (MinMax(2,3) - demand(2,2))/Sbase;
 trace(Ythree{3}*W) >= (MinMax(3,3) - demand(3,2))/Sbase;
 
 trace(Mnode{1}*W) <= (MinMax(1,6))^2;                   %equation 3(m)
 trace(Mnode{1}*W) >= (MinMax(2,6))^2;
 trace(Mnode{2}*W) <= (MinMax(3,6))^2;
 trace(Mnode{2}*W) >= (MinMax(1,5))^2;
 trace(Mnode{3}*W) <= (MinMax(1,5))^2;
 trace(Mnode{3}*W) >= (MinMax(1,5))^2;
 
 trace(Nnode{1}*W) == 0;                                 %equation 4(b)
 trace(Nnode{2}*W) == 0;
 trace(Nnode{3}*W) == 0;

cvx_end
Rango = rank(W)

The subject error is expected in this one as well.

6

Thanks, i am just trying to emulate the results of earlier published works so to have a baseline for my work and get handson on cvx.

7

Then I recommend you start diagnosing the source(s) of the infeasibility. The process of doing so will provide you valuable experience in hands-on optimization.

8

Thank you Mark. Are there any other resources for debugging issues in the code? I am referring to is The CVX Users’ Guide
Release 2.1

9
1 Like