CUDA 编程学习手册（五）

原文：annas-archive.org/md5/f6da79e769f988319eb178273ecbf55b

译者：飞龙

协议：CC BY-NC-SA 4.0

使用 OpenACC 进行 GPU 编程

每个处理器架构都提供了不同的编写代码以在处理器上运行的方法。CUDA 也不例外；它也提供了不同的编码方法。近年来变得非常流行的一种方法是使用 OpenACC，它基本上是基于指令的编程。

OpenACC 基本上是一个将异构计算作为一等公民的标准。该标准基本上规定了有两种处理器，即主机和设备/加速器，这与 CUDA 编程模型所述的概念非常相似。

对于希望获得最佳性能的程序员来说，使用诸如 C、C++、Fortran 和 Python 等语言的 CUDA 编程是表达并行性的首选方式。编程语言要求程序员从头开始重新创建他们的顺序程序，同时保持他们的关键操作的串行和并行版本。并行编程语言中创建的并行程序往往只适用于非常少数的平台。

编译器指令将编程语言的灵活性与库的易用性相结合。程序员使用高级指令对代码进行注释，编译器可以使用这些指令来并行化代码，或者可以安全地忽略。这意味着带有编译器指令的代码可以编译为许多不同的并行平台，并且无需维护代码的串行和并行版本。此外，有时需要快速测试和原型化应用程序以在 GPU 上运行。一个这样的例子是将天气代码等代码库（拥有数百万行代码）转换为在 GPU 上运行；使用流行的语言将需要大量的工作。在这种情况下，OpenACC 成为一个合乎逻辑的选择。在 OpenACC 中，开发人员以指令的形式向编译器提供提示。编译器接受这些提示并生成特定于架构的加速器代码。

OpenACC 标准还为代码的开发人员提供了供应商中立性。带有 OpenACC 指令的单一源代码可以重新编译为不同的设备。例如，PGI 编译器目前支持 OpenACC 后端，如 Intel CPU 多核、NVIDIA GPU、Intel Xeon Phi 和 FPGA/ASIC 架构。这对于希望编写供应商中立代码的开发人员来说是一个非常有吸引力的提议。高性能计算中的关键应用程序，如 Vienna Ab-initio Simulation Package（VASP）（分子动力学/量子化学）、Weather Research and Forecasting（WRF）和 ANSYS Fluent（CFD）利用 OpenACC 编程模型来针对 NVIDIA GPU。

总结 OpenACC 的关键要点：

• 当异构计算被视为新的编程模型时，OpenACC 标准得以发展。

• OpenACC 在各种加速器上提供性能可移植性。

• OpenACC 并不是 CUDA 编程语言的替代品。当选择的处理器是 NVIDIA 时，OpenACC 编译器在后台生成 CUDA 代码。

近年来，OpenMP 标准也开始纳入异构计算 API。但迄今为止，还没有支持不同处理器架构的编译器，因此我们选择在本书中坚持使用 OpenACC。

本章将涵盖以下主题：

• OpenACC 指令

• OpenACC 中的异步编程

• 额外重要的指令和子句

技术要求

本章需要一台带有现代 NVIDIA GPU（Pascal 架构或更高）的 Linux/Windows PC。

如介绍中所述，OpenACC 是一个标准，这个标准由不同的编译器实现，如 GCC、PGI 和 CRAY 编译器。我们将在本章中使用的编译器是 PGI。PGI 编译器在 Fortran 社区中非常受欢迎，并且一直在实现 OpenACC 最新规范方面处于领先地位，并且提供了一个可以从 PGI 网站免费下载的社区版。好处是在社区版和付费版本的 PGI 编译器之间在功能上基本没有变化。在本章中，您需要下载 PGI 社区版。

本章的代码也可以在 GitHub 上找到：github.com/PacktPublishing/Learn-CUDA-Programming。

示例代码示例是使用 PGI 社区版的 19.4 版本开发和测试的。但建议您使用最新的 PGI 版本。

使用 OpenACC 在 GPU 上合并图像

为了理解 OpenACC 概念，我们选择了一个简单的计算机视觉算法来合并两个图像。在这段代码中，我们试图合并两个图像，如下所示：

前面的图像演示了一个计算机视觉算法，用于合并两个图像。

我们将在本章后面更多地讨论代码结构。首先，根据以下步骤配置环境：

1. 准备您的 GPU 应用程序。例如，我们将使用一个用于合并两个图像的核算法。此代码可以在09_openacc/中找到。

2. 使用pgc++编译器编译您的应用程序：

$ pgc++ -c -acc -ta=tesla:pinned scrImagePgmPpmPackage.cpp
$ pgc++ -c -acc -ta=tesla:pinned -Minfo=accel image_merging.cpp
$ pgc++ -o merging.out -acc -ta=tesla:pinned -Minfo=accel scrImagePgmPpmPackage.o image_merging.o
$ ./merging.out

前面的命令将创建一个名为blurring.out的二进制文件。正如您可能已经观察到的，我们正在使用pgc++编译器来编译我们的代码。此外，我们向我们的代码传递了一些参数。让我们更详细地了解它们：

• -acc：此标志告诉编译器解析代码中提供的 OpenACC 指令。

• -ta：代表应该为设备代码生成的目标架构。请注意，-ta=tesla表示我们的目标是 NVIDIA GPU。其他目标的一些示例包括-ta=multi-core，它将多核作为设备目标，-ta=radeaon，它将 AMD GPU 作为目标，还有一些其他目标。此外，我们可以添加特定于设备的标志；例如，我们为分配所有 CPU 内存作为固定（不可分页）的 GPU 添加了一个固定标志。

• -Minfo：此选项告诉编译器为我们提供有关编译器采取的步骤的更多信息，使我们的代码并行化。通过说-Minfo-accel，我们要求编译器为我们提供与加速器区域相关的更多信息。我们可以将标志更改为-Minfo=all，以提供非加速器区域的详细信息。以下输出显示了向我们的代码添加Minfo标志的部分输出：

.... < More compiler output above>
merge_parallel_pragma(unsigned char *, unsigned char *, unsigned char *, long, long):
    30, Generating copyin(in1[:w*h])
    Generating copyout(out[:w*h])
    Generating copyin(in2[:w*h])
    Accelerator kernel generated
    Generating Tesla code
    30, #pragma acc loop gang /* blockIdx.x */
    32, #pragma acc loop vector(128) /* threadIdx.x */
    32, Loop is parallelizable
... < More compile output below >

要理解这个编译输出，我们需要了解 OpenACC pragma，我们将在下一节中进行。稍后我们将重新访问这个编译输出。可以使用pgc++ --help找到其他可用标志的更多详细信息。

运行二进制文件后的示例输出如下：

$ ./merging.out
Reading image width height and width [1536][2048]
Time taken for serial merge: 0.0028 seconds
Time taken for OpenACC merge(data+kernel): 0.0010 seconds
Time taken for OpenACC merge(kernel only) with Blocking: 0.0002 seconds
 Time taken for OpenACC merge(data _kernel) with blocking: 0.0014 seconds
Time taken for OpenACC merge (data+kernel)with Pipeline Async: 0.0008 seconds

前面的输出显示我们正在读取一个大小为 1536*2048 的图像。代码有一个串行实现和三个使用 OpenACC pragma 的并行实现。每个实现的时间在前面的输出中显示。最后一个使用 pipeline 方法的实现显示了最佳时间：0.0008 秒。我们将采取增量方法，并在接下来的部分详细介绍每个实现。

该算法的串行实现非常简单，如下面的代码片段所示：

void merge_serial(unsigned char *in1, unsigned char*in2, unsigned char *out, long w, long h)
{
    long x, y;
     for(y = 0; y < h; y++) {
         for(x = 0; x < w; x++) {
             out[y * w + x] = (in1[y * w + x]+in2[y * w + x])/2;
         }
     }
}

代码没有什么特别之处；基本上，它接受两个输入图像数据（in1和in2），执行平均操作以合并两个输入，最后存储输出。对于我们来说，关键的是循环是尴尬并行的，适合于 GPU 等架构。如上面的代码输出所示，串行实现花费了0.0028秒。请注意，计时可能会因运行代码的系统而略有不同。

在下一节中，我们将向您介绍 OpenACC 指令，以便将示例代码转换为在 GPU 上运行所需的指令。

OpenACC 指令

在本节中，我们将尝试理解 OpenACC pragma 的语法，并为合并操作实现基本的并行和数据指令。OpenACC pragma 的基本语法如下：

#pragma acc <directive> <clauses> 
!$acc parallel [clause [[,] clause]…] 

上述命令解释如下：

• 在 C/C++中的#pragma被称为“编译器提示”。这些与程序员注释非常相似；但是，编译器实际上会读取我们的 pragma。如果编译器不理解 pragma，它可以忽略它，而不是抛出语法错误。

• acc是我们 pragma 的一个补充。它指定这是一个 OpenACC pragma。任何非 OpenACC 编译器都会忽略此 pragma。

• 指令是 OpenACC 中的一个命令，它告诉编译器执行某些操作。目前，我们只会使用允许编译器并行化我们的代码的指令。

• 子句是对我们的指令的补充/修改。这些包括但不限于优化。

在本节中，我们将介绍三个指令：parallel，loop和data。我们将展示它们各自的用法，并最终将它们应用到我们的合并算法中。

并行和循环指令

并行指令是最直接的指令。它将标记代码的一个区域进行并行化（通常只涉及并行化一个for循环），如下面的代码所示：

#pragma acc parallel loop 
for (int i = 0; i < N; i++ ) {  
    //loop code 
}

我们还可以定义一个并行区域。并行区域可以有多个循环（尽管这通常不推荐！）。并行区域是指最外层花括号内的所有内容，如下面的代码片段所示：

#pragma acc parallel
{
    #pragma acc loop
    for (int i = 0; i < N; i++ )
    {
        < loop code >
    }
}

包含循环非常重要；否则，您将无法正确地并行化循环。并行指令告诉编译器冗余地并行化代码，如下所示：

循环指令明确告诉编译器我们希望并行化循环，如下面的屏幕截图所示：

循环指令有两个主要用途：

• 标记单个循环进行并行化

• 允许我们明确定义循环的优化/修改

我们将在本章后面讨论循环优化，以及 gang 和 vector；目前，我们将专注于并行化方面。循环指令要正常工作，必须包含在并行指令内：

#pragma acc parallel loop
for (int i = 0; i < N; i++ )
{
    //loop code 
}

使用并行指令时，必须包含循环指令才能使代码正常运行。我们还可以使用循环指令来并行化多维循环嵌套。在下面的代码片段中，我们看到了一个嵌套循环，并且我们明确为第二个循环提到了循环子句：

#pragma acc parallel loop
for (int i = 0; i < N; i++ )
{
    #pragma acc loop
    for( int j = 0; j < M; j++ )
    {
        //loop code
    }
}

请注意，在上面的代码片段中，我们没有在内部循环中再次放置并行子句，因为我们已经在从外部循环开始的范围中提到了它。

数据指令

OpenACC 并行模型规定我们有一个主机，运行我们的顺序代码（通常是 CPU）。然后我们有我们的设备，这是某种并行硬件。主机和设备通常（虽然并非总是）有单独的内存，程序员可以使用 OpenACC 在两个内存之间移动数据。

正如在第一章中讨论的，GPU 和 CPU 架构在根本上是不同的。GPU 作为吞吐量架构，具有大量计算单元和高速内存带宽。另一方面，CPU 是一种减少延迟的架构，具有大型缓存层次结构，并且提供大容量的主存储器。需要操作的任何数据都需要首先复制到 GPU 内存。（请注意，即使在统一内存的情况下，数据也会在后台以页面的形式由驱动程序复制。）

如下图所示，两种架构（CPU 和 GPU）之间的数据传输通过 I/O 总线进行：

在 OpenACC 中使用 GPU 作为目标架构的目标是仅将并行代码卸载到 GPU 上，而顺序代码将继续在 CPU 上运行。OpenACC 标准允许程序员通过使用 OpenACC 数据指令和数据子句 显式定义数据管理。数据子句允许程序员在主机和设备（或在我们的情况下，CPU 和 GPU）之间指定数据传输。

**隐式数据管理：**我们可以将数据传输留给编译器，如下例所示：

int *A = (int*) malloc(N * sizeof(int));

#pragma acc parallel loop
for( int i = 0; i < N; i++ )
{
    A[i] = 0;
}

在前面的代码中，编译器将理解需要从 GPU 复制A向量，并为开发人员生成隐式传输。

**显式数据管理：**最好使用显式数据传输来获得对传输更多控制，如下面的代码中使用复制数据子句所示：

int *a = (int*) malloc(N * sizeof(int));
#pragma acc parallel loop copy(a[0:N])
for( int i = 0; i < N; i++ )
{
     a[i] = 0;
}

在前面的代码片段中，我们使用了复制数据子句。下图解释了运行时到达复制数据指令时执行的步骤：

我们将通过合并代码的详细步骤来解释这些步骤，其中我们将应用数据子句。

其他可用的数据子句如下所列：

数据子句	描述	关键用法
copy(list)

• 在设备上分配内存

• 在进入区域时，从主机复制数据到设备

• 在退出区域时，将数据复制到主机

这是默认的输入数据结构，被修改后从函数返回
copyin(list)

• 在设备上分配内存

• 在进入区域时，从主机复制数据到设备

作为子例程的输入的向量
copyout(list)

• 在设备上分配内存

• 在退出区域时，将数据复制到主机

不覆盖输入数据结构的结果
create(list)

• 仅在设备上分配内存

• 不进行复制

临时数组

为了最大化性能，程序员应避免所有不必要的数据传输，因此显式内存管理优于隐式数据管理。

**数组形状：**数组形状是指定数组大小的方式。如果不指定形状，编译器将尝试假定大小。这在 Fortran 中效果很好，因为 Fortran 跟踪数组的大小；然而，在 C/C++中可能不起作用。数组形状也是从数组复制数据的唯一方式（例如，如果只需要复制数组的一半，这可能提高性能，减少不必要的复制），如下面的代码片段所示：

