Incorrect results for multiple PSD cones

[archenal]

Hi，

Thanks for this insightful work.

I am trying to use factorwidth.m to reformulate my SDPs. I found that an original infeasible problem will be reported as feasible by mosek if factorwidth.m is used.
Example_SDPfw.zip

One guess for this issue is that factorwidth.m works imperfectly with multiple PSD cones, as I noticed that all the examples you gave only have 1 cone. I did notice you have some code dealing with the multiple PSD cones, so this guess might not be correct.

Thanks

[Jarmill]

Archenal,

Thanks for your interest in our work.

I'm not sure this is necessarily a bug. Factor width/structured subset approximations return an upper bound on the SDP objective. A dual infeasible objective means that the optimum objective value is -infinity. An upper bound of this may be -infinity (tight), or may be zero (returned by factor-width as you've observed). With a block-size of 5 or 10 (equivalent to the original problem), a dual infeasible certificate is returned. For any other block-size (like 1-3), a cost of zero is returned since b = 0 (and the returned x = 0). With a block-size of 4 the dual infeasible certificate is returned, not sure why that is the case.

Try it out by changing the BlockSize variable

load('Example_SDPfw.mat')

s = svd(full(A));
n_null = nnz(s <= 1e-8);

BlockSize = 4;
pars.fid = 0;

[x0, y0, info0] = sedumi(A, b, c, K, pars);
cost0 = c'*x0;
disp(info0)

% Same process through factorwidth
% opts=struct;
% opts.bfw = 1;
% opts.block = BlockSize;
% 
% [Af, bf, cf, Kf] = factorwidth(A, b, c, K, opts);
% [xf, yf, infof] = sedumi(An, bn, cn, Kn, pars);
% costf = cf'*xf;
%disp(infof)

%Same process through decomposed_subset
[An, bn, cn, Kn, info_dec] = decomposed_subset(A, b, c, K, BlockSize);
[xn, yn, infon] = sedumi(An, bn, cn, Kn, pars);
xnr = decomposed_recover(xn, info_dec);
costn = c'*xnr;
disp(infon)


%eigenvalues of returned solution are in eX
count = 0;
X = cell(length(K.s), 1);
for i = 1:length(K.s)
    Ksi = K.s(i);
    xi = xnr(count + (1:Ksi^2));
    X{i} = reshape(xi, Ksi, Ksi);
    count = count + Ksi^2;
end

eX = cellfun(@eig, X, 'UniformOutput', false);

One issue with this problem is that A is rank deficient, and has two zero singular values. It looks like there is no avenue for facial reduction here, SieveSDP came up dry.

Hopefully this answers your question.

Best,

Jared

[zhengy09]

Hi Archenal,

Jared has given a detailed explanation. I didn't see an error in our code factorwidth.m either. The problem you have is dual infeasible but primal feasible. In fact, the primal problem has an unbounded primal objective, which is -infinity. So, when using the factor-width for inner approximation, it will return an upper bound for the primal problem, which can be -infinity or a bounded value ( it can be zero in some tests).

If you change the primal cost function a little bit, you may see the problem become both primal and dual feasible. Here is the code:

load Example_SDPfw
opts.bfw = 1;
opts.nop = 10;
c = rand(size(A,1),1); % change the cost function, won't change primal feasibility
[x,y,info] = sedumi(A,b,c,K);
[Anew, bnew, cnew, Knew, infofw] = factorwidth(A,b,c,K,opts);
% reformulating the SDP
[xn,yn,infon] = sedumi(Anew,bnew,cnew,Knew);

Let's know if you have further questions.
Best,
Yang

[archenal]

Hi Jared and Yang,

Thanks for your kind reply and the helpful explanations.

In fact, the problem I had should be primal infeasible, and I used YALMIP to convert it into sedumi format. Based on your reply, I realized that YALMIP automatically dualized my problem, that was why the results returned by the solver confused me. Now, I use YALMIP to convert my problem in the primal format, and factorwidth.m works perfectly.

Again, thank you for developing and providing these useful tools.