How to prevent matrices decaying into vectors?

New to Julia
1

I have a function f that can take matrices of various sizes as inputs. When I try to call it with matrices of size d \times 1, produced by m[:, k] where k is an integer, they get converted to vectors of size d. How to prevent this? This must have been discussed many times but I could not find any solutions. Below is an MWE:

function f(x1, x2)
    d, n1 = size(x1)
    d, n2 = size(x2)

    out = zeros(d, n1 * n2)
end

function test()
    x1 = randn(2, 3)
    x2 = randn(2, 5)

    # Works fine
    y = f(x1, x2[:, 1:2])

    # Does not work
    y = f(x1, x2[:, 1])

    display(y)
end

test()
2

Like this:

julia> A = randn(2,2)
2Ă—2 Array{Float64,2}:
  2.66127  0.568908
 -1.40935  0.118752

julia> A[:,1]
2-element Array{Float64,1}:
  2.6612671809730895
 -1.409351022092897 

julia> A[:,1:1]
2Ă—1 Array{Float64,2}:
  2.6612671809730895
 -1.409351022092897
3 Likes
3

Many routines in Base accept either a matrix or a vector, where the latter is treated as a matrix with one column. This is a pretty convenient pattern for linear-algebra code.

The way to do this is to simply pass dimension arguments to size, for example changing your code above to:

d, n1 = size(x1,1), size(x1,2)

The size function returns 1 for dimensions after the last one, so size(x1,2) will be 1 if x1 is a vector (1d array). Note if x1 is a vector, you are also allowed to index it as x1[i,1], i.e. you can pass 1 for “extra” dimensions, to treat it as a matrix with one column.

You can even declare this explicitly in your argument signature as function f(x1::AbstractVecOrMat, x2::AbstractVecOrMat), where AbstractVecOrMat{T} is a built-in type alias for Union{AbstractVector{T},AbstractMatrix{T}}.

4

I just tried this and was a little surprised that we throw an error:

julia> v = rand(5);

julia> d1, d2 = size(v)
ERROR: BoundsError: attempt to access (5,)
  at index [2]

We could potentially allow it and just keep returning 1 as long as someone’s asking.

2 Likes
5

Specifically for size? How would that generalize?

1 Like
6

I have the same question as @tkoolen : size(v) is an expression that has a value regardless of the LHS, how could it know how long to run?

I consider automatic expansion by 1s a wart in the linear algebra design, but I am fine with it on demand. Eg expandedsize(v, n) could

  1. expand with 1s whenever ndims(v) <= n,
  2. error otherwise.

This would still allow combining code for vectors and matrices.

7

Yeah, it’s probably a bad idea now that I think about it a bit more. It would require returning a special Dims object that acts like an infinite tuple.

8

We need to be able to dispatch on the number of return values.

9

As I understand it, multiple return values are basically returned as a Tuple - do you mean dispatching on the return type?

10

It was mostly a joke (Number of function output arguments - MATLAB nargout - MathWorks Nordic).

2 Likes