2.1.2 double不够准

C/C++有两种浮点数据类型：float型和double型，它们是进行浮点运算的基础。每个程序设计员在学习C/C++时都被告知：float型变量的取值范围是-3.4e-38～3.4e+38，有效数字包含小数点后6位；double型变量的取值范围是1.7e-308～1.7e+308，有效数字包含小数点后15位。但是，在算法设计和问题求解时，往往会忘记浮点数据类型在有效数字方面的限制，导致错误的结果，如下面的两段程序。

程序2-3运行后得到如下结果：

程序2-4运行后得到如下结果：

为什么输出1.2345却得到结果1.2344999999999999？为什么1.234f+1.2345f≠2.4685f？计算机浮点数的表示精度是产生上述现象的本质原因。运用float和double进行浮点运算，大家应该记住：“精确是偶然的、误差是必然的”。在计算机中，float型和double型数据都是通过二进制数字序列进行编码和表示，其中float型是32位，double型是64位。从组合数学角度分析，float型最多能构成232个离散的实数（对应32位0-1序列的个数），但是实数是连续的，任意两个不相等的实数之间存在无穷多个实数。显然，从float型和double型的表示能力来说，它们不能精确表示所有实数，只能近似表示相应的实数。下面详细阐述C/C++中浮点数的表示方法。

1．浮点数的存储格式

在计算机系统的发展过程中，曾经提出过多种表示实数的方法，但是到目前为止使用最广泛的是浮点表示法。相对于定点数而言，浮点数利用指数使小数点的位置可以根据需要而上下浮动，从而可以灵活地表示更大范围的实数。

浮点数表示法利用科学计数法来表达实数。通常，将浮点数表示为±d.dd…d×βe，其中d.dd…d称为有效数字（Significant），它具有p个数字（称p位有效数字精度），β为基数（Base），e为指数（Exponent），±表示实数的正负。更精确地，±d₀.d₁d₂…d_p−1×βe表示以下数字，即

显然，对实数的浮点表示仅作如上规定是不够的，因为同一实数的浮点表示还不是唯一的。例如，1.0×102、0.1×103和0.01×104都可以表示100.0。为了达到表示单一性的目的，有必要对其进一步规范。规定有效数字的最高位（即前导有效位）必须非零，即0<d₀<β。符合该标准的数称为规格化数（Normalized Numbers），否则称为非规格化数（Denormalized Numbers）。

电子电气工程师协会（Institute of Electrical and Electronics Engineers，IEEE）1985年制定了IEEE 754二进制浮点运算规范，它是浮点运算部件事实上的工业标准。一个实数R在IEEE 754标准中可以用R=(−1)s×M×2E的形式表示，说明如下：

☢ 符号s（Sign）决定实数是正数还是负数，s为0表示该实数是正数，为1则表示负数。注意，数值0的符号位属于特殊情况，有特定的表示方法。

☢ 有效数字M（Significant）是二进制小数，M的取值范围为1≤M≤2或0≤M<1。

☢ 指数E（Exponent）是2的幂，它的作用是对浮点数加权。

C/C++中的浮点数本质上对应一种数据结构，它定义了浮点数在计算机中的内部表示，规定了构成浮点数的字段及其布局、算术解释。IEEE定义的754浮点数的数据位划分为3个字段，如图2-1所示，各参数值的编码规则如下：

☢ 一个单独的符号位直接编码符号s。

☢ k位的偏置指数e，记为e=e_k−1…e₁e₀，它是编码指数E的移码表示：E=e−(2k−1−1)。

☢ n位的小数f，记为f=f_n−1…f₁f₀，编码有效数字M，以原码表示。

根据偏置指数e的值，被编码的浮点数可分为3种类型。

（1）规格化数

当有效数字M在范围1≤M<2中且指数e的位模式e_k−1…e₁e₀既不全是0也不全是1时，浮点格式所表示的数都属于规格化数。这种情况下，小数f（0≤f<1）的二进制表示为0.f_n−1…f₁f₀。有效数字M=1+f，即M=1.f_n−1…f₁f₀（其中，小数点左侧的数值位称为前导有效位）。我们总是能调整指数E，使得有效数字M在范围1≤M<2中，这样有效数字的前导有效位总是1，因此该位不需显式表示，只需通过指数隐式给出。

