Why does Julia support `::Vector * ::Matrix`?

For example,

julia> [1,2] * [2;;]
2×1 Matrix{Int64}:
 2
 4

This reshapes the vector as a n x 1 matrix and multiplies it with the 1 x 1 matrix. This seems to go against the general idea of treating Vectors and Matrices as separate entities.

In short, it’s because it’s convenient to have vectors participate in the algebra of matrices as though they are Nx1 columns and because — with transposed “row vectors” — there’s a comprehensive and consistent set of behaviors that allow it to do so sensibly. See Julia#4774 for the in-depth design discussion.

6 Likes

More generally, we treat AbstractArrays as having an infinite amount of length-1 trailing dimensions:

julia> size([1,2], 1)
2

julia> size([1,2], 2)
1

julia> size([1,2], 3)
1

julia> size([1,2], 4)
1

julia> [1, 2][2, 1, 1, 1, 1, 1, 1, 1]
2
1 Like

To be honest, I’m a bit surprised by vector-matrix multiply being accepted with “compatible inferred dimensions.” I suppose the reason it’s a thing is that it’s mostly harmless – i.e., there isn’t something else we would do with it.

Although it does give rise to an inconsistency with respect to vector-vector multiplication:

julia> [1;] * [1;;]
1×1 Matrix{Int64}:
 1

julia> [1;] * [1;] # why not return a 1-vector?
ERROR: MethodError: no method matching *(::Vector{Int64}, ::Vector{Int64})

An inferred-dimension interpretation (like justifies vector-matrix) would have accepted this vector-vector multiplication and returned a 1-vector. A justification against is that x * y was probably a typo of x' * y or x * y', so this is a nice backstop.

If I were Czar I would probably have left vector-matrix multiplication as an error at any size. Alas, to make it so now would be breaking (even if it’s likely uncommon).

3 Likes

Right, it’s not exactly like they’re Nx1 matrices. Again, I really encourage you to read through Julia#4774 or at least watch Jiahao’s :4774: JuliaCon presentation for more on the rationale here.

For Vector * Matrix specifically, see Julia#20389.

5 Likes

I’d revolt for sure. Being able to write v1 * v2' to quickly create a rank-1 matrix is incredibly useful and nice. Having to turn that into reshape(v1, length(v1), 1) * v2' is so much worse, especially when v_1 \, v_2^\dagger is already pretty common mathematical notation for the outer product of two vectors.

I’d say that the reason there’s no method for this is probably just that it’s not useful since it only comes up if the vectors are both of length 1, whereas multiplying a vector by a covector is suuuper common.

2 Likes

Fear not – I would revolt, too. Even I never imagined to remove x * y' for outer products. But y' is a Adjoint{<:Any,<:AbstractVector}, which we can (and do) specialize beyond what we do for general AbstractMatrix. There isn’t an incompatibility here.

But like I said: there’s little harm in the current behavior, no possibility to change it in v1.x, and still almost no reason to change it beyond that.

But then if you did collect(v2’) it’d just suddenly stop working. It seems quite weird and pointless to try and omit this method to me.