《OpenACC并行编程实战》—— 3.5 计算构件parallel

本节书摘来自华章出版社《OpenACC并行编程实战》一 书中的第3章，第3.5节，作者何沧平，更多章节内容可以访问云栖社区“华章计算机”公众号查看。

3.5　计算构件parallel

OpenACC中的计算构件有两个，一个是前面介绍的kernels构件，一个就是这里要介绍的parallel构件。两个计算构件的作用都是将循环并行化，但有一些重要区别。本节将结合一些例子详细对比介绍。
parallel这个基本构件开启加速器设备上的并行执行。

在C和C++中，parallel导语的语法是：

#pragma acc parallel [子语列表] 换行
  结构块

在Fortran中的语法是：

!$acc parallel [子语列表]
  结构块
!$acc end parallel 

这里的子语是下列中的一个:

async[( 整数表达式 )]
wait[(整数表达式列表 )]
num_gangs(整数表达式)
num_workers(整数表达式 )
vector_length(整数表达式 )
device_type(设备类型列表)
if(条件)
reduction(操作符:变量列表)
copy(变量列表)
copyin(变量列表)
copyout(变量列表)
create(变量列表)
present(变量列表)
deviceptr(变量列表)
private(变量列表)
firstprivate(变量列表)
default(none|present)

遇到一个加速器parallel构件的时候，程序就创建一个或多个gang来执行这个加速器并行区域。gang的数量、每个gang中的worker数量和每个worker中的向量通道数量，在这个并行区域的存续期内保持为常数。每个gang都开始以gang冗余模式执行结构块。这意味着，在并行区域之内，带loop导语且以gang层次工作分摊的循环之外的代码，将被所有的gang冗余地执行。
在每个gang中都有一个worker开始执行本构件的结构块中的代码。

如果没有使用async子句，加速器parallel区域出口处将会有一个隐式障碍，此时本地线程不再向前执行，直到所有的gang都已到达parallel区域的出口。
假设计算构件里面引用一个变量，该变量既没有预置数据属性也没有出现在计算构件的数据子语之中，还没有出现在data构件、可见的declare导语之中，如果构件上没有default（none）子语，那么编译器将隐式地确定该变量的数据属性。构件上没有default（present）子语的时候，如果这个parallel构件内引用的数组或者聚合数据类型变量既没有出现在本构件的数据子语中，也没有出现在包围本构件的data构件中，那么它们将被视为出现在parallel构件的copy子语中。如果构件上有一个default（present）子语，那么对那些没有预置数据类型的所有数组和聚合数据类型的变量，编译器会将它们看做好像出在一个present子语之中。如果一个标量变量在parallel构件中被引用，并且没有出现在该parallel构件的数据子语中，也没有包含在任何data构件中，那么它等同于出现在parallel构件的firstprivate子语中。

使用限制：

• 一个程序不能通过分支转入或跳出parallel构件。
• 一个程序决不能依赖于子语的求值顺序或求值的副作用。
• 只有async、wait、num_gangs、num_workers和vector_length这些子语可以跟在一个device_type子语后面。
• 至多出现一个if子语。Fortran中，条件的运算结果必须是一个标量逻辑值；C和C++中，条件的运算结果必须是一个标量整数值。
• 至多只能出现一个default子语，而且必须带有一个值，要么是none要么是present。

数据子语copy、copyin、copyout、create、present、deviceptr、firstprivate和private子语将在第4章中讲述。

3.5.1　gang单独模式

与kernels导语一样，parallel导语也会自动探测结构块代码内的并行性，存在数据依赖时拒绝并行化。与kernels导语不同的地方是，在没有loop导语的情况，parallel构件只有一个gang来并行化。例3.30使用parallel导语来开启2个并行区域，第11～13行的并行区域做向量相加，第14～16行的并行区域平移数组元素。Fortran 版例3.31操作相同。

例3.30的编译反馈（此处略去了数据管理信息）如下：

         12, #pragma acc loop vector(128) /* threadIdx.x */
     12, Loop is parallelizable
     15, Loop carried dependence of b prevents parallelization
         Loop carried backward dependence of b prevents vectorization

