changes to cuda notebooks

Author: HrisShterev

cuda/04_Speed-comparison.ipynb cuda/03_Speed-comparison.ipynbcuda/04_Speed-comparison.ipynb renamed to cuda/03_Speed-comparison.ipynb

@@ -0,0 +1,139 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "b",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "#undef __noinline__"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "a",
+   "metadata": {},
+   "source": [
+    "# Thread Cooperation - The Grid-Stride Loop\n",
+    "\n",
+    "In the previous notebooks each GPU thread handled exactly one array element.\n",
+    "That only works when you launch exactly as many threads as you have data points.\n",
+    "\n",
+    "A more flexible pattern is the **grid-stride loop**: launch a fixed number of threads\n",
+    "(typically chosen to saturate the GPU) and have each thread walk through the array\n",
+    "in steps equal to the total number of threads in the grid.\n",
+    "\n",
+    "This means:\n",
+    "- You never need to know N at launch time\n",
+    "- You can reuse the same launch config for any array size\n",
+    "- Each thread does a balanced share of the work"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "ec24e43e",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "#include <cstdio>\n",
+    "#include <cmath>\n",
+    "\n",
+    "static const int N   = 1 << 17;   // 131 072 elements\n",
+    "static const float ALPHA = 2.5f;     // scalar for SAXPY\n",
+    "\n",
+    "__global__ void saxpy(float alpha, const float *x, float *y, int n) {\n",
+    "    int start  = threadIdx.x + blockIdx.x * blockDim.x;\n",
+    "    int step   = blockDim.x  * gridDim.x;          // total threads in the grid\n",
+    "\n",
+    "    for (int i = start; i < n; i += step)\n",
+    "        y[i] = alpha * x[i] + y[i];\n",
+    "}"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "9d46b9e7",
+   "metadata": {},
+   "source": [
+    "We launch with 64 blocks × 256 threads = 16 384 threads, but the array has 131 072 elements.\n",
+    "Each thread therefore handles roughly 8 elements via the loop. The result is verified\n",
+    "element-by-element on the CPU."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "id": "d080db69",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Verified - 131072 elements correct, max error = 0.00e+00\n"
+     ]
+    }
+   ],
+   "source": [
+    "// Allocate, fill, launch, verifi\n",
+    "float *h_x = (float*)malloc(N * sizeof(float));\n",
+    "float *h_y = (float*)malloc(N * sizeof(float));\n",
+    "float *h_ref = (float*)malloc(N * sizeof(float));\n",
+    "\n",
+    "for (int i = 0; i < N; i++) {\n",
+    "    h_x[i]   = sinf((float)i);\n",
+    "    h_y[i]   = cosf((float)i);\n",
+    "    h_ref[i] = ALPHA * h_x[i] + h_y[i];   // CPU reference\n",
+    "}\n",
+    "\n",
+    "float *d_x, *d_y;\n",
+    "cudaMalloc(&d_x, N * sizeof(float));\n",
+    "cudaMalloc(&d_y, N * sizeof(float));\n",
+    "cudaMemcpy(d_x, h_x, N * sizeof(float), cudaMemcpyHostToDevice);\n",
+    "cudaMemcpy(d_y, h_y, N * sizeof(float), cudaMemcpyHostToDevice);\n",
+    "\n",
+    "// 256 threads/block, 64 blocks — 16 384 total threads for 131 072 elements\n",
+    "saxpy<<<64, 256>>>(ALPHA, d_x, d_y, N);\n",
+    "\n",
+    "cudaMemcpy(h_y, d_y, N * sizeof(float), cudaMemcpyDeviceToHost);\n",
+    "\n",
+    "// Verify (allow small floating-point tolerance)\n",
+    "int errors = 0;\n",
+    "float max_err = 0.f;\n",
+    "for (int i = 0; i < N; i++) {\n",
+    "    float err = fabsf(h_y[i] - h_ref[i]);\n",
+    "    if (err > 1e-4f) errors++;\n",
+    "    if (err > max_err) max_err = err;\n",
+    "}\n",
+    "\n",
+    "if (errors == 0)\n",
+    "    printf(\"Verified - %d elements correct, max error = %.2e\\n\", N, max_err);\n",
+    "else\n",
+    "    printf(\"FAIL - %d mismatches\\n\", errors);\n",
+    "\n",
+    "cudaFree(d_x);  cudaFree(d_y);\n",
+    "free(h_x);  free(h_y);  free(h_ref);"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "C++23 CUDA",
+   "language": "cpp",
+   "name": "xcpp23-cuda"
+  },
+  "language_info": {
+   "codemirror_mode": "text/x-c++src",
+   "file_extension": ".cpp",
+   "mimetype": "text/x-c++src",
+   "name": "CUDA",
+   "nbconvert_exporter": "",
+   "pygments_lexer": "",
+   "version": "cxx23"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

cuda/03_The_Julia_Set.ipynb

@@ -74,13 +74,13 @@
    "source": [
     "__global__ void dot_product(float *a, float *b, float *partial) {\n",
     "    __shared__ float tile[THREADS];   // one float per thread, on-chip\n",
+    "    float running = 0.0f;\n",
     "\n",
     "    int tid    = threadIdx.x + blockIdx.x * blockDim.x;\n",
     "    int stride = blockDim.x * gridDim.x;\n",
     "    int lane   = threadIdx.x;\n",
     "\n",
     "    // Phase 1: each thread accumulates its share of element-wise products\n",
-    "    float running = 0.0f;\n",
     "    for (int i = tid; i < N; i += stride)\n",
     "        running += a[i] * b[i];\n",
     "    tile[lane] = running;\n",