[翻译] 溢出原理及如何防止溢出
本文首发于 算法社区 dspstack.com，转载请注明出处，谢谢。
感恩节快乐!也许吃太多火鸡在你脑中还记忆犹新。如果是这样，这将是讨论溢出的好时机。
在浮点运算的世界中，溢出是可能的，但不是特别常见。当数字变得太大时，就会溢出；IEEE双精度浮点数支持2^1024以下的范围，如果超过这个范围，就会出现问题：
for k in [10, 100, 1000, 1020, 1023, 1023.9, 1023.9999, 1024]:
try:
print "2^%.4f = %g" % (k, 2.0 ** k)
except OverflowError, e:
print "2^%.4f ---> %s" % (k,e)
2^10.0000 = 1024
2^100.0000 = 1.26765e+30
2^1000.0000 = 1.07151e+301
2^1020.0000 = 1.12356e+307
2^1023.0000 = 8.98847e+307
2^1023.9000 = 1.67731e+308
2^1023.9999 = 1.79757e+308
2^1024.0000 ---> (34, 'Result too large')
2^1024 是一个非常非常大的数字。你可能永远不会需要像2^1024这样的数字，除非你在用组合学。为了方便比较，这里有一些物理量的数值区间：
宇宙中的质子数估计为 1080≈2266
宇宙大小与普朗克长度之比约为 46×10^9 光年/ 1.62×10^-35 m≈2^204
宇宙年龄与普朗克时间之比约为138 亿年 / 5.39×10^-44 秒 ≈2^202
因此，双精度数字应该是相当安全的。
单精度浮点数能到略小于2^128，在巨大的物理计算中容易溢出。更可能的是，由于精度而不是数值区间，不会使用它们。
在嵌入式系统世界中，我们通常没有奢侈去使用浮点计算。它无论在 CPU 执行时间，还是在钱上都需要额外的开销。我们只能使用整数算术。（实际上我们中的一些人喜欢整数运算！）
所以本文的其余部分是关于整数算术中的溢出。接下来的事情可能会让你吃惊。
溢出源
常见的疑点是：加法和乘法
什么会导致整数数学溢出？常见的疑点为加法和乘法：
import numpy as np
a = np.int16(32767)
b = np.int16(2)
print "32767 + 2 = %s" % (a+b)
print "32767 * 2 = %s" % (a*b)
32767 + 2 = -32767
32767 * 2 = -2
-c:3: RuntimeWarning: overflow encountered in short_scalars
-c:4: RuntimeWarning: overflow encountered in short_scalars
每个整数格式都有一个可表示的范围：
有符号 16 位整数支持范围 [ -32768, 32767 ]
无符号 16 位整数支持范围 [ 0, 65535 ]
有符号 32 位整数支持范围 [ -2147483648, 21474833647 ]
无符号 32 位整数支持范围 [ 0, 4294967295]
如果超出这些范围，即使是暂时的，你也需要非常小心。大多数环境都能很好地处理溢出，并提供“舍入”或“模”行为（在有符号 16 位环境中为32767+1=-32768），其中，超出范围的进位将消失，剩下的是与精确结果对应的低阶位。
如果你使用 C 语言编程，并且有个满足 C99 标准的编译器，你可能会惊讶的发现，无符号整数数学在溢出条件下保证具有模行为，但是有符号整数数学的溢出行为是未定义的。
对于想成为编程语言律师的人来说，C99 的相关章节是：
6.5 节 第 5 项：
如果在表达式求值期间出现异常条件(也就是说，如果结果在数学上没有定义或其值不在其类型的可表示值范围内)，那么该行为是未定义的。
6.2.5 小节 第 9 项：
涉及到无符号操作数的计算是永远不会溢出的，因为不能由无符号整数类型所表示的结果值会通过取模来进行简化，而取模的分母是比结果值类型所能表示的最大值更大的数值。
(实践中，现代编译器，比如 gcc 和 clang，都倾向于对基本的算术计算使用取模行为。但是要小心，因为有时候编译器优化会错误地处理比较或条件表达式，即使基本算术是“正确的”。 这种行为的症状很难识别。LLVM 网站对此有一些解释，同时也有对其他未定义行为的缺陷的说明。 如果想要安全，使用 gcc 和 clang 都支持的 -fwrapv 编译器标志，它保证了有符号运算的所有方面都具有取模行为。)
C 标准这个小缺陷的历史原因是，在过去，一些处理器使用 1 补数来表示有符号数，这种情况的有符号和无符号的算术运算有所不同。如今，2 补数是标准， 当对结果使用基本字长时，有符号和无符号的加法、减法和乘法的实现是相同的。 但是无论好坏，C 语言的其中一个主要准则是允许特定机器的行为来维护有效代码，所以 编写标准时不会保证有符号运算所产生的溢出结果是何种形式。另一方面，像 Java 语言具有独立于机器的算术定义。 Java 语言的任何算术操作无论是在哪种处理器上运行都会给你相同的答案。这个保证的代价是，在某些处理器上，需要额外的指令。
顺便说一下，如果您使用 C 进行整数运算，请确保 #include <stdint.h>，并使用它包含的类型定义 typedefs，如 int16_t、uint16_t、int32_t 和 uint32_t。这些都是可移植的，而短整型、整型或长整型的位数可能会随着处理器架构的不同而变化。
如果您使用 MATLAB 的定点特性，请注意整数溢出的默认行为是在整数范围的极限处饱和。这避免了一些溢出问题，但是如果您习惯于使用概括语义的语言，那么它可能不会给您所期望的结果。
防止溢出
仅仅因为我们可以在 gcc 或者 clang 编译器使用 -fwrapv 而保证截断行为（wraparound behavior），这在应用程序并不会是正确的行为。
如果我控制着一个值，并且我想使输出增加，并不断给过程变量加 1 ，那么我会得到 32765，32766，32767，-32768，-32767，等等。这样就产生了一个跳跃式不连续，这很糟糕。防止这种情况发生的与机器无关的唯一方法是避免溢出。具体方法为向上转型到更大的类型大小，检查溢出，并使结果饱和：
#include <stdint.h>
int16_t add_and_saturate(int16_t x, int16_t y)
{
int32_t z = (int32_t)x + y;
if (z > INT16_MAX)
{
z = INT16_MAX;
}
else if (z < INT16_MIN)
{
z = INT16_MIN;
}
return (int16_t)z;
}
您也可以这样做，即保持在 16 位类型的范围内，但它会变得有点棘手，使用更多的操作，并且我不能 100% 确信我的代码是正确的：
#include <stdint.h>
int16_t add_and_saturate(int16_t x, int16_t y)
{
int16_t z = x+y;
if ((y > 0) && (x > (INT16_MAX - y)))
{
z = INT16_MAX;
}
else if ((y < 0) && (x < (INT16_MIN - y)))
{
z = INT16_MIN;
}
return z;
}
处理乘法运算也是这样，并且实际上只使用类型扩展方法。
减法?
下面的 C 代码有个 bug，你能明白为什么吗？
int16_t calcPI(int16_t command, int16_t measurement)
{
int16_t error = command - measurement;
/* other stuff */
}
问题是减法。 如果 command = 32767 并且 measurement = -32768，那么误差的“实际”值是 65535 。如果 command = -32768 and measurement = 32767, 那么误差的“实际”值是 -65535。 刚和加法一样,在减法中两个 k-位数字操作的结果是 (k+1)-位数字。有几个方法可以避免不正确的结果。
首先，我们可以使用32-位中间计算值：
int32_t error = (int32_t)command - measurement;
这有个缺陷，就是加法操作的越多，输出值的范围就越大。但是只要中间计算值不陷入溢出，就没事。
其次，我们可以把误差值规范到 int16_t 的区间，所以上下限为-32768 到 +32767。这有个缺陷，就是极小边缘值和极大边缘值没什么区别。 (比如: command = 32767 和 measurement = 0 得出 error = 32767，但是 measurement = -32768 也是这个结果。) 实际上， 我们会想在整个区间保持线性关系（线性增长或线性缩小）。
第三，我们可以修改这个计算操作，使得计算结果在上下限范围内：
int16_t error = ((int32_t)command - measurement) / 2；
这个减少了净增益，但是避免了溢出和饱和。
除法?
除法是算术运算中最难看的东西。幸运的是，溢出情况并不太复杂。
如果您使用的是 C ，那么这里的除法意味着取一个n位数字并除以另一个n位数字，那么只有一个例子可以导致溢出，稍后我们将讨论这个问题。（它不被零除。如果在不首先排除d=0的情况下对n/d进行除法，那么你会得到任何计算结果。）
一些处理器具有内置的除法功能，编译器通过直接映射到处理器的除法特性所拥有的“内置”或“固有”功能来支持这些运算。在这里，您可以做一些事情，比如用一个 32 位整数除以一个16位整数，得到一个 16 位的结果，但是您需要确保避免用一个会产生一个不符合16位结果的商的数字来做除法。例如：180000/2=90000，90000 超出了 16 位数字（有符号或无符号）的界限。
移位?
不要忘了移位操作符 << 和 >>。这些不会引起溢出，但是在 C 语言，不要忘了讨厌的未定义行为和具体实现的行为。以下是从 C99 标准中与之相关6.5.7 节的 3-5 段：
3 对每个操作数执行整数提升。结果类型是左操作数提升后的类型。如果右操作数的值是负数或者大于等于提升后的左操作数宽度，那么操作行为是未定义的。
4 E1 << E2 的结果是 E1 左移 E2 个位；空出的位用 0 补充。如果 E1 是非负类型，那么结果是 E1 × 2^E2，是对结果类型所能表示的最大值进行一个或多个取模的值。如果 E1 是有符号的并且是非负值，那么 E1 × 2^E2 就是结果类型可表示的，并且该值就是结果值；否则这个行为未定义。
5 E1 >> E2 的结果是 E1 右移 E2 个位。如果 E1 是无符号类型或者 E1 是有符号类型并且是非负值，那么结果值就是 E1 / 2^E2 的整数部分。如果 E1 是有符号类型并且是复数，那么结果值未定义。
就是这样！如果一个有符号整数且其值为负数的 E1，进行右移位，那么这结果实现是已经定义的。“明智的”做法是算术右移，将 1 位放入最重要的位置来处理符号扩展。 但是 C 标准允许编译器产生一个逻辑右移，并把 0 位放在最重要的位上，而这个很可能不是你想要的结果。
最近的一些语言，像 Java 和 Javascript 提供了两个右移操作符：常规右移操作符 >>用来计算算术右移的，而 >>> 操作符用来计算逻辑右移来产生一个无符号二进制整数，并将 0 移动到最重要的位置。### 就这样吗? 噢!
不是，不仅仅是这样！还有其他什么吗？是的，还有自增和自减操作符，但是使用 -fwrapv 标志，它们会按预期使用截断语义（wraparound semantics）。
其余的溢出项可能属于所谓的病理例子。这些都涉及到宏常量 INT_MIN 和 INT_MAX 之间的不对称，这会导致不必要的混叠。
我们之前讨论过除法，并且当 -32768 除以 -1 时整数除法会溢出。
# Let's divide -32768 / -1
a=np.int16(-32768)
b=np.int16(-1)
a/b
-32768
噢! 我们应该得到 +32768，但是这不符合16-位有符号整数的取值范围，所以得到别名为 int16_t 位的同等值为 -32768。
相同的事情会发生在一元减法上：
# Let's negate -32768
-a
-32768
然后在同一条直线上还有定点乘法的缺陷。假设您正在做 Q15（小数点后有 15 位） 数学，其中整数表示形式为 \(\frac{k}{2^{15}}\)。C代码是这样的：
int16_t a = ...;
int16_t b = ...;
int16_t c = ((int32_t)a * b) >> 15;
这个能正常运算而不会溢出吗？是的，除了一种情况。\(-32768 \times -32768 = 2^{30}\)，并且如果我们向右移动 15 位，得到 \(2^{15} = 32768\)。但是对一个有符号 16-位整数而言，32768 是 -32768 的别名。哟！
我们该怎么处理这个？
Q15 乘法
Q15 乘法比较难。如果您绝对确定永远不会计算 -32768×-32768，那么继续使用通常的 C 代码。其中一个例子就是 PI 反馈控制，其中增益总是非负的。或者如果你知道其中一个数字的范围不包括 -32768，那么没问题。
或者，你可以做一个额外的右移位，像这样：
int16_t a = ...;
int16_t b = ...;
int16_t c = ((int32_t)a * b) >> 16;
在嵌入式系统中，右移 16 的操作通常比右移 15 的指令周期快，因为它通常映射到汇编指令的“抓取高字”操作，您在存储内存时通常就可以免费获取该操作。如果 a 和 b 是 Q15 定点整数，则该操作表示 \(a \times b \times \frac{1}{2}\)，或者如果其中一个值是 Q15，另一个值是 Q16，那么该操作表示 \(a \times b\)。
除此之外，你需要某种方法使中间值饱和，这样就可以安全地右移位：
int16_t a = ...;
int16_t b = ...;
int16_t c = sat30((int32_t)a * b) >> 15;
其中，sat30(x) 将结果限制在 \(-2^{30} + k_1\) 到 \(2^{30} - 1 - k_2\)，其中 k1 和 k2 是 0 到 32767 的任意数。（为什么这么模棱两可？因为最重要的 15 位是没关系的，并且某些处理器可能会有小机制允许你选择不同的 k 值使代码的执行速度更快。比如 \(2^{30} - 1\) 不可能用一条指令执行完毕，而 \(2^{30} - 32768\)= 32767 << 15 可能可以用一条指令来执行。）
一元减法
一元减法的情况有些类似。我们有几个选项：
如果我们完全确认输入不会是 -32768，那么我们就可以不用管它。
否则我们必须对 -32768 进行判断，然后转成 +32767。
或者，在 C 语言我们可以使用位补码操作符，因为对于有符号整数x = -x-1：~-32768 = 32767, ~0 = -1, ~-1 = 0, ~32767 = -32768。所有输入被 1 进行调整在某些应用可能不能接受，但是在其他应用这样做没问题，这只是表示一个小的偏差，并且补码操作往往是一个非常快的单指令操作。
相同的思想应用于绝对值，代替
#define ABS(x) ((x) < 0 ? (-x) : (x))
我们可以使用位补码操作符：
#define ABS_Q15(x) ((x) < 0 ? (~x) : (x))
要么我们对 -32768 特殊例子做判断。
只是一个提醒：#define 预定义宏对编译器是不可见的（因为它们出现在预编译期中），并且使用边际效应的参数是不安全的，比如上述定义的宏定义 ABS(++x) 会增加三次。
