# Most efficient way to find cholesky decomposition of slices of a 3D array in KernelAbstractions

**URL:** <https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292>\
**Category:** GPU\
**Tags:** question\
**Created:** [August 16, 2024, 8:33pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292 "2024-08-16T20:33:41Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![Jarrod-Angove](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarrod-angove/32/42663_2.png) [@Jarrod-Angove](https://discourse.julialang.org/u/Jarrod-Angove)\
**Post date:** [August 16, 2024, 8:33pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292/1 "2024-08-16T20:33:41Z")

</div>

I am using [KernelAbstractions.jl](https://github.com/JuliaGPU/KernelAbstractions.jl) to write a machine learning package that is device-agnostic. At one point in my algorithm, I end up with a NxNxn tensor, let’s call it `T`, on the device. I need to create a new vector, let’s call it `L`, that contains the lower triangular component of the cholesky decomposition of the slice `i`: `L[i] = cholesky(T[:,:,i]).L.data`. Is there a way to implement this in `KernelAbstractions.jl` that doesn’t involve me implementing my own cholesky kernel? For clarity, here is a sequential example of the full problem:

```julia
logp = Vector{Float32}(undef, n)
for k in 1:n
    L = cholesky(T[:,:,k])
    a = L.U \ (L.L \ y)
    logp[k] = -0.5 * y' * a - tr(L.L) - N/2 log(2pi)
end
return logp

```

Thank you!

---

<div class="post-metadata">

**Author:** ![amontoison](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/amontoison/32/218741_2.png) [@amontoison](https://discourse.julialang.org/u/amontoison)\
**Post date:** [August 17, 2024, 6:06pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292/2 "2024-08-17T18:06:18Z")

</div>

We have a dense method `cholesky` for NVIDIA and AMD arrays.  
It should work if you just rely on the multiple dispatch of `cholesky` and `ldiv!`.

---

<div class="post-metadata">

**Author:** ![biona001](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/biona001/32/16497_2.png) [@biona001](https://discourse.julialang.org/u/biona001)\
**Post date:** [August 17, 2024, 8:55pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292/3 "2024-08-17T20:55:58Z")

</div>

Btw, last time I checked `L.L` creates a copy. You can avoid this by `UpperTriangular(L.factors)'` (if `L.uplo=='U'`) or `LowerTriangular(L.factors)` (if L.uplo==‘L’)

---

<div class="post-metadata">

**Author:** ![Jarrod-Angove](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarrod-angove/32/42663_2.png) [@Jarrod-Angove](https://discourse.julialang.org/u/Jarrod-Angove)\
**Post date:** [August 21, 2024, 11:00pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292/4 "2024-08-21T23:00:59Z")

</div>

Sorry, I should have made my question a little more clear. I have already been running my code on GPU using multiple dispatch, I’m just concerned about the latency from spinning up n kernels. In my case, n is on the order of 10,000, while N is only on the order of 100, so I’m worried that launching the kernels is creating a lot of overhead. I guess the better question would have been this:

Is there some way to fuse the kernels for all n cholesky calls? In the KernelAbstractions.jl docs it just has an empty section for fusing kernels, so I was wondering if there was some undocumented way to go about this.

Thanks for your help!

---

<div class="post-metadata">

**Author:** ![Jarrod-Angove](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarrod-angove/32/42663_2.png) [@Jarrod-Angove](https://discourse.julialang.org/u/Jarrod-Angove)\
**Post date:** [August 21, 2024, 11:01pm UTC](https://discourse.julialang.org/t/most-efficient-way-to-find-cholesky-decomposition-of-slices-of-a-3d-array-in-kernelabstractions/118292/5 "2024-08-21T23:01:30Z")

</div>

This is a good suggestion, thanks! I’ll check the performance improvement on this