如果 C 是行优先顺序,为什么 ARM 内部代码采用列优先顺序?
If C is row-major order, why does ARM intrinsic code assume column-major order?
我不确定问这个问题的最佳地点在哪里,但我目前正在研究使用 ARM 内在函数并遵循本指南:https://developer.arm.com/documentation/102467/0100/Matrix-multiplication-example
但是,此处编写的代码假定数组是按列优先顺序存储的。我一直认为 C 数组是按行优先存储的。他们为什么会这样认为?
编辑:
例如,如果不是这样:
void matrix_multiply_c(float32_t *A, float32_t *B, float32_t *C, uint32_t n, uint32_t m, uint32_t k) {
for (int i_idx=0; i_idx < n; i_idx++) {
for (int j_idx=0; j_idx < m; j_idx++) {
for (int k_idx=0; k_idx < k; k_idx++) {
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
}
}
}
}
他们这样做了:
void matrix_multiply_c(float32_t *A, float32_t *B, float32_t *C, uint32_t n, uint32_t m, uint32_t k) {
for (int i_idx=0; i_idx < n; i_idx++) {
for (int k_idx=0; k_idx < k; k_idx++) {
for (int j_idx=0; j_idx < m; j_idx++) {
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
}
}
}
}
由于按 C[0]、C[1]、C[2]、C[3] 的顺序而不是按 C[ 的顺序访问 C 的空间局部性,代码会 运行 更快0]、C[2]、C[1]、C[3](其中 C[0]、C[1]、C[2]、C[3] 在内存中是连续的)。
C 本质上不是行优先或列优先。
写a[i][j]
时,由您决定i
是行索引还是列索引。
虽然先写行索引(使数组成为行优先)是一种常见的约定,但没有什么可以阻止您做相反的事情。
此外,请记住 A × B = C
等同于 Bt × At = Ct
(t
表示转置矩阵),并且读取行主矩阵就好像它是列主矩阵一样(或反之亦然)转置它,这意味着如果你想保持你的矩阵行优先,你可以反转操作数的顺序。
您没有使用像 C[i][j]
这样的 C 二维数组,所以这不是 C 如何存储任何东西的问题,而是二维索引如何手动完成的代码,使用 n * idx_1 + idx_2
,您可以选择在内循环还是外循环中循环。
但是两个非转置矩阵的 matmul 的困难部分是你需要为两个输入矩阵做出相反的选择: 一个朴素的 matmul 必须跨越远距离输入矩阵之一的元素,所以它本质上是搞砸了。这就是为什么仔细的缓存阻塞/循环平铺对于矩阵乘法很重要的主要原因。 (O(n^3) 处理 O(n^2) 数据 - 每次将其放入 L1d 缓存,and/or 放入寄存器时,您都希望充分利用它。)
循环交换 可以 加快速度以利用最内层循环中的空间局部性,如果你做得对的话。
参见 What Every Programmer Should Know About Memory? 中的缓存阻塞 matmul 示例,它在内部几个循环中遍历所有 3 个输入中的连续内存,选择在 3 个中的任何一个中都没有缩放的索引矩阵作为内部矩阵。看起来像这样:
for (j_idx)
for (k_idx)
for (i_idx)
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
请注意,B[k * j_idx + k_idx]
在循环内部循环中是不变的,并且您正在对连续内存(这很容易进行 SIMD 向量化)执行简单的 dst[0..n] += const * src[0..n]
操作,尽管您是仍在为每个 FMA 执行 2 次加载 + 1 次存储,因此这不会最大化您的 FP 吞吐量。
与缓存访问模式分开,这也避免了长依赖链进入单个累加器(C 元素)。但这对于优化的实现来说并不是真正的问题:您当然可以使用多个累加器。由于舍入误差,FP 数学不是严格关联的,但多个累加器是 并且可能比连续添加行 x 列点积的每个元素具有更小的 FP 舍入误差。
与按标准简单 C 循环的顺序添加相比,它会有 不同的 结果,但通常 更接近 确切答案。
您建议的循环顺序 i、k、j 可能是最差的。
您正在内部循环中跨越 3 个矩阵中的 2 个的远距离元素,包括对 C[] 的不连续访问,这与您在上一段中所说的相反。
将 j
作为最内层循环,您将在第一个外部迭代中访问 C[0]
、C[n]
、C[2n]
等。 B[]
也是如此,这真的很糟糕。
互换 i
和 j
循环 会让你在中间循环中连续访问 C[]
而不是跨步,而且仍然一个的行,另一个的列,在最内层的循环中。所以这将是一个严格的改进:是的,你是对的,这个天真的例子比它需要的更糟糕。
但关键问题是在 内部 循环中跨步访问某些内容:这是性能灾难;这就是为什么仔细的缓存阻塞/循环平铺对于矩阵乘法很重要的主要原因。唯一从未与比例因子一起使用的索引是 i
.
我不确定问这个问题的最佳地点在哪里,但我目前正在研究使用 ARM 内在函数并遵循本指南:https://developer.arm.com/documentation/102467/0100/Matrix-multiplication-example
但是,此处编写的代码假定数组是按列优先顺序存储的。我一直认为 C 数组是按行优先存储的。他们为什么会这样认为?
编辑: 例如,如果不是这样:
void matrix_multiply_c(float32_t *A, float32_t *B, float32_t *C, uint32_t n, uint32_t m, uint32_t k) {
for (int i_idx=0; i_idx < n; i_idx++) {
for (int j_idx=0; j_idx < m; j_idx++) {
for (int k_idx=0; k_idx < k; k_idx++) {
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
}
}
}
}
他们这样做了:
void matrix_multiply_c(float32_t *A, float32_t *B, float32_t *C, uint32_t n, uint32_t m, uint32_t k) {
for (int i_idx=0; i_idx < n; i_idx++) {
for (int k_idx=0; k_idx < k; k_idx++) {
for (int j_idx=0; j_idx < m; j_idx++) {
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
}
}
}
}
由于按 C[0]、C[1]、C[2]、C[3] 的顺序而不是按 C[ 的顺序访问 C 的空间局部性,代码会 运行 更快0]、C[2]、C[1]、C[3](其中 C[0]、C[1]、C[2]、C[3] 在内存中是连续的)。
C 本质上不是行优先或列优先。
写a[i][j]
时,由您决定i
是行索引还是列索引。
虽然先写行索引(使数组成为行优先)是一种常见的约定,但没有什么可以阻止您做相反的事情。
此外,请记住 A × B = C
等同于 Bt × At = Ct
(t
表示转置矩阵),并且读取行主矩阵就好像它是列主矩阵一样(或反之亦然)转置它,这意味着如果你想保持你的矩阵行优先,你可以反转操作数的顺序。
您没有使用像 C[i][j]
这样的 C 二维数组,所以这不是 C 如何存储任何东西的问题,而是二维索引如何手动完成的代码,使用 n * idx_1 + idx_2
,您可以选择在内循环还是外循环中循环。
但是两个非转置矩阵的 matmul 的困难部分是你需要为两个输入矩阵做出相反的选择: 一个朴素的 matmul 必须跨越远距离输入矩阵之一的元素,所以它本质上是搞砸了。这就是为什么仔细的缓存阻塞/循环平铺对于矩阵乘法很重要的主要原因。 (O(n^3) 处理 O(n^2) 数据 - 每次将其放入 L1d 缓存,and/or 放入寄存器时,您都希望充分利用它。)
循环交换 可以 加快速度以利用最内层循环中的空间局部性,如果你做得对的话。
参见 What Every Programmer Should Know About Memory? 中的缓存阻塞 matmul 示例,它在内部几个循环中遍历所有 3 个输入中的连续内存,选择在 3 个中的任何一个中都没有缩放的索引矩阵作为内部矩阵。看起来像这样:
for (j_idx)
for (k_idx)
for (i_idx)
C[n*j_idx + i_idx] += A[n*k_idx + i_idx]*B[k*j_idx + k_idx];
请注意,B[k * j_idx + k_idx]
在循环内部循环中是不变的,并且您正在对连续内存(这很容易进行 SIMD 向量化)执行简单的 dst[0..n] += const * src[0..n]
操作,尽管您是仍在为每个 FMA 执行 2 次加载 + 1 次存储,因此这不会最大化您的 FP 吞吐量。
与缓存访问模式分开,这也避免了长依赖链进入单个累加器(C 元素)。但这对于优化的实现来说并不是真正的问题:您当然可以使用多个累加器。由于舍入误差,FP 数学不是严格关联的,但多个累加器是
您建议的循环顺序 i、k、j 可能是最差的。
您正在内部循环中跨越 3 个矩阵中的 2 个的远距离元素,包括对 C[] 的不连续访问,这与您在上一段中所说的相反。
将 j
作为最内层循环,您将在第一个外部迭代中访问 C[0]
、C[n]
、C[2n]
等。 B[]
也是如此,这真的很糟糕。
互换 i
和 j
循环 会让你在中间循环中连续访问 C[]
而不是跨步,而且仍然一个的行,另一个的列,在最内层的循环中。所以这将是一个严格的改进:是的,你是对的,这个天真的例子比它需要的更糟糕。
但关键问题是在 内部 循环中跨步访问某些内容:这是性能灾难;这就是为什么仔细的缓存阻塞/循环平铺对于矩阵乘法很重要的主要原因。唯一从未与比例因子一起使用的索引是 i
.