intel MacBook 和 M1 之间的 np.float32 浮点差异

Question

我最近将我的 Intel MacBook Pro 13" 升级为配备 M1 Pro 的 MacBook Pro 14"。一直在努力让我的软件重新编译和工作。幸运的是，除了一些晦涩的 Fortran 代码和 Python 中的浮点问题之外，没有什么大问题。关于 python/numpy 我有以下问题。

我有一个很大的代码库，但为了简单起见，我将使用这个简单的函数，将飞行高度转换为压力来显示问题。

def fl2pres(FL):
    P0=101325
    T0=288.15
    T1=216.65
    g=9.80665
    R=287.0528742
    GAMMA=0.0065
    P11=P0*np.exp(-g/GAMMA/R*np.log(T0/T1))

    h=FL*30.48

    return np.where(h<=11000, \
        P0*np.exp(-g/GAMMA/R*np.log((T0/(T0-GAMMA*h) ))),\
            P11*np.exp(-g/R/T1*(h-11000)) )

当我在 M1 Pro 上运行代码时，我得到：

In [2]: fl2pres(np.float64([400, 200]))
Out[3]: array([18753.90334892, 46563.239766  ])

和;

In [3]: fl2pres(np.float32([400, 200]))
Out[3]: array([18753.90234375, 46563.25080916])

在我的旧款 Intel MacBook Pro 上执行同样的操作，我得到：

In [2]: fl2pres(np.float64([400, 200]))
Out[2]: array([18753.90334892, 46563.239766  ])

和;

In [3]: fl2pres(np.float32([400, 200]))
Out[3]: array([18753.904296888, 46563.24778944])

float64 计算匹配，但 float32 不匹配。我们在代码中大量使用 float32 来优化内存。我知道由于架构差异，可能会发生这种浮点错误，但想知道是否可以进行简单的修复，因为目前某些单元测试失败了。我可以将架构包含在这些测试中，但我希望有一个更简单的解决方案？

将所有输入转换为 float64 使我的单元测试通过，从而解决了这个问题，但由于我们有相当多的大型数组和数据帧，因此对内存的影响是不必要的。

两台笔记本电脑都运行通过 homebrew 安装的 python 3.9.10、pandas 1.4.1 和 numpy 1.22.3（安装用于映射加速和 blas）。

编辑 我更改了打印中间值的函数以查看发生变化的位置：

def fl2pres(FL):
    P0=101325
    T0=288.15
    T1=216.65
    g=9.80665
    R=287.0528742
    GAMMA=0.0065
    P11=P0*np.exp(-g/GAMMA/R*np.log(T0/T1))

    h=FL*30.48
    A = np.log((T0/(T0-GAMMA*h)))
    B = np.exp(-g/GAMMA/R*A)
    C = np.exp(-g/R/T1*(h-11000))
    print(f"P11:{P11}, h:{h}, A:{A}, B:{B}, C:{C}")
    return np.where(h<=11000, P0*B, P11*C)

使用与上述 float32 情况相同的输入运行此函数，我在 M1 Pro 上得到：

P11:22632.040591374975, h:[12192.  6096.], A:[0.32161594 0.14793371], B:[0.1844504  0.45954345], C:[0.82864394 2.16691503]
array([18753.90334892, 46563.239766  ])

关于英特尔：

P11:22632.040591374975, h:[12192.  6096.], A:[0.32161596 0.14793368], B:[0.18445034 0.45954353], C:[0.828644 2.166915]
array([18753.90429688, 46563.24778944])

Answer 1

根据我在 numpy 的 GitHub 上创建的问题：

您遇到的差异似乎都在一个“ULP”（最后一个单位）内，也许是2？对于特殊的数学函数，如 exp 或 sin，不幸的是，小错误是预料之中的，并且可能取决于系统（硬件和操作系统/数学库）。

可能会产生稍大影响的一件事可能是在较新的机器上使用 NumPy 的 SVML（即仅在英特尔机器上）。可以在构建时使用 NPY_DISABLE_SVML=1 作为环境变量来禁用它，但我认为您不能在不构建 NumPy 的情况下禁用它的使用。（然而，目前来看，M1 机器很可能是精度较低的机器，或者它们都大致相同，只是不同）

我还没有尝试使用编译 numpy NPY_DISABLE_SVML=1，我现在的计划是使用一个可以在我的所有平台上运行的 docker 容器，并使用单个“真相”进行测试。