cholesky on Hermitian matrix leads to scalar indexing error

[mfherbst]

trafficstars

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]

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