Cvxpy: multiply doesn't broadcast correctly

[dave31415]

`def test_cvxpy_multiply():
# fails in cvxpy 1.0.6
# python 2.7.15
# ValueError: Cannot broadcast dimensions (4,) (4, 1)
x = np.ones(4)
y = cvxpy.Variable((4, 1))
cvxpy.multiply(x, y)

def test_cvxpy_multiply_2():
# passes in cvxpy 1.0.6, python 2.7.15
x = np.ones((4, 1))
y = cvxpy.Variable((4, 1))
cvxpy.multiply(x, y)

def test_cvxpy_multiply_3():
# passes in cvxpy 1.0.6, python 2.7.15
x = np.ones((4, 2))
y = cvxpy.Variable((4, 2))
cvxpy.multiply(x, y)`

[SteveDiamond]

This kind of broadcasting is intentionally not supported, because I think it adds more confusion than functionality.

[lumbric]

@SteveDiamond should this be made more explicit in the documentation? I can prepare a PR if you want, but tbh I haven't fully understand yet when numpy broadcasting rules are not followed precisly.

The documentation says at the moment:

The semantics for how shapes behave under arithmetic operations are the same as for NumPy ndarrays (except some broadcasting is banned).

The examples above can be easily fixed by changing shape of x or by use of np.newaxis. But there are cases where I don't see an easy replacement for complex broadcasting rules other than (slow?) loops.

Example use case:
I'm currently struggling with an error of the form ValueError: Cannot broadcast dimensions (589, 1) (1, 589). I have a 1d boolean Variable x of length 589 and need a constraint for each i,j in range(589) if x[i] != 0 and x[j] != 0. So I thought I could do something like y >= x[np.newaxis, :] + x[:, np.newaxis] - 1 (for a non-negative constant array y of shape (589, 589)). But cvxpy's broadcasting rules forbid this.
_Edit:_ I've managed to avoid broadcasting by using y >= cp.atoms.affine.vstack.vstack([x] * len(x)).T + cp.atoms.affine.vstack.vstack([x] * len(x)) - 1, but this vectorized form turns out to be slower by a factor of two in comparison to the form [y[i] >= x[i] + x - 1 for i in range(len(x))] where one axis is replaced by a loop comprehension. I'll keep testing this a bit and summarize results here.

[lumbric]

Just had the same troubles again and found my own comment.

@SteveDiamond should we add this in the documentation? I'd be happy to help here, but haven't fully understood the precise differences to numpy broadcasting yet.

[SteveDiamond]

The problem is everyone I've taught in person has found the broadcasting behavior: (n,1) + (n,) = (n,n) incredibly confusing. I would be fine with adding all other broadcasting.

So it's fine to add (n,1) + (1,n) = (n,n).

[SteveDiamond]

Another comment is my experience implementing indexing is that matching numpy behavior exactly is way way more difficult than you would think. But it would be ok to add in broadcasting rules one by one.

[SteveDiamond]

I added broadcasting multiply( (n,1) , (1,m) ) to master. Getting broadcasting working for addition is a little more complicated, but the basic principle is to replicate using np.ones((589, 1)) @ x[None, :] + x[:, None] @ np.ones((1, 589)).

[SteveDiamond]

I added broadcasting for addition on master, so I'm going to close this issue.