mirror of
https://github.com/andreinechaev/nvcc4jupyter.git
synced 2026-06-14 03:00:47 +05:30
Cosmin Ștefan Ciocan
579 lines
29 KiB
Plaintext
579 lines
29 KiB
Plaintext
{
|
|
"cells": [
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"# Introduction\n",
|
|
"\n",
|
|
"This notebook is an adaptation of [this session](https://www.olcf.ornl.gov/calendar/introduction-to-cuda-c/) (presentation and assignments) from [the CUDA training series](https://www.olcf.ornl.gov/cuda-training-series/) provided to the Oak Ridge National Laboratory by [Bob Crovella](https://developer.nvidia.com/blog/author/bob-crovella/), who is on the Solution Architecture team at NVIDIA. While not meant as a replacement to the course, this notebook goes over the main points and acts as a way to quickly put them into practice. Remember that the theoretical part may be a bit overwhelming at first, so make sure you run the examples and play around with them to get a better understanding of the concepts."
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"# Setup\n",
|
|
"- Make sure you run the notebook on either **Kaggle**, **Colab** or anywhere there is an NVIDIA GPU readily available.\n",
|
|
"- We need to install and load the [nvcc4jupyter](https://github.com/andreinechaev/nvcc4jupyter) extension, which enables us to run CUDA C++ code directly from the notebook cells"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"!pip install nvcc4jupyter"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%load_ext nvcc4jupyter"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"# Introduction to CUDA C++\n",
|
|
"\n",
|
|
"This course covers what CUDA is, the CUDA programming model, CUDA syntax, the basics of error handling and running simple CUDA programs with little care for performance, as that will be covered in the next courses. It contains:\n",
|
|
"- The theoretical basics of CUDA\n",
|
|
"- Learning CUDA syntax\n",
|
|
"- Vector addition example\n",
|
|
"- Writing and launching CUDA C++ kernels\n",
|
|
"- Managing GPU Memory\n",
|
|
"- Matrix multiplication assignment"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## CUDA basics\n",
|
|
"\n",
|
|
"CUDA is a **parallel computing platform and API** created by NVIDIA which allows developers to use NVIDIA cards for **accelerating repetitive computational tasks**. The most well known applications of this are **video games**, where each pixel is processed in parallel with many others, and in running **artificial neural networks**. Both applications use lots of matrix multiplications and other mathematical operations which can be run in parallel. We will call this accelerator card **\"the device\"** (also known as GPU).\n",
|
|
"\n",
|
|
"The device is a separate hardware component from the CPU (which we will call **\"the host\"**), but they are similar in many ways. Both have memory and compute resources to manage, which is the task of their operating system. Both have the concept of **cores**, which are \"independent\" processors which execute instructions. We will see that device cores are less independent than host cores, and they also are less complex, meaning you can fit a lot more of them in the same space. Having more cores means the GPU is capable of processing the same workload a lot faster if it can be parallelized. The host operating system is much more complex as well, having to manage the interactions with the network, disk, etc.\n",
|
|
"\n",
|
|
"<div style=\"text-align:center\"><img src=\"https://nvcc4jupyter.s3.eu-central-1.amazonaws.com/notebooks/introduction-to-cuda-cpp/host_device_communication.png\"/></div>\n",
|
|
"\n",
|
|
"The figure above shows the basic interactions between device and host. The **host is what decides what needs to be done (the kernel code) and on what data**, and the **device will execute that code on all of the data with as much parallelism as it can**. This interaction makes sense only if it is faster to send the data to the device, process it with its numerous cores, and have it sent back to the host, than to directly process it on the host. This is a matter of **exposing enough parallelism** in your problem to justify sending the data. In this session we will cover only [embarrassingly parallel](https://en.wikipedia.org/wiki/Embarrassingly_parallel) problems such as **vector addition** and **matrix multiplication** which should easily benefit from device acceleration.\n",
|
|
"\n",
|
|
"<div style=\"text-align:center\"><img src=\"https://nvcc4jupyter.s3.eu-central-1.amazonaws.com/notebooks/introduction-to-cuda-cpp/porting_to_cuda.png\"/></div>\n",
|
|
"\n",
|
|
"In most applications, there will be sequences of code that **cannot be parallelized** or **require I/O** (as shown in blue in the figure above) and which must be run on the much fewer in number, but faster and smarter cores of the host (CPU). The developer has to identify those functions which are **computationally intensive and parallelizable** (shown in green).\n",
|
|
"\n",
|
|
"## Cuda programming model\n",
|
|
"\n",
|
|
"### 3-step processing flow\n",
|
|
"\n",
|
|
"Using a GPU to accelerate a computation can be reduced into **3 main steps** (with unified memory the story is a bit more complex, but that is a topic of discussion for future sessions):\n",
|
|
"\n",
|
|
"1. Copy input data from host to device.\n",
|
|
"2. Send kernel code from host to device and process the data from step 1.\n",
|
|
"3. Copy output data from device to host.\n",
|
|
"\n",
|
|
"Communication is done through a high performance bus such as [PCIe](https://en.wikipedia.org/wiki/PCI_Express) or [NVLink](https://en.wikipedia.org/wiki/NVLink), who ensure data transfer at high bandwidth and as low latency as possible.\n",
|
|
"\n",
|
|
"It is important to keep in mind the fact that host and device memory are separate entities, so a pointer to host memory will not be usable in device code and vice versa. \n",
|
|
"\n",
|
|
"### Host function definition\n",
|
|
"\n",
|
|
"CUDA code is compiled by [nvcc](https://docs.nvidia.com/cuda/cuda-compiler-driver-nvcc/index.html), which is a **compiler driver** that **calls multiple compilers under the hood** to perform various compilation stages. It will separate host and device functions, host functions being compiled by standard host compilers such as [gcc](https://gcc.gnu.org/), and device functions by the NVIDIA compiler. Now we will take a look at how to write syntactically correct CUDA C++ code. The first step is writing a simple function to run on the device:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"/* \n",
|
|
" * the \"__global__\" keyword signals to \"nvcc\" that this code must be compiled to\n",
|
|
" * run on the device, not the host; other than that, this is very much the same\n",
|
|
" * way you would define a function in C/C++\n",
|
|
" */\n",
|
|
"__global__ void my_kernel(void) {\n",
|
|
" // FUNCTION CODE\n",
|
|
"}\n",
|
|
"```\n",
|
|
"\n",
|
|
"Now for a more complex example:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"/*\n",
|
|
" * this is the signature of the vector addition function we will encounter in \n",
|
|
" * the examples; we're passing the pointers for the two input vectors, A and \n",
|
|
" * B, one for the output vector C, and the length of those vectors so we know\n",
|
|
" * where to stop;\n",
|
|
" */\n",
|
|
"__global__ void vector_add_device(const float *A, const float *B, float *C, int vector_size){\n",
|
|
" // FUNCTION CODE\n",
|
|
"}\n",
|
|
"```\n",
|
|
"\n",
|
|
"### Host function calls\n",
|
|
"\n",
|
|
"Calling a device function is a topic on its own. Compared to C/C++, device functions are called with the **kernel launch** syntax:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"vector_add_device<<1, 1>>(A, B, C, vector_size);\n",
|
|
"```\n",
|
|
"\n",
|
|
"For reference, here is how it would look for C/C++:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"vector_add_host(A, B, C, vector_size);\n",
|
|
"```\n",
|
|
"\n",
|
|
"The kernel launch configuration is one of the ways we tell the device how to parallelize our function when we want to use multiple workers at the same time. Just as a quick example to get some intuition, this is how you would tell it to run **N workers in parallel** to perform vector addition:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"vector_add_device<<N, 1>>(A, B, C, vector_size);\n",
|
|
"```\n",
|
|
"\n",
|
|
"### Block and thread indices\n",
|
|
"\n",
|
|
"Up until now we have only looked at function signatures, but we have not actually completed the `vector_add_device` function with its code. Let's take a look without worrying about understanding everything from the start:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"__global__ void vector_add_device(const float *A, const float *B, float *C, int vector_size){\n",
|
|
" int idx = threadIdx.x + blockIdx.x * blockDim.x;\n",
|
|
" if (idx < vector_size) {\n",
|
|
" C[idx] = A[idx] + B[idx];\n",
|
|
" }\n",
|
|
"}\n",
|
|
"```\n",
|
|
"\n",
|
|
"It is time to have a brief explanation on what **threads** and **blocks** are. A thread is the **smallest unit of work** in CUDA. Take a quick look at the code above and you will notice there is no for-loop to process all elements of the two input vectors. The function only adds and saves the result of the elements at the index `idx`. What you see here is the view of one of the many threads that exist in a device function call. Each thread computes an index in the array, adds the elements at that index and saves the result in the output array at that index. This is an essential part of writing device functions: figuring out how to divide work among threads by giving each of them an index to process. Of course, we are simplifying a bit.\n",
|
|
"\n",
|
|
"Each thread is assigned into a block of multiple threads (at most 1024 threads per block). This is where understanding the hardware starts being important and we will cover the details in the following courses. For now, you only need to know that threads that run in the same block can more easily communicate and synchronize with one another because they are executed on the same Streaming Multiprocessor.\n",
|
|
"\n",
|
|
"<div style=\"text-align:center\"><img src=\"https://nvcc4jupyter.s3.eu-central-1.amazonaws.com/notebooks/introduction-to-cuda-cpp/kernel_execution_on_gpu.png\"/></div>\n",
|
|
"\n",
|
|
"Now take a look at this line:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"int idx = threadIdx.x + blockIdx.x * blockDim.x;\n",
|
|
"```\n",
|
|
"\n",
|
|
"- `threadIdx`: Is a built-in structure (it is not defined anywhere in your code, but it is provided at runtime) that contains the index of the thread inside its own block. This index is not unique by itself as there are threads in other blocks that will share the same index.\n",
|
|
"- `blockIdx`: Is a built-in structure that contains the index of the block that contains this thread. Is constant for all threads in the same block.\n",
|
|
"- `blockDim`: Is a built-in structure that contains the number of threads per block. Is constant for all threads, no matter the block they belong to. When multiplied with `blockIdx`, it yields the index of the first thread in the block with index `blockIdx`. When also adding `threadIdx` you get the global index of the given thread. This example matches the figure above which represents a 1-dimensional grid of 1-dimensional blocks. Multiple dimensions are possible, but we only use one, hence accessing the first value of all 3 structures: `threadIdx.x`, `blockIdx.x`, and `blockDim.x`.\n",
|
|
"\n",
|
|
"\n",
|
|
"At this point we know enough to run our hello world example:"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%%cuda\n",
|
|
"#include <stdio.h>\n",
|
|
"\n",
|
|
"__global__ void hello(){\n",
|
|
" int global_idx = threadIdx.x + blockIdx.x * blockDim.x;\n",
|
|
" printf(\n",
|
|
" \"Hello from thread with global index %u (threadIdx.x: %u, blockIdx.x: %u, blockDim.x: %u)\\n\", \n",
|
|
" global_idx, threadIdx.x, blockIdx.x, blockDim.x\n",
|
|
" );\n",
|
|
"}\n",
|
|
"\n",
|
|
"int main() {\n",
|
|
" int n_blocks = 2;\n",
|
|
" int n_threads_per_block = 3;\n",
|
|
" hello<<<n_blocks, n_threads_per_block>>>();\n",
|
|
" \n",
|
|
" // wait for the execution of the asynchronous \"hello\" function call to finish\n",
|
|
" cudaDeviceSynchronize();\n",
|
|
"}"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"If you run the hello world example with `n_blocks = 2` and `n_threads_per_block = 3` then you will be able to see that there are 6 threads with global indices from 0 to 5. There are two threads with thread id (`threadIdx.x`) equal to 0, one of them being on the first position of the first block and the second thread being on the first position of the second block. The block id (`blockIdx.x`) is shared between threads from the same block (3 threads have block id 0 and 3 have id 1). The block size (`blockDim.x`) is equal to 3 for all threads since this is a global setting (the number of threads per block) that is set in the kernel launch configuration.\n",
|
|
"\n",
|
|
"<div style=\"text-align:center\"><img src=\"https://nvcc4jupyter.s3.eu-central-1.amazonaws.com/notebooks/introduction-to-cuda-cpp/block_thread_indexing.png\"/></div>\n",
|
|
"\n",
|
|
"If you change `n_blocks = 4` and `n_threads_per_block = 8` then you will see how the globally unique index is computed for the example in the figure above. \n",
|
|
"\n",
|
|
"### Memory management\n",
|
|
"\n",
|
|
"To be able to transition from our simple hello world program to vector addition we need to see how to handle memory. The most important thing to remember is that **device memory and host memory are two different entities**. This means you cannot use pointers for one memory space in the other. This also means you need to have **two copies of each input and output array**, one in host memory, and one in device memory. Usually, memory handling follows this flow: \n",
|
|
"\n",
|
|
"1. Allocate host memory (input and output arrays)\n",
|
|
" - `malloc()`\n",
|
|
"2. Allocate device memory (input and output arrays)\n",
|
|
" - `cudaMalloc()`\n",
|
|
"3. Fill host input arrays from any of a number of sources: files, network, etc.\n",
|
|
" - `memcpy()` or other means of populating the arrays\n",
|
|
"4. Copy host input arrays' content to device input arrays (this is done over the NVLink / PCIe bus)\n",
|
|
" - `cudaMemcpy()`\n",
|
|
"5. Run kernel code to turn inputs into outputs\n",
|
|
"6. Copy device output arrays' content to the host output arrays (NVLink / PCIe bus)\n",
|
|
" - `cudaMemcpy()`\n",
|
|
"7. Free device memory\n",
|
|
" - `cudaFree()`\n",
|
|
"8. Free host memory\n",
|
|
" - `free()`\n",
|
|
" \n",
|
|
"We have **equivalents** between **C/C++** and the **CUDA API**:\n",
|
|
"- Memory allocation: \n",
|
|
" - [malloc](https://en.cppreference.com/w/c/memory/malloc) / [cudaMalloc](https://docs.nvidia.com/cuda/cuda-runtime-api/group__CUDART__MEMORY.html)\n",
|
|
"- Copying data from one array to another: \n",
|
|
" - [memcpy](https://en.cppreference.com/w/c/memory/memcpy) / [cudaMemcpy](https://docs.nvidia.com/cuda/cuda-runtime-api/group__CUDART__MEMORY.html)\n",
|
|
"- Memory deallocation: \n",
|
|
" - [free](https://en.cppreference.com/w/c/memory/free) / [cudaFree](https://docs.nvidia.com/cuda/cuda-runtime-api/group__CUDART__MEMORY.html)\n",
|
|
" \n",
|
|
"Unlike **memcpy**, **cudaMemcpy** can copy and write data to/from either host or device memory. This means there are 4 types of data transfer, **host-to-host** (same as the regular memcpy), **host-to-device** (same as step 4 from above), **device-to-host** (same as step 6 from above), and **device-to-device**. The developer always needs to pay attention to what kind of pointers he provides to cudaMemcpy.\n",
|
|
"\n"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Useful functions\n",
|
|
"\n",
|
|
"We will define some functions which will be useful for both the vector addition example and the matrix multiplication assignment. You may skip reading those (especially those marked as assignment solutions). These functions will be available when adding `#include \"utils.h\"`."
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%%cuda_group_save --group shared --name \"utils.h\"\n",
|
|
"#include <math.h>\n",
|
|
"\n",
|
|
"// error checking macro\n",
|
|
"#define cudaCheckErrors(msg) \\\n",
|
|
" do { \\\n",
|
|
" cudaError_t __err = cudaGetLastError(); \\\n",
|
|
" if (__err != cudaSuccess) { \\\n",
|
|
" fprintf(stderr, \"Fatal error: %s (%s at %s:%d)\\n\", \\\n",
|
|
" msg, cudaGetErrorString(__err), \\\n",
|
|
" __FILE__, __LINE__); \\\n",
|
|
" fprintf(stderr, \"*** FAILED - ABORTING\\n\"); \\\n",
|
|
" exit(1); \\\n",
|
|
" } \\\n",
|
|
" } while (0)\n",
|
|
" \n",
|
|
"int divide_ceil(int dividend, int divisor) {\n",
|
|
" // equivalent to \"ceil(dividend / divisor)\" converted to int, but fast\n",
|
|
" return (dividend + divisor - 1) / divisor;\n",
|
|
"}\n",
|
|
"\n",
|
|
"bool almost_equal(float first, float second, float abs_tol = 0.001, float rel_tol = 0.001)\n",
|
|
"{\n",
|
|
" float diff = fabs(first - second);\n",
|
|
" if (diff <= abs_tol)\n",
|
|
" return true;\n",
|
|
"\n",
|
|
" first = fabs(first);\n",
|
|
" second = fabs(second);\n",
|
|
" float largest = (second > first) ? second : first;\n",
|
|
"\n",
|
|
" if (diff <= largest * rel_tol)\n",
|
|
" return true;\n",
|
|
" return false;\n",
|
|
"}\n",
|
|
"\n",
|
|
"void compare_arrays(float *x, float *y, int n, float abs_tol = 0.001, float rel_tol = 0.001) {\n",
|
|
" for (int i = 0; i < n; i++) {\n",
|
|
" if (!almost_equal(x[i], y[i], abs_tol, rel_tol)) {\n",
|
|
" printf(\"[ERROR] Arrays/matrices are not equal. At index %d expected \\\"%f\\\", got \\\"%f\\\"\\n\", i, x[i], y[i]);\n",
|
|
" return;\n",
|
|
" }\n",
|
|
" }\n",
|
|
" printf(\"[SUCCESS] Arrays/matrices are equal.\\n\");\n",
|
|
"}\n",
|
|
"\n",
|
|
"// this was obfuscated so you do not accidentally see the solution; it is also\n",
|
|
"// a purely C/C++ solution so not a 1:1 replacement for the assignment\n",
|
|
"void assignment_solution(float *x, float *y, float *z, int w) {\n",
|
|
" for (int i = 0; i < w; ++i) { for (int j = 0; j < w; ++j) { z[i*w+j] = 0; for (int k = 0; k < w; ++k) { z[i*w+j] += x[i*w+k] * y[k*w+j]; } } }\n",
|
|
"}"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Vector addition example\n",
|
|
"\n",
|
|
"Now we can finally implement the vector addition program which will also give us a view into how to manage memory in CUDA C++. Read the comments for detailed explanations:"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%%cuda\n",
|
|
"#include <stdio.h>\n",
|
|
"#include \"utils.h\"\n",
|
|
"\n",
|
|
"const int VECTOR_SIZE = 4096; // the size of the input and output vectors\n",
|
|
"const int BLOCK_SIZE = 256; // the number of threads per block; limited to 1024\n",
|
|
"\n",
|
|
"// vector add kernel: C[i] <- A[i] + B[i] for i in [0, 1, ..., VECTOR_SIZE - 1]\n",
|
|
"__global__ void vadd(const float *A, const float *B, float *C, int vector_size) {\n",
|
|
" // compute the globally unique index which tells the current thread \n",
|
|
" // what element it is supposed to be processing\n",
|
|
" int global_idx = threadIdx.x + blockIdx.x * blockDim.x;\n",
|
|
" \n",
|
|
" // if BLOCK_SIZE does not divide VECTOR_SIZE exactly, the last block\n",
|
|
" // will have threads whose global index will be out of bounds\n",
|
|
" if (global_idx < vector_size) {\n",
|
|
" C[global_idx] = A[global_idx] + B[global_idx];\n",
|
|
" }\n",
|
|
"}\n",
|
|
"\n",
|
|
"int main(){\n",
|
|
" // all vectors are float pointers, but some of them will contain an\n",
|
|
" // address in host space (h_A, h_B, h_C) and the others an address in\n",
|
|
" // device space (d_A, d_B, d_C)\n",
|
|
" float *h_A, *h_B, *h_C, *d_A, *d_B, *d_C;\n",
|
|
"\n",
|
|
" // allocate space for vectors in host memory\n",
|
|
" h_A = new float[VECTOR_SIZE];\n",
|
|
" h_B = new float[VECTOR_SIZE];\n",
|
|
" h_C = new float[VECTOR_SIZE];\n",
|
|
"\n",
|
|
" // assign random values to input vectors in host memory\n",
|
|
" for (int i = 0; i < VECTOR_SIZE; i++) {\n",
|
|
" h_A[i] = rand()/(float)RAND_MAX;\n",
|
|
" h_B[i] = rand()/(float)RAND_MAX;\n",
|
|
" }\n",
|
|
"\n",
|
|
" // allocate space for vectors in device memory\n",
|
|
" cudaMalloc(&d_A, VECTOR_SIZE*sizeof(float));\n",
|
|
" cudaMalloc(&d_B, VECTOR_SIZE*sizeof(float));\n",
|
|
" cudaMalloc(&d_C, VECTOR_SIZE*sizeof(float));\n",
|
|
" // check that the last CUDA API call was successful; if not, exit\n",
|
|
" cudaCheckErrors(\"cudaMalloc failure\"); \n",
|
|
"\n",
|
|
" // copy vectors A and B from host to device:\n",
|
|
" cudaMemcpy(d_A, h_A, VECTOR_SIZE*sizeof(float), cudaMemcpyHostToDevice);\n",
|
|
" cudaMemcpy(d_B, h_B, VECTOR_SIZE*sizeof(float), cudaMemcpyHostToDevice);\n",
|
|
" cudaCheckErrors(\"cudaMemcpy H2D failure\");\n",
|
|
"\n",
|
|
" // we need to have as many threads as there are elements in the vector;\n",
|
|
" // one option would be to call the function like this:\n",
|
|
" // vadd<<<VECTOR_SIZE, 1>>>(d_A, d_B, d_C, VECTOR_SIZE);\n",
|
|
" // however, we want to showcase the use of multiple threads per block,\n",
|
|
" // which is a feature that will benefit us in future problems for performance\n",
|
|
" // optimization; the number of threads is N_BLOCKS * BLOCK_SIZE, both of which\n",
|
|
" // are integers; this means we need to choose the number of blocks such that\n",
|
|
" // N_BLOCKS * BLOCK_SIZE >= VECTOR_SIZE; this is the reason for the if \n",
|
|
" // statement in the vector addition kernel function (we may have more threads\n",
|
|
" // than elements in the vector so the extra threads will do nothing useful)\n",
|
|
" int n_blocks = divide_ceil(VECTOR_SIZE, BLOCK_SIZE);\n",
|
|
"\n",
|
|
" // launch the vector adding kernel\n",
|
|
" vadd<<<n_blocks, BLOCK_SIZE>>>(d_A, d_B, d_C, VECTOR_SIZE);\n",
|
|
" cudaCheckErrors(\"kernel launch failure\");\n",
|
|
"\n",
|
|
" // wait for the kernel to finish execution\n",
|
|
" cudaDeviceSynchronize();\n",
|
|
" cudaCheckErrors(\"kernel execution failure\");\n",
|
|
"\n",
|
|
" // copy output array data from device memory to host memory\n",
|
|
" cudaMemcpy(h_C, d_C, VECTOR_SIZE*sizeof(float), cudaMemcpyDeviceToHost);\n",
|
|
" cudaCheckErrors(\"cudaMemcpy D2H failure\");\n",
|
|
" \n",
|
|
" // verify results\n",
|
|
" float *ground_truth = new float[VECTOR_SIZE];\n",
|
|
" for (int i = 0; i < VECTOR_SIZE; i++) {\n",
|
|
" ground_truth[i] = h_A[i] + h_B[i];\n",
|
|
" }\n",
|
|
" compare_arrays(ground_truth, h_C, VECTOR_SIZE);\n",
|
|
" \n",
|
|
" // print the first element of each host array to see that C[0] == A[0] + B[0]\n",
|
|
" printf(\"A[0] = %f\\n\", h_A[0]);\n",
|
|
" printf(\"B[0] = %f\\n\", h_B[0]);\n",
|
|
" printf(\"C[0] = %f\\n\", h_C[0]);\n",
|
|
" return 0;\n",
|
|
"}\n"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Assignment - Matrix Multiplication\n",
|
|
"\n",
|
|
"- **TODO_01** - For this task you need to add the indices for the dot product operation between the row of the first matrix and the column of the second matrix in order to finish the matrix multiplication code. Keep in mind that the matrices are stored in "
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%%cuda\n",
|
|
"#include <stdio.h>\n",
|
|
"#include <time.h>\n",
|
|
"#include \"utils.h\"\n",
|
|
"\n",
|
|
"const int MATRIX_SIZE = 512;\n",
|
|
"const int BLOCK_SIZE = 16;\n",
|
|
"\n",
|
|
"// matrix multiply (naive) kernel: C <- A * B\n",
|
|
"__global__ void matrix_multiply(const float *A, const float *B, float *C, int matrix_size) {\n",
|
|
" // each thread must compute one element in the output matrix, which is \n",
|
|
" // determined by the column and row index; together they form a globally \n",
|
|
" // unique pair of indices since there will be only one thread that will\n",
|
|
" // compute C[c_row][c_col]\n",
|
|
" int c_col = threadIdx.x + blockDim.x * blockIdx.x;\n",
|
|
" int c_row = threadIdx.y + blockDim.y * blockIdx.y;\n",
|
|
"\n",
|
|
" if ((c_col < matrix_size) && (c_row < matrix_size)) {\n",
|
|
" // compute in \"temp\" the dot product of row \"c_row\" from A \n",
|
|
" // and column \"c_col\" in B; save the result in C[c_row][c_col]\n",
|
|
" // assuming row-major ordering\n",
|
|
" float temp = 0;\n",
|
|
" /**************** TODO_01 *****************/\n",
|
|
" for (int i = 0; i < matrix_size; i++) {\n",
|
|
" temp += A[FIXME] * B[FIXME];\n",
|
|
" }\n",
|
|
" /******************************************/\n",
|
|
" C[c_row * matrix_size + c_col] = temp;\n",
|
|
" }\n",
|
|
"}\n",
|
|
"\n",
|
|
"int main(){\n",
|
|
" float *h_A, *h_B, *h_C, *d_A, *d_B, *d_C;\n",
|
|
" int n_elements = MATRIX_SIZE * MATRIX_SIZE;\n",
|
|
"\n",
|
|
" // variables used for computing the execution time\n",
|
|
" clock_t t0, t1, t2, t3;\n",
|
|
" double t1_diff = 0.0;\n",
|
|
" double t2_diff = 0.0;\n",
|
|
" double t3_diff = 0.0;\n",
|
|
" t0 = clock();\n",
|
|
" \n",
|
|
" h_A = new float[n_elements];\n",
|
|
" h_B = new float[n_elements];\n",
|
|
" h_C = new float[n_elements];\n",
|
|
" for (int i = 0; i < n_elements; i++){\n",
|
|
" h_A[i] = rand()/(float)RAND_MAX;\n",
|
|
" h_B[i] = rand()/(float)RAND_MAX;\n",
|
|
" }\n",
|
|
"\n",
|
|
" // Initialization timing\n",
|
|
" t1 = clock();\n",
|
|
" t1_diff = ((double)(t1-t0))/CLOCKS_PER_SEC;\n",
|
|
" printf(\"Init took %f seconds. Begin compute\\n\", t1_diff);\n",
|
|
"\n",
|
|
" // Allocate device memory and copy input data over to GPU\n",
|
|
" cudaMalloc(&d_A, n_elements*sizeof(float));\n",
|
|
" cudaMalloc(&d_B, n_elements*sizeof(float));\n",
|
|
" cudaMalloc(&d_C, n_elements*sizeof(float));\n",
|
|
" cudaCheckErrors(\"cudaMalloc failure\");\n",
|
|
" cudaMemcpy(d_A, h_A, n_elements*sizeof(float), cudaMemcpyHostToDevice);\n",
|
|
" cudaMemcpy(d_B, h_B, n_elements*sizeof(float), cudaMemcpyHostToDevice);\n",
|
|
" cudaCheckErrors(\"cudaMemcpy H2D failure\");\n",
|
|
"\n",
|
|
" // Launch kernel\n",
|
|
" dim3 block(BLOCK_SIZE, BLOCK_SIZE); // dim3 variable holds 3 dimensions\n",
|
|
" dim3 grid(divide_ceil(MATRIX_SIZE, block.x), divide_ceil(MATRIX_SIZE, block.y));\n",
|
|
" matrix_multiply<<<grid, block>>>(d_A, d_B, d_C, MATRIX_SIZE);\n",
|
|
" cudaCheckErrors(\"kernel launch failure\");\n",
|
|
"\n",
|
|
" // Copy results back to host\n",
|
|
" cudaMemcpy(h_C, d_C, n_elements*sizeof(float), cudaMemcpyDeviceToHost);\n",
|
|
"\n",
|
|
" // GPU timing\n",
|
|
" t2 = clock();\n",
|
|
" t2_diff = ((double)(t2-t1))/CLOCKS_PER_SEC;\n",
|
|
" printf (\"Done. Compute took %f seconds\\n\", t2_diff);\n",
|
|
"\n",
|
|
" cudaCheckErrors(\"kernel execution failure or cudaMemcpy H2D failure\");\n",
|
|
"\n",
|
|
" // verify results\n",
|
|
" float *ground_truth = new float[n_elements];\n",
|
|
" assignment_solution(h_A, h_B, ground_truth, MATRIX_SIZE);\n",
|
|
" t3 = clock();\n",
|
|
" t3_diff = ((double)(t3-t2))/CLOCKS_PER_SEC;\n",
|
|
" printf (\"Verifying results on CPU took %f seconds\\n\", t3_diff);\n",
|
|
" compare_arrays(ground_truth, h_C, n_elements);\n",
|
|
"\n",
|
|
" return 0;\n",
|
|
"}"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"### Solution\n",
|
|
"\n",
|
|
"We are not yet dealing with CUDA optimization problems. We just want to make this program produce the correct answer. The TODO is only a matter of mapping a 2D matrix into a 1D vector using [row-major ordering](https://en.wikipedia.org/wiki/Row-_and_column-major_order) like this:\n",
|
|
"\n",
|
|
"```cpp\n",
|
|
"/**************** TODO_01 *****************/\n",
|
|
"for (int i = 0; i < matrix_size; i++) {\n",
|
|
" temp += A[c_row * matrix_size + i] * B[i * matrix_size + c_col];\n",
|
|
"}\n",
|
|
"/******************************************/\n",
|
|
"```"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Future sessions\n",
|
|
"\n",
|
|
"- CUDA Shared Memory\n",
|
|
"- CUDA GPU architecture and basic optimizations\n",
|
|
"- Atomics, Reductions, Warp Shuffle\n",
|
|
"- Using Managed Memory\n",
|
|
"- Concurrency (streams, copy/compute overlap, multi-GPU)\n",
|
|
"- Analysis Driven Optimization\n",
|
|
"- Cooperative Groups"
|
|
]
|
|
}
|
|
],
|
|
"metadata": {
|
|
"kaggle": {
|
|
"accelerator": "nvidiaTeslaT4",
|
|
"dataSources": [],
|
|
"dockerImageVersionId": 30698,
|
|
"isGpuEnabled": true,
|
|
"isInternetEnabled": true,
|
|
"language": "python",
|
|
"sourceType": "notebook"
|
|
},
|
|
"kernelspec": {
|
|
"display_name": "Python 3",
|
|
"language": "python",
|
|
"name": "python3"
|
|
},
|
|
"language_info": {
|
|
"codemirror_mode": {
|
|
"name": "ipython",
|
|
"version": 3
|
|
},
|
|
"file_extension": ".py",
|
|
"mimetype": "text/x-python",
|
|
"name": "python",
|
|
"nbconvert_exporter": "python",
|
|
"pygments_lexer": "ipython3",
|
|
"version": "3.10.13"
|
|
}
|
|
},
|
|
"nbformat": 4,
|
|
"nbformat_minor": 4
|
|
}
|