从中可知，第12行的循环成功并行化，而第15行的循环因为有数据依赖不能并行，只能串行执行。
例3.30的运行时反馈如下：

launch CUDA kernel  file=C:\cygwin64\home\he\accbook\parallel\p1c.c function=main line=11 device=0 threadid=1 num_gangs=1 num_workers=1 vector_length=128 grid=1 block=128
launch CUDA kernel  file=C:\cygwin64\home\he\accbook\parallel\p1c.c function=main line=14 device=0 threadid=1 num_gangs=1 num_workers=1 vector_length=1 grid=1 block=1

从中可知，第11行的循环只使用了1个gang，gang中包含一个长度为128的向量。对长度为1024的数组，循环体在这1个gang上并行运行。第14行的循环只能串行执行，所以只有1个gang、1个worker，而且1个vector的长度为1。也就是说，只有1个线程来执行整个循环。

3.5.2　gang分裂模式

为了使用多个gang来并行执行parallel构件中的循环，就必须使用loop导语。这里的loop导语默认自己的关联循环是数据独立的，即各个迭代步之间没有数据依赖，可以安全地并行化。如果确实有数据依赖，那么程序员要承担错误后果。
作为对比，此处再提一次，kernels构件会要求编译器自行检测loop导语的关联循环的数据依赖性，无法确定数据独立的时候就拒绝并行化。
例3.32和例3.33中的每个parallel区域内均有一个loop导语，编译器会认为紧接着的关联循环中无数据依赖。

例3.32的编译反馈如下：

         13, #pragma acc loop gang, vector(128) /* blockIdx.x threadIdx.x */
         17, #pragma acc loop gang, vector(128) /* blockIdx.x threadIdx.x */

由此可知，例3.32的第17行与例3.30的第15行完全相同，循环中有数据依赖，例3.30中没有loop导语，编译器自动检测到不可并行化。而例3.32中采用了loop导语，隐含independent子语，由程序员保证无数据依赖。而实际上，第17～18行的循环中有数据依赖，强制并行化将导致错误结果，图3.10直观地解释了错误原因，不再重述。
例3.32的运行时反馈如下：

a[N-1] = 2046
b[2] = 1   错误，正确的结果是b[2] = 0
launch CUDA kernel  file=C:\cygwin64\home\he\accbook\parallel\p2c.c function=main line=11 device=0 threadid=1 num_gangs=8 num_workers=1 vector_length=128 grid=8 block=128
launch CUDA kernel  file=C:\cygwin64\home\he\accbook\parallel\p2c.c function=main line=15 device=0 threadid=1 num_gangs=8 num_workers=1 vector_length=128 grid=8 block=128

由此可知，两个parallel构件各生成1个CUDA内核，运行在8个gang上，实现了多个gang并行计算，不再是只有1个gang。
无论parallel构件里包含多少个循环，编译器都将该构件中的所有循环编译成仅仅一个CUDA内核。因此，无论各个循环的迭代步数是否相同，它们对应的CUDA内核都将采用相同的执行参数。如果不同循环的迭代步数差别太大，性能可能无法达到最优，请看例3.34和例3.35。

例3.34的编译反馈如下：

         14, #pragma acc loop gang, vector(128) /* blockIdx.x threadIdx.x */
         17, #pragma acc loop gang, vector(128) /* blockIdx.x threadIdx.x */

例3.34的运行时反馈如下：

launch CUDA kernel  file=C:\cygwin64\home\he\accbook\parallel\p3c.c function=main line=11 device=0 threadid=1 num_gangs=5 num_workers=1 vector_length=128 grid=5 block=128

可以看出，两个循环的并行化方式完全一样，运行时果然只启动1个CUDA内核。执行配置为：（CUDA C术语）网格中包含5个线程块，每个线程块里有128个线程；（OpenACC术语）一共启用5个gang，每个gang里只包含1个worker，worker中的向量长度为128。

3.5.3　二重循环

parallel构件不会自动并行化本区域内的循环，必须使用loop导语显式地指定将哪些循环并行化，以及如何并行化。例3.36和例3.37将两个二维数组的对应元素相加，结果保存在第三个数组中，需要用到二重循环。第10行的parallel导语开启一个并行区域，区域范围是第10～16行。第11行上的loop导语作用到第12行的for循环上，第14行的内层for循环上没有放loop导语。

