Try this code. load variables from this mat file https://uploadfiles.io/38e79 . Mosek and Gurobi report optimal value is 118229(I think it’s correct). Now uncomment %x==[0,0,0,1;1,0,0,0;0,0,1,0]; Gurobi fails and Mosek reports 110757 as optimal value!. I think it’s an infeasible point!

cvx_clear;

for sample=1:nrofsamp

for u=1:v_nodes

for v=1:v_nodes

for i=1:s_nodes

for j=1:s_nodes

li_cost(sample,u,v,i,j)=s_link_cost(u,v,i,j);

`end end end end`

end

%%

[row_s,col_s]=find((s_net_topology)~=0);

[row_v,col_v]=find((v_net_topology)~=0);

Timestart=tic;

cvx_beginvariable x(v_nodes,s_nodes) integer nonnegative;%nodes

%variable y(v_nodes,v_nodes,s_nodes,s_nodes) nonnegative;%links

variable y(nrofsamp,v_nodes,v_nodes,s_nodes,s_nodes) nonnegative;%linksminimize sum(sum(s_processing_cost.*x))+(1/nrofsamp)*sum(sum(sum(sum(sum(li_cost.*y)))));

subject tofor i=1:v_nodes

sum(x(i,:))==1;

endfor i=1:s_nodes

sum(x(:,i))<=1;

endfor i=1:s_nodes

for j=1:v_nodes

x(j,i)*v_processing_capacity(j)<=s_processing_capacity(i);

end

endfor sample=1:nrofsamp

for i=1:s_nodes

for j=1:sum(nonzeros(v_net_topology))%# of virtual links

[ph1]=find(s_net_topology(i,:)~=0);

[ph2]=find(s_net_topology(:,i)~=0);

ind1=zeros;

ind2=zeros;

if(~isempty(ph1)&&~isempty(ph2))

for t=1:length(ph1)

ind1(t)=sub2ind(size(y),sample,row_v(j),col_v(j),i,ph1(t));

end

for t=1:length(ph2)

ind2(t)=sub2ind(size(y),sample,row_v(j),col_v(j),ph2(t),i);

end`%sum(y(ind1(:)))-sum(y(ind2(:)))== v_link_capacity(row_v(j),col_v(j))*(x(row_v(j),i)-x(col_v(j),i)); sum(y(ind1(:)))-sum(y(ind2(:)))==v_link_capacity{sample}(row_v(j),col_v(j))*(x(row_v(j),i)-x(col_v(j),i)); end if (isempty(ph1)&&~isempty(ph2)) for t=1:length(ph2) ind2(t)=sub2ind(size(y),sample,row_v(j),col_v(j),ph2(t),i); end -sum(y(ind2(:)))- v_link_capacity{sample}(row_v(j),col_v(j))*(x(row_v(j),i)-x(col_v(j),i))==0; end if (isempty(ph2)&&~isempty(ph1)) for t=1:length(ph1) ind1(t)=sub2ind(size(y),sample,row_v(j),col_v(j),i,ph1(t)); end sum(y(ind1(:)))- v_link_capacity{sample}(row_v(j),col_v(j))*(x(row_v(j),i)-x(col_v(j),i))==0; end end`

end

endfor sample=1:nrofsamp

for i=1:sum(nonzeros(s_net_topology))`arrind=null(1,4); for j=1:sum(nonzeros(v_net_topology)) a=[row_v(j),col_v(j),row_s(i),col_s(i)]; arrind=vertcat(arrind,a); end [s z]=size(arrind); if(~isempty(arrind)) indices(1)=sub2ind(size(y),sample,arrind(1),arrind(2),arrind(3),arrind(4)); for row=2:s indices(row)=sub2ind(size(y),sample,arrind(row-1,1),arrind(row-1,2),arrind(row-1,3),arrind(row-1,4)); end end sum(y(indices(:)))<=s_link_capacity(row_s(i),col_s(i));`

end

end

%x==[0,0,0,1;1,0,0,0;0,0,1,0];

t1=toc(Timestart);

cvx_end