#pragma acc parallel loop copy(A[1:N-2])

这将复制A的所有元素，除了第一个和最后一个元素。

将并行、循环和数据指令应用于合并图像代码

现在让我们尝试将并行、循环和数据指令应用于合并顺序代码：

void merge_parallel_pragma(unsigned char *in1, unsigned char*in2,unsigned char *out, long w, long h)
{
    long x, y;
    #pragma acc parallel loop gang copyin(in1[:h*w],
                                          in2[:h*w]) 
                                          copyout(out[:h*w])
     for(y = 0; y < h; y++) {
        #pragma acc loop vector
        for(x = 0; x < w; x++) {
            out[y * w + x] = (in1[y * w + x]+in2[y * w + x])/2;
        }
    }
}

我们已经使用并行循环指令并行化了两个循环（高度：y和宽度：x）。此外，我们还明确地添加了数据子句来复制数据。请注意，由于in1和in2向量只是输入，它们是使用copyin()数据子句进行复制的。out向量是输出，使用copyout()数据子句进行复制。让我们试着理解这个函数的编译器输出：

merge_parallel_pragma(unsigned char *, unsigned char *, unsigned char *, long, long):
    30, Generating copyin(in1[:w*h])
        Generating copyout(out[:w*h])
        Generating copyin(in2[:w*h])
        Accelerator kernel generated
        Generating Tesla code
        30, #pragma acc loop gang /* blockIdx.x */
        32, #pragma acc loop vector(128) /* threadIdx.x */
32, Loop is parallelizable

前面的编译器输出显示，对于merge_parallel_pragma函数，编译器生成了以下操作：

• 在第 30 行，为in1和in2变量生成了copyin。在内核启动前将被复制到 GPU 的数组大小将是[0:w*h]。

• 在第 30 行，为out变量生成了copyout。在 GPU 内核启动后将被复制的数组大小将是[0:w*h]。

• 在第 30 和 32 行，生成了 Tesla 内核代码：

• 在第 30 行，外部循环使用了 gang 级并行化。

• 在第 32 行，内部循环使用了矢量级并行化

当代码在 V100 上运行时，整个内核所花费的时间为0.0010s。这基本上是串行代码的两倍快。这可能听起来并不令人印象深刻。原因是大部分时间花在了数据传输上，而不是内核计算。为了确认这一点，让我们使用nvprof：

$ nvprof ./merging.out
==26601== DoneProfiling application: ./merging.out
==26601== Profiling result:
Type Time(%) Time Calls Avg Min Max Name
GPU activities: 67.36% 609.41us 2 304.71us 286.34us 323.08us [CUDA memcpy HtoD]
27.63% 250.02us 1 250.02us 250.02us 250.02us [CUDA memcpy DtoH]
5.01% 45.344us 1 45.344us 45.344us 45.344us merge_parallel_pragma_30_gpu(unsigned char*, unsigned char*, unsigned char*, long, long)
...

正如您在前面的分析输出中所观察到的，94%的时间花在了数据传输上，而只有 5%的时间（45 微秒）花在了内核执行上。您可能会问：我怎么知道这是哪个内核？如果您仔细看 GPU 内核的名称，merge_parallel_pragma_30_gpu，PGI 编译器在merge_parallel_pragma函数的第 30 行生成了一个 CUDA 内核，这就是我们如何将其与在该行号放置的编译指示联系起来的方式。

所以我们知道问题在哪里，但解决方案呢？我们将使用的优化技术是 blocking 来隐藏这种延迟。我们将在接下来的章节中更多地介绍 blocking 技术，并使用异步子句来重叠这个传输。

OpenACC 中的异步编程

为了实现合并并行代码的更好性能，我们将利用一个叫做 blocking 的概念。Blocking 基本上意味着，我们可以创建数组的块，而不是一次性传输整个输入和输出数组，这些块可以并行传输和操作。以下图表演示了创建块并在内核执行时重叠数据传输：

前面的图表显示了不同的块被传输，这些块的内核执行可以独立于每个块。为了实现这一点，我们需要数据传输命令和内核调用被异步地触发和执行。为了实现 blocking，我们将在本节中引入更多的指令/子句：结构化/非结构化数据指令和async子句。我们将展示它们的每一个，并最终将它们应用到我们的基本 OpenACC 合并并行代码中。

结构化数据指令

OpenACC 数据指令允许程序员显式地管理设备上的数据（在我们的例子中是 GPU）。以下代码片段显示了标记结构化数据区域的示例：

< Initialize data on host (CPU) >
#pragma acc data < data clauses >
{
    //< Code >
}

设备内存分配发生在区域的开始，设备内存释放发生在区域的结束。此外，从主机到设备（CPU 到 GPU）的任何数据移动发生在区域的开始，从设备到主机（GPU 到 CPU）的任何数据移动发生在区域的结束。内存分配/释放和数据移动是由程序员包含的子句定义的。

**包含多个计算区域：**一个数据区域可以包含任意数量的并行/内核区域，如下例所示：

#pragma acc data copyin(A[0:N]) create(C[0:N])
{
    #pragma acc parallel loop
    for( int i = 0; i < N; i++ )
    {
        C[i] = A[i] + 10;
    }
    #pragma acc parallel loop
    for( int i = 0; i < N; i++ )
    {
        C[i] = C[i] / 10;
    }
}

非结构化数据指令

有两个非结构化数据指令：

• 进入数据：处理设备内存分配，并从主机复制到设备。您可以在进入数据中使用的两个子句是：

• create：这将只执行设备内存分配。

• copyin：这将执行分配以及从设备到设备的内存复制。

• 退出数据：处理设备内存释放，并从设备复制到主机。您可以在退出数据中使用的两个子句是：

• delete：这将仅执行设备内存释放。

• copyout：这将首先从设备复制内存到主机，然后执行设备内存释放。

非结构化数据指令不会将数据区域标记为您可以在代码中有多个进入数据和退出数据指令。最好将它们纯粹视为内存分配和释放。使用非结构化数据指令的最大优势是它们能够跨多个函数进行分支。您可以在一个函数中分配数据，并在另一个函数中释放它。我们可以看一个简单的例子：

#define N 1024
int* allocate(int size)
{
    int *ptr = (int*) malloc(size * sizeof(int));
    #pragma acc enter data create(ptr[0:size])
    return ptr;
} 
void deallocate(int *ptr)
{
    #pragma acc exit data delete(ptr)
    free(ptr);
}
int main()
{
    int *ptr = allocate(N);
    #pragma acc parallel loop
    for( int i = 0; i < N; i++ )
    {
        ptr[i] = 0;
    }
    deallocate(ptr);
}

上面的代码片段显示了分配发生在单独的allocate()函数中，删除发生在deallocate()中。您可以将相同的概念链接到 C++中构造函数的一部分enter data create和析构函数的一部分exit data delete。

OpenACC 中的异步编程

默认情况下，所有 OpenACC 调用都是同步的。这意味着，在每次数据传输或每次对 GPU 的内核调用之后，都会隐式添加同步。CPU 将等待直到 OpenACC 调用完成，然后开始执行下一条指令。为了使调用异步，我们可以在数据和并行指令中使用async子句，如下面的代码所示：

#pragma acc data copyin(a[:N]) async 
// performing copyin asynchronously 
#pragma acc parallel loop async 
//performing parallel loop asynchronously. 

使用async的主要好处可以总结如下：

• 如果我们想要同时执行主机和设备代码，我们可以使用async启动我们的设备代码，而在执行时我们可以返回到主机继续不相关（非设备相关）的代码。

• 我们可以排队多个设备内核启动，以便它们连续执行，这在某些情况下可以减少启动设备内核的开销。

• 我们可以在主机和设备之间同时执行数据移动和设备计算**。**这是我们将应用于我们的代码的优化，并且是async的最常见用例。

在幕后，每当我们使用async子句时，我们都会向队列添加一些工作。提交给不同队列的工作可以异步执行，而在同一队列中的工作将顺序执行（一个接着一个）。当我们使用async时，我们可以指定队列号。如果未指定队列号，则将自动使用默认值。

将非结构化数据和异步指令应用于合并图像代码

现在让我们尝试将数据指令与async子句一起应用于合并并行代码：

void merge_async_pipelined(unsigned char *in1, unsigned char*in2,unsigned char *out, long w, long h)
{
    long x, y;
    #pragma acc enter data create(in1[:w*h], in2[:h*w], out[:w*h])
    const long numBlocks = 8;
    const long rowsPerBlock = (h+(numBlocks-1))/numBlocks;
    for(long block = 0; block < numBlocks; block++) {
        long lower = block*rowsPerBlock; // Compute Lower
        long upper = MIN(h, lower+rowsPerBlock); // Compute Upper
        #pragma acc update device(in1[lower*w:(upper-lower)*w],
                                  in2[lower*w:(upper-lower)*w]) 
                                  async(block%2)
        #pragma acc parallel loop present(in1,in2, out) async(block%2)
        for(y = lower; y < upper; y++) {
            #pragma acc loop
            for(x = 0; x < w; x++) {
                out[y * w + x] = (in1[y * w + x]+in2[y * w + x])/2;
            }
        }
        #pragma acc update self(out[lower*w:(upper-lower)*w]) 
                                async(block%2)
    }
#pragma acc wait
#pragma acc exit data delete(in1, in2, out)
}

我们已经使用了数据指令和async子句来实现阻塞概念。让我们分解整体实现，这将使其更容易理解：

1. 进入数据区域：enter data create子句在 GPU 中为in1和in2变量以及out分配内存。

2. 创建块：我们决定将图像分成八个块。这些块分布在行中。外部的for循环用于此目的添加了这个原因。

3. 异步从主机传输数据到设备：acc update device基本上将数据从主机异步复制到设备，因为我们已经在其中添加了一个async子句。

4. 异步启动并行循环：async子句被添加到并行子句中，以异步启动 GPU 内核。

5. 异步从设备传输数据到主机：acc update self基本上是将数据从设备异步地复制到主机，因为我们已经在同一个地方添加了一个async子句。

6. 等待：acc wait将确保 CPU 等待，直到所有 OpenACC 启动都完成，然后在所有队列中继续前进。

7. 退出数据区域：acc exit data delete将删除在enter data子句中分配的in1和in2向量以及out。

让我们试着理解merge_async_pipelined函数的编译器输出：

merge_async_pipelined(unsigned char *, unsigned char *, 
                      unsigned char *, long, long):
     67, Generating enter data create(out[:h*w],in2[:h*w],in1[:h*w])
     74, Generating update device(in1[w*lower:w*(upper-lower)],
                                  in2[w*lower:w*(upper-lower)])
         Generating present(in1[:],out[:],in2[:])
         Accelerator kernel generated
         Generating Tesla code
         74, #pragma acc loop gang /* blockIdx.x */
         76, #pragma acc loop vector(128) /* threadIdx.x */
     76, Loop is parallelizable
     81, Generating update self(out[w*lower:w*(upper-lower)])
     84, Generating exit data delete(out[:1],in2[:1],in1[:1])

前面的编译器输出显示，对于merge_async_pipelined函数，编译器生成了以下操作：

• 在第 67 行，为in1、in2和out变量生成了data create区域。

• 在第 74 行，为in1和in2调用了update device，并且数据传输到设备被限制在上下界之间：in1[w*lower:w*(upper-lower)],in2[w*lower:w*(upper-lower)]。

• 在第 74 和 76 行，Tesla 内核代码已经生成。

• 在第 81 行，为out变量调用了update self，并且数据从设备传输被限制在上下界之间：out[w*lower:w*(upper-lower)]。

• 在第 84 行，数据区域结束，并调用delete来释放在 GPU 上分配的内存。

当代码在 V100 上运行时，整个内核所花费的时间为 0.0008 秒。为了更详细地了解这一点，让我们回到分析器。这次我们将利用 NVIDIA Visual Profiler 来可视化输出：

使用 NVIDIA Visual Profiler 输出

前面的屏幕截图显示了使用async和阻塞后的 Visual Profiler 输出。来自分析器窗口的关键消息如下：

1. 我们看到有三个流被创建和使用。这是因为我们的代码使用了async(block%2)，这意味着我们请求了最大 2 个队列。第三个队列是默认队列，在管道执行期间不被使用。

2. 我们看到主机到设备和设备到主机的传输也重叠了，因为 GPU 有两个直接内存访问（DMA）引擎，因此反向的数据传输也可以重叠。

3. 我们还看到我们的内核执行与数据传输重叠。

到目前为止，我们已经看到了帮助我们将顺序代码转换为在 GPU 上运行的图像合并的关键指令。在下一节中，我们将向您介绍更多的子句，这些子句将帮助您进一步优化您的 OpenACC 代码。

其他重要的指令和子句

在本节中，我们将介绍其他重要的广泛使用的指令，可以应用到我们的合并算法中。

Gang/vector/worker

Gang/worker/vector 定义了我们可以在 OpenACC 中实现的各种并行级别。这种并行在并行化多维循环嵌套时非常有用。OpenACC 允许我们定义一个通用的 gang/worker/vector 模型，适用于各种硬件，但我们将更多地专注于 GPU 特定的实现。下图显示了 OpenACC 并行编程模型：

这个前面的图表代表了一个单一的 gang。当我们并行化我们的for循环时，循环迭代将会被均匀地分配给多个 gang。每个 gang 将包含一定数量的线程。这些线程被组织成块。一个 worker 是一行线程。

在前面的图中，有三个 worker，这意味着有三行线程。向量指的是每行有多长。所以在前面的图中，向量是八，因为每行有八个线程。在为 GPU 编程时，默认情况下会自动应用 gang 和 vector 并行。