例3.36的编译反馈如下：

         12, #pragma acc loop gang /* blockIdx.x */
         14, #pragma acc loop vector(128) /* threadIdx.x */

编译器不但将第12行的循环并行化了，竟然将第14行的循环也并行化了，但并行线程是按一维组织的。运行时反馈信息中更详细地给出执行配置：
num_gangs=32 num_workers=1 vector_length=128 grid=32 block=128
线程网格包含32个线程块，每个线程块包含128个一维线程。作为对比，例3.9中用kernels导语对相同的二重循环并行化，将二重循环映射成了二维线程网格，需要注意这个区别。
例3.36的第12～16行的循环，还可以用其他方式加导语，例如：

  #pragma acc parallel
  #pragma acc loop
  for(i=0; i < M; i++)
  {
    #pragma acc loop
    for(j=0; j<N; j++)
      a[i][j] = b[i][j] + c[i][j];
  }

和

  #pragma acc parallel
  #pragma acc loop gang
  for(i=0; i < M; i++)
  {
    #pragma acc loop gang, vector 
    for(j=0; j<N; j++)
      a[i][j] = b[i][j] + c[i][j];
  }

它们的编译反馈和运行时反馈与例3.36一致。导语

#pragma acc loop gang, vector 

的含义是要求使用gang和vector来并行化内层循环。其实还可以要求用多少个gang、多长的vector来并行化，但建议不要这么做，交由编译器选择最好的执行配置，通常比程序员手动选的好，而且可移植性更好。
Fortran版例3.37的二重循环也可以用其他方式加导语，例如：

  !$acc parallel
  !$acc loop
  do j = 1, N
    !$acc loop
    do i = 1, M
      a(i,j) = b(i,j) + c(i,j)
    enddo
  enddo
  !$acc end parallel

3.5.4　三重循环

很多行业中用到大量的三重循环，例如天气预报中需要对地理坐标、时间点进行循环。循环次数越多，计算量越大，并行化加速潜力越大，请看例3.38和例3.39。

例3.38的编译反馈如下：

         15, #pragma acc loop gang /* blockIdx.x */
         17, #pragma acc loop vector(128) /* threadIdx.x */

例3.38的运行时反馈如下：

num_gangs=32 num_workers=1 vector_length=128 grid=32 block=128

例3.38中用三重循环实现三维数组加操作，第14行的loop导语作用到第15行的for循环上，但编译器竟然将第15、17行循环并行化，跳过了第16行的循环，而且并行线程按照一维组织。挺奇怪，不明白编译器是怎么想的。作为对比，例3.11中，kernels循环将三维循环中的两个内层循环映射成二维线程网格，一定要注意这个区别。
例3.38中第12～19行还有其他的并行化方法，效果相同：

  #pragma acc parallel
  {
    #pragma acc loop
    for(i=0; i<L; i++)
    { 
      #pragma acc loop
      for(j=0; j<M; j++)
      for(k=0; k<N; k++)
        a[i][j][k] = b[i][j][k] + c[i][j][k];
    }
  }

和

  #pragma acc parallel
  {
    #pragma acc loop
    for(i=0; i<L; i++)
    { 
      #pragma acc loop
      for(j=0; j<M; j++)
      {
        #pragma acc loop
        for(k=0; k<N; k++)
          a[i][j][k] = b[i][j][k] + c[i][j][k];
      }
    }
  }

Fortran版例3.39中三重循环的并行方法还有：

!$acc parallel
!$acc loop
do k = 1, N
  !$acc loop
  do j = 1, M
  do i = 1, L
    a(i,j,k) = b(i,j,k) + c(i,j,k);
  enddo
  enddo
enddo
!$acc end parallel

或者

!$acc parallel
!$acc loop
do k = 1, N
  !$acc loop
  do j = 1, M
    !$acc loop
    do i = 1, L
      a(i,j,k) = b(i,j,k) + c(i,j,k);
    enddo
  enddo
enddo
!$acc end parallel