cuBLAS（numba）中的非方形C顺序矩阵

Question

我正在尝试使用Anaconda的Numba软件包中的cuBLAS函数，遇到了一个问题。我需要将输入矩阵按C顺序排列，输出可以按Fortran顺序排列。

我可以运行软件包提供的示例脚本（此处）。脚本有两个函数gemm_v1和gemm_v2。在 gemm_v1 中，用户必须以Fortran顺序创建输入矩阵。在gemm_v2中，它们可以传递给GEMM的cuda实现并在设备上进行转置。我可以让这些示例与正方形矩阵一起工作。但是，我无法弄清楚如何让gemm_v2与非正方形输入矩阵一起工作。有没有办法使用按C顺序排列的非正方形输入矩阵？

注意：
理想情况下，在调用GEMM后，输入和输出矩阵都应保留在设备上，供其他计算使用（这是迭代方法的一部分）。

- user1554752

---

在调用blas时，您需要指定gemm(transa, transb, m, n, k, alpha, A:r, B:r, beta, C:w);其中transa和transb是要应用于矩阵的操作。在gemm_v1示例中，这是恒等操作，在gemm_v2示例中，它是转置操作。然后，您需要指定m、n和k。这些是A的行数(m)、A的列数/B的行数(n)和B的列数(k)。如果您将其保持在示例的语法中，则指定为方阵，因此这就是更改的地方。请确保您的矩阵形状与声明相匹配。 - Uvar

Answer 1

这个例子的问题在于，它只适用于方阵。如果矩阵不是方阵，由于维度不匹配（假设维度对于 A*B 是正确的），你就无法计算出 A^t*B^t。

我手头没有可用的 cuBLAS 安装程序，所以这有点像一次尝试，但我会非常惊讶如果 cuBLAS 的工作方式与通常的 BLAS 不同。BLAS 期望矩阵按列主序（即 Fortran 主序）排列，但也可以用于按行主序（即 C 主序）排列的矩阵。

在我看来（这可能完全错误），gemm_v2 不是处理两个 C 主序矩阵乘法的通常/最佳方法，例如因为如果一个人将两个 C 主序矩阵相乘，答案也将是一个 C 主序矩阵。

使用 gemm 计算两个 C 主序矩阵的积的技巧如下：

即使这可能已经为你所知，我还是想首先详细说明一下行主序（c-memory-layout）和列主序（fortran-memory-layout），以便充分解释我的答案。

因此，如果我们有一个 2x3（即 2 行 3 列）矩阵 A，并将其存储在一些连续的内存中，我们得到：

row-major-order(A) = A11, A12, A13, A21, A22, A23
col-major-order(A) = A11, A21, A12, A22, A13, A33

这意味着，如果我们得到一个连续的内存块，它按行主序表示矩阵，我们将其解释为按列主序表示矩阵，那么我们将得到完全不同的矩阵！

然而，如果我们看一下转置矩阵A^t，我们可以很容易地看出：

row-major-order(A) = col-major-order(A^t)
col-major-order(A) = row-major-order(A^t)

这意味着，如果我们想要以行主序（row-major-order）的方式得到矩阵C作为结果，blas例程应该将转置后的矩阵C以列主序（column-major-order）的方式（毕竟我们无法改变）写入这个内存。然而，C^t=(AB)^t=B^t*A^t 且 B^t 和 A^t 是在列主序中重新解释的原始矩阵。

现在，设A是一个n x k的矩阵，B是一个k x m的矩阵，则调用gemm例程应该按如下方式进行：

gemm('N', 'N', m, n, k, 1.0, B, m, A, k, 0.0, C, m)

请注意：

1. 我们不必转置矩阵 A 和 B，因为通过将 C-order 重新解释为 Fortran-order 可以处理它。
2. 我们必须交换矩阵 A 和 B 的位置，以便得到以 Fortran-order 形式表示的 C^t 作为结果。
3. 由于重新解释为 C-order，所以得到的矩阵 C 是按 C-order 排列的（从 Fortran-order 转换）并且没有 ^t。

- ead