由于 OpenACC 是一个开放标准并且面向多种硬件，它提供了通用构造。但是这个构造如何映射到特定的目标设备呢？答案很简单；这取决于架构和编译器，因此提供了性能可移植性。如果我们要映射当前 PGI 编译器如何将这个概念映射到 CUDA（NVIDIA GPU），那么它将如下所示：

• OpenACC gang 映射到 CUDA 块。

• worker 本质上映射到 CUDA 线程束。

• OpenACC 向量映射到threadIdx.x和（X 维度）。

• OpenACC worker 映射到threadIdx.y（Y 维度）。

再次强调，这是 PGI 编译器如何映射 OpenACC 构造的方式。其他编译器可能会以不同的方式进行映射。特别是对于 NVIDIA GPU，gang worker vector 将定义我们的 GPU 线程的组织。通过添加以下子句，开发人员可以告诉编译器在给定的循环上使用哪些并行级别：

• gang: 标记用于 gang 并行的循环。

• worker: 标记用于工作并行的循环。

• vector: 标记用于向量并行的循环。

以下代码片段有三个循环，并且每个循环的并行性都已经明确定义：外循环为gang，中间循环为worker循环，最内层循环为vector循环：

#pragma acc parallel loop gang
for( i = 0; i < size; i++ )
    #pragma acc loop worker
    for( j = 0; j < size; j++ )
        #pragma acc loop vector
        for( k = 0; k < size; k++ )
          c[i][j] += a[i][k] * b[k][j];

**调整 gangs、workers 和 vectors：**编译器将为您选择一定数量的 gangs 和 workers 以及向量长度，但您可以使用以下子句进行更改：

• num_gangs(N): 为并行区域生成N个 gangs

• num_workers(M): 为并行区域生成M个 workers。

• vector_length(Q): 为并行区域使用向量长度Q

在以下代码片段的示例中，我们将 gangs 的数量设置为2，workers 的数量设置为2，向量长度设置为32：

#pragma acc parallel num_gangs(2) \
  num_workers(2) vector_length(32)
{
  #pragma acc loop gang worker
  for(int x = 0; x < 4; x++){
    #pragma acc loop vector
    for(int y = 0; y < 32; y++){
      array[x][y]++;
    }
  }
}

在代码中设置 gangs 的数量很少是一个好主意——让编译器决定。大多数情况下，您可以通过调整向量长度有效地调整循环嵌套。此外，在 GPU 上很少使用 worker 循环。

托管内存

OpenACC 提供了一个选项，允许编译器处理内存管理。通过自己管理内存，我们将能够获得更好的性能；但是，允许编译器使用托管内存非常简单。我们不需要对我们的代码进行任何更改，就可以让托管内存正常工作。

为了使用托管内存，我们可以像这样将托管标志传递给pgc++编译器：

$ pgc++ -c -acc -ta=tesla:managed scrImagePgmPpmPackage.cpp
$ pgc++ -c -acc -ta=tesla:managed -Minfo=accel image_merging.cpp
$ pgc++ -o merging.out -acc -ta=tesla:managed -Minfo=accel scrImagePgmPpmPackage.o image_merging.o
$ ./blurring.out

添加了托管子句后，编译器基本上会忽略数据子句，并且托管内存用于在 CPU 和 GPU 之间传输数据。请注意，托管内存仅用于堆数据，而不是栈/静态数据。我们在上一章介绍的统一内存概念将保持不变。

内核指令

内核指令允许程序员退一步，完全依赖编译器。使用内核指令的一些示例代码如下：

#pragma acc kernels 
for (int i = 0; i < N; i++ ) 
{ 
    //< loop code > 
}

就像并行指令示例中一样，我们正在并行化一个循环。请记住，使用并行指令时，必须始终与循环指令配对；否则，代码将无法正确并行化。内核指令不遵循相同的规则；在一些编译器中，添加循环指令可能会限制编译器优化代码的能力。

内核指令是并行指令的完全相反。这意味着编译器做出了很多假设，甚至可能覆盖程序员并行化代码的决定。此外，默认情况下，编译器将尝试优化循环。编译器通常很擅长优化循环，并且有时甚至可以以程序员无法描述的方式优化循环。然而，通常程序员可以通过自己优化循环来获得更好的性能。

如果您遇到编译器拒绝并行化循环的情况，您可以覆盖编译器的决定。（但请记住，通过覆盖编译器的决定，您要对并行化代码造成的任何错误负责！）在这段代码中，我们使用独立子句来向编译器保证我们认为该循环是可以并行化的：

#pragma acc kernels loop independent
for (int i = 0; i < N; i++ )
{
    //< loop code >
}

Kernel 指令最明显的优势之一是它能够同时并行化许多循环。例如，在下面的代码段中，我们能够通过利用内核区域同时有效地并行化两个循环：

#pragma acc kernels
{
    for (int i = 0; i < N; i++ )
    {
        //< loop code >
    } 
... some other sequential code
    for (int j = 0; j < M; j++ )
    {
        //< loop code >
    }
}

Collapse 子句

collapse 子句允许我们将多维循环嵌套转换为单一维度循环。这个过程对于增加循环的整体长度（通常增加并行性）和通常有助于内存局部性。让我们看一下语法：

#pragma acc parallel loop collapse( 3 )
for(int i = 0; i < N; i++)
{
    for(int j = 0; j < M; j++)
    {
        for(int k = 0; k < Q; k++)
        {
            < loop code >
        }
    }
}

该代码将三维循环嵌套合并为单一维度循环。

Tile 子句

tile 子句允许我们将多维循环分解为瓦片或块。这通常对于增加某些代码的内存局部性很有用。让我们看一下语法：

#pragma acc parallel loop tile( 32, 32 )
for(int i = 0; i < N; i++)
{
    for(int j = 0; j < M; j++)
    {
        < loop code >
    }
}

前面的代码将我们的循环迭代分成 32 x 32 个瓦片（或块），然后并行执行这些块。

CUDA 互操作性

正如本章前面提到的，OpenACC 并不是 CUDA 语言的替代品；事实上，开发人员可以开始利用 OpenACC 将热点部分移植到 GPU 上。他们可以开始仅集成 CUDA 内核以用于最关键的功能。有几种方法可以将 OpenACC/CUDA 转换为可互操作的代码。我们将在本节中介绍其中一些。

DevicePtr 子句

这个子句可以用来映射使用cudaMalloc分配的 CUDA 设备指针，并将其传递给 OpenACC。以下代码片段展示了deviceptr子句的使用：

double *cuda_allocate(int size) {
    double *ptr;
    cudaMalloc((void**) &ptr, size * sizeof(double));
    return ptr;
}
int main() {
    double *cuda_ptr = cuda_allocate(100); 
    // Allocated on the device, but not the host!

    #pragma acc parallel loop deviceptr(cuda_ptr)
    for(int i = 0; i < 100; i++) {
        cuda_ptr[i] = 0.0;
    }
}

通常，OpenACC 运行时期望得到一个主机指针，然后将其转换为一些相关的设备指针。deviceptr子句是一种告诉 OpenACC 运行时一个给定指针不应该被转换，因为它已经是一个设备指针的方法。

Routine 指令

最后要讨论的话题是在 OpenACC 并行和内核区域内使用 CUDA 设备函数。这些函数是编译为由 GPU 内核或 OpenACC 区域调用的。为了在我们的 OpenACC 循环中使用 CUDA __device__函数，我们还可以使用 routine 指令：

//In CUDA code
extern "C" __device__
int cuda_func(int x) {
        return x*x;
}

//In OpenACC Code
#pragma acc routine seq
extern int cuda_func(int);

...

int main() {
    A = (int*) malloc(100 * sizeof(int));
    #pragma acc parallel loop copyout(A[:100])
    for(int i = 0; i < 100; i++) {
        A[i] = cuda_func(i);
    }
}

请注意，本章提供了一种实际利用 OpenACC 的方法，不涵盖整个标准 API。有关广泛的 API 信息，请参阅www.openacc.org/.

总结

在本章中，我们为您提供了一种利用 GPU 的替代方法。使用 OpenACC 的基于指令的编程方法对于传统应用程序非常受欢迎，对于新应用程序也提供了一种非常简单和可移植的方法。使用这种方法，您可以看到编译器变得更加先进。用户对指令的反馈已经被使用，通过利用指令可以为不同的架构生成最佳的并行代码。

我们介绍了提供指示/提示给编译器的并行指令。我们还利用数据指令来控制数据传输，而不是依赖于托管内存。通过使用异步子句，我们还尝试通过重叠内核和数据传输来优化我们的应用程序。我们探讨了将 OpenACC 构造映射到 CUDA 层次结构，以及 OpenACC 和 CUDA C/C++代码之间的互操作性。

在下一章中，我们将开始将我们对 CUDA 的知识应用于深度学习。

使用 CUDA 加速深度学习

深度学习是一种可以根据人工神经网络解释数据的机器学习方法。具体来说，我们提供机器可以理解的数据，并构建学习数据表示的神经网络模型。我们可以使用这种技术构建识别语音、从图像中分类对象、理解文本、翻译语言、转换数据域等模型。基本的神经网络包括全连接层（FCL）、卷积神经网络（CNN）和循环神经网络（RNN）。这些架构在数据分类、区域理解和顺序关系方面显示出强大的准确性。

深度学习需要大量计算，以便广泛应用。然而，通过使用 GPU 计算能力，我们可以显著减少训练时间，从而解决了这个问题。这是因为神经网络的基本架构是基于矩阵运算的，而 GPU 是一个针对此进行了优化的硬件平台。具体来说，深度学习的创新是通过 NVIDIA CUDA 加速来解决的，因为深度学习中的许多算法可以加速。

在本章中，我们将简要回顾神经网络操作，并讨论如何在 GPU 上加速这些操作。作为实践，我们将使用 cuDNN 和 cuBLAS CUDA 库实现一个卷积网络。cuDNN 库是 NVIDIA 的 CUDA 库，专门优化了深度学习操作。我们将在三个部分中介绍其实现。我们还将介绍 GPU 如何优化所需的操作。然后，我们将通过比较长短期记忆（LSTM）网络的性能来介绍使用 cuDNN 库的有效性。然后，我们将介绍使用NVIDIA 工具扩展（NVTX）进行深度学习的性能分析。这可以测量 GPU 上的网络操作，以便我们可以分析时间线上的操作并了解其性能。

在本章中，我们将涵盖以下主题：

• 使用 CUBLAS 加速全连接层

• 使用 cuDNN 的逐元素层

• cuDNN/CUDA 中的 Softmax 和损失函数

• 使用 cuDNN 的卷积神经网络

• 使用 CUDA 的循环神经网络

• 深度学习框架的性能分析

技术要求

本章需要安装 cuDNN 库和 CUDA 工具包。我们还需要 CUDA 启用的 GPU。本章将介绍深度学习的基础知识和性能，因此不需要新的 GPU 功能。换句话说，如果您已经涵盖了前几章的大部分内容，您将拥有一个适当的 GPU 来使用。

要安装 cuDNN 库，您需要从developer.nvidia.com/cudnn下载软件包。您需要登录 NVIDIA 开发者网站才能访问下载页面。如果您还没有帐户，您需要注册一个 NVIDIA 开发者帐户。确保 cuDNN 与您安装的 CUDA 版本编译一致。

使用 cuBLAS 加速全连接层

全连接层是深度学习的基本架构。让我们回顾一下它的操作，并看看 CUDA 如何加速神经网络的前向和反向传播过程。然后，我们将把它们应用到 GPU 上。

神经网络操作

神经网络的基本操作是在输入数据和参数之间执行点操作。我们称之为感知。在深度学习中，神经网络以分层方式连接多个感知。我们称这些为前馈神经网络。以下图表显示了一个感知和基本神经网络：

感知器的基本操作是使用输入数据和适当的权重创建点积。然后，它使用激活函数进行非线性操作，例如 sigmoid 或整流线性单元（ReLU）。在前馈神经网络中，操作只是一个仿射变换，然后是激活函数的应用。一个向量将被馈送到神经网络作为输入，并与两层中每个节点之间的权重参数相乘。

为了训练神经网络，我们进行前向传播、损失计算和梯度反向传播，然后使用更新参数。让我们简要介绍一下它们。然后，我们将使用 cuBLAS 和其他 CUDA 操作来匹配每个步骤。

前向操作可以用以下方程表示：

这里， 是给定输入向量的预测结果， 是权重参数矩阵， 是激活函数。正如我们所看到的，全连接层中的基本操作是矩阵运算。因此，我们需要对输入和激活函数实现矩阵乘法运算。因为我们进行分类任务，所以我们使用 softmax 函数来规范化输出，并在下一层获得概率分布结果。

为了获得真实值之间的损失，我们对标签应用 one-hot 编码，并通过从每个元素获得熵来获得交叉熵损失，如下所示：

我们可以通过每个交叉熵损失的总和来获得总损失值。然后，我们可以从前述方程中获得梯度。这看起来像一个复杂的操作，但可以简化如下：

现在，我们将梯度传播到前一层，这被称为反向传播。在这个任务中，我们使用链式法则来获得每个权重和偏差参数的梯度。然后，我们可以更新权重参数集和偏差。例如，我们可以通过以下方程获得权重和偏差的梯度：

我们可以通过以下方程获得梯度传播到前一层：

这里， 是激活函数的梯度。因此，我们需要从第二层获得 用于第一层。然后，可以通过以下方程获得第一层的权重和偏差的梯度：

现在，我们可以根据梯度下降规则更新权重和偏差，如下所示：

,

这里， 是迭代步骤。

激活函数的梯度可能不同，其类型也可能不同。这个激活层的实现将在下一节中介绍。激活函数的导数可以用以下方程表示：

,

因此，神经网络操作是一组线性代数操作，并且可以使用 cuBLAS 库进行覆盖。实现的代码可以在01_ann中找到。我们将在实现全连接层、实现层操作和实现 softmax 层部分介绍这些实现细节。

神经网络层的设计

