Vector of vectors non-allocating

General Usage
,
1

How do I make the function that follows non-allocating?

function genvec2(N)
    result = [zeros(Int64,2) for i=1:N] 
    for i=1:N
        result[i] = [2*i, 2i+1]  
    end 
    return result
end
2

In the case that you can fix N,

using StaticArrays

function genvec2()
           @SVector [ ( 2*i, 2*i+1 ) for i=1:10 ]
end

I’m not sure otherwise…

2 Likes
3

Like the above poster said, with a dynamic N it’s impossible to make this non-allocating in Julia (and, in general, for any language). However, there are many ways to make this happen with fewer allocations.

One way is to not generate all length=2 vectors in the results initialization, since you simply throw them away when you assign them in the loop. This will reduce your number of allocations by almost half. Alternatively, if you write result[i] .= (2*i, 2i+1) it will take the arrays you were allocating with zeros and overwrite them with those two values without generating a new (redundant) vector.

It’s also possible to reduce the number of allocations here to just one. I include two others here that manage this:

function genvec3(N)
    result = Vector{Vector{Int64}}(undef, N) # allocate one Vector
    for i = 1:N
        result[i] = [2i, 2i+1] # allocate one Vector per loop
    end 
    return result
end

function genvec4(N)
    result = eachcol(zeros(Int64, 2, N)) # allocate one Matrix but read it like many Vectors
    for i = 1:N
        result[i] .= (2i, 2i+1) # fill the sliced matrix with values
    end 
    return result
end

using StaticArrays
function genvec5(N)
    result = map(i -> @SVector[2i,2i+1], 1:N) # allocate one Vector
    return result
end

julia> @time genvec3(3) # 1+N allocations
  0.000003 seconds (4 allocations: 320 bytes)
3-element Vector{Vector{Int64}}:
 [2, 3]
 [4, 5]
 [6, 7]

julia> @time genvec4(3)
  0.000006 seconds (1 allocation: 112 bytes)
3-element ColumnSlices{Matrix{Int64}, Tuple{Base.OneTo{Int64}}, SubArray{Int64, 1, Matrix{Int64}, Tuple{Base.Slice{Base.OneTo{Int64}}, Int64}, true}}:
 [2, 3]
 [4, 5]
 [6, 7]

julia> @time genvec5(3)
  0.000003 seconds (1 allocation: 112 bytes)
3-element Vector{SVector{2, Int64}}:
 [2, 3]
 [4, 5]
 [6, 7]

The runtimes on all of these are too short for an actual comparison (try with bigger N), but note that the latter two allocate less. Note the return types of these functions are all slightly different, which might matter to you.

Personally, I prefer the version using StaticArrays.SVector.

2 Likes
4

If you need many of these vectors, it might make sense to allocate them outside of this function (in the place where you actually use it) and reuse that memory. In that case you would implement a mutating version of genvec2:

function genvec2!(vec)
    for i in eachindex(vec)
        vec[i] .= (2i, 2i+1)
    end
    return vec # optionally return
end

You probably have a bit more code that you want to optimize. If you want help with that, it is much easier for us when you post your actual code. Maybe simplify a bit if necessary but most importantly see that it can be executed so we can profile ourselves :slight_smile:

2 Likes
5

Many thanks! Will go back to the drawing board and get back to you.

6

I am still looking into the issue. I explored various options. My notebook on this issue is at Notebook Vector of Subvectors. Input is most welcome.

7

The last solution doesn’t need map. I think the rather straight-forward is just as fast:

function genvec8(N) 
    result = [@SVector[2*i, 2*i+1] for i=1:N ]
    return result
end

or equivalently (closer to your original code)

function genvec9(N)
    elemtype = typeof(@SVector[1,2])
    result = Vector{elemtype}(undef,N) # vector of SVector
    for i = eachindex(result)
        result[i] = @SVector[2i, 2i+1] # fill the slots (no broadcast!)
    end
    return result
end

To comment a bit on the other versions you have there:

  • genvec2 and genvec3 both use a Vector of Vector and thus need to allocate 1 outer Vector and N small Vectors.
  • genvec4 works for me (are you on an older Julia perhaps?) and needs only 1 allocation because it then uses views into the memory. However each view is an indirection and so it is not as fast as it could be.
  • genvec5 is really the very same thing as genvec3. In some sense a SVector is just a fancy tuple, so genvec5 and genvec3 identical.
  • genvec6 doesn’t work for me. But conceptually it would use a SVector of SVectors which is theoretically fine, but would cause a lot of trouble for the compiler as it now has to compile a new method for every length of the outer vector. That’s why StaticArrays.jl recommends not using it beyond ~100 elements or so.
  • genvec7 (and my genvec8 and genvec9) now do the most sensible thing: They use a Vector of SVector. That is best of both worlds. We can have a varying amount of small vectors without troubling the compiler but operations on the small vectors is still fast because Julia knows their length. In fact Julia can optimize away even the allocation of the small vectors because they all have the same known size and are thus stored inline inside the outer array (this is not possible with Vector because a Vector can grow in size thus needs to live in the heap separately). Btw: This approach is conceptually very close to genvec4. The data of the inner vectors is packed together in a single large array in memory. The main difference is just the way of accessing it. Here SVector makes a large difference compared to views, because again views can have varying sizes and thus Julia cannot optimize them away and so there is an additional layer of indirection.
2 Likes
8

@abraemer : sincere thanks!

It the actual use case we are considering to store the geometry of a block as a set of N sub-blocks. Each block will have 8 vertices. Each vertex will be stored as a 3-vector of x, y and z coordinates.

Can static vectors of length 8 contain static vectors of size 3 of type Float64?

9

I don’t see why not. Perhaps you want a MVector on the outside if you wish to change the contents at some point.

1 Like
10

@abraemer : Question: in line 2 of your genvec9() the call to @SVector has as arguments [1,2]. Why is this? Thx!

11

@SVector is a macro that helps constructing SVectors. You can find its documentation here or look at what it expands to in the REPL with:

julia> @macroexpand @SVector [1,2]
:((SVector){2}((tuple)(1, 2)))

which means that this is just a shortcut to writing: SVector{2}((1,2)) but you could alternatively also write any of the form listed in the documentation under “Create an SVector…”.

1 Like