使用 FORTRAN 排序的 Numpy 数组到 ctypes

Question

是否有一种高性能方法可以将 numpy 数组转换为 FORTRAN 有序 ctypes 数组，理想情况下不需要副本，并且不会触发与步幅相关的问题？

\n

import numpy as np\n\n# Sample data\nn = 10000\nA = np.zeros((n,n), dtype=np.int8)\nA[0,1] = 1\n\ndef slow_conversion(A):\n    return np.ctypeslib.as_ctypes(np.ascontiguousarray(A.T))\n\nassert slow_conversion(A)[1][0] == 1\n

\n

仅调用as_ctypes的性能分析：

\n

%%timeit\nnp.ctypeslib.as_ctypes(A)\n

\n

3.35 \xc2\xb5s \xc2\xb1 每个循环 10.5 ns（意味着 \xc2\xb1 标准偏差 7 次运行，每次 100000 次循环）

\n

所提供的（慢速）转换的性能分析

\n

%%timeit\nslow_conversion(A)\n

\n

206 ms \xc2\xb1 每个循环 10.4 ms（平均 \xc2\xb1 标准偏差 7 次运行，每次 1 次循环）

\n

理想的结果是获得与调用类似的性能as_ctypes。

\n

Answer 1

默认情况下，Numpy 创建 C 序数组（因为它是用 C 编写的），即行优先数组。进行转置A.T会创建一个数组视图，其中步长是相反的（即没有复制）。话虽如此，np.ascontiguousarray进行复制是因为数组现在不再连续，并且复制成本很高。这就是为什么slow_conversion慢的原因。yourarray.flags['F_CONTIGUOUS']请注意，可以使用和 通过检查来测试连续性yourarray.strides。另请注意yourarray.flags，并yourarray.__array_interface__提供有关数组是否已被复制的信息以及有关步幅的信息。

np.asfortranarray返回有关文档的内存中按 Fortran 顺序排列的数组。如果需要，它可以执行复制。事实上，np.asfortranarray(A)复制而不复制np.asfortranarray(A.T)。您可以检查该函数的C 代码以获取有关此行为的更多信息。由于两者都被视为 FORTRAN 连续，因此在这种情况下最好使用np.asfortranarray(A.T)不进行任何复制。

关于 ctypes，它处理以行优先顺序存储的 C 数组，而不是以列优先顺序存储的 FORTRAN 数组。此外，与 FORTRAN 数组不同，C 数组本身不支持步幅。这意味着行基本上是存储在内存中的连续数据的内存视图。由于slow_conversion(A)[1][0] == 1要求为 true，这意味着返回的对象第二行的第一项应该为 1，因此这些列必须连续存储在内存中。问题是初始数组不是 FORTRAN 连续的而是 C 连续的，因此需要转置。转置非常昂贵（尽管 Numpy 实现不是最优的）。

假设您不想支付复制/转置的开销，则需要放宽该问题。有几个可能的选择：

• 例如使用 Numpy 直接创建 FORTRAN 有序数组 np.zeros((n,n), dtype=np.int8, order='F')。这将创建一个具有转置步长的 C 数组，以便其行为类似于 FORTRAN 数组，其中对列进行的计算速度很快（请记住，Numpy 是用 C 编写的，因此行优先有序数组是参考）。这样，ctypes 中的第一行实际上是一列。请注意，出于性能考虑，在混合 C 和 FORTRAN 有序数组时应非常小心，因为非连续访问速度要慢得多。
• 使用跨步 FORTRAN 数组。这一解决方案基本上意味着基本的基于列的计算将慢得多，并且需要编写基于行的计算，这在 FORTRAN 中非常不寻常。您可以使用 提取指向 C 连续数组的指针A.ctypes.data_as(POINTER(c_double))，使用 提取步长A.strides，使用 提取形状A.shape。话虽这么说，这个解决方案似乎并不是真正的便携/标准。标准方法似乎是在 FORTRAN 中使用 C 绑定。我对此不是很熟悉，但您可以在此答案中找到完整的示例。

最后一个解决方案是使用快速转置算法手动转置数据。就地转置比异地转置更快（通常是 2 倍），但这需要一个平方数组，并且它会改变稍后不应使用的输入数组（除非可以在转置数组上进行操作） 。快速转置不能直接使用 Numpy 编写。一种解决方案是在 Numba 中执行（这是一个次优的异地解决方案，但已经相当快了），或者在 C 中或直接在 FORTRAN 中执行（在所有情况下都使用包装函数）。这应该比 Numpy 快得多，但仍然比基本的 ctypes 包装慢得多。

Answer 2

如果您的起始数组必须是 C 序的：

\n

正如您所指出的，np.ctypeslib.as_ctypes(...)速度很快。

\n

你的计算的瓶颈是np.ascontiguousarray(A.T)- 这相当于np.asfortranarray(A)同样慢。

\n\n

这让我相信仅使用 numpy 无法使这变得更快。我的意思是，由于已经存在一个完整的专用函数来执行此操作（np.asfortranarray） - 我们假设它已经过优化以具有最佳性能。

\n

如果您可以从 F 有序数组开始：

\n

这样的话，问题就解决了！您可以获得与调用一样的性能，as_ctypes因为数组在内存中已经按照所需的顺序排列：

\n

n = 10000\nA = np.zeros((n, n), dtype=np.int8, order=\'F\')\nA[0, 1] = 1\n\ndef conversion(A):\n    return np.ctypeslib.as_ctypes(np.ascontiguousarray(A.T))\n\nassert conversion(A)[1][0] == 1\n

\n

%%timeit\nconversion(A)\n\n# 4.76 \xc2\xb5s \xc2\xb1 68.4 ns per loop (mean \xc2\xb1 std. dev. of 7 runs, 100,000 loops each)\n

\n

感谢@J\xc3\xa9r\xc3\xb4meRichard 和@StephanSchlecht 激发了这个想法。阅读@J\xc3\xa9r\xc3\xb4meRichard\ 的答案以获得更多见解！

\n