需要特别指出的是，指数E要加上一个偏置值（Bias），转换成无符号的偏置指数e，也就是说，指数E要以移码的形式存放在计算机中。e、E和Bias三者的对应关系为e=E+Bias，其中Bias=2k−1−1。

（2）非规格化数

当指数e的位模式e_k−1…e₁e₀全为0（即e=0）时，浮点格式表示的数是非规格化数。这种情况下，E=1−Bias，有效数字M=f=0.f_n−1…f₁f₀，有效数字的前导有效位为0。

非规格化数的引入有两个目的：一是它提供了一种表示数值0的方法，二是它可用来表示那些非常接近于0.0的数。

（3）特殊数

当指数e的位模式e_k−1…e₁e₀全为1，小数f的位模式f_n−1…f₁f₀全为0（即f=0）时，该浮点格式所表示的值表示无穷，s为0时是+∞，s为1时是−∞。

当指数e的位模式e_k−1…e₁e₀全为1，小数f的位模式f_n−1…f₁f₀不全为0（即f_n−1,…,f₁,f₀至少有一个非零，即f≠0）时，该浮点格式表示的值被称为NaN（Not a Number）。

① 单精度浮点数（float）存储格式。

对于float型，IEEE 754标准规定用32位二进制编码表示，具体如下：

☢ 最高位D₃₁，保存符号位s，为0表示正数，为1表示负数。

☢ D₃₀…D₂₃，共8位，移码方式（指数值加上偏移量127）保存指数部分，称为阶码。

☢ D₂₂…D₀，共23位，保存系数部分，称为尾数，对于规范化二进制数，整数位的前导“1”不保存（隐含），直接保存小数部分f₂₂…f₁f₀。

并有以下规定：

☢ 若1≤E≤254，则R=(−1)s×(1+0.M)×2E−127，此时R为规范化数。

☢ 若E=0且M=0，则R=0。

☢ 若E=0且M≠0，则R=(−1)s×(0.M)×2−126，此时R为非规范化数。

☢ 若E=255且M≠0，则为非数值的编码组合（用NaN表示）。

☢ 若E=255且M=0，则表示R溢出（用+INF、-INF分别表示正、负无穷大）。

② 双精度浮点数（double）存储格式。

对于double型，IEEE 754标准规定用64位二进制编码表示，具体如下：

☢ 最高位D₆₃，保存符号位s，为0表示正数，为1表示负数。

☢ D₆₂…D₅₂，共11位，移码方式（指数值加上偏移量1023）保存指数部分，称为阶码。

☢ D₅₁…D₀，共52位，保存系数部分，称为尾数，对于规范化二进制数，整数位的前导“1”不保存（隐含），直接保存小数部分f₅₁…f₁f₀。

并有以下规定：

☢ 若1≤E≤2046，则R=(−1)s×(1+0.M)×2E−1023，此时R为规范化数。

☢ 若E=0且M=0，则R=0。

☢ 若E=0且M≠0，则R=(−1)s×(0.M)×2−1022，此时R为非规范化数。

☢ 若E=2047且M≠0，则为非数值的编码组合（用NaN表示）。

☢ 若E=2047且M=0，则表示R溢出（用+INF、-INF分别表示正、负无穷大）。

2．浮点数的有效数字

理论上，有限位数的二进制浮点数必能转换为有限位数的十进制数，但有限位数的十进制浮点数转换为二进制数不能保证是有限位数，且多数情况下不是有限位数。本质上，十进制浮点数在计算机内存中保存的并不是具体的值，而是一个计算其近似值的如下公式，即

R=(−1)s×M×2E

s、M和E的值如浮点数的存储格式所定义。因此自然有一个如何看待和精确使用数值的问题，即有效数字位数的问题。尤其对于各种数值计算，必须掌握其有效数字位数。下面从相对误差（误差绝对值与数值绝对值的比值）的角度给出一般性的结论。

对于具有严格的n位有效数字的十进制数D的最大绝对误差（用λ表示）为末位数字的半个单位，其最大相对误差为λ/D，显然该值不超过0.5×10−n。反之，可得出相对误差不大于等于0.5×10−n的十进制数必能提供至少n位有效数字。因此，从相对误差的大小就可确定：二进制数与其能提供的十进制数的有效数字位数的关系。