在编写代码之前，让我们来看看如何将操作打包成一个层配置：

1. 首先，我们执行前向操作。

2. 然后，我们执行反向操作。

3. 然后我们从梯度中得到一个权重更新。

4. 最后，输出层将获得损失。

这样，层可以配置如下：

它具有标准化的输入和输出，以及两种类型的输入，取决于工作流程。左侧数据路径将被命名为输入，而右侧将被命名为输出。数据分为两个阶段（前向和后向）。我们将使用 blob 来管理参数和输入/输出数据。blob 是跨层处理的数据的包装器，并帮助管理内存空间。我们将使用这种设计来简化网络的配置。每个层都将有每个 blob 的描述符和前向/后向处理操作。

现在，让我们创建一个层类，它将是所有层的基类。以下代码显示了class公共函数的堆叠。而且，你可以在01_ann/src/ directory的layer.h和layer.cu中找到它的实现。这不仅有前向和后向操作，还有权重更新控制和损失计算：

class Layer
{
public:
    Layer();
    ~Layer();

    std::string get_name() { return name_; }

    virtual Blob<float> *forward(Blob<float> *input) = 0;
    virtual Blob<float> *backward(Blob<float> *grad_input) = 0;

    virtual float get_loss(Blob<float> *target);
    virtual int   get_accuracy(Blob<float> *target);

    void set_cuda_context(CudaContext *context) { cuda_ = context; }

    /* weights update control */
    void freeze() { freeze_ = true; }
    void unfreeze() { freeze_ = false;}
    void set_load_pretrain() { load_pretrain_ = true; }
    void set_gradient_stop() { gradient_stop_ = true; }

为了支持这些操作，层类维护了几个 cuDNN 描述符、blob 指针和权重更新控制器。当我们涵盖网络实现时，详细的实现将会被涵盖：

protected:
    std::string name_;

    // Tensor descriptor for the input/output tensor
    cudnnTensorDescriptor_t input_desc_;
    cudnnTensorDescriptor_t output_desc_;
    // filter and bias descriptor for weights and biases
    cudnnFilterDescriptor_t filter_desc_;
    cudnnTensorDescriptor_t bias_desc_;

    // output memory
    Blob<float> *input_ = nullptr;       /* x */
    Blob<float> *output_ = nullptr;      /* y */
    Blob<float> *grad_input_ = nullptr;  /* dx */
    Blob<float> *grad_output_ = nullptr; /* dy */

    // master weights & bias
    bool freeze_ = false;               /* control parameter updates */
    Blob<float> *weights_ = nullptr;      /* w */
    Blob<float> *biases_  = nullptr;      /* b */
    Blob<float> *grad_weights_ = nullptr; /* dw */
    Blob<float> *grad_biases_  = nullptr; /* db */

    int batch_size_ = 0; // mini-batch size

    // cuda handle container
    CudaContext *cuda_ = nullptr;

    // initialize weights along with the input size
    void init_weight_bias(unsigned int seed = 0);
    void update_weights_biases(float learning_rate);

    // pretrain parameters
    bool load_pretrain_ = false;
    int load_parameter();
    int save_parameter();

    // gradient stop tagging
    bool gradient_stop_ = false;

    friend class Network;
}

这个层类将在其他部分的深度学习网络实现中使用。因此，它具有用于 cuDNN 操作的cudnnTensorDescriptor_t变量，以及get_loss()和get_accuracy()函数。

张量和参数容器

在我们的实现中，我们将使用一个名为Blob的数据容器。它的名称是从 Caffe 借来的。这使我们能够存储张量或网络参数以及其维度大小信息和内存点。我们将使用这个来连接每一层。这有助于每一层根据输入张量的大小信息初始化其权重。此外，每一层都可以根据Blob的信息验证其结果。

这个 blob 将需要神经网络中的维度大小信息，如下一行代码所示。然后，它的构造函数将根据大小信息创建一个主机端缓冲区：

Blob<T>(int n, int c, int h, int w)

Blob还可以处理主机和设备上的内存，并帮助我们访问这些内存。Blob具有以下内存访问辅助函数：

// get specified memory pointer
ftype *ptr() { return h_ptr_; }

// get cuda memory
ftype *cuda() 
{ 
    if (d_ptr_ == nullptr) 
        cudaMalloc((void**)&d_ptr_, sizeof(ftype) * len());
    return d_ptr_;
}

// transfer data between memory
ftype *to(DeviceType target) { 
    ftype *ptr = nullptr;
    if (target == host)
    {
        cudaMemcpy(h_ptr_, cuda(), sizeof(ftype) * len(), 
                   cudaMemcpyDeviceToHost);
        ptr = h_ptr_;
    }
    else // DeviceType::cuda
    {
        cudaMemcpy(cuda(), h_ptr_, sizeof(ftype) * len(), 
                   cudaMemcpyHostToDevice);
        ptr = d_ptr_;
    }
    return ptr;
}

正如我们之前讨论的，Blob可以存储张量，我们还需要提供张量形状信息，作为 cuDNN API 所需的描述符。因此，Blob可以使用以下代码创建和设置张量描述符：

/* Tensor Control */
bool is_tensor_ = false;
cudnnTensorDescriptor_t tensor_desc_;
cudnnTensorDescriptor_t tensor()
{
    if (is_tensor_)
        return tensor_desc_;

    cudnnCreateTensorDescriptor(&tensor_desc_);
    cudnnSetTensor4dDescriptor(tensor_desc_, 
                                CUDNN_TENSOR_NCHW, CUDNN_DATA_FLOAT,
                                n_, c_, h_, w_);
    is_tensor_ = true;
    return tensor_desc_;
}

现在，让我们使用Blob来实现一个全连接层。

实现一个全连接层

在这一部分，我们将使用 cuBLAS 编写一个全连接网络。对于这个层，我们将创建一个从Layer类派生出来的Dense类。类构造函数将接收默认的层配置信息，如下所示：

Dense::Dense(std::string name, int output_size)
{
    name_ = name;
    output_size_ = output_size;
}

但这还不足以配置整个层。缺失的信息将从输入中提供，因为输入大小将由前一层确定。现在，让我们来看看前向传播。

实现前向传播

在前向传播中，我们可以将前向过程分为两个步骤，如下所示：

由于权重大小不必受批量大小的影响，我们只考虑输入权重和输出权重的数量。另一方面，数据馈送 blob，如输入和输出，受批量大小的影响。因此，我们的 GEMM 操作与过滤器和输入数据可以设计如下：

隐藏的输出将与偏置值相加。输入数据不仅限于数据加载器中的数据。当我们堆叠层时，上一层的输出将成为当前层的输入数据。前向操作可以实现如下：

Blob<float> *Dense::forward(Blob<float> *input) {
  .. { blob initialization } ..

  // output = weights^T * input (without biases)
  cublasSgemm(cuda_->cublas(),
        CUBLAS_OP_T, CUBLAS_OP_N, output_size_, 
        batch_size_, input_size_,
        &cuda_->one, weights_->cuda(), input_size_,
        input_->cuda(), input_size_,
        &cuda_->zero, output_->cuda(), output_size_);

  // output += biases * one_vec^T
  cublasSgemm(cuda_->cublas(), 
        CUBLAS_OP_N, CUBLAS_OP_N, output_size_, batch_size_, 1,
        &cuda_->one, biases_->cuda(), output_size_, one_vec, 1, 
        &cuda_->one, output_->cuda(), output_size_);
  return output_;
}

在第一次迭代中，每个层都需要初始化其权重和偏置。例如，这个Dense层可以初始化其权重、偏置和输出张量元素。我们可以将这个初始化任务分为两个阶段。第一个是权重和偏置，如下所示：

// initialize weights and biases
if (weights_ == nullptr)
{
    // setup parameter size information
    input_size_ = input->c() * input->h() * input->w();

    // initialize weight, bias, and output
    weights_ = new Blob<float>(1, 1, input_size_, output_size_);
    biases_ = new Blob<float>(1, 1, output_size_);
}

接下来的阶段是关于更新输入信息和初始化输出 blob。当它是新的或需要重新配置时，我们需要做以下工作。在这个任务中，我们还需要创建一个填满我们批量大小的向量。这将用于偏置的添加：

// initilaize input and output
if (input_ == nullptr || batch_size_ != input->n())
{
  input_ = input;
  batch_size_ = input->n();

  if (output_ == nullptr)
    output_ = new Blob<float>(batch_size_, output_size_);
  else
    output_->reset(batch_size_, output_size_);

  output_->tensor();

  if (d_one_vec != nullptr)
    cudaFree(d_one_vec);
  checkCudaErrors(cudaMalloc((void**)&d_one_vec, sizeof(float) * batch_size_));
  init_one_vec<<< (batch_size_+BLOCK_DIM_1D-1)/BLOCK_DIM_1D, BLOCK_DIM_1D >>>(d_one_vec, batch_size_);

  if (!freeze_)
    init_weight_bias();
}

这个初始化任务不仅触发了第一次迭代，还触发了批量大小的变化。在训练阶段不需要检查批量大小，但在测试阶段会很有用。这是因为训练和推断阶段的批量大小是不同的。在这种情况下，我们需要根据新的批量大小创建一个输出 blob。输出张量的大小是由通道大小确定的。以下代码创建了一个大小为（batch_size_，output_size_，1，1）的 blob：

output_ = new Blob<float>(batch_size_, output_size_);

这将创建扁平化张量。然后，我们将馈送这些张量，这要求它们在通道中对齐。这种对齐在 softmax 层中是特别需要的。我们将在 softmax 层的实现中进行讨论。

在这个阶段的另一个重要任务是初始化权重和偏置。在我们的实现中，我们将使用 ReLU 作为激活函数。我们将使用正常的初始化器（arxiv.org/abs/1502.01852）技术使网络可训练。根据前述论文的指导，所需的权重值可以用以下方程生成：

是来自上一层的输入数量。因此，我们可以在更新输入张量信息后初始化参数。此外，偏置值将被初始化为0。以下代码显示了这一实现：

void Layer::init_weight_bias(unsigned int seed)
{
    // Create random network
    std::random_device rd;
    std::mt19937 gen(seed == 0 ? rd() : static_cast<unsigned int>
                                        (seed));

    // He normal distribution
    float range = sqrt(6.f / input_->size());
    std::uniform_real_distribution<> dis(-range, range);

    for (int i = 0; i < weights_->len(); i++)
        weights_->ptr()[i] = static_cast<float>(dis(gen));
    for (int i = 0; i < biases_->len(); i++)
        biases_->ptr()[i] = 0.f;

    // copy initialized value to the device
    weights_->to(DeviceType::cuda);
    biases_->to(DeviceType::cuda);
}

现在，让我们来讨论反向传播。

实现反向传播

正如我们之前讨论的，来自下一层的梯度被传播到这一层。基于传播的梯度，我们需要获得权重、偏置和数据（输入梯度）的三个梯度。我们需要创建可以存储它们的 blob。它们的大小不取决于批量大小，所以我们只需要确保创建它们。以下代码显示了我们如何为此目的创建 blob：

if (grad_weights_ == nullptr) {
  grad_output_ = grad_output;
  grad_weights_ = new Blob<float>(weights_->shape());
  grad_biases_ = new Blob<float>(biases_->shape());
  grad_input_ = new Blob<float>(input_->shape());
}

在上述代码中，grad_output_表示从下一层传播的输出数据的梯度，grad_input_表示将传播到上一层的输入数据的梯度。因此，我们不需要创建grad_output_ blob。如果您觉得这些命名约定令人困惑，也许更容易理解grad_input_为，grad_input_为。

以下代码显示了我们如何实现这一点：

Blob<float> *Dense::backward(Blob<float> *grad_output) {
  .. { blob initialization } ..

  // db = (dy) * one_vec
  cublasSgemv(cuda_->cublas(),
    CUBLAS_OP_N,
    output_size_, batch_size_,
    &cuda_->one,
    grad_output_->cuda(), output_size_,
    one_vec, 1,
    &cuda_->zero,
    grad_biases_->cuda(), 1); 

  // dw = x * (dy)^T
  cublasSgemm(cuda_->cublas(),
    CUBLAS_OP_N, CUBLAS_OP_T,
    input_size_, output_size_, batch_size_,
    &cuda_->one,
    input_->cuda(), input_size_,
    grad_output_->cuda(), output_size_,
    &cuda_->zero,
    grad_weights_->cuda(), input_size_);

  // dx = W * dy
  if (!gradients_stop_)
    cublasSgemm(cuda_->cublas(),
      CUBLAS_OP_N, CUBLAS_OP_N,
      input_size_, batch_size_, output_size_,
      &cuda_->one,
      weights_->cuda(), input_size_,
      grad_output_->cuda(), output_size_,
      &cuda_->zero, 
      grad_input_->cuda(), input_size_);

  return grad_input_;
}

如果这一层是模型中的第一层，我们也可以跳过计算输入数据的梯度，因为我们不需要对其进行任何操作。

当我们想要更新权重时，将会更新权重和偏置值。在本节中，我们将使用随机梯度下降（SGD）来实现这一点。这个操作也可以在其他层中使用。在这里，我们将把这个函数放在Layer类中。权重更新也可以使用cublas函数来完成，如下所示：

void Layer::update_weights_biases(float learning_rate)
{
  float eps = -1.f * learning_rate;
  if (weights_ != nullptr && grad_weights_ != nullptr) {
    // w = w + eps * dw
    cublasSaxpy(cuda_->cublas(),
      weights_->len(),
      &eps,
      grad_weights_->cuda(), 1,
      weights_->cuda(), 1);
  }

  if (biases_ != nullptr && grad_biases_ != nullptr)
  {
    // b = b + eps * db
    cublasSaxpy(cuda_->cublas(),
      biases_->b(),
      &eps,
      grad_biases_->cuda(), 1,
      biases_->cuda(), 1);
  }
}

