cuBLAS dsyrk 比 dgemm 慢

cuBLAS dsyrk slower than dgemm

我正在尝试使用 cuBLAS 在 GPU 上计算 C = A*A',我发现 rank-k 更新 cublasDsyrk 运行ning 比一般矩阵慢大约 5 倍-矩阵乘法例程cublasDgemm.

这让我感到惊讶;我认为 syrk 会更快,因为它是一段更专业的代码。这是一个不合理的期望吗?我做错了吗?

为代码计时

最终我正在编写要编译为 MATLAB 的 MEX 文件的 CUDA 代码,因此很抱歉没有提供完整的工作示例(会有很多与 MATLAB 对象争论的无关代码)。

我知道这可能不是最好的方法,但我正在使用 clock() 来计算代码需要多长时间 运行:

// Start of main function
clock_t tic = clock();
clock_t toc;

/* ---- snip ---- */

cudaDeviceSynchronize();

toc = clock();
printf("%8d (%7.3f ms) Allocated memory on GPU for output matrix\n",
        toc-tic,1000*(double)(toc-tic)/CLOCKS_PER_SEC);

// Compute the upper triangle of C = alpha*A*A' + beta*C
stat = cublasDsyrk(handle, CUBLAS_FILL_MODE_UPPER, CUBLAS_OP_N, 
        M, N, &alpha, A, M, &beta, C, M);

toc = clock();
printf("%8d (%7.3f ms) cublasDsyrk launched\n",
        toc-tic,1000*(double)(toc-tic)/CLOCKS_PER_SEC);

cudaDeviceSynchronize();

toc = clock();
printf("%8d (%7.3f ms) cublasDsyrk completed\n",
        toc-tic,1000*(double)(toc-tic)/CLOCKS_PER_SEC);

/* ----- snip ----- */

运行时间

运行在 [12 x 500,000] 随机矩阵(列主存储)上的输出:

  911 (  0.911 ms) Loaded inputs, initialized cuBLAS context
 1111 (  1.111 ms) Allocated memory on GPU for output matrix
 1352 (  1.352 ms) cublasDsyrk launched
85269 ( 85.269 ms) cublasDsyrk completed
85374 ( 85.374 ms) Launched fillLowerTriangle kernel
85399 ( 85.399 ms) kernel completed
85721 ( 85.721 ms) Finished and cleaned up

syrk 调用替换为

stat = cublasDgemm(handle, CUBLAS_OP_N, CUBLAS_OP_T, M, M, N, 
        &alpha, A, M, A, M, &beta, C, M);

整个过程 运行 速度更快:

  664 (  0.664 ms) Loaded inputs, initialized cuBLAS context
  796 (  0.796 ms) Allocated memory on GPU for output matrix
  941 (  0.941 ms) cublasDgemm launched
16787 ( 16.787 ms) cublasDgemm completed
16837 ( 16.837 ms) Launched fillLowerTriangle kernel
16859 ( 16.859 ms) kernel completed
17263 ( 17.263 ms) Finished and cleaned up

我尝试了一些其他大小的矩阵;有趣的是,当矩阵的行数很少时,速度差异似乎最为明显。在 100 行时,gemm 仅快 2 倍,在 1000 行时它稍慢(这是我一直期望的)。

其他详情

我正在使用 CUDA Toolkit 7.5,GPU 设备是 NVIDIA Grid K520(Kepler,计算能力 3.0)。我 运行正在使用 Amazon EC2 g2.x2large 实例。

[n x 500,000] for n=12,100,1000 都是 非常 宽的矩阵。在这些极端情况下,gemm()syrk() 可能无法达到其最佳性能,其中 syrk() 几乎是 gemm() 的两倍(因为结果矩阵是对称的,所以你可以节省一半的计算量)。

另外一个考虑是CUDA gemm()/syrk() 通常将矩阵划分为固定大小的子矩阵作为基本计算单元以实现高性能。 dgemm()子矩阵的大小最大可达32x64,如下图link.

http://www.netlib.org/lapack/lawnspdf/lawn267.pdf

如果您的大小(12 或 100)既不比子矩阵大也不是它的倍数,则性能通常会下降很多。