On Sat, May 30, 2015 at 3:23 PM, Charles R Harris <charlesr.harris@gmail.com> wrote:

The problem arises when multiplying a stack of matrices times a vector. PEP465 defines this as appending a '1' to the dimensions of the vector and doing the defined stacked matrix multiply, then removing the last dimension from the result. Note that in the middle step we have a stack of matrices and after removing the last dimension we will still have a stack of matrices. What we want is a stack of vectors, but we can't have those with our conventions. This makes the result somewhat unexpected. How should we resolve this?

I'm afraid I don't quite understand the issue. Maybe a more specific example of the shapes you have in mind would help? Here's my attempt.

Suppose we have two arrays:

a with shape (i, j, k)

b with shape (k,)

Following the logic you describe from PEP465, for a @ b we have shapes transform like so:

(i, j, k,) @ (k, 1) -> (i, j, 1) -> (i, j)

This makes sense to me as a stack of vectors, as long as you are imagining the original stack of matrices as along the first dimension. Which I'll note is the default behavior for the new np.stack (https://github.com/numpy/numpy/pull/5605).