正如你所看到的，我们可以使用学习率更新权重和偏差。当然，你也可以改变eps操作以应用其他优化算法。

层终止

在 C/C++编程中，程序员应该覆盖如何在终止类实例时返回所使用的资源。根据我们的设计，如果层具有权重参数并且可以从梯度中更新它们，该层最多会创建六个 blob。以下代码显示了终止 blob 的层终止代码，这些 blob 是在内部创建的：

Layer::~Layer()
{
  if (output_ != nullptr) delete output_;
  if (grad_input_ != nullptr) delete grad_input_;

  if (weights_ != nullptr) delete weights_;
  if (biases_ != nullptr) delete biases_;
  if (grad_weights_ != nullptr) delete grad_weights_;
  if (grad_biases_ != nullptr) delete grad_biases_;
}

输入 blob 或张量描述符将由其他层或 blob 终止处理。层类是其他层的基类。因此，我们可以专注于终止自定义创建的资源，因为当我们终止任何派生层时，这个终止代码将一起被调用。

尽管我们已经设计了网络和层，但我们还应该开发一些额外的层来完成网络。例如，我们没有实现激活、softmax 和损失计算层。我们将在接下来的部分中介绍这些层。

使用 cuDNN 的激活层

神经网络层中有许多逐元素操作。激活函数是这些操作之一。cuDNN 库提供了六种激活函数：sigmoid、ReLU、tanh、clipped ReLU、ELU 和 identity。在 cuDNN 库中，cudnnActivationForward()执行前向操作，cudnnActivationBackward()执行后向操作。

让我们看一下cuddnnActivationForward()函数的接口，如下所示：

cudnnStatus_t cudnnActivationForward( cudnnHandle_t handle,
    cudnnActivationDescriptor_t activationDesc,
    const void *alpha, const cudnnTensorDescriptor_t xDesc, 
    const void *x, const void *beta,  
    const cudnnTensorDescriptor_t yDesc, void *y)

使用cudnnActivationDescriptor_t，我们可以确定激活函数的类型。Alpha 和 beta 是标量值，用于确定要添加的输入速率。xDesc和yDesc保存张量的形状信息。它们可以使用cudnnCreateTensorDescriptor()创建。

当你看cudnnActivationBackward()函数时，dy是来自下一层的梯度输入，dx是输出到上一层的梯度。在这种情况下，y变成了输入。这样，dyDesc提供了梯度输入形状信息，而dxDesc提供了梯度输出形状信息：

cudnnStatus_t cudnnActivationBackward( cudnnHandle_t handle,
    cudnnActivationDescriptor_t activationDesc,
    const void *alpha, const cudnnTensorDescriptor_t yDesc,  
    const void *y,
    const cudnnTensorDescriptor_t dyDesc, const void *dy,
    const cudnnTensorDescriptor_t xDesc,  const void *x,
    const void *beta,  const cudnnTensorDescriptor_t dxDesc, void *dx)

一般来说，我们可以期望层之间的张量形状不会改变。因此，我们可以对x和dx使用相同的张量描述符。这与使用y和dy是一样的。

现在，让我们使用 cuDNN API 实现启用 cuDNN 的激活函数。要使用 cuDNN API，我们需要提供一个张量描述符来指定输入和输出张量的维度给 cuDNN 函数。我们还需要指定激活操作。

层配置和初始化

虽然我们的示例实现没有使用层接口，但我们需要将我们的示例集成到层接口中。在我们的层设计中，激活层可以这样实现：

class Activation: public Layer
{
public:
  Activation(std::string name, cudnnActivationMode_t mode, 
             float coef = 0.f);
  ~Activation();

  Blob<float> *forward(Blob<float> *input);
  Blob<float> *backward(Blob<float> *grad_input);

private:
  cudnnActivationDescriptor_t act_desc_;
  cudnnActivationMode_t mode_;
  float coef_;
};

在初始化步骤中，我们需要创建几个张量描述符和一个激活描述符。cuDNN 库要求开发人员提供与 API 对应的张量大小或任何其他操作句柄：

Activation::Activation(std::string name, cudnnActivationMode_t mode, float coef)
{
  name_ = name;
  mode_ = mode;
  coef_ = coef;

  cudnnCreateActivationDescriptor(&act_desc_);
  cudnnSetActivationDescriptor(act_desc_, mode, CUDNN_PROPAGATE_NAN, coef);
}

在 cuDNN 中，我们使用激活描述符来指定激活函数操作。我们使用cudnnSetActivationDescriptor()函数来实现这一点。然后，它可以确定cudnnActivationForward/Backward()函数的操作。我们将在下一节中介绍这一点。然而，在这之前，我们需要实现类析构函数，以便它销毁激活描述符，如下所示：

cudnnDestroyActivationDescriptor(activation_desc);

现在，让我们介绍激活层的前向和后向操作。

实现层操作

这也被称为警告操作。这个层不需要我们处理权重和偏差，因此比密集层更容易实现。

实现前向传播

在第一次迭代中，我们需要初始化输入描述符、输出描述符和输出 blob。当批处理大小改变时，我们将更新输出 blob。然而，我们不需要初始化权重和偏差，因为它们没有。以下代码显示了它的实现：

if (input_ == nullptr || batch_size_ != input->n())
{
  input_ = input;
  input_desc_ = input->tensor();
  batch_size_ = input->n();

  if (output_ == nullptr)
    output_ = new Blob<float>(input->shape());
  else
    output_->reset(input->shape());

  output_desc_ = output_->tensor();
}

初始化后，我们使用 cuDNN 中的cudnnActivationForward()函数进行激活过程，如下所示：

cudnnActivationForward(cudnnHandle, act_desc_, 
    &one, input_desc_, d_input, &zero, output_desc_, d_output);

这个激活函数的操作是在我们初始化这个层时确定的，正如我们之前讨论的。

实现反向传播

下一步是实现反向传播。我们将重用我们已经拥有的输入/输出张量描述符。现在，我们必须初始化我们希望反向传播的梯度：

if (grad_input_ != grad_output_)
{
  grad_output_ = grad_output;
  grad_input_ = new Blob<float>(input_->shape());
  grad_input_->reset(input_->shape()); 
}

初始化后，我们可以调用cudnnActivationBackward()函数，如下所示：

cudnnActivationBackward(cudnnHandle, activation_desc, 
    &one, output_desc_, output_->cuda(), output_desc_, 
    d_grad_output, input_desc_, input_->cuda(),
    &zero, input_desc_, grad_input_->cuda());

请注意，我们重用了在前向传递中创建的输入张量描述符和输出张量描述符。我们之所以能够这样做，是因为激活操作不会改变张量的大小。我们可以通过在激活反向传播中使用 cuDNN API 来简化我们的实现。

cudnnActivationBackward()函数的输出是d_grad_input。正如我们在前一节中描述的，这个梯度将传递给下一层。

现在，我们将实现 softmax 层，并将我们的层实现集成为一个网络。然后，我们将讨论图像分类任务中全连接层的准确性。

cuDNN/CUDA 中的 softmax 和损失函数

对于 MNIST 数据集分类，我们将使用 softmax 分类器。softmax 函数对输入进行归一化，并生成概率的概率分布。softmax 操作可以表示如下：

cuDNN 的 softmax 前向函数支持此操作，以及通道和所有实例。之前，我们将密集层的输出与通道对齐。因此，我们将沿着通道应用 softmax 操作。

为了确认我们的训练有效完成，我们需要计算损失函数。由于 softmax 损失函数用于获取跨概率的损失，所以 softmax 损失函数被称为交叉熵损失。损失函数如下：

我们需要获得这个 softmax 损失的梯度以更新神经网络。幸运的是，softmax 损失的梯度在求导后很简单，如下所示：

对于前向操作，我们将使用 cuDNN 函数来获取 softmax 的输出。为了获得梯度，拥有自定义操作更直观和简单。

实现 softmax 层

现在，让我们看看如何使用 cuDNN 和 CUDA 代码来实现 softmax 层。

实现前向传播

我们可以使用 cuDNN 库中的cudnnSoftmaxForward()来获得 softmax 成本函数的输出：

cudnnSoftmaxForward(cudnnHandle, CUDNN_SOFTMAX_ACCURATE, 
      CUDNN_SOFTMAX_MODE_CHANNEL,
      &one,  input_desc,  d_input, &zero, output_desc, d_output);

在这种情况下使用的最重要的参数设置之一是CUDNN_SOFTMAX_MODE_CHANNEL。此选项使得在输入张量描述符信息后面进行通道级别的 softmax 操作。通过这样做，我们可以提供已经通过密集层的小批量输入按通道对齐的张量。

实现反向传播

softmax 层的反向传递与其他层的实现不同。这个操作将输入数据的标签作为输入，并获得适当的梯度。正如我们之前讨论的，softmax 损失的梯度可以使用以下方程获得：

我们可以使用cublasSaxpy()来实现这个操作，如下所示：

// set grad_input_ as predict
cudaMemcpyAsync(grad_input_->cuda(), output_->cuda(), 
                output_->buf_size(), cudaMemcpyDeviceToDevice));
// set grad_input_ = predict - target 
cublasSaxpy(cuda_->cublas(), target->len(), &cuda_->minus_one,
            target->cuda(), 1, grad_input_->cuda(), 1));

在前面的代码中，目标 blob 包含了 one-hot 编码的目标向量，因此将负目标向量添加到预测值中会产生适当的梯度。之后，我们需要在传播到前一层之前对批次梯度进行归一化，如下所示：

int grad_output_size = target->n() * target->c() * target->h() * target->w();
float scale = 1.0f / static_cast<float>(target->n());
cublasSscal(cuda_->cublas(), grad_output_size, &scale, grad_input_->cuda(), 1);

由于这引入了加权和的均值，我们可以期望每个批次的梯度被归一化。

实现损失函数

计算 softmax 的损失值是可选的。这意味着它的值在训练和推断中不被考虑。然而，我们可以将其用作训练的指标。

如我们之前讨论的，softmax 损失函数应该实现以下方程：

我们可以通过一个核函数从每个样本的输出中获得损失并累积它们，如下所示：

__global__ void
softmax_loss_kernel(float *reduced_loss, float *predict, 
                    float *target, int size)
{
  int batch_idx = blockDim.x * blockIdx.x + threadIdx.x;

  extern __shared__ float s_data[];
  float loss = 0.f;

  // each thread calculate entropy for each data 
  // and accumulate to shared memory
  if (batch_idx > 0)
    return;

  for (int c = 0; c < num_outputs; c++)
    loss += target[batch_idx * num_outputs + c] * \
                logf(predict[batch_idx * num_outputs + c]);
                workspace[batch_idx] = -loss;

  // Then, we do reduction the result to calculate loss 
  // Using 1 thread block
  if (blockIdx.x > 0) return;

  // Cumulate workspace data
  s_data[threadIdx.x] = 0.f;
  for (int i = 0; i < batch_size; i += blockDim.x)
    s_data[threadIdx.x] += workspace[threadIdx.x + i];

  __syncthreads();

  // Reduction
  for (unsigned int stride = blockDim.x / 2; stride > 0; stride >>= 1)
  {
    if (threadIdx.x + stride < batch_size)
      s_data[threadIdx.x] += s_data[threadIdx.x + stride];
    __syncthreads();
  }

  if (threadIdx.x == 0)
    reduced_loss[blockIdx.x] = s_data[0];
}

这个操作使用并行归约，在第三章 CUDA 线程编程中介绍过，用于获取一个批次中的累积损失值。由于我们只会使用这个减少的损失值来确认训练，所以我们只会监视它的输出而不是取平均值。

现在，让我们将我们实现的所有层与 MNIST 数据集加载器集成在一起。

MNIST 数据加载器

整个过程中一个重要的部分是为特定数据集创建一个数据加载器。在这个实验室中，我们将使用包含 60,000 个样本的 MNIST 数据集。在初始化时，我们告诉数据加载器它应该加载训练集还是测试集。之后，数据加载器将加载数据集中的一些魔术数字，以及所有样本和它们的标签。加载的数据将被存储在向量中，并使用相同的随机种子进行洗牌。由于数据加载器构建和洗牌样本向量，训练循环或测试循环可能会在每次迭代时获得随机化的输入数据。完整的实现代码可以在本书的 GitHub 存储库中的src/mnist.cpp文件中找到。

管理和创建模型

当我们有多个层时，我们需要一个可以管理这些层的对象，进行神经网络操作，即前向/后向传播和权重更新。在这个实验室中，我们将有一个层的数组，并迭代数组进行前向处理。例如，前向操作可以用以下代码执行：

Blob<float> *Network::forward(Blob<float> *input) {
  output_ = input;
  for (auto layer : layers_)
    output_ = layer->forward(output_);

  return output_;
}

反向传播也可以通过以相反顺序迭代数组来完成：

void Network::backward(Blob<float> *target) {
  Blob<float> *gradient = target;
  // back propagation.. update weights internally.....
  for (auto layer = layers_.rbegin(); layer != layers_.rend(); layer++) {
    // getting back propagation status with gradient size
    gradient = (*layer)->backward(gradient);
  }
}

如您所见，我们在向量中管理层，并具有每个层的操作。将新层添加到网络中甚至更简单，如下面的代码所示：

void Network::add_layer(Layer *layer) {
  layers_.push_back(layer);
}

通过使用Network类，我们可以使用各种模型管理函数，如参数更新，层注册，层初始化等。此外，我们可以构建一个像现代深度学习框架一样的神经网络。例如，我们可以创建一个模型如下：

// step 1\. loading dataset
MNIST data_loader = MNIST("./dataset");
// create training dataset loader and shuffling the data
data_loader.train(batch_size, true);  

// step 2\. model initialization
Network model;
model.add_layer(new Dense("dense1", 500));  // 1st layer
model.add_layer(new Dense("dense2", 10));   // 2nd layer
model.cuda();     // set cuda context for each layer

我们还可以有以下训练循环：

