如何提高与 Python 接口的 Fortran 代码的性能？

Question

我正在研究将 Fortran 代码连接到 Python 中的可能性。 我知道 f2py，但是，由于我没有成功地将它与外部库（例如 lapack）一起使用，所以我转而使用 ctypes。 我一切正常，但我注意到这种行为：对于相同的算法（两个矩阵元素的简单乘法），Fortran 代码比 numpy 慢，而我预计它会更快。 我尝试了两台不同的 PC，一台是 Windows 10，另一台是 Windows 11。 我使用通过 winlibs.com 版本 gcc-12.2.0-llvm-15.0.7-mingw-w64ucrt-10.0.0-r4 安装的 MinGW-64 的 gfortran，我的 Python 版本是 3.9.13

带有 timeit 的测量时间是（对于矩阵 1000x1000）

使用 Numpy：每个循环 2.77 ms ± 110 µs（7 次运行的平均值 ± 标准偏差，每次 100 次循环）

使用 Fortran：每个循环 12.4 ms ± 314 µs（7 次运行的平均值 ± 标准偏差，每次 100 次循环）

Fortran 代码在文件“libfortran3.f90”中

subroutine prodotto(a, b, n) bind(c, name='prodotto')
  ! declarations
  use iso_c_binding, only: c_double, c_int

  integer(c_int), intent(in) :: n
  real(c_double), dimension(n,n), intent(in) :: b
  real(c_double), dimension(n,n), intent(inout) :: a
  
  do i = 1,n
    do j = 1,n
        a(i,j) = a(i,j)*b(i,j)
    end do
  end do
  
end subroutine prodotto

它是用

gfortran -O3 -funroll-loops -ffast-math -fPIC -shared -o libfortran3.so libfortran3.f90

编译的。 请注意，我对输入和输出使用了相同的矩阵，因此我不必为结果传递另一个矩阵（因为我预计这会使问题恶化）。

Python代码（在笔记本中运行）是

from __future__ import print_function
from ctypes import CDLL, POINTER, c_int, c_double
import numpy as np
import time
import sys

fortran = CDLL('./libfortran3.so')
fortran.prodotto.argtypes = [ POINTER(c_double), 
                                 POINTER(c_double), 
                                 POINTER(c_int)]
fortran.prodotto.restype = None
#
N = 10
A = np.random.rand(N,N)
B = np.random.rand(N,N)
#
print('With Numpy')
%timeit A*B
#
print('With Fortran')
A = np.asfortranarray(A, dtype=c_double)
B = np.asfortranarray(B, dtype=c_double)
Act = A.ctypes.data_as(POINTER(c_double))
Bct = B.ctypes.data_as(POINTER(c_double))
%timeit fortran.prodotto( Act, Bct, c_int(N) )

当然，我的 Fortran 代码没有 Numpy 优化得那么好，但我想知道是否可以修改某些内容以获得更好的性能。

Answer 1

如评论中所述，您在 Fortran 代码中以错误的方式循环，因为 Fortran 是列优先的。话虽这么说，在这种情况下，因为它是一个元素一个元素地乘法，你可以用一个简单的

替换嵌套的 do 循环

a = a * b

其次，

N=10

很短，调用外部函数的额外开销掩盖了实际代码的性能。

设置

N=100

并用数组乘法替换手动循环，我得到

With Numpy
5.00679704999493
With Fortran
3.699547360010911