higher dimensional variables / variable indexing in cvxpy

[arkha]

Hello,

I need to use variables with three or more indices. To simplify modelling, it would be thus useful if I could create optimization variables with three or more dimensions; as far as I can see, however, in cvxpy it is only possible to create up to two dimensional arrays.

To handle the extra dimensions using 2 dimensional arrays, I thus have to squeeze the extra dimensions in the two I have. Then, to perform summations only over certain dimensions, I am thinking of creating lists of tuples that will select only certain elements of the optimization variable.

Something along the lines of:

x = cvx.Variable(2,2)
interesting_indices = [(0,0),(1,1)]
constraints = [sum([x[i] for i in interesting_indices]) <= 0.5]

to, e.g., perform the summation only over certain selected entries - which seems to work.

Is there a better way to handle higher dimensions / variables indexing in cvxpy? 

For instance, it would be useful if we could plug Variable objects into numpy structures, as in

x = np.empty((3,3,3), dtype=object)
x[:,:,:] = cvx.Variable()

to create a three dimensional array, and then work with x as any other numpy array (also in terms of indexing).

Thank you!

[Steven Diamond]

You can create a map of tuple to matrix variable. Like this:

# x is 3D array with dimensions p x q x r.

x = {}

for i in range(r):

  x[i] = Variable(p, q)

For 4D arrays the keys would be tuples with two elements, etc.

Would you mind telling me more about how you're using the ND arrays? I may add ND arrays to cvxpy at some point, but I want to understand the use cases better first.

[arkha]

Thank you for the quick reply.

I am trying to migrate some code from Yalmip for a scheduling application. Variables there could be assignments of a certain task i \in I, to a certain machine j \in J in a given time step t \in T.

In Yalmip in order to, e.g., ensure that each machine is assigned to one task only at each time step, I could simply do

x_ijt = sdpvar(3,3,3)
for j = 1:J
  for t = 1:T
    sum(x_ijt(:,j,t)) <= 1
  end
end

What would be a good way of implementing something like this in cvxpy?

By the way, I am planning to use cvxpy only for the modelling part, then extract the model in standard form and plug it into Gurobi or CPLEX. I am not sure ECOS BNB can handle the sizes I am looking into - I will try though.

[Steven Diamond]

cvxpy supports gurobi. I just added discussion of how to use gurobi as a solver to the website.

So the point of packing things in a 3D array is to do vectorized operations? You could do the operations in a nonvectorized way with a map in Python:

x_ijt = {}

for i in range(I):

  for j in range(J):

    for t in range(T):

      x_ijt[(i,j,t)] = Variable()


for j in range(J):
  for t in range(T):
    sum([x_ijt[(i,j,t)] for i in range(I)]) <= 1






[arkha]

Steven,

Since you are managing the website: I noticed that the instructions to set solvers never pass arguments as strings. On my installation this is however necessary.

For instance, on http://www.cvxpy.org/en/latest/tutorial/advanced/index.html you have:

prob.solve(solver=ECOS)

while it should be (at least for me):

prob.solve(solver='ECOS')

This is in the sections "Solve method options", "Viewing solver output", "Setting solver options" and "Getting the standard form". 

It's very minor, just thought I'd let you know.