GPU memory from Python: reading a scene object where it already is¶

A recipe that works on the GPU should not be handed a host copy. Every uLib object that holds bulk numbers — a VoxImage, a collection, an EventTable — has a PreferredDevice property, and when its numbers are in GPU memory it offers a __cuda_array_interface__: the handshake numba, CuPy and PyTorch all read. Pass the object itself and nothing crosses the bus.

image.SetPreferredDevice("Device")        # the voxels move to the GPU

from numba import cuda
values = cuda.as_cuda_array(image)        # a view of VRAM, not a copy
print(values.shape, values.dtype, values.strides)   # (32768,) float32 (8,)

The C++ side of this — DeviceMemory, DeviceStorable, IDeviceArray, the placement rules — is Data on the GPU.

Two things to know first¶

No CUDA toolkit is needed to keep data on the GPU. Allocating device memory and copying to and from it is the driver's job, and uLib loads the driver at run time. Compiling a kernel is what needs the toolkit, and that is numba's or CuPy's business. So on a machine with a GPU but no toolkit, cuda.as_cuda_array(image) gives you a real device array that you can copy from and hand to a library that brings its own kernels — while @cuda.jit will not compile. uLib.device.kernels_available() tells the two apart; uLib.device.available() is about the memory alone.

With the toolkit installed (condaenv.yml), a kernel reads and writes the voxels in place:

from numba import cuda

image.SetPreferredDevice("Device")
values = cuda.as_cuda_array(image)       # a view of VRAM

@cuda.jit
def scale(values, factor):
    i = cuda.grid(1)
    if i < values.size:
        values[i] *= factor

scale[(values.size + 127) // 128, 128](values, 3.0)
print(image.GetValue(7))                 # 3x what it was: the host read
                                         # brings the block home, change and all

If kernels_available() is false where you expected true, it is almost always libcudart rather than libNVVM: numba asks the runtime which architectures NVVM accepts and reports none when it cannot find it, which reads like "your GPU is too old". A toolkit outside the environment needs CUDA_HOME=<prefix> and <prefix>/targets/x86_64-linux/lib on LD_LIBRARY_PATH.

A host read brings the data home. image.GetValue(i) copies the block back to host memory; the object's next update puts it back where it is preferred. So Device pays off when what consumes the data works on the device too, and costs two copies a round when it does not. Ask uLib.device.location(image) at any moment:

Device: Tesla P40
Host (no NVIDIA driver library: libcuda.so.1)
Host (these voxels cannot be copied to a device)

The property¶

image.GetPreferredDevice()        # "Host" or "Device"
image.SetPreferredDevice("Device")  # -> True if the voxels actually got there
image.IsDeviceResident()          # where they are right now
image.CanUseDevice()              # whether Device could be honoured at all
image.GetDataLocation()           # the sentence above
image.ReleaseDevice()             # bring them home and give the memory back

SetPreferredDevice returns where the data ended up, not whether the preference was recorded — it is always recorded. Device is honoured for data whose bytes are the whole of it (voxels, floats, lines, MuonScatter; not a GeantEvent, which holds a path and a name) on a machine whose driver gives out memory. Otherwise the data stays on the host and says so.

In a .gcx the preference is saved; the location is not, because it is a fact about one machine. A document written before the property existed still loads — the field is behind the image's class version.

`__cuda_array_interface__`, and what the stride is for¶

image.__cuda_array_interface__
# {'shape': (32768,), 'typestr': '<f4', 'data': (139576886689792, False),
#  'strides': (8,), 'version': 3, 'stream': None}

A voxel is a Value and a Count, and only the Value is wanted — so the view exposes one float every 8 bytes rather than a repacked array. numba and CuPy both handle that. strides is None where the data really is packed (an EventTable's doubles, a collection of floats).

The attribute is None, not a broken dictionary, when the object cannot be given a device array — no GPU, or data that cannot be byte-copied. A consumer that finds the attribute finds something usable behind it.

What each object offers:

Object	Exposed	`typestr`
`VoxImage`	each voxel's `Value`, strided	`<f4`
`MultiVoxImage`	the active channel and component	`<f4`
`EventTable`	the values row-major, packed	`<f8`
`MuonScatterCollection`	the events' raw bytes	`\\|u1`

`uLib.device`¶

The helpers, for when you do not want to write the fallbacks yourself:

from uLib import device

device.available()            # can data be kept on the GPU at all
device.unavailable_reason()   # why not, in words
device.name()                 # "Tesla P40"
device.memory_info()          # (free, total) bytes
device.held()                 # (bytes, allocations) uLib holds right now
device.describe()             # one line for a log or a panel
device.kernels_available()    # can a kernel be *compiled* here

device.location(obj)          # where this object's data is
device.is_device_resident(obj)
device.array_info(obj)        # pointer, count, itemsize, stride, typestr
device.as_array(obj)          # numba or CuPy array, whichever is installed
device.as_array(obj, copy_to_host=True)   # numpy, in one bulk copy

as_array returns None for an object that offers no device array, and raises RuntimeError only when the data is on the device and neither numba nor CuPy is installed to take it. copy_to_host=True is the pragmatic path on a machine with no kernel compiler: one bulk copy still beats reading 32768 voxels one call at a time.

Every one of these works through a recipe's object proxy as well as on a raw bound object.

In a recipe¶

ctx.<name> hands back the scene object, and __cuda_array_interface__ reaches through the proxy, so a recipe reads the voxels where they are:

import numpy as np

SIDE = 32

def run(ctx):
    image = ctx['phantom']
    print('voxels are', image.GetDataLocation())

    try:
        import cupy
        block = cupy.asarray(image).reshape((SIDE, SIDE, SIDE))
        profile = cupy.asnumpy(block.sum(axis=(0, 1)))   # summed on the GPU
        how = 'on the GPU with CuPy, reading the voxels in place'
    except ImportError:
        from numba import cuda
        view = cuda.as_cuda_array(image)                 # still no host copy
        profile = view.copy_to_host().reshape(
            (SIDE, SIDE, SIDE)).sum(axis=(0, 1))         # one bulk copy back
        how = 'numba took the device view; summed on the host'

    print('summed', how)
    table = ctx['profile']
    table.Columns = 'slice,total'
    table.Clear()
    for index, total in enumerate(profile):
        table.AddRow([float(index), float(total)])
    ctx.refresh()

Set the recipe's own PreferredDevice to Device as well, and the scene reads as what it is: a GPU stage between two objects that keep their data there.

A whole scene to open¶

gcompose assets/examples/math/02_images_on_the_gpu.gcx

Built by src/Math/examples/02_images_on_the_gpu.cpp, which is also run by ctest: an image of 32³ voxels with PreferredDevice = Device, an EventTable beside it with the same setting, and the recipe above. It prints where everything is at each step, and on a machine with no GPU it prints that instead and computes the same numbers on the host.