N-dimensional variables

[SteveDiamond]

Currently cvxpy only supports scalar, vector, and matrix variables. Additional dimensions can be mimicked using a dict:

# Create a 3D k-by-m-by-n variable.
x = {}
for i in range(k):
  x[i] = Variable(m,n)

[marziacescon]

Hey Steven!

I'm trying to translate an optimization problem I've written in matlab cvx into python cvxpy.
My optimization variable x is a column vector which appears in the objective function is the format x^TRx, where R is a square matrix. Now, I encounter problems in transposing the variable x. Specifically, I can't get it transposed. The code I'm using is the following:

import numpy as np
import cvxpy as cvx

create optimization variables

	x = cvx.Variable(nf_1) # a vector
            xTRx = np.dot(np.dot(np.transpose(x),R_ins),p_ins)

I'm aware of the fact that 1D array transposed returns itself, but I haven't been able to find a way around it.
Any help is appreciated.

[SteveDiamond]

Hi Marzia,

You can't ever use numpy functions on cvxpy objects. You have to use the cvxpy functions:

	x = cvx.Variable(nf_1) # a vector
            xTRx = cvx.quad_form(x, R_ins)

[rileyjmurray]

@SteveDiamond it seems that adding support for this is pretty straightforward in 1.0.

expression.py doesn't actually impose the at-most-2d restriction; that's done in atom.py and leaf.py (note that restrictions imposed on Leaf objects are inherited by Variables, Constants, and Parameters).

When I removed these restrictions (on my fork) N-d Expressions worked with every atom I tested (all elementwise atoms, some affine atoms, some non-affine atoms). N-d Variables also behaved correctly, in the ways that I tested them.

[SteveDiamond]

That's great! I guess we did a better job than I thought with the refactor. We will need to add checks to certain operations like matrix multiplication, but it could be enough to simply ban N>2 in all relevant locations. We may not need any special case handling for N>2.

[rileyjmurray]

Update: we should check to make sure that @ with cvxpy Expressions of constant value behaves in the same way as @ with numpy ndarrays of higher dimensions. Reason being: @ and np.dot behave differently for higher-dimensional arrays, and the original implementation of matmul in cvxpy might boil down to numpy's dot function, rather than the actual matmul function. (I haven't looked into this, but it seems like something that could easily happen.)

[vchinde]

Hi Steve, I am new to cvxpy and python, I have created a three dimensional matrix using your suggestion.

Create a 3D k-by-m-by-n variable.

x = {}
for i in range(k):
x[i] = Variable(m,n)

I have problem slicing the 3D matrix. Like how to extract x(1, :, 1) (matlab notation) using python commands

[SteveDiamond]

You can't slice it. That's one advantage of truly supporting ND matrices.

[SteveDiamond]

Sorry, in your specific example you would do x[0][:,0]. But you can't slice the dictionary dimension.

[vchinde]

I am running into errors if I do what you suggested

import cvxpy

N_sch=2
N_perhour=2

x = {}

for i in range(N_perhour):
x[i] = cvxpy.Variable(N_sch,19)

print(x[0][:,0])
Traceback (most recent call last):

File "", line 1, in
print(x[0][:,0])

File "C:\Users\vchinde\AppData\Local\Continuum\anaconda2\lib\site-packages\cvxpy\expressions\expression.py", line 315, in getitem
return cvxtypes.index()(self, key)

File "C:\Users\vchinde\AppData\Local\Continuum\anaconda2\lib\site-packages\cvxpy\atoms\affine\index.py", line 48, in init
self.key = ku.validate_key(key, expr.shape)

File "C:\Users\vchinde\AppData\Local\Continuum\anaconda2\lib\site-packages\cvxpy\utilities\key_utils.py", line 49, in validate_key
raise IndexError("Too many indices for expression.")

IndexError: Too many indices for expression.

[SteveDiamond]

