【发布时间】：2023-12-18 16:14:01
【问题描述】：
我的程序执行 NxN 矩阵乘法，其中两个矩阵的元素都使用 for 循环初始化为值 (0, 1, 2, ... N)。两个矩阵元素都是浮点类型。不存在内存分配问题。矩阵大小输入为 4 的倍数，例如：4x4 或 8x8 等。通过顺序计算验证答案。在 64x64 的矩阵大小下一切正常。仅当矩阵大小超过 64（例如：68 x 68）时，才会观察到顺序版本和 SSE 版本之间的差异。
SSE sn-p如图（大小=68）：
void matrix_mult_sse(int size, float *mat1_in, float *mat2_in, float *ans_out) {
__m128 a_line, b_line, r_line;
int i, j, k;
for (k = 0; k < size * size; k += size) {
for (i = 0; i < size; i += 4) {
j = 0;
b_line = _mm_load_ps(&mat2_in[i]);
a_line = _mm_set1_ps(mat1_in[j + k]);
r_line = _mm_mul_ps(a_line, b_line);
for (j = 1; j < size; j++) {
b_line = _mm_load_ps(&mat2_in[j * size + i]);
a_line = _mm_set1_ps(mat1_in[j + k]);
r_line = _mm_add_ps(_mm_mul_ps(a_line, b_line), r_line);
}
_mm_store_ps(&ans_out[i + k], r_line);
}
}
}
有了这个，答案在元素 3673 处有所不同，我得到的乘法答案如下
标量：576030144.000000 & SSE：576030208.000000
我还用相同的初始化和设置以及 N = 68 用 Java 编写了一个类似的程序，对于元素 3673，我得到的答案是 576030210.000000
现在有三个不同的答案，我不知道如何继续。为什么会出现这种差异，我们如何消除这种差异？
【问题讨论】：
-
只要结果精确到大约 6 位有效数字，就没有什么可担心的——毕竟这是单精度浮点数。
-
但它不是失去精度的小数部分。它是整数部分。
-
无论小数点在哪里，您都会得到大约 6 个有效数字 - 这就是浮点的工作原理。我建议在编写浮点代码之前阅读docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html。
-
@Jagruth.P，代码是在 32 位还是 64 位模式下编译的？如果它是 32 位代码，那么编译器可能正在使用 x87 代码，该代码将使用 80 位进行内部计算，然后返回浮点数。您可以查看程序集或使用 SSE（例如使用
mm_add_ss）编写标量代码，以确保您使用的是相同的硬件。 -
@Z boson - 我使用的是在 Oracle Virtual Box 上运行的 32 位 ubuntu 映像（不过我的笔记本电脑是 64 位的）。我不确定虚拟 Ubuntu 会以何种方式影响这一点，但奇怪的是，当我在正确引导分区的 Ubuntu 上运行相同的东西时，代码运行顺利，没有错误，我什至得到了相同的顺序和 SSE 部分（意味着没有精度差异观察到）。
标签： c sse precision matrix-multiplication rounding-error