JuliaLinearAlgebra / Octavian.jl

Multi-threaded BLAS-like library that provides pure Julia matrix multiplication
https://julialinearalgebra.github.io/Octavian.jl/stable/
Other
222 stars 17 forks source link

nested task error: ReadOnlyMemoryError() when using `matmul!` #191

Open JinyanTeng opened 1 month ago

JinyanTeng commented 1 month ago

When using 30 threads to run the command matmul!(xtx, X', X), I encountered an error. However, when using 20 threads, no error was observed. It is important to note that this error is not consistently reproducible and may vary for different matrices X. In this case, the matrix X has dimensions of 271 by 24.

ERROR: LoadError: TaskFailedException

nested task error: ReadOnlyMemoryError()
Stacktrace:
  [1] macro expansion
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/llvm_intrin/memory_addr.jl:1482 [inlined]
  [2] __vstore!
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/llvm_intrin/memory_addr.jl:1482 [inlined]
  [3] _vstore!
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/strided_pointers/stridedpointers.jl:148 [inlined]
  [4] macro expansion
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/vecunroll/memory.jl:872 [inlined]
  [5] _vstore_unroll!
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/vecunroll/memory.jl:1174 [inlined]
  [6] _vstore!
    @ ~/.julia/packages/VectorizationBase/6AO0m/src/vecunroll/memory.jl:1803 [inlined]
  [7] macro expansion
    @ ~/.julia/packages/LoopVectorization/QgYWB/src/reconstruct_loopset.jl:1107 [inlined]
  [8] _turbo_!
    @ ~/.julia/packages/LoopVectorization/QgYWB/src/reconstruct_loopset.jl:1107 [inlined]
  [9] macro expansion
    @ ~/.julia/packages/LoopVectorization/QgYWB/src/condense_loopset.jl:1179 [inlined]
 [10] packamul!
    @ ~/.julia/packages/Octavian/iOEo8/src/macrokernels.jl:62 [inlined]
 [11] packaloopmul!(C::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, A::LayoutPointers.StridedPointer{Float64, 2, 2, 0, (2, 1), Tuple{Int64, Static.StaticInt{8}}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, B::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, α::Static.StaticInt{1}, β::Static.StaticInt{0}, M::Int64, K::Int64, N::Int64)
    @ Octavian ~/.julia/packages/Octavian/iOEo8/src/macrokernels.jl:209
 [12] matmul_st_only_pack_A!(C::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, A::LayoutPointers.StridedPointer{Float64, 2, 2, 0, (2, 1), Tuple{Int64, Static.StaticInt{8}}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, B::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, α::Static.StaticInt{1}, β::Static.StaticInt{0}, M::Int64, K::Int64, N::Int64, ::Static.StaticFloat64{0.0007423708195588264}, ::Static.StaticFloat64{0.7757548987718677}, ::Static.StaticFloat64{0.7936663315339363}, ::Static.StaticFloat64{0.7144577794375783})
    @ Octavian ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:48
 [13] matmul_st_pack_dispatcher!(pC::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, pA::LayoutPointers.StridedPointer{Float64, 2, 2, 0, (2, 1), Tuple{Int64, Static.StaticInt{8}}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, pB::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, α::Static.StaticInt{1}, β::Static.StaticInt{0}, M::Int64, K::Int64, N::Int64)
    @ Octavian ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:375
 [14] __matmul!(C::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, A::LayoutPointers.StridedPointer{Float64, 2, 2, 0, (2, 1), Tuple{Int64, Static.StaticInt{8}}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, B::LayoutPointers.StridedPointer{Float64, 2, 1, 0, (1, 2), Tuple{Static.StaticInt{8}, Int64}, Tuple{Static.StaticInt{0}, Static.StaticInt{0}}}, α::Static.StaticInt{1}, β::Static.StaticInt{0}, M::Int64, K::Int64, N::Int64, nthread::Nothing, W₁::Static.StaticFloat64{0.0007423708195588264}, W₂::Static.StaticFloat64{0.7757548987718677}, R₁::Static.StaticFloat64{0.7936663315339363}, R₂::Static.StaticFloat64{0.7144577794375783})
    @ Octavian ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:721
 [15] __matmul!
    @ ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:710 [inlined]
 [16] _matmul!(C::Matrix{Float64}, A::LinearAlgebra.Adjoint{Float64, Matrix{Float64}}, B::Matrix{Float64}, α::Static.StaticInt{1}, β::Static.StaticInt{0}, nthread::Nothing, MKN::Nothing)
    @ Octavian ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:576
 [17] matmul!
    @ ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:520 [inlined]
 [18] matmul!
    @ ~/.julia/packages/Octavian/iOEo8/src/matmul.jl:437 [inlined]
chriselrod commented 1 month ago

When using 30 threads to run the command matmul!(xtx, X', X), I encountered an error. However, when using 20 threads, no error was observed.

How many physical cores does your system have? Octavian should limit itself to that many. For a 24x271 * 271x24 matrix, Octavian will probably only use a single thread, so that's probably a red haring.

__vstore! showing the problem is a little surprising.

Can you

julia> using Octavian

julia> X = rand(271, 24);

julia> xtxb = fill(-12345.6789, 100, 100);

julia> xtx2 = @view(xtxb[40:63, 40:63]);

julia> Octavian.matmul!(xtx2, X', X);

julia> xtx2 ≈ X' * X
true

julia> all(==(first(xtxb)), @view(xtxb[1:39,:]))
true

julia> all(==(first(xtxb)), @view(xtxb[64:end,:]))
true

julia> all(==(first(xtxb)), @view(xtxb[40:63,1:39]))
true

julia> all(==(first(xtxb)), @view(xtxb[40:63,64:end]))
true 

Can you see if any of these checks fail, implying there was an out of bound access?

I tried on both an AVX512 and an AVX2 system, getting all true on both.

What is your

versioninfo()

?

chriselrod commented 1 month ago

I suspect the illegal write isn't to xtx, but to the temporary packed matrix X'.

What do you get for

julia> Octavian.first_cache_size(Val(eltype(X)))
static(65536)

julia> Octavian.first_cache_size(Val(eltype(X))) / length(X)
10.076260762607626

?

JinyanTeng commented 1 month ago

When using 30 threads to run the command matmul!(xtx, X', X), I encountered an error. However, when using 20 threads, no error was observed.

How many physical cores does your system have? Octavian should limit itself to that many. For a 24x271 * 271x24 matrix, Octavian will probably only use a single thread, so that's probably a red haring.

__vstore! showing the problem is a little surprising.

Can you

julia> using Octavian

julia> X = rand(271, 24);

julia> xtxb = fill(-12345.6789, 100, 100);

julia> xtx2 = @view(xtxb[40:63, 40:63]);

julia> Octavian.matmul!(xtx2, X', X);

julia> xtx2 ≈ X' * X
true

julia> all(==(first(xtxb)), @view(xtxb[1:39,:]))
true

julia> all(==(first(xtxb)), @view(xtxb[64:end,:]))
true

julia> all(==(first(xtxb)), @view(xtxb[40:63,1:39]))
true

julia> all(==(first(xtxb)), @view(xtxb[40:63,64:end]))
true 

Can you see if any of these checks fail, implying there was an out of bound access?

I tried on both an AVX512 and an AVX2 system, getting all true on both.

What is your

versioninfo()

?

cbcfea83f2465132e34d01f01e42635

0c6fa78c81a36547b78dd577e56311b

FYI.

chriselrod commented 1 month ago

Do you think you could get an rr trace of the error? That is, start Julia with --bug-report=rr and upload the trace?