<jupyter_start><jupyter_text>Learning ObjectivesIn this tutorial, the goal is to:* Parallelize the single-GPU code using CUDA Memcpy and streams* Understand intra-node topology and underlying technologies like GPUDirect P2P and their implication on program performance Multi-GPU ProgrammingIn this section we first cover the principle behind decomposing data among the GPUs, known as domain decomposition. Then, we understand and implement the baseline multi-GPU code using `cudaSetDevice` and `cudaMemcpy` functions.  Domain DecompositionBefore we begin, we define two important terms:* **Latency:** The amount of time it takes to take a unit of data from point A to point B. For example, if 4B of data can be transferred from point A to B in 4 $\mu$s, that is the latency of transfer.* **Bandwidth:** The amount of data that can be transferred from point A to point B in a unit of time. For example, if the width of the bus is 64KB and latency of transfer between point A and B is 4 $\mu$s, the bandwidth is 64KB * (1/4$\mu$s) = 1.6 GB/s.To parallelize our application to multi-GPUs, we first review the different methods of domain decomposition available to us for splitting the data among the GPUs, thereby distributing the work. Broadly, we can divide data into either stripes or tiles.* **Stripes**: They minimize the number of neighbours, require communication among less neighbours, and are optimal for latency bound communication.* **Tiles**: They minimize surface area/ volume ratio of the grid, require communicating less data, and are optimal for bandwidth bound communication.![domain_decomposition](../../images/domain_decomposition.png)When we divide the global grid between GPUs, only the boundaries of each GPU-local grid need to be communicated with the neighboring GPUs, as they need the updated grid-point values for the next iteration. Therefore, we use horizontal stripes (as C/ C++ are row-major) in our tutorials for domain decomposition, enabling data parallelism. Halo ExchangeWe term the exchange of top and bottom rows after each iterations the "halo exchange". Review the image below and notice that we update the topmost and bottomost rows of the grid to implement the periodic boundary condition. Recall that the left and right columns of the grid constitute Dirichlet boundary conditions (that is, constant value).![halo_exchange](../../images/halo_exchange.png) CUDA concepts: Part 1 Setting the GPUTo verify that our system has multiple GPUs in each node, run the command below:<jupyter_code>!nvidia-smi<jupyter_output><empty_output><jupyter_text>`nvidia-smi` utility shows the available GPU on a node, but inside a CUDA program, the number of GPU in the node can be obtained using the `cudaGetDeviceCount(int *count)` function. To perform any task, like running a CUDA kernel, copy operation, etc. on a particular GPU, we use the `cudaSetDevice(int device)` function. Copying between GPUsThe `cudaMemcpy` function supports GPU to GPU copy using the `cudaMemcpyDeviceToDevice` flag and the source and destination memory addresses should reside in GPU devices. For example, if we want to copy 1000 floats from the array `arr_gpu_0` allocated on GPU 0 to the array `arr_gpu_1`, the function call is:```ccudaMemcpy(arr_gpu_1, arr_gpu_0, 1000 * sizeof(float), cudaMemcpyDeviceToDevice);```Recall that CUDA kernel calls made from the host are non-blocking (asynchronous) by default. That is, the control may return back to the host thread before the device kernel finishes execution. To perform the halo exchange, we need to perform copy operations between each GPU and its neighbours. However, for large copy sizes, `cudaMemcpy` is blocking with respect to the host. Thus, we cannot use the following code snippet:```cfor (int i = 0; i < 2; i++) {    // Set current device    cudaSetDevice(i);    // Define row number of top and bottom neighbours, etc.    TopNeighbour = ...; BotNeighbour = ...; // and so-on    // Launch device kernel on GPU i    jacobi_kernel>>(...);    // Halo exchange    cudaMemcpy(grid_rows[TopNeighbour], grid_rows[myTop], size, cudaMemcpyDeviceToDevice);    cudaMemcpy(grid_rows[BotNeighbour], grid_rows[myBot], size, cudaMemcpyDeviceToDevice);    // Norm check, swapping current and previous grid arrays, etc.} // Serializes operations with respect to the host```This code results in serialized execution as shown in diagram below:![memcpy_serialized](../../images/memcpy_serialized.png) Asynchronous operationsInstead of `cudaMemcpy`, we can use the `cudaMemcpyAsync` function which is asynchronous with respect to the host. This allows the host to launch device kernels and copy operations concurrently, enabling parallel execution across GPUs. The correct code snippet is as follows:```cfor (int i = 0; i < 2; i++) {    // Set current device    cudaSetDevice(i);    // Launch device kernel on GPU i    jacobi_kernel>>(...);}for (int i = 0; i < 2; i++) {    // Define row number of top and bottom neighbours, etc.    TopNeighbour = ...; BotNeighbour = ...; // and so-on    // Halo exchange, notice the use of Async function    cudaMemcpyAsync(grid_rows[TopNeighbour], grid_rows[myTop], size, cudaMemcpyDeviceToDevice);    cudaMemcpyAsync(grid_rows[BotNeighbour], grid_rows[myBot], size, cudaMemcpyDeviceToDevice);    // Norm check, swapping current and previous grid arrays, etc.} // Parallel execution across multiple GPUs```And the execution time of the application is reduced:![memcpyasync_parallel](../../images/memcpyasync_parallel.png) Implementation exercise: Part 1Now, let's parallelize our code across multiple GPUs by using `cudaSetDevice` and `cudaMemcpyAsync` operations. Open the [jacobi_memcpy.cu](../../source_code/cuda/jacobi_memcpy.cu) file.Alternatively, you can navigate to `CFD/English/C/source_code/cuda/` directory in Jupyter's file browser in the left pane. Then, click to open the `jacobi_memcpy.cu` file. Understand the flow of the program from within the `main` function. Review the following pre-Jacobi-computation steps:1. Computation of the memory chunk size to be allocated on each GPU stored in the `chunk_size` integer array.2. Allocation of memory on each GPU: Notice the use of array pointers like `a_new`, `l2_norm_d`, `iy_start`, etc. that point to device arrays allocated on GPU pointed to by `dev_id` variable.3. Initialization of Dirichlet boundary conditions on left and right boundaries.4. Share of initial top and bottom local grid-point values between neighbours.Now, within the iterative Jacobi loop (the `while` loop), implement the following marked as `TODO: Part 1-`:1. Set current GPU and call device kernel with correct device arrays in function arguments.2. Asynchronously copy GPU-local L2 norm back to CPU and implement top and bottom halo exchanges.3. Synchronize the devices at the end of each iteration using `cudaDeviceSynchronize` function.Review the topic above on Asynchronous Operations if in doubt. We will be using separate `for` loops for launching device kernels and initiating copy operations.After implementing these, let's compile the code:<jupyter_code>!cd ../../source_code/cuda && make clean && make jacobi_memcpy<jupyter_output><empty_output><jupyter_text>Ensure there are no compiler warnings or errors. Validate the implementation by running the binary:<jupyter_code>!cd ../../source_code/cuda && ./jacobi_memcpy<jupyter_output><empty_output><jupyter_text>The last couple of lines of the output will give the number and IDs of GPUs used, execution timings, speedup, and efficiency metrics. Review Metrics of Interest section in [single GPU overview](../single_gpu/single_gpu_overview.ipynb) tutorial for more information). We tested the code on a DGX-1 system with 8 Tesla V100 16GB GPUs, and we got the following output:```bashNum GPUs: 8. Using GPU ID: 0, 1, 2, 3, 4, 5, 6, 7, 16384x16384: 1 GPU:   4.4485 s, 8 GPUs:   1.0951 s, speedup:     4.06, efficiency:    50.78 ```Notice that we got a speed-up of $4.06\times$ using 8 GPUs and a corresponding efficiency of $50.78\%$. The numbers will vary depending on number of available GPUs in your system, the communication topology, GPU type, etc. ProfilingNow, profile the execution with `nsys`:<jupyter_code>!cd ../../source_code/cuda/ && nsys profile --trace=cuda,nvtx --stats=true -o jacobi_memcpy_report --force-overwrite true ./jacobi_memcpy<jupyter_output><empty_output><jupyter_text>[Download the profiler report here for visualization](../../source_code/cuda/jacobi_memcpy_report.qdrep). In the profiler timeline, the first few seconds denote the single-GPU code running on one of the GPUs. This version is executed so we can compare the multi-GPU version with it. Let's analyze the multi-GPU timeline.![jacobi_memcpy_report_overview](../../images/jacobi_memcpy_report_overview.png)The next iteration of the device kernel is not run till all inter-GPU copy operations are complete because we need to synchronize all GPUs at the end of each iteration. The total time taken by the Jacobi Solver loop (`jacobi_solve` NVTX annotatation) is visible and is 1.278 seconds. Also, notice the we have labelled halo exchanges as Device-to-Host (DtoH) and Host-to-Device) copies. Now, right click on `CUDA HW` tab and select `Show in Events View` option. ![jacobi_memcpy_report_events](../../images/jacobi_memcpy_report_events.png)The "Source Memory Kind" and "Destination Memory Kind" of the selected DtoH operation are both "Device". However the copy operation is marked as "Memcpy DtoH". By default, the device-to-device copy operation uses a temporary CPU buffer internally. Let us understand more about this CPU buffer and how we can eliminate it to improve performance. CUDA concepts: Part 2 Host Staging of Copy OperationsUsing `cudaMemcpyAsync` instead of `cudaMemcpy` allows us to issue copy and compute operations on multiple GPUs concurrently. The path taken by the data in both the cases is denoted by the red arrow as follows:![memcpy_host_staging](../../images/memcpy_host_staging.png)That is, in the GPU-to-GPU memory copy, the data traverses from GPU 0 the PCIe bus to the CPU, where it is staged in a buffer before being copied to GPU 1. This is called "host staging" and it decreases the bandwidth while increasing the latency of the operation. If we eliminate host staging, we can usually improve the performance of our application. Peer-to-Peer Memory AccessP2P allows devices to address each other's memory from within device kernels and eliminates host staging by transferring data either through the PCIe switch or through NVLink as denoted by the red arrow below. ![memcpy_p2p_overview](../../images/memcpy_p2p_overview.png)Peer-to-Peer (P2P) memory access requires GPUs to share a Unified Virtual Address Space (UVA). UVA means that a single address space is used for the host and all modern NVIDIA GPU devices (specifically, those with compute capibility of 2.0 or higher).This P2P memory access feature is supported between two devices if `cudaDeviceCanAccessPeer()` returns true for these two devices. P2P must be enabled between two devices by calling `cudaDeviceEnablePeerAccess()` as illustrated in the following code sample:```ccudaSetDevice(currDevice);int canAccessPeer = 0;cudaDeviceCanAccessPeer(&canAccessPeer, currDevice, PeerDevice);if (canAccessPeer) {    cudaDeviceEnablePeerAccess(PeerDevice, 0);}```Note that this enables a unidirectional P2P access where `currDevice` can perform memory access to `PeerDevice`. If we want `PeerDevice` to be able to access `currDevice` via P2P, then we need to use the code accordingly.First, let's check if P2P is supported between the GPUs:<jupyter_code>!nvidia-smi topo -p2p r<jupyter_output><empty_output><jupyter_text>The `topo` sub-command requests information on the GPU communication topology, `-p2p` flag requests P2P status, and `r` asks whether P2P reads are supported. Change `r` to `w` to check whether writes are supported. We share our output on a DGX-1 system with 8 Tesla V100s, focusing on the capabilities of GPU 0:![nvidia_smi_p2p_gpu0](../../images/nvidia_smi_p2p_gpu0.png)This means GPU 0 can communicate via P2P with GPUs 1 through 4. For GPUs 5 through 7, it must use host staging.To check whether P2P via NVLink is supported, run the command below:<jupyter_code>!nvidia-smi topo -p2p n<jupyter_output><empty_output><jupyter_text>In our DGX-1 system, the result is similar as before. Even if P2P via NVLink is not supported on your system, as long as `-p2p r` and `-p2p w` are supported between GPUs, P2P capability is available. Implementation Exercise: Part 2Now, let us improve our program performance by enabling P2P access between GPUs, wherever possible. The `jacobi_memcpy.cu` code accepts a runtime argument `-p2p` which should enable P2P access between GPUs. Modify the code by searching for `TODO: Part 2` and enabling GPU `devices[dev_id]` to access peer GPUs `devices[top]` and `devices[bottom]`, whenever possible. Notice that the code snippet is within a `for` loop which sets and iterates over each GPU, which is why bidirectional P2P will be enabled. Take help from the code sample in the previous section.Now, let's compile the code again:<jupyter_code>!cd ../../source_code/cuda && make clean && make jacobi_memcpy<jupyter_output><empty_output><jupyter_text>Ensure there are no compiler warnings or errors. Validate the implementation by running the binary with P2P enabled:<jupyter_code>!cd ../../source_code/cuda && ./jacobi_memcpy -p2p<jupyter_output><empty_output><jupyter_text>The output we got on our DGX-1 system is:```bashNum GPUs: 8. Using GPU ID: 0, 1, 2, 3, 4, 5, 6, 7, 16384x16384: 1 GPU:   4.4487 s, 8 GPUs:   0.8798 s, speedup:     5.06, efficiency:    63.21 ```Notice that the efficiency increased by about $8\%$ to $63.21\%$ compared to our baseline implementation. You can run the baseline again by removing the `-p2p` flag. Note that if P2P is not supported on your system, you will likely not experience any performance gain. ProfilingLet us profile the execution with `nsys`:<jupyter_code>!cd ../../source_code/cuda/ && nsys profile --trace=cuda,nvtx --stats=true -o jacobi_memcpy_p2p_report --force-overwrite true ./jacobi_memcpy -p2p<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.0754477506698632, 'fraction_numerical': 0.014031871386264279, 'mean_word_length': 4.924394319131161, 'pattern_counts': {'":': 0, '<': 40, '<?xml version=': 0, '>': 41, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '542015', 'license': 'Apache-2.0', 'n_tokens_mistral': 4045, 'n_tokens_neox': 3761, 'n_words': 1877, 'path': 'hpc/multi_gpu_nways/labs/CFD/English/C/jupyter_notebook/cuda/memcpy.ipynb', 'repo_name': 'AditiRM/gpubootcamp', 'size': '14182'}