Allocation goes up with integration time when using CUDA for mesolve

[Lightup1]

Bug Description

The allocation goes up with integration time when using CUDA for mesolve.

Code to Reproduce the Bug

using QuantumToolbox
using CUDA
CUDA.allowscalar(false) # Avoid unexpected scalar indexing
using BenchmarkTools 
##
N = 20 # cutoff of the Hilbert space dimension
ω = 1.0 # frequency of the harmonic oscillator
γ = 0.1 # damping rate

a_gpu = cu(destroy(N)) # The only difference in the code is the cu() function

H_gpu = ω * a_gpu' * a_gpu

ψ0_gpu = cu(fock(N, 3))

c_ops = [sqrt(γ) * a_gpu]
e_ops = [a_gpu' * a_gpu]

tlist = [0,100] # time list
sol = mesolve(H_gpu, ψ0_gpu, tlist, c_ops, e_ops = e_ops)
@benchmark mesolve($H_gpu, $ψ0_gpu, $tlist, $c_ops, e_ops = $e_ops,progress_bar=$Val(false))

Code Output

BenchmarkTools.Trial: 34 samples with 1 evaluation per sample.
 Range (min … max):  117.296 ms … 212.421 ms  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     144.434 ms               ┊ GC (median):    0.00%
 Time  (mean ± σ):   149.778 ms ±  18.603 ms  ┊ GC (mean ± σ):  0.00% ± 0.00%

                █▄▄▁  █              ▁                           
  ▆▁▆▁▁▁▁▆▁▆▁▆▁▆████▆▁█▆▆▁▆▁▁▁▁▁▆▁▆▆▁█▁▆▁▁▆▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▆ ▁
  117 ms           Histogram: frequency by time          212 ms <

 Memory estimate: 2.93 MiB, allocs estimate: 125233.

Expected Behaviour

The allocation should not go up with integration time when there is no additional saving.

Your Environment

julia> QuantumToolbox.about()

 QuantumToolbox.jl: Quantum Toolbox in Julia
≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡
Copyright © QuTiP team 2022 and later.
Current admin team:
    Alberto Mercurio and Yi-Te Huang

Package information:
====================================
Julia              Ver. 1.11.6
QuantumToolbox     Ver. 0.34.0
SciMLOperators     Ver. 1.5.0
LinearSolve        Ver. 3.26.0
OrdinaryDiffEqCore Ver. 1.30.0

System information:
====================================
OS       : Linux (x86_64-linux-gnu)
CPU      : 12 × AMD Ryzen 5 2600X Six-Core Processor
Memory   : 7.726 GB
WORD_SIZE: 64
LIBM     : libopenlibm
LLVM     : libLLVM-16.0.6 (ORCJIT, znver1)
BLAS     : libopenblas64_.so (ilp64)
Threads  : 1 (on 12 virtual cores)

julia> CUDA.versioninfo()
CUDA toolchain: 
- runtime 13.0, local installation
- driver 580.97.0 for 13.0
- compiler 13.0

CUDA libraries: 
- CUBLAS: 13.0.0
- CURAND: 10.4.0
- CUFFT: 12.0.0
- CUSOLVER: 12.0.3
- CUSPARSE: 12.6.2
- CUPTI: 2025.3.0 (API 130000.0.0)
- NVML: 13.0.0+580.76.4

Julia packages: 
- CUDA: 5.8.3
- CUDA_Driver_jll: 13.0.0+0
- CUDA_Compiler_jll: 0.2.0+1
- CUDA_Runtime_jll: 0.19.0+0
- CUDA_Runtime_Discovery: 1.0.0

Toolchain:
- Julia: 1.11.6
- LLVM: 16.0.6

Preferences:
- CUDA_Runtime_jll.local: true

1 device:
  0: NVIDIA GeForce RTX 2060 (sm_75, 1.769 GiB / 6.000 GiB available)

Additional Context

No response

[albertomercurio]

It is a "problem" of CUDA.jl, as the mul! function for CUDA arrays allocates

using LinearAlgebra
using SparseArrays
using CUDA
using CUDA.CUSPARSE
using BenchmarkTools

# %%

N = 1000
A = CuSparseMatrixCSC(sprand(N, N, 0.01))
x = CuVector(rand(N))
y = similar(x)

mul!(y, A, x)

@benchmark mul!($y, $A, $x)

BenchmarkTools.Trial: 10000 samples with 4 evaluations per sample.
 Range (min … max):  7.496 μs …  29.064 μs  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     7.940 μs               ┊ GC (median):    0.00%
 Time  (mean ± σ):   8.181 μs ± 733.870 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

     ▁▇█▇▃                                                     
  ▂▃▄██████▆▅▅▆▆▆▆▅▄▄▄▄▄▄▅▅▄▄▄▄▄▃▃▃▃▃▂▂▂▂▂▂▂▂▂▂▂▂▂▂▁▂▁▁▁▁▁▂▁▂ ▃
  7.5 μs          Histogram: frequency by time        10.5 μs <

 Memory estimate: 1.69 KiB, allocs estimate: 80.

However, it is a very small quantity, probably originating from the call of CUDA internal routines. Also, it doesn't scale with the array size.

Long time ago I opened this issue

• Memory allocations during in-place sparse matrix-vector multiplication JuliaGPU/CUDA.jl#1982

[ytdHuang]

@Lightup1

any further questions or comments ?

If not, maybe we can close this issue.

[Lightup1]

Thanks.