如何获得多线程的 dot() 函数？

Question

当我在我的 64 位 Windows 具有 8 个逻辑 CPU 核心的 8.1 PC 上的 Julia v0.3.7 REPL 中执行此命令时：

blas_set_num_threads(CPU_CORES)
const v=ones(Float64,100000)
@time for k=1:1000000;s=dot(v,v);end

我在任务管理器或进程资源管理器的 CPU 仪表中观察到只有 12.5% CPU（1 个逻辑 CPU 核心）被使用。我还在 Windows 7 和 Windows 8.1 上观察到 Julia v0.3.5 的相同情况。我还通过在命令行上使用 "Julia -p 8" 启动来观察到相同的行为。返回到运行在没有“-p 8”命令行选项的情况下使用 Julia REPL，我尝试了这个测试：

blas_set_num_threads(CPU_CORES)
@time peakflops(10000)

在这种情况下，CPU 仪表显示 100% CPU 使用率。

由于 dot() 和 peakflops() 都使用 BLAS（在我的例子中是 OpenBLAS），我希望两者都使用 blas_set_num_threads() 指定的线程数。但是，实际上只有后一个函数起作用。 dot() 的行为是由于错误造成的，可能是在 OpenBLAS 中？

我试图通过使用矩阵乘法函数来解决 Julia 的缺点。但是，我在 GB 大小的二维数组的子向量上迭代 dot() 操作，其中子向量使用连续内存。矩阵乘法迫使我转置每个向量，从而创建一个副本。这在内部循环中是昂贵的。因此，我的选择似乎是要么学习如何使用 Julia 的并行处理 commands/macros，要么回到 Python（Intel 的 MKL BLAS 按预期运行 ddot()）。由于 dot() 是我尝试编写的例程中消耗 99% CPU 的函数，我希望 OpenBLAS 能在 Julia 中为我提供一个用户友好的最佳解决方案，而不必担心引擎盖下的所有并行处理复杂性。但也许它并没有那么糟糕......

我可以使用一些创建多线程 dot() 函数。一些示例代码将是最好的帮助。当所有线程运行在同一台机器上时，我需要使用 SharedArray 吗？如果是这样，从 Array 转换为 SharedArray 是否会创建一个副本？由于我的数组非常大，我不想同时在内存中保存两个副本。因为我的向量长度约为 100,000，并且来自一个数组的向量以不可预测的顺序使用，所以对我来说最好的多线程解决方案是 dot() 函数，它在可用核心之间拆分任务并对来自的结果求和每个核心。如何在 Julia 中像 BLAS 一样高效地做到这一点？

我在这里看到了 Tim Holy 的回答：但是，他的示例（我认为）在一个核心上执行了整个 dot() 功能，并且没有回答我的其他问题。

编辑 1：我尝试使用“-p 8”命令行选项运行ning Julia，并将上面示例中的 dot() 替换为此处的 innersimd()： http://docs.julialang.org/en/release-0.3/manual/performance-tips/ 仍然只使用 1 个核心。我修改了 innersimd() 以将其参数键入 ::Array{Float64, 1} 然后 ::SharedArray{Float64, 1}，但这仍然使用 1 个核心。 :(

编辑 2：我测试了 Julia 的矩阵乘法（BLAS 的 gemm!() 函数）：

blas_set_num_threads(CPU_CORES)
const A=ones(Float64,(4,100000))
const B=ones(Float64,(100000,4))
@time for k=1:100000;s=A*B;end

即使 Julia 在没有“-p”命令行选项的情况下启动，Julia 也只使用一个核心。

编辑 3：这是 Python:

的可比测试代码

import numpy as np
from scipy.linalg.blas import ddot
from timeit import default_timer as timer
v = np.ones(100000)
start = timer()
for k in range(1000000):
    s = ddot(v,v)
exec_time=(timer() - start)
print
print("Execution took", str(round(exec_time, 3)), "seconds")

我在 64 位 Anaconda3 v2.1.0 和 Win 上得到了相同的结果Python：7.5 秒。与我问题中的第一个示例相比，使用 Julia 0.3.7 和 OpenBLAS，执行时间为 28 秒。这使得 Python 比 Julia 快 4 倍，OpenBLAS 的单线程实现 ddot().

解释了这一点

编辑 4：我在 Python 中做了一些测试，在我的内部循环 (N=100000) 中有效地执行了 (4xN)*(Nx2) 矩阵乘法，我发现它的速度要慢两倍多。我怀疑这是因为缓存需求现在大了 8 倍，所以缓存命中率更差。为了保持 Julia 中的缓存需求与 Python 相同，Julia 的问题是：如何将我的 100000 长度向量分成 4 个块并并行执行 OpenBLAS 的单线程 ddot()（并求和结果）？长度可能不是 4 的倍数。由于 OpenBLAS 是单线程的，我可以运行 Julia 使用“-p 8”命令行参数，并在同一会话中多线程其他自定义函数。

编辑 5：运行 Julia v0.3.7 与“-p 8”命令行参数，我观察到 OpenBLAS gemm!() 函数 does 运行作为多线程（对于某些矩阵维度）：

blas_set_num_threads(CPU_CORES)
const a = rand(10000, 10000)
@time a * a

Answer 1

原因是 OpenBLAS 的 ddot 似乎不是多线程的。我认为他们主要为xgemm等BLAS-3函数实现多线程，所以如果可能的话，最好的办法是用矩阵乘法来写你的问题。

如果您能举例说明您如何尝试矩阵乘法，将会有所帮助。也许可以避免换位。 Python 你是怎么做到的？这一切都只是 BLAS，所以 Julia 应该能够像 Python 一样好。

运行带有 -p 8 标志的 Julia 会启动九个 Julia 副本并关闭 BLAS 中的多线程，因此这可能会使您的情况变得更糟。

编辑 1：我在单线程和多线程模式下使用 OpenBLAS 和 MKL 尝试了您的示例，如果有任何差异，单线程更快。不过，我不确定 OpenBLAS 是否为您的示例使用了多个线程。可能是对于瘦问题很难实现多线程加速，因此他们只为更胖的问题打开多线程。如果您的示例在同一台机器上 Python when 运行明显更快？如果是这样，请提供 Python 示例，以便我可以与我的机器上的示例进行比较。

编辑 2：我可以确认 MKL 中 ddot 的加速。准确地说，编辑 1 中的评论是关于 dgemm。正在积极开展使 Julia 多线程化的工作，当这种情况发生时，线程化的 Julia 实现可能比 OpenBLAS 的 ddot 更快。您现在的选择是

如果可能，请重写您的问题，使 OpenBLAS 的 dgemm 使用多线程，即使用更胖的矩阵。在我的机器上，多线程似乎开始 n=16。如果这是可能的，我认为它不会比 MKL 慢。通常，OpenBLAS 的 dgemm 与 MKL 的比较好。
请 OpenBLAS 开发人员将 ddot 设为多线程。如果能给项目提供一个版本就更好了，但是高性能的BLAS代码很难写。
购买 MKL 并开始销售 Revolution Julia，包括 MKL 以支付许可证费用。

编辑 3：所有比较都是在具有 Nahalem 架构和 80 个内核（但我最多只使用了 16 个）的 Ubuntu 机器上用 Julia 进行的。我使用了相同开发版本的 Julia 的两种不同构建：一种与来自源代码的 OpenBLAS 构建链接，另一种与 MKL 链接。我还在最近的 Haswell 系统上尝试了 Julia 和 OpenBLAS，以确保在该架构上没有为 ddot 实现线程。

如何获得多线程的 dot() 函数？

How to get a multi-threaded dot() function?

multithreading

julia