与Mac和Linux上的exp函数略有不同

The following C program produces different results on my Mac and on Linux. I'm suprised because I assumed that the implementation of libm is somehow standardized

以下C程序在我的Mac和Linux上产生不同的结果。我很惊讶,因为我认为libm的实现在某种程度上是标准化的

#include<math.h>
#include<stdio.h>

int main()
{
  double x18=-6.899495205106946e+01;
  double x19 = exp(-x18);
  printf("x19     = %.15e\n", x19);
  printf("x19 hex = %llx\n", *((unsigned long long *)(&x19)));
}

The output on Mac is

Mac上的输出是

x19     = 9.207186811339878e+29
x19 hex = 46273e0149095886

and on Linux

在Linux上

x19     = 9.207186811339876e+29
x19 hex = 46273e0149095885

Both were compiled without any optimization flags as follows:

两者都编译时没有任何优化标志,如下所示:

gcc -lm ....

I know that I never should compare floats to be excatly the same.

我知道我永远都不应该比较漂浮物。

This issue came up during debugging, regrettably the algorithm using this calculation proofs to be numerically unstable and this slight difference leads to significant deviations in the final result. But this is a different problem.

在调试期间出现了这个问题,遗憾的是使用这种计算证明的算法在数值上是不稳定的,这种微小的差异导致最终结果的显着偏差。但这是一个不同的问题。

I'm just surprised that such basic operations as exp are not standardized as I can expect for the basic algebraic operations specified by IEEE 754.

我很惊讶exp这样的基本操作没有标准化,因为我可以期待IEEE 754规定的基本代数运算。

Are there any assumptions about precision I can rely on for different implementations of libm for different machines or for different versions ?

对于不同机器或不同版本的不同libm实现,是否有任何关于精度的假设?

Because of the discussion below I used mpmath to compute the value with higher than machine precision and I get with two more figures the result 9.2071868113398768244, so for both of my results the last figure is already wrong. The result on linux can be explained by down rounding this value, the Mac result is also off if the computer uses up rounding.

由于下面的讨论,我使用mpmath计算高于机器精度的值,我得到两个数字结果9.2071868113398768244,所以对于我的两个结果,最后一个数字已经错了。 linux上的结果可以通过向下舍入此值来解释,如果计算机使用舍入,则Mac结果也会关闭。

2 个解决方案

#1

The C99 Specification states (other version should be similar):

C99规范声明(其他版本应该类似):

J.3 Implementation-defined behavior

J.3实现定义的行为

1 A conforming implementation is required to document its choice of behavior in each of the areas listed in this subclause. The following are implementation-defined:

1需要符合规定的实施方案,以记录本子条款中列出的每个领域的行为选择。以下是实现定义:

...

J.3.6 Floating point

J.3.6浮点

— The accuracy of the floating-point operations and of the library functions in <math.h> and <complex.h> that return floating-point results (5.2.4.2.2).

- 返回浮点结果的和中浮点运算和库函数的准确性(5.2.4.2.2)。

Meaning GNU libm and BSD libm are free to have different levels of accuracy. Likely what is happening, is that the BSD implemention on OSX rounds to the nearest (unit in the last place) ULP, and the GNU implementation truncates down to the next ULP.

意思是GNU libm和BSD libm可以*地具有不同的准确度。可能正在发生的事情是,OSX上的BSD实现向最近的(最后一个单元)ULP舍入,并且GNU实现截断到下一个ULP。

#2

IEEE-754 behavior is specified at the binary level. Using a Linux, I get identical values for Python's native math library, mpmath, and MPFR (via gmpy2). However, conversion to decimal varies between the three methods.

IEEE-754行为在二进制级别指定。使用Linux,我获得Python的本机数学库,mpmath和MPFR(通过gmpy2)的相同值。但是,转换为十进制在三种方法之间有所不同。

>>> import mpmath, gmpy2
>>> import mpmath, gmpy2, math
>>> x18=68.99495205106946
>>> x19=math.exp(x18)
>>> mp18=mpmath.mpf("68.99495205106946")
>>> mp19=mpmath.exp(mp18)
>>> gp18=gmpy2.mpfr("68.99495205106946")
>>> gp19=gmpy2.exp(gp18)
>>> x18 == mp18
True
>>> x18 == gp18
True
>>> x19 == mp19
True
>>> x19 == gp19
True
>>> print(x18, mp18, gp18)
68.99495205106946 68.9949520510695 68.994952051069461
>>> print(x19, mp19, gp19)
9.207186811339876e+29 9.20718681133988e+29 9.2071868113398761e+29

After conversion to Python's arbitrary precision integer form, all three results also show as exact.

转换为Python的任意精度整数形式后,所有三个结果也显示为精确。

>>> hex(int(x19))
'0xb9f00a484ac42800000000000'
>>> hex(int(mp19))
'0xb9f00a484ac42800000000000'
>>> hex(int(gp19))
'0xb9f00a484ac42800000000000'

So (at least one) Linux math library, mpmath, and gmpy2.mpfr agree.

所以(至少有一个)Linux数学库,mpmath和gmpy2.mpfr同意。

Disclaimer: I maintain gmpy2 and have contributed to mpmath in the past.

免责声明:我维持gmpy2并且过去曾为mpmath做出过贡献。

#1