// get data sample's shared buffer
Blob<float> *train_data   = data_loader.get_data();   
// get target's shared buffer
Blob<float> *train_target = data_loader.get_target(); 
// load data and targets with the batch size
data_loader.get_batch();    
tp_count = 0;  step = 0;
while (step < num_steps)
{
  // transfer loaded data to the GPU
  train_data->to(cuda);
  train_target->to(cuda);

  model.forward(train_data);    // forward
  model.backward(train_target); // backward
  learning_rate *= 1.f / (1.f + lr_decay * step);
  model.update(learning_rate);  // update

  step = data_loader.next(true); // load next data

  ... monitoring logic ...
}

对于测试阶段，我们为测试数据集创建另一个数据集加载器，并只进行前向传播的迭代。以下代码显示了它的实现：

test_data_loader.test(batch_size_test);                   // create test dataset loader
Blob<float> *test_data = test_data_loader.get_data();     // get sample data shared buffer
Blob<float> *test_target = test_data_loader.get_target(); // get target shared buffer
test_data_loader.get_batch();    // load samples and targets with the batch size
tp_count = 0; step = 0;
while (step < num_steps_test) {
  // transfer loaded data to the GPU
  test_data->to(cuda);
  test_target->to(cuda);

  model.forward(test_data);  // forward
  tp_count += model.get_accuracy(test_target);

  step = test_data_loader.next(); // load next data
}
float accuracy = 100.f * tp_count / num_steps_test / batch_size_test;

在测试阶段，我们将在完成对测试数据集中所有样本的测试后获得准确率。现在，我们需要在测试循环之后获得准确率。

使用 MNIST 数据集进行网络训练

现在，让我们运行我们实现的代码并查看其结果。对于训练阶段，我们将迭代 2,400 步，批量大小为 256。MNIST 数据集在训练集中有 60,000 个样本。2,400 步意味着我们将进行大约 10 个 epochs 的迭代。样本代码可以用以下命令编译：

$ nvcc -run -m64 -std=c++11 -I/usr/local/cuda/samples/common/inc -gencode arch=compute_70,code=sm_70 -lcublas -lcudnn -lnvToolsExt -o train ./train.cpp ./src/layer.cu ./src/loss.cu ./src/mnist.cpp ./src/network.cpp

以下截图显示了我们实现的训练和测试输出：

在训练迭代中，网络从训练数据集中获得了 92%的准确率。然而，测试准确率只有 77%，这与训练结果相比是一个相对较低的分数。推断显示训练和推断之间的准确率差距很大可能有很多原因。一个可能的原因是全连接层没有考虑到前面截图中显示的区域信息。在深度学习中，我们使用卷积层来使网络学习空间信息。

现在，让我们使用 cuDNN 实现卷积层，将其添加到网络中，并比较模型的性能。

使用 cuDNN 的卷积神经网络

cuDNN 库为卷积操作提供了优化的性能。通过创建一个卷积层，我们将覆盖 API 的配置，用于前向和后向操作。

卷积网络层对输入数据进行卷积处理。当你想要构建一个了解区域信息的神经网络时，这种网络架构是很有用的。回想一下，在第七章中的卷积实现，CUDA 中的并行编程模式，它需要相当大的内存带宽，并需要进一步优化以获得最佳性能。然而，使用 cuDNN 库，我们也可以获得最佳性能，因为我们不必重新发明轮子。

卷积层的实现与全连接层的实现类似。然而，由于 cuDNN 库的存在，有两个不同之处：我们不必像以前那样完全实现那么多细节，我们需要为操作分配一个工作空间大小。对于每个卷积操作——前向、反向滤波器和反向输入——都需要额外的内存空间，取决于它们的算法。算法可以根据给定的输入/输出/滤波器张量维度而变化。详细的 API 调用将在稍后处理。

与其他层一样，它有三个工作阶段。对于推理阶段，我们将调用cudnnConvolutionForward()和cudnnAddTensor()。对于反向阶段，我们将调用cudnnConvolutionBackwardData()、cudnnConvolutionBackwardFilter()和cudnnConvolutionBackwardBias()。最后，对于更新阶段，我们可以重用全连接层的代码。该层的配置概述如下：

实现前向传播

在深度学习神经网络中，通常会与卷积网络一起使用池化层。池化层只是根据简单的规则选择输入数据进行输出。以下图示显示了最大池化的例子：

使用 cuDNN 库，我们将实现这两个卷积操作。

卷积层

与全连接层类似，这个卷积层有权重和偏置参数。在全连接层中，我们使用了 cuBLAS，它不需要 cuDNN 相关的描述符。然而，我们将使用 cuDNN 卷积函数，因此需要使用滤波器描述符和卷积操作描述符。以下代码显示了在构建层时应该初始化的资源：

Conv2D::Conv2D(std::string name,
        int out_channels, kernel_size, stride, padding, dilation):
        out_channels_(out_channels), kernel_size_(kernel_size),
        stride_(stride), padding_(padding), dilation_(dilation) {
  name_ = name;
  cudnnCreateFilterDescriptor(&filter_desc_);
  cudnnCreateConvolutionDescriptor(&conv_desc_);
  cudnnSetConvolution2dDescriptor(conv_desc_,
    padding_, padding_, stride_, stride_, dilation_,dilation_,
    CUDNN_CROSS_CORRELATION, CUDNN_DATA_FLOAT);
}

由于我们在模型构建时提供了卷积操作信息，我们可以指定卷积描述符。然而，滤波器的操作可以在推断时指定，因为我们可以在那时学习输入张量的大小。现在，让我们实现卷积层的前向传递。

正如我们之前讨论的，我们可以用输入张量大小初始化卷积层。这个输入张量大小会影响输出张量的大小。以下代码显示了前向传递中的参数初始化步骤：

// initialize weights and bias
if (weights_ == nullptr) {
  // initialize containers handles
  cudnnSetFilter4dDescriptor(filter_desc_, 
    CUDNN_DATA_FLOAT, CUDNN_TENSOR_NCHW,
    out_channels_, input->c(), kernel_size_, kernel_size_);

  weights_ = new Blob<float>(out_channels_, input->c(), kernel_size_, kernel_size_);
  biases_ = new Blob<float>(1, out_channels_); // bias size
  bias_desc_ = biases_->tensor();
}

然后，我们需要更新输入资源，初始化输出 blob，创建 cuDNN 工作空间，并初始化权重参数，如下所示：

// initilaize input and output
if (input_ == nullptr || batch_size_ != input->n()) {
  // initialize input
  input_ = input;
  input_desc_ = input->tensor();
  batch_size_ = input->n();

  // getting output tensor size
  cudnnGetConvolution2dForwardOutputDim(
    conv_desc_, input_desc_, filter_desc_,
    &output_size_[0], &output_size_[1], 
    &output_size_[2], &output_size_[3]);

  // initialize output blob
  if (output_ == nullptr)
    output_ = new Blob<float>(output_size_);
  else
    output_->reset(output_size_);
  output_desc_ = output_->tensor();

  // initialize weights
  if (!freeze_)
    init_weight_bias();

  // initialize workspace for cudnn
  set_workspace();
}

为了获得输出张量大小，我们使用cudnnGetConvolution2dForwardOutputDim()函数。该函数根据输入张量大小、卷积操作和滤波器大小输出维度大小信息。然后，我们重用了在全连接层中使用的相同参数初始化代码。

要调用 cuDNN 的卷积 API，我们需要提供其工作算法和工作空间内存。我们这样做是因为 cuDNN 根据卷积大小选择最佳卷积算法，并且需要立即进行测量。确定算法后，cuDNN 可以确定工作空间大小。卷积层需要进行前向传播的卷积操作、输入数据的梯度和权重的梯度。我们需要分别处理每个算法，但我们可以分配一个工作空间，因为工作空间专门用于每个卷积操作。

因此，我们创建的工作空间需要具有每个卷积算法所需的最大大小。以下代码显示了我们如何使用它们并管理工作空间：

Conv2d::set_workspace() {
  size_t temp_size = 0;

  // fwd
  cudnnGetConvolutionForwardAlgorithm(cuda_->cudnn(),
    input_desc_, filter_desc_, conv_desc_, output_desc_,
    CUDNN_CONVOLUTION_FWD_PREFER_FASTEST, 0, &conv_fwd_algo_);
  cudnnGetConvolutionForwardWorkspaceSize(cuda_->cudnn(),
    input_desc_, filter_desc_, conv_desc_, output_desc_, 
    conv_fwd_algo_, &temp_size);
  workspace_size = std::max(workspace_size, temp_size);

  // bwd - data
  cudnnGetConvolutionBackwardDataAlgorithm(cuda_->cudnn(), 
    filter_desc_, output_desc_, conv_desc_, input_desc_, 
    CUDNN_CONVOLUTION_BWD_DATA_PREFER_FASTEST, 0, 
    &conv_bwd_data_algo_);
  cudnnGetConvolutionBackwardDataWorkspaceSize(cuda_->cudnn(),
    filter_desc_, output_desc_, conv_desc_, input_desc_, 
    conv_bwd_data_algo_, &temp_size);
  workspace_size = std::max(workspace_size, temp_size);

  // bwd - filter
  cudnnGetConvolutionBackwardFilterAlgorithm(cuda_->cudnn(),
    input_desc_, output_desc_, conv_desc_, filter_desc_,
    CUDNN_CONVOLUTION_BWD_FILTER_PREFER_FASTEST, 0, 
    &conv_bwd_filter_algo_);
  cudnnGetConvolutionBackwardFilterWorkspaceSize(cuda_->cudnn(),
    input_desc_, output_desc_, conv_desc_, filter_desc_, 
    conv_bwd_filter_algo_, &temp_size);
  workspace_size = std::max(workspace_size, temp_size);

  if (workspace_size > 0) {
    if (d_workspace != nullptr)
      cudaFree(d_workspace);
    cudaMalloc((void**)&d_workspace, workspace_size);
  }
}

每个卷积算法都使用单独的类型进行指定，即cudnnConvolutionFwdAlgo_t、cudnnConvolutionBwdDataAlgo_t和cudnnConvolutionBwdFilterAlgo_t。我们可以通过将它们声明为类成员变量来使用它们，即conv_fwd_algo_、conv_bwd_data_algo_和conv_bwd_filter_algo_。

现在，在初始化后，我们编写前向处理代码。我们使用滤波器进行卷积并添加偏差。以下代码显示了 cuDNN 卷积前向实现：

cudnnConvolutionForward(cuda_->cudnn(), &cuda_->one, input_desc_, input_->cuda(), \
    filter_desc_, weights_->cuda(), conv_desc_, conv_fwd_algo_, d_workspace, workspace_size, \
    &cuda_->zero, output_desc_, output_->cuda());
cudnnAddTensor(cuda_->cudnn(), &cuda_->one, bias_desc_, biases_->cuda(), \
    &cuda_->one, output_desc_, output_->cuda());

卷积的结果将使用输出 blob 传递到下一层。

实现反向传播

在反向传播中，我们应该计算偏差的梯度、权重的梯度和输入数据的梯度。为此，我们需要在第一次迭代中创建 blob 以便我们可以存储它们。它们的大小不取决于批处理大小，所以我们只需要确保它们被创建。初始化步骤可以实现如下：

// initialize grad_output back-propagation space
if (grad_weights_ == nullptr) {
  grad_output_  = grad_output;
  grad_weights_ = new Blob<float>(weights_->shape());
  grad_biases_  = new Blob<float>(1, biases_->c());
  grad_input_   = new Blob<float>(input_->shape());
}

然后，我们调用 cuDNN 反向卷积 API，如下所示：