It should be x[i] = cvxpy.Variable((N_sch,19)). Writing cvxpy.Variable(N_sch,19) actually makes a variables with dimensions (N,) and raises an error in the current version, but you may not have that version on windows.

[vchinde]

Thanks Steve. I got it

[vchinde]

Hi Steve, I have the following question

m_z=cvxpy.Variable(shape=(24, 19))
for i_sch in xrange(24):
m_i=m_z[i_sch,:]
constr +=[Simconst['m_min'] <= m_i]

The dimension of Simconst['m_min'] is (19L, 1L). When I run the code I get the following error "ValueError: Cannot broadcast dimensions (19, 1) (19,)". How to overcome this error? I did tried using m_i.T but has the same error

[SteveDiamond]

Write Simconst['m_min'][:,0] <= m_i

[rotcx]

cvxpy.Variable(dim1,dim2)

some of the (dim1_index,dim2_index) combinations are meaning less

how can we "do not define variables on those meaningless combinations" instead of "define variables on the whole combination, then use constraints to fix the meaningless variables "
?

[rotcx]

Another question:
how can I vectorize constraints instead of using for-cycling clauses? Tedious when the problem is large scale.

[SteveDiamond]

The simplest thing is to make a variable vector x representing all the meaningful combinations, then expand that to a sparse matrix using cvxpy.reshape(A@x, (dim1, dim2)), where A is a dim1*dim2 by T matrix, where T is the number of meaningful combinations. You can figure out what the entries of A should be. That exercise will hopefully teach you how to vectorize.

[carlosgmartin]

Will multidimensional variables ever be added to cvxpy?

[SteveDiamond]

Not any time soon.

[rotcx]

had to turn to or-tools Steven Diamond <notifications@github.com> 于2020年2月13日周四 上午10:26写道：

…

Not any time soon. — You are receiving this because you commented. Reply to this email directly, view it on GitHub <https://github.com/cvxgrp/cvxpy/issues/198?email_source=notifications&email_token=ABCNQ6OXXSWIGRNT5MLOMY3RCSVVNA5CNFSM4BJCB5JKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOELTEWGY#issuecomment-585517851>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABCNQ6PLLAASOVH3NEW35ILRCSVVNANCNFSM4BJCB5JA> .

-- =============================== Xuan cx_2016@foxmail.com

[rileyjmurray]

@chenxuanhanhao thanks for letting us know. I'm curious: what's your application where the high-level logic offered by python makes it unrealistic to work only with 1d or 2d arrays?

[teichert]

Thanks for cvxpy! With the uses that cvxpylayers is likely to see, I'm sure higher order tensors are going to be increasingly common (e.g. anything you wanted to do with a matrix before, you might want to do batched). In my current use case, I want to find a tensor that is close to another tensor (in the sense of kl_div) while respecting a number of constraints that deal with marginalizing out various subsets of the dimensions. (number of dimensions can be between 2 and 20.)

Allowing cp.sum to work on higher-order ndarrays would satisfy my current use case (even if I am only allowed to sum over one dimension at a time). Alternatively, if just cp.transpose and cp.reshape allowed higher order ndarray, I'm sure that would facilitate easier work-arounds in many cases.

[SteveDiamond]

You can transpose or reshape outside of cvxpy, in pytorch/tensorflow.

[teichert]

You can transpose or reshape outside of cvxpy, in pytorch/tensorflow.

Thanks; that's right, and I can do fancy indexing in order to sum out "dimensions" even if my representation in cvxpy is flat. In my case, the constraints I needed to add were easy to express as summations on subsets of dimensions, but difficult to write in terms of the flat representation. So, what I ended up doing was to make my own sum function that took in a flattened ndarray (i.e. transposing and flattening having already been applied outside of cvxpy), a shape tuple, and a set of dimensions to sum over and then (within my sum function) I did the indexing (via np.ravel_multi_index, np.ix_, ndarray.transpose, and ndarray.reshape) to achieve the equivalent of multi-dimensional summations on a flat cvxpy array.

