OpenACC Fortran 循环中的 Matmul答案

【问题标题】：Matmul in OpenACC Fortran loopOpenACC Fortran 循环中的 Matmul
【发布时间】：2021-04-24 16:37:35
【问题描述】：

使用 PGI 编译器通过 OpenACC 加速 Fortran 代码，我在加速循环中遇到了 matmul 调用问题。

在简化示例中，我将单位矩阵应用于两个向量，因此输入和输出值应该相同：

program test
        implicit none
        integer :: a(3, 3)
        integer :: v1(3, 2), v2(3, 2)
        integer :: i

        a = reshape([1, 0, 0, 0, 1, 0, 0, 0, 1], [3, 3])
        v1 = reshape([1, 2, 3, 4, 5, 6], [3, 2])

        print *, v1

        !$acc kernels copyin(a, v1) copyout(v2)
        !$acc loop independent
        do i = 1, 2
                v2(:, i) = matmul(a, v1(:, i))
        enddo
        !$acc end kernels

        print *, v2
endprogram

使用 PGI 编译器版本 20.9 进行编译时，我得到了以下信息：

test:
     12, Generating copyin(a(:,:),v1(:,:)) [if not already present]
         Generating implicit copyout(z_a_0(:)) [if not already present]
         Generating copyout(v2(:,:)) [if not already present]
     14, Loop is parallelizable
         Generating Tesla code
         14, !$acc loop gang ! blockidx%x
         15, !$acc loop vector(32) ! threadidx%x
     15, Loop is parallelizable

运行代码给出以下值：

1 2 3 4 5 6
4 5 6 4 5 6

第二行应该和第一行一样，这是顺序执行的情况。代码有什么问题？

【问题讨论】：

看起来两个线程都被送入v1的第二列

标签： fortran openacc pgi-accelerator

【解决方案1】：

@Mat Colgrove 解释了错误行为的原因。我发现的解决方法是显式编写矩阵向量乘法：

program test
        implicit none
        integer :: a(3, 3)
        integer :: v1(3, 2), v2(3, 2)
        integer :: i, j, k

        a = reshape([1, 0, 0, 0, 1, 0, 0, 0, 1], [3, 3])
        v1 = reshape([1, 2, 3, 4, 5, 6], [3, 2])

        print *, v1

        !$acc kernels copyin(a, v1) copyout(v2)
        !$acc loop independent
        do i = 1, 2
                !$acc loop seq
                do k = 1, 3
                        v2(k, i) = 0
                        !$acc loop seq
                        do j = 1, 3
                                v2(k, i) = v2(k, i) + a(j, k) * v1(j, i)
                        enddo
                enddo
        enddo
        !$acc end kernels

        print *, v2
endprogram

【讨论】：

【解决方案2】：

看起来是编译器问题。问题线是：

Generating implicit copyout(z_a_0(:))

“z_a_0”是正在创建的编译器临时数组，用于保存调用 matmul 的中间结果。它的声明被提升出循环，然后作为共享数组复制回来。由于它是共享的，因此会导致竞争条件。

我已提交问题报告 (TPR #29482) 并将其发送给我们的工程师进行进一步调查。

【讨论】：