1.4 Your first useful kernel: vector add

A complete vector-add example that walks the full CUDA workflow: memory allocation, data transfer, writing the kernel, error checking, and timing, then comparing CPU vs GPU.

Contents

• 1. From CPU to GPU: a full development flow
• 2. Implementing the vector-add kernel
• 3. Complete host-side code
• 4. CUDA error checking
• 5. Timing: CPU vs GPU
• 6. Full runnable code
• 7. Performance analysis and discussion
• 8. Advanced exercises
• Summary

---

1. From CPU to GPU: a full development flow

Moving a compute job from the CPU to the GPU is like handing work from one capable employee to a team of a thousand: you prepare the materials, assign the work, wait for it to finish, and collect the results. In CUDA that maps to:

graph LR
    A["1. Allocate GPU memory"] --> B["2. Copy CPU→GPU"]
    B --> C["3. Launch kernel"]
    C --> D["4. Copy GPU→CPU"]
    D --> E["5. Free GPU memory"]

Each step has a matching CUDA API:

---

2. Implementing the vector-add kernel

Vector add is the simplest parallel problem: every element is independent, so it maps cleanly to a GPU. Given two length-N vectors A and B, compute:

$$C[i] = A[i] + B[i], \quad i = 0, 1, \ldots, N-1$$

2.1 Kernel function

__global__ void vectorAdd(const float* A, const float* B, float* C, int N) {
    // Compute the global thread index
    int idx = blockIdx.x * blockDim.x + threadIdx.x;

    // Bounds check: avoid out-of-range access
    if (idx < N) {
        C[idx] = A[idx] + B[idx];
    }
}

Code notes:

• const float*: A and B are read-only inputs
• int idx: each thread computes a unique global index
• if (idx < N): the last block may have spare threads; you must check bounds
• C[idx] = A[idx] + B[idx]: each thread independently adds one element

2.2 Why the bounds check is needed

Suppose N = 1000 and block size = 256:

• Blocks needed = $\lceil 1000/256 \rceil = 4$
• Total threads = 4 × 256 = 1024
• The extra 24 threads have indices 1000~1023, past the end of the array

Without if (idx < N), those threads touch illegal addresses and you get undefined behavior or a crash.

---

3. Complete host-side code

3.1 Memory management

#include <cstdio>
#include <cstdlib>
#include <cuda_runtime.h>

int main() {
    int N = 1 << 20;  // about 1 million elements
    size_t size = N * sizeof(float);

    // 1. Allocate host memory
    float* h_A = (float*)malloc(size);
    float* h_B = (float*)malloc(size);
    float* h_C = (float*)malloc(size);

    // Initialize input data
    for (int i = 0; i < N; i++) {
        h_A[i] = (float)i;
        h_B[i] = (float)(i * 2);
    }

    // 2. Allocate device memory
    float *d_A, *d_B, *d_C;
    cudaMalloc(&d_A, size);
    cudaMalloc(&d_B, size);
    cudaMalloc(&d_C, size);

    // 3. Copy inputs from host to device
    cudaMemcpy(d_A, h_A, size, cudaMemcpyHostToDevice);
    cudaMemcpy(d_B, h_B, size, cudaMemcpyHostToDevice);

    // 4. Configure and launch the kernel
    int blockSize = 256;
    int gridSize = (N + blockSize - 1) / blockSize;
    vectorAdd<<<gridSize, blockSize>>>(d_A, d_B, d_C, N);

    // 5. Copy results from device back to host
    cudaMemcpy(h_C, d_C, size, cudaMemcpyDeviceToHost);

    // 6. Verify results
    for (int i = 0; i < N; i++) {
        if (fabsf(h_C[i] - (h_A[i] + h_B[i])) > 1e-5f) {
            printf("Verification FAILED at index %d!\n", i);
            break;
        }
    }
    printf("Verification PASSED!\n");

    // 7. Free memory
    cudaFree(d_A);
    cudaFree(d_B);
    cudaFree(d_C);
    free(h_A);
    free(h_B);
    free(h_C);

    return 0;
}

3.2 Build and run

$ nvcc vectorAdd.cu -o vectorAdd
$ ./vectorAdd
Verification PASSED!

---

4. CUDA error checking

CUDA API calls can fail (out of memory, bad arguments, and so on). Production code should check the return value of every call.

