Quantum memory management within CUDA-Q kernels does not work

[renezander90]

Required prerequisites

• Consult the security policy. If reporting a security vulnerability, do not report the bug using this form. Use the process described in the policy to report the issue.
• Make sure you've read the documentation. Your issue may be addressed there.
• Search the issue tracker to verify that this hasn't already been reported. +1 or comment there if it has.
• If possible, make a PR with a failing test to give us a starting point to work on!

Describe the bug

When sequentially executing a sub-kernel that allocates a local cudaq.qvector, the quantum memory management fails to actually free and reuse the allocated qubits after the sub-kernel returns. This causes total qubit usage to scale linearly with the number of calls, eventually causing a crash due to qubit exhaustion/OOM.

Crucially, inspecting the MLIR confirms that the frontend correctly inserts the quake.dealloc instruction. The issue appears to be that the execution runtime or simulator backend is ignoring this instruction or failing to return the qubits to the available pool.

Please note that many quantum algorithms will also require memory management within a single kernel (e.g., for allocating ancillary qubits for an mcx gate, and deallocating them afterwards https://www.iccs-meeting.org/archive/iccs2022/papers/133530169.pdf). But, this might rather be a feature request.

Steps to reproduce the bug

import cudaq

@cudaq.kernel
def test_kernel() -> int:
    q = cudaq.qvector(10)
    h(q[0])
    return mz(q[0])

@cudaq.kernel
def main_kernel() -> int:

    a = test_kernel()
    b = test_kernel()
    c = test_kernel()
    return a + b + c

# print(test_kernel)
cudaq.run(main_kernel, shots_count=10)

This will cause a crash due to excessive memory usage. (It will work for 1 or 2 calls of test_kernel).

Expected behavior

The quantum resources go out of scope at the end of each sub-kernel. The program contains a quake.dealloc statement:

    %0 = quake.alloca !quake.veq<10>
    %1 = quake.extract_ref %0[0] : (!quake.veq<10>) -> !quake.ref
    quake.h %1 : (!quake.ref) -> ()
    %2 = quake.extract_ref %0[0] : (!quake.veq<10>) -> !quake.ref
    %measOut = quake.mz %2 : (!quake.ref) -> !quake.measure
    %3 = quake.discriminate %measOut : (!quake.measure) -> i1
    %4 = cc.cast unsigned %3 : (i1) -> i64
    quake.dealloc %0 : !quake.veq<10>
    return %4 : i64

Expected behavior would be that the maximum number of qubits does not exceed 10. But, simulation time and memory usage grows drastically with the number of calls of test_kernelinside the main_kernel.

Is this a regression? If it is, put the last known working version (or commit) here.

Not a regression

Environment

• CUDA-Q version: 0.14.0
• Python version: 3.11.14
• Operating system: macOs 15.7.4

Suggestions

No response

[schweitzpgi]

Yes, this looks like a bug with the simulator losing track of the number of active qubits. It will need further investigation.

[schweitzpgi]

I suspect that the issue is that all of the device-only kernel calls are inlined into the main kernel and fused into a single larger kernel. (So in the case of the example, a single kernel with 30 qubits.) So it's not that they "don't work" exactly, just that the code is inadvertently trying to use way, way too many qubits.

The programming model, presently, is for the programmer to specify all the qubits needs in the entry-point kernel (i.e. main_kernel) and then pass them as arguments to the device-only kernel(s) (i.e. test_kernel). Using this style will workaround folding all the qubits defined in all the device-only kernels into one gigantic collection of qubits.

import cudaq

@cudaq.kernel
def test_kernel(q: cudaq.qview) -> int:
    h(q[0])
    return mz(q[0])

@cudaq.kernel
def main_kernel() -> int:
    qs = cudaq.qvector(10)
    a = test_kernel(qs)
    b = test_kernel(qs)
    c = test_kernel(qs)
    return a + b + c

print(test_kernel)
cudaq.run(main_kernel, shots_count=10)

[renezander90]

Thanks for investigating!

Please note that this behavior is different than allocating/deallocating new quilts for each test_kernel:

import cudaq

@cudaq.kernel
def test_kernel(q: cudaq.qview) -> int:
    x(q[0])
    return mz(q[0])

@cudaq.kernel
def main_kernel() -> int:
    qs = cudaq.qvector(10)
    a = test_kernel(qs)
    b = test_kernel(qs)
    return a + b

print(test_kernel)
cudaq.run(main_kernel, shots_count=10)

[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]

The result is always 1, because the second kernel acts on the same quantum state (flips the measured qubit).

The programming model, presently, is for the programmer to specify all the qubits needs in the entry-point kernel

Then I assume that the following would also not be possible currently:

1.) allocate 10 qubits q0
2.) allocate 5 more ancilla qubits q1, perform some computation with q0 + q1 after which q1 are safely uncomputed, deallocate qubits q1
3.) Repeat step 2. multiple times

In this case, the total qubit count would grow with each repetition of step 2.), and cannot stay 15?

[schweitzpgi]

Yes. It's likely to experience the same bug at the moment.

[schweitzpgi]

Needs verification that #4455 fixed the issue. (It added the reported reproducer to the regression tests.)