[Xeunon]

@SteveDiamond it seems that adding support for this is pretty straightforward in 1.0.

expression.py doesn't actually impose the at-most-2d restriction; that's done in atom.py and leaf.py (note that restrictions imposed on Leaf objects are inherited by Variables, Constants, and Parameters).

When I removed these restrictions (on my fork) N-d Expressions worked with every atom I tested (all elementwise atoms, some affine atoms, some non-affine atoms). N-d Variables also behaved correctly, in the ways that I tested them.

Would you mind explaining how exactly you removed the restrictions?
I'm working on a production planning optimization problem and the nature of its variables are three dimensional, i.e. Product, Machine, Period. I have found a workaround but having the ability to define 3D variables would sure make the abstraction much simpler and I wouldn't mind missing out on a couple of functions that might not work on higher dimensions.

[EmilianoG-byte]

@teichert would you mind sharing your code? I think might have a similar use case where I need to reshape and transpose a CVXPY variable (i.e. reshape into a 4-rank tensor, transpose 2 dimensions, then reshape back to a 2-rank tensor).

[rileyjmurray]

Cristian— CVXPY’s partial_transpose atom might be useful for you here.

…

On Sat, Jul 8, 2023 at 7:59 PM Cristian Emiliano Godinez Ramirez < ***@***.***> wrote: @teichert <https://github.com/teichert> would you mind sharing your code? I think might have a similar use case where I need to reshape and transpose a CVXPY variable (i.e. reshape into a 4-rank tensor, transpose 2 dimensions, then reshape back to a 2-rank tensor). — Reply to this email directly, view it on GitHub <#198 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACRLIFAK5AXJTYTNCMS4AH3XPHX4VANCNFSM4BJCB5JA> . You are receiving this because you commented.Message ID: ***@***.***>

[EmilianoG-byte]

thank you for your fast response, @rileyjmurray! I've seen the function, but I should have mention that my two legs correspond to "two different subsystems" as:

choi = np.reshape(superop, [dim] * 4)
choi = choi.swapaxes(1, 2).reshape([dim ** 2, dim ** 2])

and as I understand partial_transpose could do a swap between 0<->2 or 1<->3 only. Any ideas on how I could achieve mine: 1<->2 ?

[rileyjmurray]

@EmilianoG-byte thanks for the clarification. Can you open a separate GitHub issue requesting this functionality?

[EmilianoG-byte]

@rileyjmurray Just opened the issue :). Here: #2181

[rileyjmurray]

Thanks so much! I can’t promise quick progress on it, but that will help a lot in keeping track of your request.

…

On Wed, Jul 12, 2023 at 1:55 PM Cristian Emiliano Godinez Ramirez < ***@***.***> wrote: @rileyjmurray <https://github.com/rileyjmurray> Just opened the issue :) — Reply to this email directly, view it on GitHub <#198 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACRLIFALIYZFLJLGTPWDMY3XP3QJNANCNFSM4BJCB5JA> . You are receiving this because you were mentioned.Message ID: ***@***.***>

[teichert]

@EmilianoG-byte sorry for the slow reply---I wasn't able to easily find the code, so I'm guessing that I was just tinkering when I did this. But it looks like you all are working on a better solution anyway.

[EmilianoG-byte]

Hi @teichert! There is an open issue regarding this general transpose of two dimensions, but if you have already some sample code that can handle the partial transpose of multi dimensional arrays, I would be interested in taking a look at it :)

[teichert]

@EmilianoG-byte okay; I did end up finding my code fwiw. I'm guessing it isn't relevant to your needs, but I'll add it here just in case. There might be more transposing capabilities in cvxpy now, but, as I recall, at the time of writing, I wanted to add constraints on ndarray variables regarding summing out various subsets of the dimensions, and the way I ended up doing it was by just having a flat variable representation and forming the ndsums based on that (using np.transpose before and after if necessary). No guarantees that this code makes any sense! :)