4.1 Macro

#define CUDA_CHECK(call)                                                      \
    do {                                                                       \
        cudaError_t err = call;                                                \
        if (err != cudaSuccess) {                                              \
            fprintf(stderr, "CUDA error at %s:%d - %s\n",                     \
                    __FILE__, __LINE__, cudaGetErrorString(err));              \
            exit(EXIT_FAILURE);                                                \
        }                                                                      \
    } while (0)

4.2 How to use it

CUDA_CHECK(cudaMalloc(&d_A, size));
CUDA_CHECK(cudaMemcpy(d_A, h_A, size, cudaMemcpyHostToDevice));

// Kernel launch needs a two-step check
vectorAdd<<<gridSize, blockSize>>>(d_A, d_B, d_C, N);
CUDA_CHECK(cudaGetLastError());          // check launch configuration
CUDA_CHECK(cudaDeviceSynchronize());     // check errors during execution

4.3 Common error codes

Tip: During development, put a check after every CUDA API call. For release builds you can gate the checks with a compile-time macro so you drop the overhead.

---

5. Timing: CPU vs GPU

5.1 CPU baseline

void vectorAddCPU(const float* A, const float* B, float* C, int N) {
    for (int i = 0; i < N; i++) {
        C[i] = A[i] + B[i];
    }
}

5.2 Timing with CUDA events

CUDA events are the standard way to time work on the GPU, with microsecond-level precision:

cudaEvent_t start, stop;
cudaEventCreate(&start);
cudaEventCreate(&stop);

// Start timing
cudaEventRecord(start);

// Work to time
vectorAdd<<<gridSize, blockSize>>>(d_A, d_B, d_C, N);

// Stop timing
cudaEventRecord(stop);
cudaEventSynchronize(stop);

// Elapsed time
float milliseconds = 0;
cudaEventElapsedTime(&milliseconds, start, stop);
printf("GPU Kernel time: %.4f ms\n", milliseconds);

// Cleanup
cudaEventDestroy(start);
cudaEventDestroy(stop);

5.3 End-to-end timing including transfers

cudaEventRecord(start);

cudaMemcpy(d_A, h_A, size, cudaMemcpyHostToDevice);
cudaMemcpy(d_B, h_B, size, cudaMemcpyHostToDevice);
vectorAdd<<<gridSize, blockSize>>>(d_A, d_B, d_C, N);
cudaMemcpy(h_C, d_C, size, cudaMemcpyDeviceToHost);

cudaEventRecord(stop);
cudaEventSynchronize(stop);

float totalTime = 0;
cudaEventElapsedTime(&totalTime, start, stop);
printf("GPU total (with memcpy): %.4f ms\n", totalTime);

5.4 Sample performance numbers

On an RTX 4090, with N = $2^{24}$ (about 16 million elements), typical results look like:

Note: Vector add is a memory-bound op (one add per element, but you read two floats and write one). Once you include transfers, PCIe bandwidth becomes the bottleneck. The GPU wins most clearly when the data is already in device memory.

---

6. Full runnable code

Putting the pieces above into one program you can compile and run:

#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <chrono>
#include <cuda_runtime.h>

// Error-check macro
#define CUDA_CHECK(call)                                                      \
    do {                                                                       \
        cudaError_t err = call;                                                \
        if (err != cudaSuccess) {                                              \
            fprintf(stderr, "CUDA error at %s:%d - %s\n",                     \
                    __FILE__, __LINE__, cudaGetErrorString(err));              \
            exit(EXIT_FAILURE);                                                \
        }                                                                      \
    } while (0)

// GPU kernel
__global__ void vectorAdd(const float* A, const float* B, float* C, int N) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < N) {
        C[idx] = A[idx] + B[idx];
    }
}

// CPU baseline
void vectorAddCPU(const float* A, const float* B, float* C, int N) {
    for (int i = 0; i < N; i++) {
        C[i] = A[i] + B[i];
    }
}

