JuliaSIMD / JuliaSIMD/LoopVectorization.jl

Custom A*X .+ b returning wrong results

Aperta
#379 1 commento 0 reazioni 0 assegnatari Vedi su GitHub
Lingua principale
Julia
Stelle
789
Fork
73
Metriche di merge delle PR
Nessuna PR unita negli ultimi 30g

Descrizione

First of all, great package! I was trying to implement a slight alteration to the matrix multiplication example from the `Readme.md`:
```julia
function custom_gemm!(C::Matrix{T}, A::Matrix{T}, B::Matrix{T}, b::Vector{T}) where { T }
@turbo for m ∈ axes(A,1), n ∈ axes(B,2)
Cmn = b[m]
for k ∈ axes(A,2)
Cmn += A[m,k] * B[k,n]
end
C[m,n] = Cmn
end
end
```
where I only changed the original line `Cmn = zero(eltype(C))` to `Cmn = b[m]`. Basically I want this function to return `A*X .+ b`. Although this change seems trivial in my eyes, the implementation gives:
```julia
julia> using LoopVectorization

julia> Y, A, X, b = randn(2,2), randn(2,2), randn(2,2), ones(2)

julia> A*X .+ b
2×2 Matrix{Float64}:
1.44203 0.534377
0.931009 1.14601

julia> custom_gemm!(Y, A, X, b)
2×2 Matrix{Float64}:
1.0 1.0
1.0 1.0
```
It seems that the inner `for`-loop is ignored and `Cmn` is only set to `b[m]`. If I slightly change the above function the code works again:
```julia
function custom_gemm2!(C::Matrix{T}, A::Matrix{T}, B::Matrix{T}, b::Vector{T}) where { T }
@turbo for m ∈ axes(A,1), n ∈ axes(B,2)
Cmn = zero(T)
Cmn += b[m]
for k ∈ axes(A,2)
Cmn += A[m,k] * B[k,n]
end
C[m,n] = Cmn
end
end

julia> custom_gemm2!(Y, A, X, b)
2×2 Matrix{Float64}:
1.44203 0.534377
0.931009 1.14601
```
I am not sure what is causing this issue, but I can imagine that a lot of people might accidentally run into this.
```julia
julia> versioninfo()
Julia Version 1.7.0
Commit 3bf9d17731 (2021-11-30 12:12 UTC)
Platform Info:
OS: Windows (x86_64-w64-mingw32)
CPU: Intel(R) Core(TM) i7-9750H CPU @ 2.60GHz
WORD_SIZE: 64
LIBM: libopenlibm
LLVM: libLLVM-12.0.1 (ORCJIT, skylake)
Environment:
JULIA_NUM_THREADS = 12
```
`LoopVectorization v0.12.101`

Guida per i contributori

Nessuna guida per i contributori indicizzata per questo repository

Direzione di ricerca

Start with the matrix multiplication example in Readme.md and run the reported Julia reproduction using LoopVectorization v0.12.101. Compare custom_gemm! with custom_gemm2! and the reference A*X .+ b result. The issue is resolved when the direct b[m] initialization produces the same matrix as the reference computation.

Scritto dal modello di indicizzazione a partire dal testo della issue.

Valutazione

Stack tecnologico
julia
Ambito
performance
Tipo di issue
Bug
Difficoltà
4/5
Tempo stimato
3-5 giorni
Stato di attività
Ferma
Chiarezza
Specificata chiaramente
Idoneità per principianti
45/100

Ricevi le nuove issue nella tua casella

Un breve riepilogo di issue GitHub adatte ai principianti.