Fortran多维子数组性能

Question

在Fortran90中操作和分配多维数组中的子数组时，我偶然发现了一个有趣的性能问题。

Fortran90引入了操作数组子部分的能力，我看到一些地方建议使用这种"切片"方法来执行数组操作，而不是循环。例如，如果我必须添加两个数组，a和b，大小为10，那么最好这样写:

c(1:10) = a(1:10) + b(1:10)

或

c = a + b

不是

do i = 1, 10
    c(i) = a(i) + b(i)
end do

我对简单的一维和二维数组尝试了这种方法，发现使用"切片"表示法更快。然而，当在多维数组中分配这样的结果时，事情开始变得有点有趣了。

首先，我必须为我的相当粗糙的成绩评估练习道歉。我甚至不确定我所采用的方法是否是计时和测试代码的正确方法，但是我对测试的定性结果相当有信心。

program main
    implicit none
    integer, parameter :: mSize = 10000
    integer :: i, j
    integer :: pCnt, nCnt, cntRt, cntMx
    integer, dimension(mSize, mSize) :: a, b
    integer, dimension(mSize, mSize, 3) :: c
    pCnt = 0
    call SYSTEM_CLOCK(nCnt, cntRt, cntMx)
    print *, "First call: ", nCnt-pCnt
    pCnt = nCnt
    do j = 1, mSize
        do i = 1, mSize
            a(i, j) = i*j
            b(i, j) = i+j
        end do
    end do
    call SYSTEM_CLOCK(nCnt, cntRt, cntMx)
    print *, "Created Matrices: ", nCnt-pCnt
    pCnt = nCnt
    ! OPERATIONS BY SLICING NOTATION
    !c(1:mSize, 1:mSize, 1) = a + b
    !c(1:mSize, 1:mSize, 2) = a - b
    !c(1:mSize, 1:mSize, 3) = a * b
    ! OPERATIONS WITH LOOP
    do j = 1, mSize
        do i = 1, mSize
            c(i, j, 1) = a(i, j) + b(i, j)
            c(i, j, 2) = a(i, j) - b(i, j)
            c(i, j, 3) = a(i, j) * b(i, j)
        end do
    end do
    call SYSTEM_CLOCK(nCnt, cntRt, cntMx)
    print *, "Added Matrices: ", nCnt-pCnt
    pCnt = nCnt
end program main

可以看到，我有两种方法对两个大的2D数组进行操作并将其赋值给一个3D数组。我非常喜欢使用切片表示法，因为它帮助我编写更短、更优雅的代码。但是在观察到我的代码是多么的缓慢之后，我不得不重新检查切片符号在循环内计算的能力。

我使用GNU Fortran 4.8.4 for Ubuntu 14.04运行了上面的代码，有和没有-O3标志

1. 无-O3标志

   。切片符号

        5 Runs - 843, 842, 842, 841, 859
        Average - 845.4
   

   b。循环计算

        5 Runs - 1713, 1713, 1723, 1711, 1713
        Average - 1714.6
   

2. 带-O3标志

   。切片符号

        5 Runs - 545, 545, 544, 544, 548
        Average - 545.2
   

   b。循环计算

        5 Runs - 479, 477, 475, 472, 472
        Average - 475
   

我发现非常有趣的是，没有-O3标志，切片符号继续比循环执行得更好。然而，使用-O3标志会使这个优势完全消失。相反，在这种情况下使用数组切片表示法是有害的。

事实上，对于我相当大的3D并行计算代码，这是一个重要的瓶颈。我强烈怀疑，在将低维数组赋值到高维数组的过程中，数组临时的形成是罪魁祸首。但是为什么优化标志在这种情况下不能优化赋值呢?

此外，我觉得责备-O3标志不是一件值得尊敬的事情。那么，数组临时变量真的是罪魁祸首吗?我还遗漏了什么吗?任何见解都将对加快我的代码非常有帮助。谢谢!

Answer 1

在进行任何性能比较时，您必须将苹果与苹果进行比较，将橙子与橙子进行比较。我的意思是你并不是在比较同一件事。即使它们产生相同的结果，它们也是完全不同的。

这里起作用的是内存管理，考虑操作期间的缓存故障。如果你按照haraldkl的建议把循环版本变成3个不同的循环，你肯定会得到类似的性能。

发生的情况是，当你在同一个循环中组合3个赋值时，右侧有很多缓存重用，因为所有3个赋值在右侧共享相同的变量。对于循环版本，a或b的每个元素只被加载到缓存和寄存器中一次，而对于数组操作版本，a或b的每个元素被加载3次。这就是区别所在。数组的大小越大，差异越大，因为您将获得更多的缓存错误和更多的元素重新加载到寄存器中。

Answer 2

我不知道编译器真正做什么所以不是真正的答案，但太多的文本评论…我怀疑编译器将数组表示法扩展为这样的内容:

do j = 1, mSize
    do i = 1, mSize
        c(i, j, 1) = a(i, j) + b(i, j)
    end do
end do
do j = 1, mSize
    do i = 1, mSize
        c(i, j, 2) = a(i, j) - b(i, j)
    end do
end do
do j = 1, mSize
    do i = 1, mSize
        c(i, j, 3) = a(i, j) * b(i, j)
    end do
end do

当然，如果这样写的话，编译器可能仍然会折叠这些循环，所以你可能需要更多地迷惑他，例如，在循环之间向屏幕上写一些c。