Evaluate MvNormal pdf on gpu array

General Usage
1

Is there a way to evaluate MvNormal densities on GPU data?

using Distributions, CUDA

julia> logpdf(MvNormal((rand(2)), (diagm(rand(2)))), rand(2))
-1.251537609209612

julia> logpdf(MvNormal(cu(rand(2)), cu(diagm(rand(2)))), cu(rand(2)))
ERROR: Scalar indexing is disallowed.
Invocation of getindex resulted in scalar indexing of a GPU array.
This is typically caused by calling an iterating implementation of a method.
Such implementations *do not* execute on the GPU, but very slowly on the CPU,
and therefore should be avoided.

If you want to allow scalar iteration, use `allowscalar` or `@allowscalar`
to enable scalar iteration globally or for the operations in question.
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] errorscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:155
  [3] _assertscalar(op::String, behavior::GPUArraysCore.ScalarIndexing)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:128
  [4] assertscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:116
  [5] getindex
    @ ~/.julia/packages/GPUArrays/qt4ax/src/host/indexing.jl:50 [inlined]
  [6] scalar_getindex
    @ ~/.julia/packages/GPUArrays/qt4ax/src/host/indexing.jl:36 [inlined]
  [7] _getindex
    @ ~/.julia/packages/GPUArrays/qt4ax/src/host/indexing.jl:19 [inlined]
  [8] getindex
    @ ~/.julia/packages/GPUArrays/qt4ax/src/host/indexing.jl:17 [inlined]
  [9] logdet(C::Cholesky{Float32, CuArray{Float32, 2, CUDA.DeviceMemory}})
    @ LinearAlgebra ~/.julia/juliaup/julia-1.11.2+0.x64.linux.gnu/share/julia/stdlib/v1.11/LinearAlgebra/src/cholesky.jl:673
 [10] logdet
    @ ~/.julia/packages/PDMats/cAM9h/src/pdmat.jl:89 [inlined]
 [11] logdetcov
    @ ~/.julia/packages/Distributions/Qeo2J/src/multivariate/mvnormal.jl:263 [inlined]
 [12] mvnormal_c0
    @ ~/.julia/packages/Distributions/Qeo2J/src/multivariate/mvnormal.jl:101 [inlined]
 [13] _logpdf
    @ ~/.julia/packages/Distributions/Qeo2J/src/multivariate/mvnormal.jl:143 [inlined]
 [14] logpdf(d::MvNormal{Float32, PDMats.PDMat{…}, CuArray{…}}, x::Vector{Float64})
    @ Distributions ~/.julia/packages/Distributions/Qeo2J/src/common.jl:263
 [15] top-level scope


  [052768ef] CUDA v5.5.2
  [31c24e10] Distributions v0.25.113
Julia Version 1.11.2
Commit 5e9a32e7af2 (2024-12-01 20:02 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: 128 × AMD EPYC 7713 64-Core Processor
  WORD_SIZE: 64
  LLVM: libLLVM-16.0.6 (ORCJIT, znver3)
Threads: 1 default, 0 interactive, 1 GC (on 128 virtual cores)
Environment:
  LD_LIBRARY_PATH = /usr/local/nvidia/lib:/usr/local/nvidia/lib64
  JULIA_EDITOR = code
  JULIA_NUM_THREADS =