Robust Grid Stride Loops

Welcome back! So far, we’ve used a "one-thread-per-element" approach, which assumes our grid is always large enough to cover the data. But what happens if the dataset is massive or the GPU hardware limits our grid size?

To build professional, "battle-ready" CUDA code, we use the Grid Stride Loop. This pattern allows threads to process multiple elements by "stepping" through the array. By the end of this lesson, your kernels will be hardware-agnostic—capable of handling any data size correctly, regardless of how many threads you launch. Let’s upgrade our indexing for maximum robustness!

Now we add the key idea of this lesson. We compute the same starting index, but then we loop forward by a stride of gridDim.x * blockDim.x.

__global__ void vectorAdd(const float* a, const float* b, float* c, int n) {
    int i = blockDim.x * blockIdx.x + threadIdx.x;
    int stride = blockDim.x * gridDim.x;

    for (; i < n; i += stride) {
        c[i] = a[i] + b[i];
    }
}

Why this is more robust:

• Even if the grid has fewer threads than n, threads “come back around” to handle later elements.
• The loop condition i < n becomes our safety check, so every write is still in bounds.
• We can tune blocksPerGrid for performance later without breaking correctness.

Next, the program creates a large problem size and prepares host and device memory. Notice how we compute bytes once, since CUDA allocation and copies use byte counts.

int main() {
    const int N = 50000;
    size_t bytes = N * sizeof(float);

    std::vector<float> h_a(N, 1.0f), h_b(N, 2.0f), h_c(N, 0.0f);

    float *d_a, *d_b, *d_c;
    CUDA_CHECK(cudaMalloc(&d_a, bytes));
    CUDA_CHECK(cudaMalloc(&d_b, bytes));
    CUDA_CHECK(cudaMalloc(&d_c, bytes));

Here we allocate:

• h_a, h_b as inputs, filled with constants.
• h_c as output, initialized to zero.
• d_a, d_b, d_c as device arrays with matching sizes.

Now we launch the kernel, check for launch errors, and synchronize. The same launch works whether we use the baseline kernel or the grid stride version.

    vectorAdd<<<blocksPerGrid, threadsPerBlock>>>(d_a, d_b, d_c, N);
    CUDA_CHECK(cudaGetLastError());
    CUDA_CHECK(cudaDeviceSynchronize());

Two important checks happen here:

• cudaGetLastError() catches configuration and launch issues.
• cudaDeviceSynchronize() ensures the kernel is finished before we copy results back.

Finally, we copy d_c back, verify every element is 3.0f, print success or failure, and free device memory.

    CUDA_CHECK(cudaMemcpy(h_c.data(), d_c, bytes, cudaMemcpyDeviceToHost));

    bool success = true;
    for (int i = 0; i < N; i++) {
        if (std::fabs(h_c[i] - 3.0f) > 1e-5f) {
            success = false;
            break;
        }
    }

    std::cout << "Final Result: " << (success ? "SUCCESS" : "FAILURE") << std::endl;

    CUDA_CHECK(cudaFree(d_a));
    CUDA_CHECK(cudaFree(d_b));
    CUDA_CHECK(cudaFree(d_c));
    return success ? 0 : 1;
}

This verification step performs real work for us: it confirms that our indexing logic and our launch setup actually produced correct results for all 50,000 elements.

When everything is configured correctly, the program prints the launch size and the final verification result.

Launching 50000 operations using 196 blocks.
Final Result: SUCCESS

This tells us two things: we computed a grid large enough to cover N, and the CPU-side check confirmed that every output element matched the expected sum.

We now have a stronger mental model for 1D CUDA work: global index gives each thread a starting point, and a grid stride loop lets that thread safely handle more elements by stepping forward with a fixed stride. This pattern makes kernels correct across many sizes, and it keeps working even when we change the grid size for tuning.

Next, we will move into practice problems where we apply grid stride loops and launch math ourselves, and you will see how quickly this becomes a dependable habit.