import cvxpy as cp
import numpy as np

def cvx_flat_sum(flattened_ndarray, dim=None, shape=None):
    """
    args:
        a: a cvxpy variable that is a flattened ndarray
        dim: a tuple of dimensions to sum over; if None, sum over all dimensions
        shape: the shape of the ndarray that flat represents
    """
    assert len(flattened_ndarray.shape) == 1
    if dim is None:
        return cp.sum(flattened_ndarray)
    if not hasattr(dim, '__len__'):
        dim = (dim,)
    dim = tuple(dim)
    set_to_drop = set(dim)
    size_to_keep = 1
    size_to_drop = 1
    dim_to_keep = []
    for d, size in enumerate(shape):
        if d in set_to_drop:
            size_to_drop *= size
        else:
            dim_to_keep.append(d)
            size_to_keep *= size
    dim_to_keep = tuple(dim_to_keep)
    ixs = (np.ravel_multi_index(np.ix_(*[range(s) for s in shape]), shape)
             .transpose(dim_to_keep + dim)
             .reshape(size_to_keep, size_to_drop))
    tmp = flattened_ndarray[ixs]
    out = cp.sum(tmp, axis=1)
    return out

def test_simple_first_dim():
    original = np.array([[1, 2],
                         [3, 4]])
    flat = original.flatten()
    shape = original.shape
    assert np.allclose(cvx_flat_sum(flat, dim=0, shape=shape).value, np.array([4, 6]))

def test_simple_last_dim():
    original = np.array([[1, 2],
                         [3, 4]])
    flat = original.flatten()
    shape = original.shape
    assert np.allclose(cvx_flat_sum(flat, dim=1, shape=shape).value, np.array([3, 7]))

def check_large_generalized_sum(dim):
    original = np.arange(3*4*5*6).reshape((3,4,5,6))
    flat = original.flatten()
    shape = original.shape
    flat_sum = cvx_flat_sum(flat, dim=dim, shape=shape).value
    expected_flat_sum = original.sum(axis=dim).flatten()
    assert np.allclose(expected_flat_sum, flat_sum)

def test_sum_all_dims():
    check_large_generalized_sum(dim=())

def test_sum_single_internal_dim():
    check_large_generalized_sum(dim=2)

def test_sum_multiple_dims():
    check_large_generalized_sum(dim=(3, 1))

def test_simple_cvx_shape():
    flat = cp.Variable(4)
    shape = (2, 2)
    assert cvx_flat_sum(flat, dim=0, shape=shape).shape == (2,)

def test_large_cvx_shape():
    flat = cp.Variable(2*3*4*5)
    shape = (2, 3, 4, 5)
    assert cvx_flat_sum(flat, dim=(0, 2), shape=shape).shape == (3 * 5,)

[yasirroni]

This issue seems to loss its focus from the title of N-dimensional variables
#198

Can we back to that or make new issue to track that?

Re-upping @Xeunon question here.

@SteveDiamond it seems that adding support for this is pretty straightforward in 1.0.
expression.py doesn't actually impose the at-most-2d restriction; that's done in atom.py and leaf.py (note that restrictions imposed on Leaf objects are inherited by Variables, Constants, and Parameters).
When I removed these restrictions (on my fork) N-d Expressions worked with every atom I tested (all elementwise atoms, some affine atoms, some non-affine atoms). N-d Variables also behaved correctly, in the ways that I tested them.

Would you mind explaining how exactly you removed the restrictions? I'm working on a production planning optimization problem and the nature of its variables are three dimensional, i.e. Product, Machine, Period. I have found a workaround but having the ability to define 3D variables would sure make the abstraction much simpler and I wouldn't mind missing out on a couple of functions that might not work on higher dimensions

We might be able to help to check and track list of functions that works in higher dimension and not, also make tests for that.

CC: @SteveDiamond