Chapter 1: Thinking in Threads: The CUDA Model with Numba
Write your first real CUDA kernels in Python with Numba. Understand the thread/block/grid launch hierarchy, how each thread finds its own global index, why every kernel needs a bounds check, and how grid-stride loops make a kernel work for any input size.
- Topics
- 4
- Videos
- 4
- Min read
- 12
Chapter Overview
CUDA is the language GPUs speak, but you do not need C++ to start speaking it. Numba compiles a subset of Python decorated with @cuda.jit straight to GPU machine code (PTX), so you can write genuine CUDA kernels — with threadIdx, blockIdx, shared memory, and atomics — in plain Python.
The mental shift this chapter asks of you is the one that matters most: on a CPU you write a loop that visits each element one at a time; on a GPU you write the body of the loop once and launch thousands of threads that each run that body on a different element, all at the same time.
This chapter builds the foundation everything else rests on:
- Why Numba: real CUDA kernels without leaving Python
- Threads, Blocks & Grids: the launch hierarchy you configure on every kernel
- The Global Index: how a thread figures out which element it owns
- Bounds Checks & Grid-Stride Loops: making a kernel correct for any size, with any launch configuration
Chapter Roadmap
Click any topic to jump in
Why Numba for CUDA
@cuda.jit compiles a Python subset to PTX — real CUDA kernels without a C++ toolchain.
Threads, Blocks & Grids
kernel[blocks, threads](...) configures the launch hierarchy; threadIdx/blockIdx/blockDim/gridDim locate a thread.
The Global Index
i = blockIdx.x*blockDim.x + threadIdx.x (or cuda.grid(1)) maps each thread to one element.
Bounds Checks & Grid-Stride Loops
if i < n guards the ragged tail; grid-stride loops decouple launch size from data size.
Topics
4 topics in this chapter
Programming an NVIDIA GPU has traditionally meant writing CUDA C/C++, installing the nvcc toolchain, and leaving Python behind to manage kernels, memory, and build steps by hand. That barrier keeps many numerical Python programmers away from the hardware their workloads would benefit from most. Numba removes it: it is a just-in-time compiler that turns a decorated subset of Python directly into GPU machine code, so you can write, launch, and debug real kernels without ever opening a C++ file. This first topic sets up why that matters before you write any parallel code. You will see that the CUDA programming model — threads, blocks, shared memory, atomics — is identical whether you express it in C++ or in Python, so nothing you learn here is throwaway. What changes is only the surface syntax and the speed of iteration. The rest of the chapter builds on this foundation: once you trust that @cuda.jit emits genuine device code, the launch hierarchy, the global index, and bounds checking all follow as the mechanics of telling thousands of threads what to do.
Definition
Numba is a just-in-time compiler that translates a restricted subset of Python and NumPy into optimized machine code. Its cuda submodule targets NVIDIA GPUs, compiling functions decorated with @cuda.jit into PTX device assembly so genuine CUDA kernels can be authored, launched, and debugged entirely from Python.
@cuda.jit: Python In, PTX Out
Numba is a just-in-time compiler. When you decorate a function with @cuda.jit, Numba does not run it as Python — the first time you launch it, Numba compiles that function to PTX (NVIDIA's GPU assembly) and runs it on the device. Inside the kernel you get the real CUDA vocabulary: cuda.threadIdx.x, cuda.blockIdx.x, cuda.shared.array, cuda.syncthreads(), cuda.atomic.add.
The trade-offs versus CUDA C++ are worth knowing. You give up some low-level control and the absolute peak of hand-tuned C++, but you gain the ability to write, launch, and debug a kernel without leaving Python or installing nvcc. The programming model is identical, so everything you learn here transfers directly to CUDA C++.
A kernel is a function executed times in parallel, once per thread. If launching the kernel costs a fixed overhead and each thread does work , the wall-clock time is roughly where is the number of threads running concurrently. Numba's one-time compile adds to on the first call only; amortized over many launches it disappears.
What actually happens the first time you call a @cuda.jit kernel, and why is the second call faster?
Related Problems on PixelBank
A CPU attacks a large array by looping over it, visiting one element after another on a handful of cores. A GPU attacks the same array by launching an army of lightweight threads that each handle a slice of the work simultaneously. To make that army manageable, CUDA organizes it as a two-level hierarchy, and this topic is where you learn to configure it. Threads are grouped into blocks, and blocks are grouped into a grid, and you choose the shape of both on every single kernel launch. That choice is not cosmetic: threads within one block run on the same Streaming Multiprocessor and can cooperate through fast shared memory and barriers, while threads in different blocks cannot synchronize during a launch at all. Getting a feel for the hierarchy now — how many threads a block holds, how many blocks cover your data, which built-ins report a thread's coordinates — is what makes the next topics tractable. Everything downstream, from computing a global index to tiling a matrix multiply, is ultimately a decision about how to map data onto this grid of blocks of threads.
Definition
A CUDA launch is structured as a grid of blocks of threads. A thread runs one instance of the kernel; a block groups threads that share an SM, fast shared memory, and barrier synchronization; the grid is the full set of blocks. The launch syntax kernel[blocks_per_grid, threads_per_block](args) fixes both dimensions.
The Launch Hierarchy
When you launch a kernel you write kernel[blocks_per_grid, threads_per_block](args). That single line configures the whole hierarchy:
- A thread is the smallest unit; it runs one copy of your kernel body.
- A block is a group of threads (commonly 128 or 256) that run on the same Streaming Multiprocessor and can cooperate through shared memory and barriers.
- The grid is all the blocks for this launch.
Inside the kernel, four built-ins tell a thread where it sits: cuda.threadIdx.x (its index within its block), cuda.blockIdx.x (its block's index within the grid), cuda.blockDim.x (threads per block), and cuda.gridDim.x (blocks per grid). Threads in the same block share fast on-chip memory; threads in different blocks cannot synchronize with each other during a launch.
Choosing threads_per_block = T and needing to cover n elements, the number of blocks is . Writing this with integer math, blocks = (n + T - 1) // T, guarantees at least n threads are launched (usually a few extra, which the bounds check handles).
You launch kernel[4, 128](...). How many threads run in total, and what is the range of cuda.blockIdx.x and cuda.threadIdx.x?
Here is the puzzle at the center of almost every CUDA kernel: you launch thousands of threads, every one of them runs the exact same function body, yet each must operate on a different piece of data. Nothing in the source code distinguishes one thread from another — so how does thread number 5,000 know it should touch element 5,000 and not trample element 0? The answer is the global index, and computing it is the single most common move in GPU programming. A thread knows two local facts: which block it belongs to and its position inside that block. Multiplying the block index by the block size skips past every earlier block, and adding the in-block position lands on a unique global slot. This topic teaches that arithmetic and the cuda.grid(1) helper that packages it. Once a thread can name its own element, a one-thread-per-element kernel becomes trivial — read input i, compute, write output i — and, crucially, there is no loop inside the kernel. The grid itself is the loop, unrolled across the hardware, which is the mental inversion that makes GPUs fast.
Definition
The global thread index is the unique position a thread occupies across the entire grid, computed as blockIdx.x * blockDim.x + threadIdx.x. It maps each thread to one element of a flat array. Numba's cuda.grid(1) returns exactly this value, and cuda.grid(2) returns a 2D coordinate pair.
From Block Coordinates to a Global Index
To map a thread to a unique element of a 1D array, combine its block index and thread index: multiply the block index by the block size to skip past all earlier blocks, then add the position within the block. In Numba you can write this by hand as i = cuda.blockIdx.x * cuda.blockDim.x + cuda.threadIdx.x, or use the convenience helper i = cuda.grid(1), which returns exactly that value.
With the global index in hand, a one-thread-per-element kernel is trivial: thread i reads input i, computes, and writes output i. There is no loop in the kernel — the grid is the loop.
cuda.grid(1) is pure index arithmetic: it partitions into contiguous chunks of size blockDim.x, one chunk per block. Block owns indices , and within it thread owns — a bijection between (block, thread) pairs and global positions.
Write the body of a kernel that doubles every element of an array x into out, using the global index.
Related Problems on PixelBank
The first kernel you write will probably work on a friendly input size and then corrupt memory on a real one, and this topic is the vaccine. Two small habits separate a toy kernel from a correct one. The first is the bounds check. Because the number of blocks is always rounded up to cover the data, the final block nearly always contains threads whose global index runs off the end of the array — threads with no valid element to touch. Let them write anyway and they scribble past the buffer, an out-of-bounds write that may silently overwrite whatever sits next in memory. A single if i < n: guard masks them off. The second habit is the grid-stride loop, which lets a fixed, modest number of threads process an array of any size by having each thread step forward by the total thread count until it falls off the end. Together these decouple correctness from the launch configuration: the same kernel runs safely whether the input has one element, a thousand, or ten million, and whether you launched a handful of threads or a flood of them.
Definition
A bounds check is the guard if i < n: that stops threads in the ragged final block from accessing memory past the array's end. A grid-stride loop has each thread iterate from its global index in steps of the total thread count, so any launch size correctly and disjointly covers any data size.
In this topic
- 1The Bounds Check
- 2The Grid-Stride Loop
The Bounds Check
Because you launch blocks, the last block almost always contains threads whose global index is — there is no element for them to touch. If those threads write anyway, they corrupt memory past the end of the array (an out-of-bounds write). The fix is one line: if i < n: before any access. Threads past the end simply do nothing.
This is why a non-multiple size like n = 1000 with blockDim = 256 needs the guard: the fourth block runs threads 768–1023, but only 768–999 are valid.
Make the doubling kernel safe for any n.
The Grid-Stride Loop
Sometimes you want to launch a fixed, modest number of threads and have each one process several elements — for example to reuse a tuned launch configuration across very different input sizes. The grid-stride loop does this: each thread starts at its global index and steps forward by the total number of threads (cuda.gridDim.x * cuda.blockDim.x) until it runs off the end.
This decouples the launch size from the data size. The same kernel correctly processes an array smaller than the thread count (most threads do one or zero iterations) or far larger (each thread loops many times).
If you launch total threads over elements, each thread runs the loop body about times. The stride guarantees the sets of indices touched by different threads are disjoint and together cover — a clean partition with no overlaps and no gaps.
Rewrite doubling as a grid-stride loop so it works even when there are fewer threads than elements.
Theory Exercise
Problem:
A teammate launches a one-thread-per-element kernel with kernel[(n + 255) // 256, 256](x, out) and no if guard inside. It works for n = 4096 but intermittently corrupts memory for n = 5000. Explain why, and give two different ways to fix it.
Hints:
- Is 5000 a multiple of 256?
- How many threads does the launch actually start, and how many elements exist?
- What do the extra threads in the last block do without a guard?