Using permute
and repmat
, most constraint loops in CVX can be eliminated. For example, the following:
variables X(N) Y(M) Z(N);
for n=1:N
for m=1:M
if n~=m && n^3>m^2
X(n)>=Y(m)+Z(n); % (1)
end
end
end
The first step is to expand. (1) actually means there are M*N constraints, so X
is expanded to repmat(X,1,M)
, Y
is expanded to repmat(Y,1,N)
, Z(n) is expanded to repmat(Z,1,M)
, according to their original dimensions. Now, repmat(X,1,M)
is a N-by-M matrix, repmat(Y,1,N)
is M-by-N, repmat(Z,1,M)
is N-by-M. Note that some matrixs are M-by-N whlie others are N-by-M.
Before doing operations between them (in this case, add them), their dimensions need to be permuted to a uniform order. We choose the order N-by-M so that only repmat(Y,1,N)
need to be permuted by using permute(repmat(Y,1,N),[2 1])
. Then the whole loop can be easily replaced by the following
repmat(X,1,M) >= permute(repmat(Y,1,N),[2 1]) + repmat(Z,1,M); %(2)
if the condition “if n~=m && n^3>m^2” is removed.
To handle the condition “if n~=m && n^3>m^2”, we create a N-by-M filter matrix F, in which F(n,m)=1 if the condition is satisfied, otherwise 0. First of all, create a matrix for n using repmat((1:N)',1,M)
, and for m using permute(repmat((1:M)',1,N),[2 1])
, permute is used here to adjust order like above. Then obviously condition1 “n~=m” can be expressd as repmat((1:N)',1,M)~=permute(repmat((1:M)',1,N),[2 1])
, condition2 “n^3>m^2” can be repmat((1:N)',1,M).^3>permute(repmat((1:M)',1,N),[2 1]).^2
, the product of them is F
F=(repmat((1:N)',1,M)~=permute(repmat((1:M)',1,N),[2 1])).*(repmat((1:N)',1,M).^3>permute(repmat((1:M)',1,N),[2 1]).^2);
since it is a “and” relationship. To filter (2), we multiply it by F. So finally, the whole loop can be expressed as
F.*(repmat(X,1,M) - permute(repmat(Y,1,N),[2 1]) - repmat(Z,1,M))>=0;
With these three steps, I think most CVX constraint loops can be vectorized, so a lot of time can be saved especially when variables are too many. Obviously what’s described here can be easily generalized to other examples.