对于能转换为float型规范化二进制编码的十进制数值，其二进制编码的尾数为23位（包括一位隐含“1”），按“最近舍入方式”，则其保存的二进制浮点数的最大绝对误差为0.5×2−23×2E−127，尾数值为最小值1时的相对误差可取得最大值，为0.5×2−23=2−24。尾数位码值为全1时（对应的尾数为2−2−23，约为2）的相对误差可取得最小值，约为2−25；总之，float型规范化数保存的二进制浮点数的相对误差不大于等于2−24成立。在不考虑数据参与运算导致的误差传递等情况下（浮点数的机内运算采取了增加附加位的做法，一般能保证运算过程不降低数据精度），将保存的float型二进制数完整转换为十进制数（必定是有限位数）时，其带来的十进制数的相对误差不大于等于2-24（约0.596×10−7）。对于float型数据能提供相对误差不超过2-24（约0.596×10−7）的十进制数，由前述的推论可得到保守的结论：其十进制数的有效数字不少于6位。另外，为了充分使用float型的有效数字，向计算机提供的十进制实数尽量为8位有效数字（超出8位的数字在转换时不起作用），以使计算机能更精确地转换，使误差最小。

对于double型浮点数，其二进制编码的尾数为53位（包括一位隐含“l”），保存的二进制相对误差不超过2-54（约为0.555×10−16）。十进制数的16位有效数字的最大相对误差不超过0.5×10−16，因此，double型数据能提供的十进制数的有效数字不少于l5位（这也是默认情况下C/C++提供l5位有效数字的原因）。另外，使用double型数据时，向计算机提供的十进制实数尽量为17位有效数字（超出17位的数字在转换时不起作用），以使double型的保存误差最小。

3．高精度浮点数处理实例

既然C/C++的浮点数本质上表示的是一个近似数，在问题求解时，如果对于数据的精度要求特别高，如对有效数字位数要求超过了20位，就需要设计和编写自己的高精度浮点数算法。尤其是在数值计算和科学计算中，高精度浮点数算法是必不可少的。下面看一个简单的例子。

【例2-2】 复利计算问题

问题描述：从投资的角度来看，以复利计算的投资报酬效果是相当惊人的，许多人都知道复利计算的公式：F=P×(1+R)N，其中F表示本利和，P表示本金，R表示利率，N表示期数。虽然复利公式并不难懂，但若是期数很多，本金很大，计算相当麻烦。

输入：本金P（0≤P≤1010），利率R（0<R<1），最长4位有效数字，期数N，1<N<40。

输出：投资N期后的本利和。

复利计算问题的有效数字最长可达4×40 = 160位，显然，它的有效数字长度超出了double的长度范围，而且本利和也可能超出long的数据范围，因此需要编写高精度的浮点算法。幂运算可以转换为乘法运算，此问题本质上就是一个高精度浮点数乘法问题。下面给出高精度浮点数乘法运算的设计思想和原理。

首先，建立如下高精度浮点数的数据结构。

该定义类似IEEE 754定义的浮点数存储格式，但是superFloat是以10为底的科学计数法表示。

然后，设计无精度损失的高精度浮点数乘法。为了不损失计算精度，可以将高精度浮点数乘法转换为大整数加法实现，基本算法流程如下：

① 将数符较多的乘数表示为X，在算法中作为加数，数符较少的乘数表示为Y，在算法中控制加法次数。

② 如果Y含有小数部分，则将Y转变为纯整数Y_D，并记录小数点的右移位数I。

③ 初始化返回值T为0，取得Y的位数Width，设计数器Count为0。

④ 取Y右侧第Count+1位，以此数字为次数加X次，再左移Count位，加到T中。

⑤ 把Count加1。

⑥ 若Count等于Width，则转下一步，否则转第4步。

⑦ 将T中的小数点左移I位。

⑧ 返回T，得到乘法结果。

本算法的特点是加法次数少，若Y的宽度为W，最多进行9×W次加法和W次移位即可。

大整数加法addBigInt的定义见程序2-2。

memset的函数原型是void * memset(void *ptr, int value, size_t num)，在内存块ptr的前num字节中填充给定值value，它是对较大的结构体或数组进行清零操作最快的一种方法。注意：memset是一种按位赋值的操作。

memcpy的函数原型是void *memcpy(void *dest, const void *src, size_t n)，它按位复制src所指的内存内容前n字节到dest所指的内存中。