int main() {
    const int N = 1 << 24;  // about 16 million elements
    const size_t size = N * sizeof(float);
    printf("Vector size: %d elements (%.1f MB)\n", N, size / (1024.0 * 1024.0));

    // ========== Allocate host memory ==========
    float* h_A = (float*)malloc(size);
    float* h_B = (float*)malloc(size);
    float* h_C_cpu = (float*)malloc(size);
    float* h_C_gpu = (float*)malloc(size);

    // Initialize data
    for (int i = 0; i < N; i++) {
        h_A[i] = sinf(i * 0.001f);
        h_B[i] = cosf(i * 0.001f);
    }

    // ========== CPU compute ==========
    auto cpu_start = std::chrono::high_resolution_clock::now();
    vectorAddCPU(h_A, h_B, h_C_cpu, N);
    auto cpu_end = std::chrono::high_resolution_clock::now();
    float cpu_ms = std::chrono::duration<float, std::milli>(cpu_end - cpu_start).count();
    printf("CPU time: %.4f ms\n", cpu_ms);

    // ========== GPU compute ==========
    float *d_A, *d_B, *d_C;
    CUDA_CHECK(cudaMalloc(&d_A, size));
    CUDA_CHECK(cudaMalloc(&d_B, size));
    CUDA_CHECK(cudaMalloc(&d_C, size));

    // Create CUDA event timers
    cudaEvent_t start, stop, kernel_start, kernel_stop;
    cudaEventCreate(&start);
    cudaEventCreate(&stop);
    cudaEventCreate(&kernel_start);
    cudaEventCreate(&kernel_stop);

    // End-to-end timing (including transfers)
    cudaEventRecord(start);

    CUDA_CHECK(cudaMemcpy(d_A, h_A, size, cudaMemcpyHostToDevice));
    CUDA_CHECK(cudaMemcpy(d_B, h_B, size, cudaMemcpyHostToDevice));

    int blockSize = 256;
    int gridSize = (N + blockSize - 1) / blockSize;

    // Kernel-only timing
    cudaEventRecord(kernel_start);
    vectorAdd<<<gridSize, blockSize>>>(d_A, d_B, d_C, N);
    cudaEventRecord(kernel_stop);

    CUDA_CHECK(cudaGetLastError());
    CUDA_CHECK(cudaMemcpy(h_C_gpu, d_C, size, cudaMemcpyDeviceToHost));

    cudaEventRecord(stop);
    cudaEventSynchronize(stop);

    float gpu_total_ms = 0, kernel_ms = 0;
    cudaEventElapsedTime(&gpu_total_ms, start, stop);
    cudaEventElapsedTime(&kernel_ms, kernel_start, kernel_stop);

    printf("GPU kernel time: %.4f ms\n", kernel_ms);
    printf("GPU total time (with memcpy): %.4f ms\n", gpu_total_ms);
    printf("Speedup (kernel only): %.1fx\n", cpu_ms / kernel_ms);
    printf("Speedup (end-to-end): %.1fx\n", cpu_ms / gpu_total_ms);

    // ========== Verify results ==========
    int errors = 0;
    for (int i = 0; i < N; i++) {
        if (fabsf(h_C_gpu[i] - h_C_cpu[i]) > 1e-5f) {
            if (errors < 5) {
                printf("Mismatch at %d: CPU=%.6f, GPU=%.6f\n",
                       i, h_C_cpu[i], h_C_gpu[i]);
            }
            errors++;
        }
    }
    if (errors == 0) {
        printf("Verification PASSED!\n");
    } else {
        printf("Verification FAILED: %d mismatches\n", errors);
    }

    // ========== Effective bandwidth ==========
    // Vector add reads 2 floats and writes 1 = 12 bytes/element
    float bandwidth_gb = (3.0f * size) / (kernel_ms * 1e-3f) / 1e9f;
    printf("Effective bandwidth: %.1f GB/s\n", bandwidth_gb);

    // ========== Cleanup ==========
    cudaEventDestroy(start);
    cudaEventDestroy(stop);
    cudaEventDestroy(kernel_start);
    cudaEventDestroy(kernel_stop);
    CUDA_CHECK(cudaFree(d_A));
    CUDA_CHECK(cudaFree(d_B));
    CUDA_CHECK(cudaFree(d_C));
    free(h_A);
    free(h_B);
    free(h_C_cpu);
    free(h_C_gpu);

    return 0;
}

# Build
nvcc -O3 -arch=sm_89 vectorAdd.cu -o vectorAdd

# Run
./vectorAdd

---

7. Performance analysis and discussion

7.1 Arithmetic intensity of vector add

Arithmetic intensity = compute / bytes moved:

