Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

cholesky on Hermitian matrix leads to scalar indexing error #734

Open
mfherbst opened this issue Feb 15, 2025 · 1 comment
Open

cholesky on Hermitian matrix leads to scalar indexing error #734

mfherbst opened this issue Feb 15, 2025 · 1 comment

Comments

@mfherbst
Copy link
Contributor

Performing a cholesky factorisation on a Hermitian matrix

using AMDGPU
using LinearAlgebra

A = ROCArray(rand(ComplexF64, 8, 6))
G = A'A + I  # Positive definite matrix
cholesky(Hermitian(G))

leads to

ERROR: Scalar indexing is disallowed.

and a stacktrace

Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] errorscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/aNaXo/src/GPUArraysCore.jl:151
  [3] _assertscalar(op::String, behavior::GPUArraysCore.ScalarIndexing)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/aNaXo/src/GPUArraysCore.jl:124
  [4] assertscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/aNaXo/src/GPUArraysCore.jl:112
  [5] getindex
    @ ~/.julia/packages/GPUArrays/uiVyU/src/host/indexing.jl:50 [inlined]
  [6] scalar_getindex
    @ ~/.julia/packages/GPUArrays/uiVyU/src/host/indexing.jl:36 [inlined]
  [7] _getindex
    @ ~/.julia/packages/GPUArrays/uiVyU/src/host/indexing.jl:19 [inlined]
  [8] getindex
    @ ~/.julia/packages/GPUArrays/uiVyU/src/host/indexing.jl:17 [inlined]
  [9] getindex
    @ ~/stuffs/julia-1.11.3/share/julia/stdlib/v1.11/LinearAlgebra/src/symmetric.jl:249 [inlined]
 [10] _getindex
    @ ./abstractarray.jl:1358 [inlined]
 [11] getindex
    @ ./abstractarray.jl:1312 [inlined]
 [12] iterate
    @ ./abstractarray.jl:1209 [inlined]
 [13] iterate
    @ ./abstractarray.jl:1207 [inlined]
 [14] copyto_unaliased!(deststyle::IndexLinear, dest::ROCArray{…}, srcstyle::IndexCartesian, src::Hermitian{…})
    @ Base ./abstractarray.jl:1086
 [15] copyto!
    @ ./abstractarray.jl:1061 [inlined]
 [16] copy_similar(A::Hermitian{ComplexF64, ROCArray{ComplexF64, 2, AMDGPU.Runtime.Mem.HIPBuffer}}, ::Type{ComplexF64})
    @ LinearAlgebra ~/stuffs/julia-1.11.3/share/julia/stdlib/v1.11/LinearAlgebra/src/LinearAlgebra.jl:459
 [17] eigencopy_oftype
    @ ~/stuffs/julia-1.11.3/share/julia/stdlib/v1.11/LinearAlgebra/src/symmetriceigen.jl:4 [inlined]
 [18] cholcopy
    @ ~/stuffs/julia-1.11.3/share/julia/stdlib/v1.11/LinearAlgebra/src/cholesky.jl:182 [inlined]
 [19] cholesky
    @ ~/stuffs/julia-1.11.3/share/julia/stdlib/v1.11/LinearAlgebra/src/cholesky.jl:401 [inlined]
 [20] cholesky(A::Hermitian{ComplexF64, ROCArray{ComplexF64, 2, AMDGPU.Runtime.Mem.HIPBuffer}})

I suppose this is due to a missing method for Hermitian{T, <:ROCArray}, because a cholesky on a plain ROCArray works fine.

@mfherbst
Copy link
Contributor Author

For reference: For now I added a workaround, manually calling the respective rocSOLVER routine:
https://github.com/JuliaMolSim/DFTK.jl/blob/907fbf50d03aecda70e4a30fdb6c59f8e00c9537/ext/DFTKAMDGPUExt.jl#L10

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant