Python 楼层划分中的舍入误差
rounding errors in Python floor division
我知道浮点运算中会出现舍入错误，但有人能解释一下原因吗：
>>> 8.0 / 0.4 # as expected
20.0
>>> floor(8.0 / 0.4) # int works too
20
>>> 8.0 // 0.4 # expecting 20.0
19.0
x64 上的 Python 2 和 3 都会发生这种情况。
据我所知，这要么是错误，要么是 // 的一个非常愚蠢的规范，因为我看不出为什么最后一个表达式的计算结果应该是 19.0.
为什么 a // b 不简单定义为 floor(a / b)？
编辑：8.0 % 0.4 也计算为 0.3999999999999996。至少这是结果，因为 8.0 // 0.4 * 0.4 + 8.0 % 0.4 计算为 8.0
EDIT：这不是 Is floating point math broken? 的副本，因为我在问为什么这个特定操作会出现（可能是可以避免的）舍入错误，为什么 a // b 未定义为/等于 floor(a / b)
REMARK：我想这不起作用的更深层原因是楼层划分是不连续的，因此有一个无限的 condition number 使它成为一个病态 -提出的问题。底除法和浮点数根本上是不兼容的，你不应该在浮点数上使用 //。只需使用整数或分数即可。
那是因为python（浮点有限表示）中没有0.4，它实际上是一个像0.4000000000000001这样的浮点数，这使得除法底数为19。
>>> floor(8//0.4000000000000001)
19.0
但真正的除法 (/) returns a reasonable approximation of the division result if the arguments are floats or complex. 这就是 8.0/0.4 结果为 20 的原因。它实际上取决于参数的大小（在 C 中双参数）。 （不四舍五入到最接近的浮点数）
阅读 Guido 本人关于 pythons integer division floors 的更多信息。
另外，有关浮点数的完整信息，您可以阅读这篇文章 https://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html
对于那些有兴趣的人，以下函数是float_div，它在Cpython的源代码中对浮点数进行真正的除法任务：
float_div(PyObject *v, PyObject *w)
{
double a,b;
CONVERT_TO_DOUBLE(v, a);
CONVERT_TO_DOUBLE(w, b);
if (b == 0.0) {
PyErr_SetString(PyExc_ZeroDivisionError,
"float division by zero");
return NULL;
}
PyFPE_START_PROTECT("divide", return 0)
a = a / b;
PyFPE_END_PROTECT(a)
return PyFloat_FromDouble(a);
}
最终结果将由函数PyFloat_FromDouble:
计算得出
PyFloat_FromDouble(double fval)
{
PyFloatObject *op = free_list;
if (op != NULL) {
free_list = (PyFloatObject *) Py_TYPE(op);
numfree--;
} else {
op = (PyFloatObject*) PyObject_MALLOC(sizeof(PyFloatObject));
if (!op)
return PyErr_NoMemory();
}
/* Inline PyObject_New */
(void)PyObject_INIT(op, &PyFloat_Type);
op->ob_fval = fval;
return (PyObject *) op;
}
好的，经过一些研究，我发现了这个 issue。
似乎正在发生的事情是，正如@khelwood 所建议的那样，0.4 在内部计算为 0.40000000000000002220，除以 8.0 时会产生比 20.0 略小的结果。 / 运算符然后舍入到最接近的浮点数，即 20.0，但 // 运算符立即截断结果，产生 19.0.
这应该更快，我想它 "close to the processor"，但我仍然不是用户想要/期望的。
@jotasi 解释了背后的真正原因。
但是如果你想阻止它，你可以使用 decimal 模块，它基本上被设计用来表示十进制浮点数，与二进制浮点表示完全不同。
所以在你的情况下，你可以这样做：
>>> from decimal import *
>>> Decimal('8.0')//Decimal('0.4')
Decimal('20')
正如您和 khelwood 已经注意到的，0.4 不能精确表示为浮点数。为什么？它是五分之二（4/10 == 2/5），它没有有限的二进制分数表示形式。
试试这个：
from fractions import Fraction
Fraction('8.0') // Fraction('0.4')
# or equivalently
# Fraction(8, 1) // Fraction(2, 5)
# or
# Fraction('8/1') // Fraction('2/5')
# 20
不过
Fraction('8') // Fraction(0.4)
# 19
在这里，0.4 被解释为需要（二进制）舍入的浮点文字（因此是浮点二进制数），只有 then 转换为有理数 Fraction(3602879701896397, 9007199254740992)，几乎但不完全是 4 / 10。然后执行底除法，因为
19 * Fraction(3602879701896397, 9007199254740992) < 8.0
和
20 * Fraction(3602879701896397, 9007199254740992) > 8.0
结果是 19，不是 20。
可能也会发生同样的情况
8.0 // 0.4
即，似乎下限除法是原子确定的（但在解释的浮点文字的唯一近似浮点值上）。
那为什么
floor(8.0 / 0.4)
给出"right"结果？因为在那里，两个舍入误差相互抵消。 首先 1) 执行除法，产生略小于 20.0 的值，但不能表示为浮点数。它会四舍五入到最接近的浮点数，恰好是 20.0。只执行 then，执行 floor 操作，但现在作用于 exactly 20.0，因此不会改变任何数字更多
1) As Kyle Strand , 即确定确切的结果 then 四舍五入 isn 't什么实际上发生在低2)级别（CPython的C代码甚至CPU 说明）。但是，它可以作为确定预期 3) 结果的有用模型。
2) 在最低 4)水平，然而，这可能不是太远了。一些芯片组通过首先计算更精确（但仍然不精确，只是有更多二进制数字）内部浮点结果然后四舍五入到 IEEE 双精度来确定浮点结果。
3) "expected" 根据 Python 规范，不一定是我们的直觉。
4) 嗯，最低级别 高于 逻辑门。我们不必考虑使半导体成为可能的量子力学就可以理解这一点。
在 github (https://github.com/python/cpython/blob/966b24071af1b320a1c7646d33474eeae057c20f/Objects/floatobject.c) 上检查了 cpython 中浮动对象的半官方来源后，可以理解这里发生了什么。
对于正常除法 [=12=] 被调用（第 560 行），它在内部将 python floats 转换为 c-doubles，进行除法然后转换结果double回到了pythonfloat。如果你只是用 8.0/0.4 在 c 中这样做，你会得到：
#include "stdio.h"
#include "math.h"
int main(){
double vx = 8.0;
double wx = 0.4;
printf("%lf\n", floor(vx/wx));
printf("%d\n", (int)(floor(vx/wx)));
}
// gives:
// 20.000000
// 20
对于楼层划分，还有一些事情发生了。在内部，float_floor_div（第 654 行）被调用，然后调用 float_divmod，一个应该 return 包含 python float 的元组的函数floored 分区，以及 mod/remainder，尽管后者刚刚被 PyTuple_GET_ITEM(t, 0) 丢弃。这些值按以下方式计算（转换为 c-doubles 后）：
- 余数是使用
double mod = fmod(numerator, denominator)计算的。
- 除法时分子减去
mod得到整数值。
- 通过有效计算
floor((numerator - mod) / denominator) 计算出底除法的结果
- 之后，@Kasramvd 的回答中已经提到的检查就完成了。但这只会将
(numerator - mod) / denominator 的结果捕捉到最接近的整数值。
这给出不同结果的原因是，由于浮点运算，fmod(8.0, 0.4) 给出 0.4 而不是 0.0。因此，计算的结果实际上是 floor((8.0 - 0.4) / 0.4) = 19 并且将 (8.0 - 0.4) / 0.4) = 19 对齐到最接近的整数值并不能修复由 fmod 的 "wrong" 结果引入的错误。你也可以很容易地在 c 中检查它：
#include "stdio.h"
#include "math.h"
int main(){
double vx = 8.0;
double wx = 0.4;
double mod = fmod(vx, wx);
printf("%lf\n", mod);
double div = (vx-mod)/wx;
printf("%lf\n", div);
}
// gives:
// 0.4
// 19.000000
我猜想，他们选择了这种计算底数除法的方式来保持 (numerator//divisor)*divisor + fmod(numerator, divisor) = numerator 的有效性（如 @0x539 的回答中的 link 中所述），即使这现在结果floor(8.0/0.4) != 8.0//0.4.
有点出乎意料的行为
我知道浮点运算中会出现舍入错误，但有人能解释一下原因吗：
>>> 8.0 / 0.4 # as expected
20.0
>>> floor(8.0 / 0.4) # int works too
20
>>> 8.0 // 0.4 # expecting 20.0
19.0
x64 上的 Python 2 和 3 都会发生这种情况。
据我所知，这要么是错误，要么是 // 的一个非常愚蠢的规范，因为我看不出为什么最后一个表达式的计算结果应该是 19.0.
为什么 a // b 不简单定义为 floor(a / b)？
编辑：8.0 % 0.4 也计算为 0.3999999999999996。至少这是结果，因为 8.0 // 0.4 * 0.4 + 8.0 % 0.4 计算为 8.0
EDIT：这不是 Is floating point math broken? 的副本，因为我在问为什么这个特定操作会出现（可能是可以避免的）舍入错误，为什么 a // b 未定义为/等于 floor(a / b)
REMARK：我想这不起作用的更深层原因是楼层划分是不连续的，因此有一个无限的 condition number 使它成为一个病态 -提出的问题。底除法和浮点数根本上是不兼容的，你不应该在浮点数上使用 //。只需使用整数或分数即可。
那是因为python（浮点有限表示）中没有0.4，它实际上是一个像0.4000000000000001这样的浮点数，这使得除法底数为19。
>>> floor(8//0.4000000000000001)
19.0
但真正的除法 (/) returns a reasonable approximation of the division result if the arguments are floats or complex. 这就是 8.0/0.4 结果为 20 的原因。它实际上取决于参数的大小（在 C 中双参数）。 （不四舍五入到最接近的浮点数）
阅读 Guido 本人关于 pythons integer division floors 的更多信息。
另外，有关浮点数的完整信息，您可以阅读这篇文章 https://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html
对于那些有兴趣的人，以下函数是float_div，它在Cpython的源代码中对浮点数进行真正的除法任务：
float_div(PyObject *v, PyObject *w)
{
double a,b;
CONVERT_TO_DOUBLE(v, a);
CONVERT_TO_DOUBLE(w, b);
if (b == 0.0) {
PyErr_SetString(PyExc_ZeroDivisionError,
"float division by zero");
return NULL;
}
PyFPE_START_PROTECT("divide", return 0)
a = a / b;
PyFPE_END_PROTECT(a)
return PyFloat_FromDouble(a);
}
最终结果将由函数PyFloat_FromDouble:
PyFloat_FromDouble(double fval)
{
PyFloatObject *op = free_list;
if (op != NULL) {
free_list = (PyFloatObject *) Py_TYPE(op);
numfree--;
} else {
op = (PyFloatObject*) PyObject_MALLOC(sizeof(PyFloatObject));
if (!op)
return PyErr_NoMemory();
}
/* Inline PyObject_New */
(void)PyObject_INIT(op, &PyFloat_Type);
op->ob_fval = fval;
return (PyObject *) op;
}
好的，经过一些研究，我发现了这个 issue。
似乎正在发生的事情是，正如@khelwood 所建议的那样，0.4 在内部计算为 0.40000000000000002220，除以 8.0 时会产生比 20.0 略小的结果。 / 运算符然后舍入到最接近的浮点数，即 20.0，但 // 运算符立即截断结果，产生 19.0.
这应该更快，我想它 "close to the processor"，但我仍然不是用户想要/期望的。
@jotasi 解释了背后的真正原因。
但是如果你想阻止它，你可以使用 decimal 模块，它基本上被设计用来表示十进制浮点数，与二进制浮点表示完全不同。
所以在你的情况下，你可以这样做：
>>> from decimal import *
>>> Decimal('8.0')//Decimal('0.4')
Decimal('20')
正如您和 khelwood 已经注意到的，0.4 不能精确表示为浮点数。为什么？它是五分之二（4/10 == 2/5），它没有有限的二进制分数表示形式。
试试这个：
from fractions import Fraction
Fraction('8.0') // Fraction('0.4')
# or equivalently
# Fraction(8, 1) // Fraction(2, 5)
# or
# Fraction('8/1') // Fraction('2/5')
# 20
不过
Fraction('8') // Fraction(0.4)
# 19
在这里，0.4 被解释为需要（二进制）舍入的浮点文字（因此是浮点二进制数），只有 then 转换为有理数 Fraction(3602879701896397, 9007199254740992)，几乎但不完全是 4 / 10。然后执行底除法，因为
19 * Fraction(3602879701896397, 9007199254740992) < 8.0
和
20 * Fraction(3602879701896397, 9007199254740992) > 8.0
结果是 19，不是 20。
可能也会发生同样的情况8.0 // 0.4
即，似乎下限除法是原子确定的（但在解释的浮点文字的唯一近似浮点值上）。
那为什么
floor(8.0 / 0.4)
给出"right"结果？因为在那里，两个舍入误差相互抵消。 首先 1) 执行除法，产生略小于 20.0 的值，但不能表示为浮点数。它会四舍五入到最接近的浮点数，恰好是 20.0。只执行 then，执行 floor 操作，但现在作用于 exactly 20.0，因此不会改变任何数字更多
1) As Kyle Strand
2) 在最低 4)水平，然而，这可能不是太远了。一些芯片组通过首先计算更精确（但仍然不精确，只是有更多二进制数字）内部浮点结果然后四舍五入到 IEEE 双精度来确定浮点结果。
3) "expected" 根据 Python 规范，不一定是我们的直觉。
4) 嗯，最低级别 高于 逻辑门。我们不必考虑使半导体成为可能的量子力学就可以理解这一点。
在 github (https://github.com/python/cpython/blob/966b24071af1b320a1c7646d33474eeae057c20f/Objects/floatobject.c) 上检查了 cpython 中浮动对象的半官方来源后，可以理解这里发生了什么。
对于正常除法 [=12=] 被调用（第 560 行），它在内部将 python floats 转换为 c-doubles，进行除法然后转换结果double回到了pythonfloat。如果你只是用 8.0/0.4 在 c 中这样做，你会得到：
#include "stdio.h"
#include "math.h"
int main(){
double vx = 8.0;
double wx = 0.4;
printf("%lf\n", floor(vx/wx));
printf("%d\n", (int)(floor(vx/wx)));
}
// gives:
// 20.000000
// 20
对于楼层划分，还有一些事情发生了。在内部，float_floor_div（第 654 行）被调用，然后调用 float_divmod，一个应该 return 包含 python float 的元组的函数floored 分区，以及 mod/remainder，尽管后者刚刚被 PyTuple_GET_ITEM(t, 0) 丢弃。这些值按以下方式计算（转换为 c-doubles 后）：
- 余数是使用
double mod = fmod(numerator, denominator)计算的。 - 除法时分子减去
mod得到整数值。 - 通过有效计算
floor((numerator - mod) / denominator)计算出底除法的结果
- 之后，@Kasramvd 的回答中已经提到的检查就完成了。但这只会将
(numerator - mod) / denominator的结果捕捉到最接近的整数值。
这给出不同结果的原因是，由于浮点运算，fmod(8.0, 0.4) 给出 0.4 而不是 0.0。因此，计算的结果实际上是 floor((8.0 - 0.4) / 0.4) = 19 并且将 (8.0 - 0.4) / 0.4) = 19 对齐到最接近的整数值并不能修复由 fmod 的 "wrong" 结果引入的错误。你也可以很容易地在 c 中检查它：
#include "stdio.h"
#include "math.h"
int main(){
double vx = 8.0;
double wx = 0.4;
double mod = fmod(vx, wx);
printf("%lf\n", mod);
double div = (vx-mod)/wx;
printf("%lf\n", div);
}
// gives:
// 0.4
// 19.000000
我猜想，他们选择了这种计算底数除法的方式来保持 (numerator//divisor)*divisor + fmod(numerator, divisor) = numerator 的有效性（如 @0x539 的回答中的 link 中所述），即使这现在结果floor(8.0/0.4) != 8.0//0.4.