$$\text{AI} = \frac{\text{FLOPs}}{\text{Bytes}} = \frac{1 \text{ FLOP}}{12 \text{ Bytes}} \approx 0.083 \text{ FLOP/Byte}$$

That number is tiny, so vector add is a classic memory-bound kernel: performance is limited by memory bandwidth, not compute.

7.2 Bandwidth utilization

If the kernel takes 0.3 ms and moves $3 \times 2^{24} \times 4$ Bytes = 201 MB:

$$\text{Effective bandwidth} = \frac{201 \text{ MB}}{0.3 \text{ ms}} \approx 670 \text{ GB/s}$$

Against the RTX 4090’s theoretical peak of 1008 GB/s, that is about 66% utilization, which is already solid for such a simple kernel.

7.3 PCIe transfer cost

Key point: For simple elementwise ops, data movement is often the main bottleneck. In real apps, keep data on the GPU: chain several kernels over the same buffers instead of shuttling them back and forth.

7.4 When the GPU is worth it

---

8. Advanced exercises

After vector add, try these to lock in what you learned:

8.1 Exercise 1: Dot product

Implement the dot product $\text{result} = \sum_{i=0}^{N-1} A[i] \times B[i]$.

Hint: this needs a reduction. Each thread computes a partial product, then you accumulate inside the block with shared memory.

8.2 Exercise 2: SAXPY

Implement SAXPY: $Y[i] = \alpha \times X[i] + Y[i]$

This is a basic BLAS op; it adds one scalar multiply on top of vector add.

8.3 Exercise 3: Grid-stride loop version

Rewrite vector add with a grid-stride loop so a fixed number of threads can cover an array of any size:

__global__ void vectorAddStride(const float* A, const float* B,
                                float* C, int N) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    int stride = blockDim.x * gridDim.x;
    for (int i = idx; i < N; i += stride) {
        C[i] = A[i] + B[i];
    }
}

// At launch you no longer need gridSize = (N + 255) / 256
// Pick any reasonable grid size
vectorAddStride<<<128, 256>>>(d_A, d_B, d_C, N);

8.4 Exercise 4: Overlap async copy with compute

Use a CUDA stream so transfers and compute can run in parallel:

cudaStream_t stream;
cudaStreamCreate(&stream);

// Async copies and compute can overlap
cudaMemcpyAsync(d_A, h_A, size, cudaMemcpyHostToDevice, stream);
cudaMemcpyAsync(d_B, h_B, size, cudaMemcpyHostToDevice, stream);
vectorAdd<<<gridSize, blockSize, 0, stream>>>(d_A, d_B, d_C, N);
cudaMemcpyAsync(h_C, d_C, size, cudaMemcpyDeviceToHost, stream);

cudaStreamSynchronize(stream);
cudaStreamDestroy(stream);

Tip: Real overlap needs host memory allocated with cudaMallocHost() (pinned memory).

---

Summary

With this vector-add walkthrough you have the core CUDA workflow:

1. Five-step standard path: cudaMalloc → cudaMemcpy(H2D) → kernel launch → cudaMemcpy(D2H) → cudaFree
2. Three kernel essentials: global index, bounds check, parallel work
3. Error checking: wrap every API call in a CUDA_CHECK macro
4. Performance measurement: CUDA events plus effective bandwidth
5. Performance intuition: for simple ops the bottleneck is transfer; the GPU wins when data stays on device or the work is compute-heavy

Vector add is simple, but it builds the full CUDA scaffolding. Later kernels (matmul, convolution, reduction, and so on) add complexity on top of this same frame.

Self-check list

• Write a full CUDA vector-add program on your own (memory management and error checks included)
• Explain why the kernel needs if (idx < N)
• Time a kernel precisely with CUDA events
• Compute a kernel’s effective memory bandwidth and compare it to peak
• Tell whether an op is compute-bound or memory-bound
• Explain why end-to-end time with transfers is much larger than kernel-only time
• Handle arbitrary sizes with a grid-stride loop
• Write robust CUDA code with the CUDA_CHECK macro

References

• NVIDIA CUDA Samples - vectorAdd
• An Even Easier Introduction to CUDA (NVIDIA Blog)
• CUDA C++ Programming Guide - Performance Guidelines
• CUDA C++ Best Practices Guide - Bandwidth