为什么 Fortran 内部函数 "spread" 通常比显式迭代慢

Why is the Fortran intrinsic function "spread" often slower than explicit iteration

我使用地球物理模型,常见情况是需要将 2D 数据与 3D 数据相乘、相加等。下面是一个例子。

module benchmarks
  implicit none
  integer, parameter :: n=500
  integer :: k
  real :: d2(n,n)
  real :: d3(n,n,n)
  contains
  ! Iteration
  subroutine benchmark_a(res)
    real, intent(out) :: res(n,n,n)
    do k = 1, size(d3,3)
      res(:,:,k) = d2*d3(:,:,k)
    end do
  end subroutine
  ! Spread
  subroutine benchmark_b(res)
    real, intent(out) :: res(n,n,n)
    res = d3*spread(d2, 3, size(d3,3))
  end subroutine
end module

program main
  use benchmarks
  real :: t, tarray(2)
  real :: res(n,n,n)
  call random_number(d2)
  call random_number(d3)
  ! Iteration
  call dtime(tarray, t)
  call benchmark_a(res)
  call dtime(tarray, t)
  write(*,*) 'Iteration', t
  ! Spread
  call dtime(tarray, t)
  call benchmark_b(res)
  call dtime(tarray, t)
  write(*,*) 'Spread', t
end program

当我 运行 具有不同维度大小 n 时,我通常发现 spread 慢得多;例如:

Spread   2.09942889
Iteration  0.458283991

有谁知道为什么 spread 方法比显式 for 循环(我认为通常要不惜一切代价避免)要慢得多?

这里的基本答案是"it isn't"。也许对于特定的编译器和特定的环境,内在函数不如显式 DO 循环优化得那么好,但它不一定是那样。我使用 ifort 19 进行了测试,即使在默认优化级别下,SPREAD 内在循环和显式循环生成了相似的代码,当我更正程序以使用结果时内在循环更快。

Iteration 0.2187500 0.1376885 Spread 9.3750000E-02 0.1376885

我还要提醒(就像我在对您的问题的评论中所做的那样)简单的基准程序通常不会衡量作者认为他们所做的事情。您的原始示例和修改后的示例都出现了最常见的错误,即从未使用过测试工作的结果,因此足够聪明的编译器可以简单地蒸发整个操作。事实上,当我使用 ifort 19 构建您的两个测试用例时,编译器完全删除了所有工作,只留下计时代码。不用说,运行速度非常快。

  implicit none
  integer, parameter :: n=500
  integer :: k
  real :: d2(n,n)
  real :: d3(n,n,n)
  contains
  ! Iteration
  subroutine benchmark_a(res)
    real, intent(out) :: res(n,n,n)
    do k = 1, size(d3,3)
      res(:,:,k) = d2*d3(:,:,k)
    end do
  end subroutine
  ! Spread
  subroutine benchmark_b(res)
    real, intent(out) :: res(n,n,n)
    res = d3*spread(d2, 3, size(d3,3))
  end subroutine
end module

program main
  use benchmarks
  real :: tstart,tend
  real :: res(n,n,n)
  call random_number(d2)
  call random_number(d3)
  ! Iteration
  call cpu_time(tstart)
  call benchmark_a(res)
  call cpu_time(tend)
  write(*,*) 'Iteration', tend-tstart, res(10,10,10)
  ! Spread
  call cpu_time(tstart)
  call benchmark_b(res)
  call cpu_time(tend)
  write(*,*) 'Spread', tend-tstart, res(10,10,10)
end program```