Threaded push! different for supposedly similar array initializations

I have been using threads to assemble sparse matrices. Now I recently wanted to change the initialization of the indices arrays of length nt=Threads.nthreads() from something like [Int[] for i=1:nt] to fill(Int[],nt), which I thought is nicer. Now I ran into the issue that the code does not seem to @sync properly anymore.

Here’s the minimal working example:

using SparseArrays

#do some calculations and push certain values + indices into arrays for each threadid.
function pusharrays!(i0,j0, is, js, xs, n) 
	@sync Threads.@threads for i=1:n
		id = Threads.threadid()
		for j=1:n
			x = randn()
			if x>3.0
				push!(is[id], i+i0)
				push!(js[id], j+j0)
				push!(xs[id], x)
			end
		end
	end
	return nothing
end

function doit_ugly(n=1000)
	nt = Threads.nthreads()

	# "ugly" array initialization
	is = [Int[] for i=1:nt] 
	js = [Int[] for i=1:nt]
	xs = [Float64[] for i=1:nt]

	pusharrays!(0,0,is,js,xs,n)
	pusharrays!(n,0,is,js,xs,n)

	return sparse(vcat(is...), vcat(js...), vcat(xs...), 2n, n)
end


function doit_fill(n=1000)
	nt = Threads.nthreads()

	# much nicer
	is = fill(Int[],nt)
	js = fill(Int[],nt)
	xs = fill(Float64[],nt)

	pusharrays!(0,0,is,js,xs,n)
	pusharrays!(n,0,is,js,xs,n)

	return sparse(vcat(is...), vcat(js...), vcat(xs...), 2n, n)  
end


doit_ugly() # ✓

doit_fill() # ⚡

The doit_fill() gives two errors between different runs:

ERROR: TaskFailedException
Stacktrace:
 [1] wait
   @ ./task.jl:322 [inlined]
 [2] threading_run(func::Function)
   @ Base.Threads ./threadingconstructs.jl:34
 [3] macro expansion
   @ ./threadingconstructs.jl:93 [inlined]
 [4] macro expansion
   @ ./task.jl:387 [inlined]
 [5] pusharrays!(i0::Int64, j0::Int64, is::Vector{Vector{Int64}}, js::Vector{Vector{Int64}}, xs::Vector{Vector{Float64}}, n::Int64)
   @ Main ./REPL[2]:3
 [6] doit_fill(n::Int64)
   @ Main ./REPL[4]:9
 [7] doit_fill()
   @ Main ./REPL[4]:2
 [8] top-level scope
   @ REPL[6]:1

    nested task error: BoundsError: attempt to access 1156-element Vector{Float64} at index [22]
    Stacktrace:
     [1] setindex!
       @ ./array.jl:839 [inlined]
     [2] push!
       @ ./array.jl:930 [inlined]
     [3] macro expansion
       @ ./REPL[2]:10 [inlined]
     [4] (::var"#3#threadsfor_fun#1"{Int64, Int64, Vector{Vector{Int64}}, Vector{Vector{Int64}}, Vector{Vector{Float64}}, Int64, UnitRange{Int64}})(onethread::Bool)
       @ Main ./threadingconstructs.jl:81
     [5] (::var"#3#threadsfor_fun#1"{Int64, Int64, Vector{Vector{Int64}}, Vector{Vector{Int64}}, Vector{Vector{Float64}}, Int64, UnitRange{Int64}})()
       @ Main ./threadingconstructs.jl:48

or

ERROR: ArgumentError: the first three arguments' lengths must match, length(I) (=40016) == length(J) (= 40816) == length(V) (= 39120)
Stacktrace:
 [1] sparse(I::Vector{Int64}, J::Vector{Int64}, V::Vector{Float64}, m::Int64, n::Int64, combine::Function)
   @ SparseArrays /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.6/SparseArrays/src/sparsematrix.jl:736
 [2] sparse(I::Vector{Int64}, J::Vector{Int64}, V::Vector{Float64}, m::Int64, n::Int64)
   @ SparseArrays /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.6/SparseArrays/src/sparsematrix.jl:956
 [3] doit_fill(n::Int64)
   @ Main ./REPL[4]:12
 [4] doit_fill()
   @ Main ./REPL[4]:2
 [5] top-level scope
   @ REPL[6]:1

I don’t see why the fill version should cause issues here, as [Int[] for i=1:nt]==fill(Int[],nt). Maybe this is already known and I don’t know the right words to search for?

julia> Y = fill(Int[],3);

julia> push!(Y[1], 9)
1-element Vector{Int64}:
 9

julia> Y
3-element Vector{Vector{Int64}}:
 [9]
 [9]
 [9]


julia> X = [Int[] for i=1:3]

julia> push!(X[1], 9)
1-element Vector{Int64}:
 9

julia> X
3-element Vector{Vector{Int64}}:
 [9]
 []
 []
1 Like

Ah, okay I didn’t know this :man_facepalming:

Is there another way to do a nicer initialization that works?

1 Like