How to model this constraint with binary variables in CVX?

MIDCP
1

I have a sparse square symmetric matrix w of size N\times N. Most of the elements in this matrix are 0.

I have a constraint as
for c=1:C
y_{ij}\ge x_{ic}+x_{jc}-1, for each (i,j) with w_{ij}>0

here, i=1,2,\cdots N and j=1,2,\cdots N

So, we need to focus on only the pairs of indices for which we have non zero elements in w.

All diagonal elements in w are 0. We only need to consider the upper right or lower left triangular matrix.

I think we also can express the constraint as

for c=1:C
for i=1:N-1
for j=i+1:N
if w_{ij}>0
y_{ij}\ge x_{ic}+x_{jc}-1,
end
end
end
end

2

If that does what you want and the constraints are generated quickly enough to suit your needs, then use it.

It can probably be vectorized, which would speed up constraint generation. @jackfsuia is the vectorization wizard.

3

The reason I am trying to generates quicker constraints is because N=500. So it is a very big and very sparse matrix.

4

Presuming all elements of w are >= 0, you can at least vectorize the inner loop.

w(i,i+1:N).*y(i,i+1:N) >= w(i,i+1:N).*(x(i,i+1:N) + x(i+1:N,c) - 1)

I don’t know whether the “constraints” having w(i,j) = 0 will cause CVX to realize they are a vacuous 0 >= 0 constraint, and therefore not provide those to the solver, nor whether the solver would be “bothered” very much f they are provided to the solver (maybe depends on the solver). Nor do I know whether CVX would be slowed down making such an assessment. So I don’t know whether this method is any good.

Perhaps @jackfsuia can figure out how to also vectorize the for loop on i, maybe using repmat and permute, and perhaps even on c. But eliminating even one level of for loop should hopefully give a big speedup.

5

write them as
for
for
for
if(wij>0)&& j>i
…
then follow the steps here How to vectorize most constraint loops in CVX

6

Would you please help me to fix this. I tried but failed. I am not so experienced with repmat and permute operators.

7

sorry, there were some typos, which I have fixed now:

(repmat(w, 1 ,1, C)>0).*(repmat((1:N-1)‘, 1,N-1, C)<permute(repmat((1:N)’,1,N-1,C),[2 1 3])).*repmat(y, 1,1, C)-permute(repmat(x, 1, 1, N), [1 3 2])-permute(repmat(x, 1 , 1, N-1),[3 1 2])+1>=0;

Please test it before using it.

8

Thanks a lot. I appreciate it very much. Would you please explain a bit. Also, I don’t see ‘<=’ in the constraint.

Why do we still have i and j in the expression. You mean we still have loops?

9

Sorry, that is plenty of typos. After fixing them, the link provided to you previously should make sense now for explanations.

10

Great. Thanks a lot. Can I use the same technique for objective function right? For example, my objective function is

\sum_{i=1}^N\sum_{j=1}^Nw_{ij}y_{ij}

Of course, again, we consider only the pairs with w_{i,j}>0 and j>i

11

Sure. lets denote the array X=w.*y. ,which is a X(N,N). What your objectivev func do is just summing on it, but with a condition on i,j. So, first use a filter array to filter it, after that, use the cvx function sum.

12

Thanks.

This is how I have done it.

cvx_begin

variable Y(N,N) binary
variable X(N,C) binary
expressions Zmat(N,N) nonegative

Zmat=W.*Y;

minimize sum((repmat(W, 1 ,1, C)>0).*(repmat((1:N-1)', N-1, C)<permute(repmat((1:N)',N-1,C),[2 1 3])).*repmat(Zmat, 1,1, C))

but it is throwing some error as below

Array dimensions must match for binary array op.

Error in minimize (line 14)
x = evalin( ‘caller’, sprintf( '%s ', varargin{:} ) );

Error in RIM_v2 (line 27)
minimize sum((repmat(W, 1 ,1, C)>0).*(repmat((1:N-1)‘, N-1, C)<permute(repmat((1:N)’,N-1,C),[2 1 3])).*repmat(Zmat, 1,1, C))

13

you need to learn how to use sum. its much simpler.

14

Sure. But there is a type/error is the Filter formulation as there is matrix dimension mismatch

(repmat(w, 1 ,1, C)>0).*(repmat((1:N-1)‘, N-1, C)<permute(repmat((1:N)’,N-1,C),[2 1 3]))

would you please check.

15

you missed a 1 in repmat((1:N-1)‘, N-1, C) and repmat((1:N)’,N-1,C). I believe it should be placed between (1:N)’ and N-1.

16

Sorry. But the problem still exists…

(repmat(W, 1 ,1, C)>0).*(repmat((1:N-1)‘,1, N-1, C)<permute(repmat((1:N)’,1,N-1,C),[2 1 3]))

Array dimensions must match for binary array op.

17

do you need (N-1,N-1,C),or(N,N-1,C)? you have both in one sentence. sorry, i am confused as to which expression you want to express? if its the following

“function is∑Ni=1∑Nj=1wijyijOf course, again, we consider only the pairs with wi,j>0 and j>i”

replace all your N-1 with N, try again

18

Great. I will definitely try.

In the meantime, I also tried to filter out the unneccessary pair.

The original problem I have is

for i=1:L
    for j=1:L
        for n=1:N
            if  W(i,j)>0 and i<j
                z(i,j)>=X(i,n)+X(j,n)-1;
            end
        end
    end
end

where W (L,L) is a known matrix and is symmetric. X(L,N) is optimization variable.

So, I tried to reduce the size of the problem by extracting all the pairs (i,j) for which W(i,j)>0. So, I removed the if condition.

Now, let’s says there are M such pairs.

veca contains the i’s

vecb contains the j’s

Now, the constraint becomes

for m=1:M
    for n=1:N
        z(m)>=X(veca(m),n)+ X(vecb(m),n)-1;
    end
end

Now, check that i still have some loop here.

is there any way to overcome that by your techniques??

19

I think the array X(veca(m),n) can be expressed by

X(veca(repmat((1:m)',1,N)).*(permute(repmat((1:N)',1,M),[2 1])-1)+veca(repmat((1:m)',1,N)))

, which is done alomost by the same technique. If this works, then other elements in your loop can be easiliy expressed following the same path.