Blob<float> *Conv2D::backward(Blob<float> *grad_output) {
  ... { initialization step } ...

  // gradients of biases
  cudnnConvolutionBackwardBias(cuda_->cudnn(),
    &cuda_->one, 
    output_desc_, grad_output->cuda(),
    &cuda_->zero, 
    bias_desc_, grad_biases_->cuda());

  // gradients of weights 
  cudnnConvolutionBackwardFilter(cuda_->cudnn(),
    &cuda_->one, 
    input_desc_, input_->cuda(), 
    output_desc_, grad_output_->cuda(),
    conv_desc_, conv_bwd_filter_algo_, d_workspace, workspace_size,
    &cuda_->zero, 
    filter_desc_, grad_weights_->cuda());

  // gradients of input data
  if (!gradient_stop_)
    cudnnConvolutionBackwardData(cuda_->cudnn(),
      &cuda_->one, 
      filter_desc_, weights_->cuda(), 
      output_desc_, grad_output->cuda(), 
      conv_desc_, conv_bwd_data_algo_, d_workspace, workspace_size,
      &cuda_->zero, 
      input_desc_, grad_input_->cuda());

然后，我们将输入数据的梯度传递给前一层以传播梯度。在更新步骤中，我们将使用基类的梯度更新代码来更新权重和偏差的梯度。在全连接层中实现反向传播时，我们已经涵盖了这一点。如果这是第一层，则我们也可以跳过计算输入数据的梯度。

使用 cuDNN 的池化层

池化层有两个特点。首先，它的输出大小与卷积层不同，cuDNN 为此提供了相应的 API。其次，它没有任何内部权重。

为了指定池化操作，我们可以使用 cuDNN 的cudnnPoolingDescriptor_t函数，并在类构造函数中创建和指定 cuDNN 的池化描述符，如下所示：

cudnnCreatePoolingDescriptor(&pool_desc_);
cudnnSetPooling2dDescriptor(pool_desc_, mode_, CUDNN_PROPAGATE_NAN,
  kernel_size_, kernel_size_, padding_, padding_, stride_, stride_);

现在，让我们实现池化层的前向和反向操作。

实现前向传播

池化层有助于减小张量的大小。因此，我们需要计算输出大小。我们可以使用cudnnGetPooling2dForwardOutputDim()函数来计算大小，就像我们在卷积层实现中所做的那样。此外，张量大小取决于批处理大小。这意味着如果批处理大小发生变化，我们需要更新张量大小。以下代码显示了我们如何初始化输入和输出 blob：

if (input_ == nullptr || batch_size_ != input->n()) {
  input_ = input;

  // resource initialize
  input_desc_ = input_->tensor();
  batch_size_ = input->n();

  // setting output
  cudnnGetPooling2dForwardOutputDim(pool_desc_, input_desc_, 
    &output_size_[0], &output_size_[1], &output_size_[2], 
    &output_size_[3]);
  if (output_ == nullptr)
    output_ = new Blob<float>(output_size_);
  else
    output_->reset(output_size_);

  output_desc_ = output_->tensor();
}

对于前向传播，我们调用cudnnPoolingForward()函数，如下所示：

Blob<float> *Pooling::forward(Blob<float> *input) {
  ... { initialization step } ...

  cudnnPoolingForward(cudnnHandle, pool_desc_, &one, 
    input_desc_, input_->cuda(),
    &zero, output_desc_, output_->cuda());
}

实现反向传播

对于反向传播步骤，我们调用cudnnPoolingBackward()函数，如下所示：

Blob<float> *Pooling::backward(Blob<float> *grad_output) {
  if (grad_input_ == nullptr)
    grad_input_ = new Blob<float>(input_->shape());

  cudnnPoolingBackward(cudnnHandle, pool_desc_,
    &one, output_desc_, output_->cuda(), 
    output_desc_, grad_output->cuda(), 
    input_desc_, input_->cuda(), 
    &zero, input_desc_, grad_input_->cuda());
}

池化层的张量形状的输入和梯度的输入是相同的，输出和梯度的输出的形状也是相同的。因此，我们可以分别重用输入和输出的张量描述符。

现在，让我们将这些集成到单个卷积层实现中。

网络配置

现在，我们将更新我们之前的网络 LeNet。网络代码可以编写如下：

Network model;
model.add_layer(new Conv2D("conv1", 20, 5));
model.add_layer(new Pooling("pool", 2, 0, 2, CUDNN_POOLING_MAX));
model.add_layer(new Conv2D("conv2", 50, 5));
model.add_layer(new Pooling("pool", 2, 0, 2, CUDNN_POOLING_MAX));
model.add_layer(new Dense("dense1", 500));
model.add_layer(new Activation("relu", CUDNN_ACTIVATION_RELU));
model.add_layer(new Dense("dense2", 10));
model.add_layer(new Softmax("softmax"));
model.cuda();

现在，我们可以开始训练和推断阶段，因为我们已经配置了我们的层，使它们彼此连接。让我们使用以下命令编译代码：

$ nvcc -run -m64 -std=c++11 -I/usr/local/cuda/samples/common/inc -gencode arch=compute_70,code=sm_70 -lcublas -lcudnn -lnvToolsExt -o train ./train.cpp ./src/layer.cu ./src/loss.cu ./src/mnist.cpp ./src/network.cpp

然后，我们可以看到训练和测试结果如下：

正如您所看到的，该网络的训练准确度和推断准确度都比仅使用全连接网络时要高。我们还可以通过查看 NVIDIA 配置文件来确认其操作，如下所示：

混合精度操作

最新的 NVIDIA GPU 支持深度学习的混合精度操作。我们不会在本书中涵盖这一点，因为它超出了我们的范围。但是，如果您希望了解更多，可以访问 NVIDIA 提供的示例，位于/usr/src/cudnn_samples_v7/conv_sample。要访问此示例，您需要从 cuDNN 网页下载示例。此示例代码显示了如何使用 cuDNN 库进行混合精度操作。

为了使 cuDNN API 与张量核心一起工作，我们需要设置数学类型，如下所示：

cudnnSetConvolutionMathType(cudnnConvDesc, CUDNN_TENSOR_OP_MATH);

然后，我们需要使用cudnnSetTensorNdDescriptor()初始化输入和输出张量的张量描述符。这为张量提供填充，以便我们获得优化的张量核心性能。

一个很好的基于 cuDNN 的实现是cudnn-training：github.com/tbennun/cudnn-training。它将 LeNet 实现为一系列 cuDNN 函数。您可以跟踪每一行，看看 CUDNN 函数是如何工作的。

如果您有兴趣使用 cuDNN 部署您的网络，请查看以下关于 GTC-CNN 推断与 cuDNN 的视频（developer.nvidia.com/gtc/2019/video/S9644/video）。这个讲座介绍了使用 cuDNN 进行 CNN 推断的有用性能优化技巧。

在深度学习训练中使用半精度需要超过 FP16 操作的利用率。我们需要在 FP16 中计算张量，同时将权重保持在 FP32 中。此外，一些操作需要 FP32。我们称之为混合精度。cuDNN 库提供了一个名为 mnistCUDNN 的混合精度推断示例。该示例显示了输入和层数据类型的转换。如果您想了解更多关于深度学习和训练中混合精度操作的信息，请阅读以下文章：devblogs.nvidia.com/video-mixed-precision-techniques-tensor-cores-deep-learning/。

现在，我们将从性能方面讨论深度学习中的其他 GPU 使用注意事项。

循环神经网络优化

RRN 允许您在深度学习中分析顺序数据。尽管该网络具有顺序依赖性，但仍有大量的优化空间。在本节中，我们将介绍其算法以及 cuDNN 如何提供优化性能。

有许多种类型的 RNN，但 cuDNN 只支持四种，即带有 ReLU 的 RNN，带有 tanh 的 RNN，LSTM 和 GRU。它们有两个输入：来自先前网络的隐藏参数和来自源的输入。根据它们的类型，它们有不同的操作。在本实验室中，我们将介绍 LSTM 操作。下图显示了 LSTM 的前向操作：

从计算的角度来看，有八个矩阵-矩阵乘法和许多逐元素操作。根据这个估计，我们可以期望 LSTM 可能是内存受限的，因为每个操作都是内存受限的。另一方面，CUDNN 提供了cudnnRNNForwardInference()和cudnnRNNFowardTraining()RNN 函数。我们将通过测量这个函数的性能和模拟 LSTM 的性能来介绍使用这个函数的好处。为了做到这一点，我们将实现一个虚拟的 LSTM 层，并将其性能与 cuDNN LSTM 函数进行比较。

为了测试目的，我们将设置超参数如下：

int mode = 2; // LSTM in CUDNN
int seq_length = 512;
int num_layers = 4;
int hidden_size = 512;
int input_size = hidden_size;
int batch_size = 32;
float dropout_rate = 0;
bool bidirectional = 0;
int persistent = 0;

序列长度或隐藏大小可能会有所不同，这取决于问题。在这个测试中，我们将使用512作为长度，在序列研究中经常使用。CUDNN API 需要更多的选项才能工作，比如 dropout 率、双向或单向以及持久 RNN。在本节中，我们只测试 vanilla LSTM。

使用 CUDNN LSTM 操作

让我们编写一些执行cudnnRNNForwardTraining()函数作为 LSTM 层的代码：

1. 我们需要初始化输入和输出内存空间。为了执行 cuDNN 的 RNN API，我们需要使用以下变量：

// hx, cx, hy, cy, dhy, dcy, dhx, and dcs can be null.
void *x;            // input
void *hx = nullptr; // input of initial hidden state
void *cx = nullptr; // input of cell state (LSTM)

void *y;            // output
void *hy = nullptr; // output of final hidden state
void *cy = nullptr; // output of final cell state (LSTM)

void *dy;            // input of gradient 
void *dhy = nullptr; // input of final hidden state
void *dcy = nullptr; // input of final cell state (LSTM)

void *dx;            // output of gradient at the input of rnn
void *dhx = nullptr; // output of gradient at the initial hidden state
void *dcx = nullptr; // output of gradient at the initial cell state

这些变量是 LSTM 的输入和输出。为了提供输入和获取输出，我们需要分配适当的内存空间。根据 LSTM 的定义，我们需要考虑输入、输出和隐藏层的长度。这些大小可以确定如下：

int input_length = seq_length * input_size * batch_size;
int output_length = seq_length * hidden_size * batch_size;
int hidden_length = hidden_size * batch_size * num_layers;

然后，我们可以为每个项目分配内存。

1. 现在，我们需要为 cuDNN RNN API 设置张量描述符。以下代码显示了我们应该设置的所需张量描述符：

cudnnTensorDescriptor_t x_desc[seq_length], y_desc[seq_length], \
                        dx_desc[seq_length], dy_desc[seq_length];
cudnnTensorDescriptor_t hx_desc, cx_desc;
cudnnTensorDescriptor_t dhx_desc, dcx_desc;
cudnnTensorDescriptor_t hy_desc, cy_desc;
cudnnTensorDescriptor_t dhy_desc, dcy_desc;

对于输入和输出描述符，我们需要初始化每个元素，即批量大小和其输入大小。其他隐藏的张量描述符是用层数、批量大小和隐藏大小进行初始化的。本节不涵盖如何编写初始化代码。但是，如果您想了解更多信息，可以查看10_deep_learning/03_rnn文件中的代码。

1. 我们还需要为 RNN 操作提供一个工作空间，就像我们为卷积操作做的那样：

void *workspace;
cudnnFilterDescriptor_t w_desc, dw_desc;
cudnnSetRNNDescriptor_v6(cudnnHandle, rnn_desc,
                         hidden_size, num_layers, dropout_desc, CUDNN_LINEAR_INPUT,
                         bidirectional ? CUDNN_BIDIRECTIONAL : CUDNN_UNIDIRECTIONAL,
                         CUDNN_LSTM, CUDNN_RNN_ALGO_STANDARD, CUDNN_DATA_FLOAT));
size_t weight_size;
cudnnGetRNNParamsSize(cudnnHandle, rnn_desc, x_desc[0], &weight_size, CUDNN_DATA_FLOAT);
cudaMalloc((void**)&workspace, weight_size);

然后，我们可以根据工作空间的大小设置滤波器描述符，如下所示：

dimW = {weight_size / sizeof(float), 1, 1}
cudnnCreateFilterDescriptor(&w_desc);
cudnnCreateFilterDescriptor(&dw_desc);
cudnnSetFilterNdDescriptor(w_desc, CUDNN_DATA_FLOAT, CUDNN_TENSOR_NCHW, 3, dimW);
cudnnSetFilterNdDescriptor(dw_desc, CUDNN_DATA_FLOAT, CUDNN_TENSOR_NCHW, 3, dimW);
cudnnRNNForwardTraining(cudnnHandle, rnn_desc, seq_length,
                x_desc, x, hx_desc, hx, cx_desc, cx,
                w_desc, w, 
                y_desc, y, hy_desc, hy, cy_desc, cy,
                workspace, workspace_size, reserved_space, 
                reserved_size);

我们可以使用cudaEvnetRecoard()和 flops 计算来衡量它们的性能。例如，前向操作可以配置为以下方程：

然后，我们将通过将批量大小从 32 增加到 256 来测试我们的实现，每次增加 32。适用的测试范围可能会有所不同，以及 GPU 的内存大小。

在本节中，我们实现了基于 LSTM 的模拟和cudnnRNNForwardTraining()调用。我们部分模拟的版本只有 GEMM 操作，这是最计算密集的操作。现在，让我们比较这些实现的性能。

实现虚拟 LSTM 操作

在我们的实现中，我们将专注于模拟 LSTM 的主要操作，而不是完全实现它。

让我们确定 LSTM 网络的超参数。一般来说，输入序列长度范围从 512 到 2,048。层数的数量是不确定的。但是，由于tanh操作，它不能太大。对于输入大小，我们将使用 512。通常情况下，批量大小在 RNN 使用方面在 32 到 256 之间。CUDNN 需要更多关于 dropout 率、双向或单向以及是否使用持久 RNN 的输入。我们现在不使用它们。我们的 LSTM 配置信息如下：

现在，我们将部分实现 LSTM 操作以测量计算强度。正如我们之前讨论的，LSTM 有两个矩阵-矩阵乘法需要计算。LSTM 操作将为输入序列的每个元素以及每个层计算。然后，操作可以配置如下：

for (int layer = 0; layer < num_layers; layer++) {
  for (int linear_layer = 0; linear_layer < 4; linear_layer++) {
    for (int sequence = 0; sequence < seq_length; sequence++) {
      cublasSgemm(cublas_handle, CUBLAS_OP_T, CUBLAS_OP_N,
      hidden_size, input_size, batch_size,
      &alpha, input_weight, input_size, x, input_size,
      &beta, h, hidden_size);
      cublasSgemm(cublas_handle, CUBLAS_OP_T, CUBLAS_OP_N,
      hidden_size, hidden_size, batch_size,
      &alpha, recurrent_weight, hidden_size,
      h, hidden_size,
      &beta, y, hidden_size);
    }
  }
}

我们可以使用更多的逐元素操作，但这只是近似计算强度，所以我们暂时不考虑它们。

比较 CUDNN 和 SGEMM LSTM 的性能

让我们比较它们的性能以及不同的批处理大小，如下所示的代码实现在main()函数中：

for (int step = 1; step <= 8; step++)
{
 batch_size = 32 * step;
 printf("Batch Size: %3d\n", batch_size);
 rnn_operation(seq_length, num_layers, hidden_size, input_size,   
   batch_size, dropout_rate, bidirectional, mode, persistent);
 cublas_operation(mode, 2ull, input_size, hidden_size, seq_length, batch_size, num_layers);
}

然后，我们可以使用以下命令编译和执行示例源代码：

$ nvcc -run -m64 -std=c++11 -I/usr/local/cuda/samples/common/inc -gencode arch=compute_70,code=sm_70 -lcublas -lcudnn -lcurand -o rnn ./rnn.cpp

以下图表显示了来自 Tesla V100 卡的 cuBLAS 和 cuDNN 的性能：

