为什么在cpu和gpu上执行pytorch matmul会得到不同的结果？

Question

我正在尝试找出 numpy/pytorch、gpu/cpu、float16/float32 数字之间的舍入差异，而我的发现让我感到困惑。

基本版是：

a = torch.rand(3, 4, dtype=torch.float32)
b = torch.rand(4, 5, dtype=torch.float32)
print(a.numpy()@b.numpy() - a@b)

结果如预期的那样全为零，但是

print((a.cuda()@b.cuda()).cpu() - a@b)

得到非零结果。 为什么 Pytorch float32 matmul 在 gpu 和 cpu 上的执行方式不同？

一个更令人困惑的实验涉及float16，如下：

a = torch.rand(3, 4, dtype=torch.float16)
b = torch.rand(4, 5, dtype=torch.float16)
print(a.numpy()@b.numpy() - a@b)
print((a.cuda()@b.cuda()).cpu() - a@b)

这两个结果都是非零的。 为什么numpy和torch处理float16数字的方式不同？我知道cpu只能做float32操作，numpy在计算前将float16转换为float32，但是torch计算也是在cpu上执行的。

你猜怎么着，print((a.cuda()@b.cuda()).cpu() - a.numpy()@b.numpy()) 得到一个全为零的结果！这对我来说纯粹是幻想......

环境如下：

• 蟒蛇：3.8.5
• 火炬：1.7.0
• numpy：1.21.2
• cuda：11.1
• GPU：GeForce RTX 3090

---

根据一些评论者的建议，我添加了以下相等测试

(a.numpy()@b.numpy() - (a@b).numpy()).any()
((a.cuda()@b.cuda()).cpu() - a@b).numpy().any()
(a.numpy()@b.numpy() - (a@b).numpy()).any()
((a.cuda()@b.cuda()).cpu() - a@b).numpy().any()
((a.cuda()@b.cuda()).cpu().numpy() - a.numpy()@b.numpy()).any()

分别直接跟上面五个打印函数，结果分别是：

False
True
True
True
False

而对于最后一个，我试了好几次，我想我可以排除运气。

Answer 1

如@talonmies 所述，差异主要是数字。 CPU/GPU 及其各自的 BLAS 库的实现方式不同，使用不同的操作/操作顺序，因此存在数值差异。

一个可能的原因是顺序操作与缩减 (https://discuss.pytorch.org/t/why-different-results-when-multiplying-in-cpu-than-in-gpu/1356/3)，例如(((a+b)+c)+d) 与 ((a+b)+(c+d)) 相比具有不同的数值属性。

This question 还谈到了可能导致数值差异的融合操作（乘加）。

我做了一点测试，发现如果我们在计算前将数据类型提升为float32，然后将其降级，可以匹配GPU在float16模式下的输出。这可能是由内部中间转换或融合操作的更好的数值稳定性引起的（torch.backends.cudnn.enabled 无关紧要）。但这并不能解决 float32 中的情况。

import torch

def test(L, M, N):
    # test (L*M) @ (M*N)
    for _ in range(5000):
        a = torch.rand(L, M, dtype=torch.float16)
        b = torch.rand(M, N, dtype=torch.float16)

        cpu_result = a@b
        gpu_result = (a.cuda()@b.cuda()).cpu()
        if (cpu_result-gpu_result).any():
            print(f'({L}x{M}) @ ({M}x{N}) failed')
            return
    else:
        print(f'({L}x{M}) @ ({M}x{N}) passed')


test(1, 1, 1)
test(1, 2, 1)
test(4, 1, 4)
test(4, 4, 4)

def test2():
    for _ in range(5000):
        a = torch.rand(1, 2, dtype=torch.float16)
        b = torch.rand(2, 1, dtype=torch.float16)

        cpu_result = a@b
        gpu_result = (a.cuda()@b.cuda()).cpu()

        half_result = a[0,0]*b[0,0] + a[0,1]*b[1,0]
        convert_result = (a[0,0].float()*b[0,0].float() + a[0,1].float()*b[1,0].float()).half()

        if ((cpu_result-half_result).any()):
            print('CPU != half')
            return
        if (gpu_result-convert_result).any():
            print('GPU != convert')
            return
    else:
        print('All passed')

test2()

输出：

(1x1) @ (1x1) passed
(1x2) @ (2x1) failed
(4x1) @ (1x4) passed
(4x4) @ (4x4) failed
All passed

你可以看出，当内部维度为1时，它通过了检查（不需要乘加/减）。