在上图中，两种实现在性能上有很大差异。cuDNN 的 LSTM 性能比使用 cuBLAS 模拟的 LSTM 要好得多。此外，LSTM 操作的性能遵循 Tesla V100 GPU 的屋顶线。另一方面，两个 SGEMM 操作并没有显示出这种性能，因为矩阵大小不够大以获得完整的性能。要从 Tesla V100 获得 10 TFlops，矩阵大小应与 1,024 的平方相似或更大。然而，正如我们所看到的，我们的矩阵大小大约是 512 的平方。

LSTM 优化在以下 NVIDIA 文章中有解释：devblogs.nvidia.com/optimizing-recurrent-neural-networks-cudnn-5。它结合了矩阵-矩阵乘法，融合逐元素操作，多个流和多层并行化。

RNN 的优化版本之一是持久 RNN（svail.github.io/persistent_rnns），由 Greg Diamos 介绍。尽管他的实现不包括 LSTM 和 GRU，但您可以了解 RNN 如何进行优化。

深度学习框架的性能分析

一般来说，我们使用 TensorFlow、PyTorch 和 MxNet 等深度学习框架开发和研究神经网络。由于这些框架，我们可以有效地开发复杂的模型。然而，当涉及性能工程时，由于性能分析工具的能力，理解框架下 GPU 操作是一个陡峭的学习曲线。例如，使用 Chrome 跟踪进行性能分析在模型简单时很有用，但在模型复杂时就不那么有用。

在第五章中，CUDA 应用程序性能分析和调试，我们介绍了NVIDIA 工具扩展（NVTX），它允许我们在 GPU 应用程序中进行自定义注释，并使用 NVIDIA Nsight Systems 查看时间轴。对于复杂的应用程序，程序员分析其性能并找到瓶颈非常有用。

在本节中，我们将介绍如何通过修改 ResNet-50 示例代码在 PyTorch 和 TensorFlow 中使用 NVTX。示例代码可以在本书的 GitHub 存储库的10_deep_learining/05_framework_profile文件夹中找到。您可以从github.com/nvidia/DeepLearningExamples获取原始源代码。

为了简化工作环境配置，我们将使用NVIDIA GPU 云（NGC）深度学习容器用于 PyTorch 和 TensorFlow。如果您需要了解 NGC 或容器的基本用法，请访问本书附录中的 NGC。

现在，让我们先从 PyTorch 开始。

对 PyTorch 模型进行性能分析

在 PyTorch 中，我们可以使用torch.cuda.nvtx.range_push("foo")和torch.cuda.nvtx.range_pop()来放置自定义标签。这保持了原始的 CUDA NVTX API，即nvtxRangePush()和nvtxRangePop()。让我们看看 NVTX 注释如何帮助我们在时间轴上理解深度学习操作。在接下来的步骤中，我们将使用05_framework_profile/pytorch/RN50v1.5文件中的 ResNet-50 示例代码：

1. 我们将在train()函数中的训练循环中放置 NVTX 注释以注释step值。该函数可以在image_classificaiton/training.py文件中找到。以下截图显示了训练循环和分别在第 234 行和第 260 行的 NVTX 注释：

在上述代码中，训练操作是在step函数中实现的，该函数由get_train_step()函数定义。因此，我们需要在该函数中放置 NVTX 注释以了解更多信息。

1. 让我们在第 164 行的get_train_step()函数中添加一些 NVTX 注释。该函数返回_step()函数，其中包括训练操作。因此，我们将在该函数中放置 NVTX 注释。训练过程包括前向和反向传播、全局归约和优化（更新权重）。以下截图显示了在第 166 行和第 171 行的前向传播的注释：

通过这种方式，我们可以在其余操作上放置其他注释。

1. 我们还可以为模型层添加 NVTX 注释。在这个例子中，ResNet-50 模型是在image_classification/resnet.py文件中实现的。以下截图显示了网络的示例注释：

正如我们所看到的，我们可以按照 ResNet 架构放置 NVTX 注释。如果我们在每个构建块中放置注释，我们可以获得更多信息。

1. 现在，让我们对模型进行分析。正如我们之前讨论的，我们将使用 NGC 深度学习容器，即 PyTorch。imagenet数据集位于/raid/datasets/imagenet/raw-data文件夹中。为了限制分析时间范围，我们将使用延迟选项（-y）和持续时间选项（-d）。以下代码显示了一个执行容器并对网络进行分析的 bash shell 脚本：

#/bin/bash

CODE_PATH="RN50v1.5"
DATASET_PATH="/raid/datasets/imagenet/raw-data/"
OUTPUT_NAME="resnet50_pyt"

# default profile
docker run --rm -ti --runtime=nvidia \
    -v $(pwd)/${CODE_PATH}:/workspace \
    -v ${DATASET_PATH}:/imagenet \
    nvcr.io/nvidia/pytorch:19.08-py3 \
       nsys profile -t cuda,nvtx,cudnn,cublas -o ${OUTPUT_NAME} 
         -f true -w true -y 60 -d 20 \
       python /workspace/main.py --arch resnet50 -b 64 
         --fp16 /imagenet

执行后，上述代码将在 RN50v1.5 目录中生成 profiled 结果，即resnet50_pyt.qdrep。

1. 最后，使用 NVIDIA Nsight Systems 打开 profiled 输出resnet50_pyt.qdrep，并查看操作。以下截图显示了带有 NVTX 注释的测量步骤：

在这里，我们可以看到反向操作所花费的时间是前向操作的两倍。此外，PyTorch 将主机线程分开用于训练循环和反向传播。从内核分析来看，耗时最长的点是逐元素的内核执行。让我们扩大前向传递以查看层的执行时间，如下截图所示：

在这里，我们可以看到第二个卷积块需要最长的时间来完成。如果这一层存在效率低下的点，我们可以进一步挖掘。如果某个操作被确定为瓶颈并需要优化，我们还可以使用 NVIDIA Nsight Compute 来分析特定的内核函数。比较主机 API 跟踪和 GPU，我们可以看到时间持续时间是不同的。这是因为主机和 GPU 操作是异步的。因此，当我们从主机测量 GPU 执行时间时，我们需要谨慎。现在，让我们看一下优化步骤，如下截图所示：

我们可以看到，从主机和 GPU 的测量执行时间中存在巨大差异。主机的测量执行时间为 25.367 毫秒，而 GPU 的时间为 4.048 毫秒。其操作主要是逐元素操作，其执行被延迟直到反向传播完成。我们还可以找到异步执行。之后，我们可以看到cudaDeviceSynchronize()函数，该函数防止当前步骤被下一步骤更新。

我们还可以通过设置环境来禁用这些异步操作，即CUDA_LAUNCH_BLOCKING=1。我们可以使用环境选项（-e）将其传递给 Nsight System 的配置选项。然后，我们可以分析应用程序的align操作与主机和内核函数。

PyTorch 在其 CUDA 对象中具有几个具有 NVTX 特色的 API。 PyTorch 文档可以在pytorch.org/docs/stable/_modules/torch/cuda/nvtx.html找到。通过直接在 PyTorch 中调用 NVTX API，将调用 CUDA NVTX API。这意味着我们可以在分析时间线中获得自定义标记的 NVTX 标记。

对 TensorFlow 模型进行分析

对 TensorFlow 图进行分析需要使用启用 NVTX 注释的 NVTX 插件。要在 TensorFlow 中使用 NVTX 注释，我们需要使用以下命令安装nvtx-plugins-tf Python 插件：

$ pip install nvtx-plugins-tf

但是，如果我们使用的是版本晚于 19.08 的 NGC TensorFlow 容器，则无需执行此操作。

TensorFlow 图形 API 是符号 API，因此它们需要特定的编程方法。 NVTX 插件为此提供了两个选项：装饰器和 Python 函数。

以下是 NVTX 装饰器的示例：

import nvtx.plugins.tf as nvtx_tf
ENABLE_NVTX=true
@nvtx_tf.ops.trace(message='Dense Block', domain_name='Forward',
        grad_domain_name='Gradient', enabled=ENABLE_NVTX, 
        trainable=True)
def dense_layer(x):
    x = tf.layers.dense(x, 1000, activation=tf.nn.relu, name='dense_1')
    x = tf.layers.dense(x, 1000, activation=tf.nn.relu, name='dense_2’) 
return x

以下是 NVTX Python 函数的示例：

import nvtx.plugins.tf as nvtx_tf
ENABLE_NVTX=true
x, nvtx_context = nvtx_tf.ops.start(x, message='Dense Block', \ 
        domain_name='Forward’, grad_domain_name='Gradient’, 
        enabled=ENABLE_NVTX, trainable=True)
x = tf.layers.dense(x, 1000, activation=tf.nn.relu, name='dense_1')
x = tf.layers.dense(x, 1000, activation=tf.nn.relu, name='dense_2’) 
x = nvtx_tf.ops.end(x, nvtx_context)

NVTX 插件提供了 NVTXHook，它允许我们对 TF 估算器和会话进行分析。例如，我们可以按以下方式使用该钩子：

from nvtx.plugins.tf.estimator import NVTXHook

nvtx_callback = NVTXHook(skip_n_steps=1, name='Train’)
training_hooks=[]
training_hooks.append(nvtx_callback)

然后，我们可以使用以下代码将其应用于任一选项：

with tf.train.MonitoredSession(hooks=training_hooks) as sess:

或者，我们可以使用以下代码：

tf.estimator.Estimator(hooks=training_hooks, ...)

现在，让我们将其应用到示例 ResNet-50 代码中并进行操作审查。示例代码可以在05_framework_profile/tensorflow/RN50v1.5文件夹中找到：

1. 让我们首先将NVTXHook应用于估算器。训练图的定义可以在runtime/runner.py文件的第 312 行找到。在构建图之前，我们将NVTXHook附加到钩子列表中，如下面的代码块所示：

1. 然后，我们将 NVTX 注释应用于模型构建函数。model_build()函数可以在model/resnet_v1_5.py文件的ResnetModel类中找到。以下代码显示了如何在model_build()函数中的conv1层上使用 Python 函数放置 NVTX 注释的示例：

在上述代码中，当使用nvtx_tf.ops.start()和nvtx_tf.ops.end()函数时，我们需要谨慎选择适当的输入和输出。只在其他层中放置 NVTX 注释。确保最终的全连接层输出是网络的输出。

我们还必须禁用用于检查可训练变量数量的代码。如果 NVTX 的trainable参数值为True，则大小会发生变化。在resnet_v1_5.py文件的第 174 行，有一段断言代码，用于检查该变量的数量。只需将其注释掉，如下所示：

1. 我们还使用 NVTX 装饰器来构建 ResNet 模块。在model/blocks目录中，我们可以在conv2d_blocks.py和resnet_bottleneck_block.py中找到conv2d和 ResNet 瓶颈块的实现。在conv2d_blocks.py文件中，我们可以装饰conv2d_block()函数以注释 NVTX 分析，如下所示：

同样，我们也可以对resnet_bottleneck_block.py文件执行相同操作：

1. 现在，让我们对模型进行性能分析。就像我们使用 PyTorch 容器一样，我们将使用 TensorFlow 的 NGC 容器。我们假设imagenet数据集的tfrecord文件位于/raid/datasets/imagenet/tfrecord目录中。以下代码显示了一个执行容器并对网络进行性能分析的 bash shell 脚本：

#/bin/bash

CODE_PATH="RN50v1.5"
DATASET_PATH="/raid/datasets/imagenet/tfrecord"
OUTPUT_NAME="resnet50_tf"

# default profile
docker run --rm -ti --runtime=nvidia \
    -v $(pwd):/result \
    -v $(pwd)/${CODE_PATH}:/workspace \
    -v ${DATASET_PATH}:/imagenet \
    nvcr.io/nvidia/tensorflow:19.08-py3 \
        nsys profile -t cuda,nvtx,cudnn,cublas -o ${OUTPUT_NAME} 
                     -f true -w true -y 40 -d 20 \
            python /workspace/main.py --mode=training_benchmark 
                                      --warmup_steps 200 \
                --num_iter 500 --iter_unit batch 
                --results_dir=results --batch_size 64

当我们执行这个函数时，我们将在RN50v1.5目录中得到resnet50_tf.qdrep文件。

1. 最后，让我们使用 NVIDIA Nsight System 审查分析输出：

在这里，我们可以确认反向传播所花费的时间是前向传播的两倍。这个示例代码与 CPU 和 GPU 不同步。因此，我们可以看到主机和 GPU 之间的时间差异更大。当我们在构建块中放置额外的注释时，我们将能够在层中看到子块的注释。

使用 NVIDIA Nsight Systems 进行性能分析在多 GPU 训练中监视所有归约操作的执行时间时提供了额外的好处。以下截图显示了一个使用两个 GPU 进行训练的 GPU 的性能分析结果：

在突出显示的行中，我们可以看到ncclAllRecude()函数，它同时调用了反向传播。通过这样做，我们不会延迟所有归约操作。这个示例代码使用 Horovod 来训练多个 GPU。如果你想了解更多，请访问 Horovod 的 GitHub 页面：github.com/horovod/horovod。你可以从这里获取文档和示例代码。

总结

在本章中，我们学习了如何使用 CUDA 库进行深度学习和性能优势。在回顾它们的用途时，我们将它们与每个步骤的深度学习机制进行匹配。由于我们可以使用的深度学习库，我们可以实现一个简单的 CNN，而不必实现算法。然后，我们使用 NVTX 注释在 PyTorch 和 TensorFlow 中对 ResNet-50 模型进行了性能分析。

对于一些深度学习工程师和研究人员来说，实现基本算法可能是不切实际的。然而，了解性能因素和基本操作可以帮助您构建高效和有效的基于深度学习的产品。如今，我们看到许多产品化的基于深度学习的服务。工程师们花费大量资源将他们训练好的模型产品化，以及训练他们的模型，以便获得尽可能低的错误率。希望您能够了解如何在深度学习应用中使用 NVTX 性能分析。利用这些知识，您可以更好地利用您的 GPU。祝你好运！