深度理解IEEE754浮点数

浮点数不是小数，计算机只能用有限的浮点数来表示无限多的小数，浮点数标准是IEEE754。

两种主要的浮点数格式，分别为单精度浮点数和双精度浮点数，即float和double，他两的格式如下：

Float:
 sign |  exponent  |     mantissa     |
 1bit |   8bits    |      23bits      |
Double:
 sign |   exponent   |           mantissa           |
 1bit |    11bits    |            52bits            |

浮点数最高位固定为符号位，0正1负。浮点数没有unsigned概念！所有浮点数都是signed！float占用32bit，double占用64bit。

浮点数名称的由来在于小数点可以浮动，但在具体存储时，需要固定一种形式。IEEE754规定，在二进制数中，通过移位，将小数点前面的值固定为1。IEEE754称这种形式的浮点数为规范浮点数（Normal Number）。

尾数用原码表示，在表示Normal Number时，第一位总是1，因而可在尾数中省略，称为隐藏位，这样使得单精度格式的23bits尾数实际表示了24bits有效数字，双精度格式的52bits的尾数实际上表示了53bit有效数字。IEEE754规定隐藏位1的位置在小数点之前。

因为指数有正有负，为了避免使用符号位，同时方便比较和排序，指数部分采用了 Biased Exponent（有偏指数，也称为移码）。IEEE754规定，\(2^{e-1}-1\)对应的值表示0（其中e表示指数部分的bit位数），小于这个值表示负数，大于这个值表示正数。因此，对于单精度浮点数而言，指数0表示为\(2^{8-1}-1=127=0b01111111\)，双精度浮点数指数0表示为\(2^{11-1}-1=1023=0b01111111111\)。

特殊浮点数

在浮点数格式中，指数部分E有特殊作用，当其取特殊值时，表示特殊浮点数，IEEE754的具体规定如下：

浮点数的分布

Normalized Number有隐藏位1，而且指数部分可变动，因此这部分的浮点数，间隔并不均匀。越是远离中心点0，前后两个浮点数的间隔就越大，因为此时的指数部分值较大。而Denormalized Number，没有隐藏位，只能用固定的最小指数来表达变动的尾数，因此这部分浮点数之间的间隔均匀。

浮点数的二进制表达

下面是一段精心调试过的，在little endian平台（x86就是）上运行的C代码，可以很方便的将一个float或double数的二进制表示打印出来：

#include <stdio.h>
#include <stdlib.h>
#include <string.h>


#define USAGE_RETURN \
    do{ \
        printf("usage: ./float_show f|d 123.456\n"); \
        return 0; \
    }while(0)


int main(int argc, char **argv) {
    float f;
    double d;
    char *p;
    int nbyte;

    if(argc != 3)
        USAGE_RETURN;

    printf("input: %s %s\n", argv[1], argv[2]);
    if(strcmp(argv[1],"f") == 0) {
        f = strtof(argv[2], NULL);
        printf("float: %.50f\n", f);
        p = (char *)&f;
        nbyte = sizeof(float);
    }
    else if(strcmp(argv[1],"d") == 0) {
        d = strtod(argv[2], NULL);
        printf("double: %.50f\n", d);
        p = (char *)&d;
        nbyte = sizeof(double);
    }
    else
        USAGE_RETURN;

    char cc[nbyte*8];  // VLA
    int k = 0;
    for(int i=nbyte-1; i>=0; --i){  // little endian
        char data = p[i];
        unsigned char pbit = 0x80;
        for(int j=0; j<8; ++j){
            if(data & pbit)
                cc[k++] = '1';
            else
                cc[k++] = '0';
            pbit >>= 1;  // logical right shift for unsigned
        }
    }

    printf("S  E  M:\n");
    if(strcmp(argv[1],"f") == 0)
        for(int i=0; i<nbyte*8; ++i) {
            if ((i==1) || (i==9)) printf("  ");
            printf("%c", cc[i]);
        }
    else if(strcmp(argv[1],"d") == 0)
        for(int i=0; i<nbyte*8; ++i) {
            if ((i==1) || (i==12)) printf("  ");
            printf("%c", cc[i]);
        }
    printf("\n");

    return 0;
}

$ gcc -Wall -Wextra float_show.c -o float_show
$ ./float_show f 0.15625
input: f 0.15625
float: 0.15625000000000000000000000000000000000000000000000
S  E  M:
0  01111100  01000000000000000000000

这段代码使用strtof和strtod来转换浮点数，这两个C接口在遇到无法转换的字符串时，会输出全0，不会报错。而且，这两个接口可以识别inf和nan（字符串）。请看下面的测试。

$ ./float_show f 99999999999999999999999999999999999999999999999999999999999999
input: f 99999999999999999999999999999999999999999999999999999999999999
float: inf
S  E  M:
0  11111111  00000000000000000000000
$ ./float_show f -99999999999999999999999999999999999999999999999999999999999999
input: f -99999999999999999999999999999999999999999999999999999999999999
float: -inf
S  E  M:
1  11111111  00000000000000000000000

前面提到了strtof和strtod可以识别inf和nan这样的字符串，我们来试一下：

$ ./float_show d inf
input: d inf
double: inf
S  E  M:
0  11111111111  0000000000000000000000000000000000000000000000000000
$ ./float_show d -inf
input: d -inf
double: -inf
S  E  M:
1  11111111111  0000000000000000000000000000000000000000000000000000

NaN也可以有正负之分，但NaN就不是数字，正负又有何意义呢？！

$ ./float_show d nan
input: d nan
double: nan
S  E  M:
0  11111111111  1000000000000000000000000000000000000000000000000000
$ ./float_show d -nan
input: d -nan
double: -nan
S  E  M:
1  11111111111  1000000000000000000000000000000000000000000000000000

inf和nan这样的字符串，可以被打印出来，可以用strtof转换，但是在代码中，并不能使用！

$ ./float_show f 0
input: f 0
float: 0.00000000000000000000000000000000000000000000000000
S  E  M:
0  00000000  00000000000000000000000
$ ./float_show f -0
input: f -0
float: -0.00000000000000000000000000000000000000000000000000
S  E  M:
1  00000000  00000000000000000000000

有了这个程序，我们就可以很容易理解，为什么0.1+0.2=0.300000004470348358154296875，请看：

$ ./float_show f 0.1
input: f 0.1
float: 0.10000000149011611938476562500000000000000000000000
S  E  M:
0  01111011  10011001100110011001101
$ ./float_show f 0.2
input: f 0.2
float: 0.20000000298023223876953125000000000000000000000000
S  E  M:
0  01111100  10011001100110011001101

浮点数并不能精确的表达输入的0.1和0.2，因此计算得到的0.3也不精确！

浮点数的最大最小值

下面这段代码（运行在x86_64这样的小端平台），打印出了浮点数的最大最小值，同时打出了紧挨着的两个浮点数：

#include <stdio.h>

int main(void) {
    float f;
    char *p = (char *)&f;
    // max
    *(p+3) = 0b01111111;
    *(p+2) = 0b01111111;
    *(p+1) = 0b11111111;
    *p     = 0b11111111;
    printf("max float: %f\n", f);
    *p     = 0b11111110;
    printf("next one:  %f\n", f);
    *p     = 0b11111101;
    printf("next one:  %f\n", f);
    // min
    *(p+3) = 0b00000000;
    *(p+2) = 0b00000000;
    *(p+1) = 0b00000000;
    *p     = 0b00000001;
    printf("min float: %.50f\n", f);
    *p     = 0b00000010;
    printf("next one:  %.50f\n", f);
    *p     = 0b00000011;
    printf("next one:  %.50f\n", f);
    return 0;
}

$ gcc t3.c -o t3
$ ./t3
max float: 340282346638528859811704183484516925440.000000
next one:  340282326356119256160033759537265639424.000000  # diff:16
next one:  340282306073709652508363335590014353408.000000
min float: 0.00000000000000000000000000000000000000000000140130
next one:  0.00000000000000000000000000000000000000000000280260
next one:  0.00000000000000000000000000000000000000000000420390

浮点数之间是有间隔的，数值越大，紧挨着的两个数的间隔越大！浮点数表达的数的范围远大于long long类型，但是数值越大，数的分布就越稀松，越不精确！

浮点数的表达范围

有了前面的知识，我们来自己动手找出浮点数的表达范围，指数全1减1，尾数全1，或正或负。同样的，下面这段代码也只能在x86这样的Little Endian平台上正确运行。

#include <stdio.h>
#include <string.h>

int main(void) {
    float f;
    char *pf = (char *)&f;
    // max float
    memset(pf, 0xFF, 4);
    *(pf+3) = 0b01111111;
    *(pf+2) = 0b01111111;
    printf("float range: %g -- %g\n", -f, f);
    double d;
    char *pd = (char *)&d;
    // max double
    memset(pd, 0xFF, 8);
    *(pd+7) = 0b01111111;
    *(pd+6) = 0b11101111;
    printf("double range: %g -- %g\n", -d, d);
    return 0;
}

float range: -3.40282e+38 -- 3.40282e+38
double range: -1.79769e+308 -- 1.79769e+308

浮点数的舍入（Rounding）

IEEE754定义了4种浮点数rouding的方式，默认采用银行家舍入法，即Round to Even，参考：各种Rounding方法详解。

由于Rounding的存在，会导致多个不同的小数对应到一个相同的浮点数上的现象。比如：

$ ./float_show f 77777776
input: f 77777776
float: 77777776.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101110
$ ./float_show f 77777777
input: f 77777777
float: 77777776.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101110
$ ./float_show f 77777778
input: f 77777778
float: 77777776.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101110
$ ./float_show f 77777779
input: f 77777779
float: 77777776.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101110
$ ./float_show f 77777780
input: f 77777780
float: 77777776.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101110
$ ./float_show f 77777781
input: f 77777781
float: 77777784.00000000000000000000000000000000000000000000000000
S  E  M:
0  10011001  00101000101100101101111

77777776到77777780，这个范围的浮点数二进制表达完全一样。而当要表达77777781时，浮点数直接跳到了77777784。

由于Rounding的存在，还导致了浮点数计算不具有可结合性（associativity）和分布性（distributivity）。请看下面这个经典case：

#include <stdio.h>

int main(){
    float a = 3.14;
    float b = 1e10;
    float c = (a+b)-b;
    float d = a+(b-b);
    printf("c = %f\n", c);
    printf("d = %f\n", d);
    float e = 1e-20;
    float f = (b*b)*e;
    float g = b*(b*e);
    printf("f = %f\n", f);
    printf("g = %f\n", g);
    return 0;
}

c = 0.000000
d = 3.140000
f = inf
g = 100000002004087734272.000000

如果只是看代码，可能我们会认为c和d应该相等，f和g应该相等，但是实际上的计算结果如上！这就是浮点数Rounding带来的不一致。

#include <stdio.h>

int main(){
    float b = 1e20;
    float c = b*(b-b);
    float d = b*b-b*b;
    printf("c = %f\n", c);
    printf("d = %f\n", d);
    return 0;
}

c = 0.000000
d = -nan

要特别注意，在并行计算的场景下，由于浮点数计算的特性，可能会导致每一次计算后的结果都有些微的不同，这些结果是否正确，需要仔细分析确认！

慎用++

在经典的《The C Programming Language》书中，看到一个case，作者编写的一个非常简单的统计输入word数量的程序，没有使用传统的unsigned int作为计数器，而是用了一个double变量，并对double变量施加++这样的操作。float和double变量是可以施加++或--这样操作的！但是非常危险....

#include <stdio.h>

int main(){
    printf("++float:\n");
    float a = 0xFFFFFFFF;
    printf("%f\n", ++a);
    printf("%f\n", ++a);
    printf("%f\n", ++a);
    printf("%f\n", ++a);
    printf("++double:\n");
    double b = 0xFFFFFFFF;
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    return 0;
}

$ gcc test.c && ./a.out
++float:
4294967296.000000
4294967296.000000
4294967296.000000
4294967296.000000
++double:
4294967296.000000
4294967297.000000
4294967298.000000
4294967299.000000

用float表达最大unsigned int时，++操作已经没有任何变化（因为有Rounding），而double还可以继续增加！这应该就是作者想要表达的point。（通过分析汇编，++float实际上加的是一个浮点数表达的1）

由于这本经典的C语言教材，写作时间比较久远，那个时候，64bit的CPU还没有出现，sizeof(long)也只等于32。因此，为了简单解决overflow的问题，double变量成为了一个选项。float不行，是因为float的精度不够，看前面的内容可以发现，至少从77777776开始，float变量的+1就开始没有变化了。但我们要清楚，就算使用double变量，也会从某一个数字开始，+1后不再变化，只是这个数字比较大，远比unsigned int要大！但是大不过unsigned long long。如下测试：

#include <stdio.h>

int main(){
    double b = 0xFFFFFFFFFFFFFFFF;
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    printf("%f\n", ++b);
    return 0;
}

$ gcc test.c && ./a.out
18446744073709551616.000000
18446744073709551616.000000
18446744073709551616.000000
18446744073709551616.000000

看起来，当有了64bit的unsigned整数类型时，再使用double变量来做counting，就不对了！就算是在没有64bit的unsigned整数类型的环境下，也不推荐这样做，很容易出错，不如使用多个unsigned int来计数，代码只是稍微复杂一点点而已！

IEEE754-2008

Name	Common Name	Base	Exponent	Mantissa
binary16	Half precision	2	5	10
binary32	Single precision	2	8	23
binary64	Double precision	2	11	52
binary128	Quadruple precision	2	15	112
binary256	Octuple precision	2	19	236

理解Binary Float Point

其实，我们日常生活中使用的浮点数计算，或者说小数计算，都是 Decimal Float Point。而计算机使用的浮点数，是 Binary Float Point。

Decimal Float Point 计算的特点是，以1/10为单位进行细分，比如0.1是对1的细分，可以有0.01，0.001等等，表示将1的长度，等分成10份，100份，1000份等等...而Binary Float Point，由于表达受限，它的计算特点是，以1/2为单位进行细分，等分成2份，4份，8份，16份等等...而小数部分的具体值，就表示取其中的多少份！

虽然在有限位数的情况下计算，不管是decimal还是binary，都有精度问题，但是这两者的精度问题的表现是不一样的。当我们需要满足日常 Decimal Float Point 的计算要求时，就要注意了：

>>> from decimal import *
>>> Decimal('0.1') + Decimal('0.1') + Decimal('0.1')
Decimal('0.3')
>>>
>>> 0.1+0.1+0.1
0.30000000000000004

Python标准库中提供的decimal模块，可以辅助我们在程序中，按照 Decimal Float Point 的方式进行计算！

不管是十进制还是二进制，都可以按照基于位置的技术系统的方法，将数字展开：

\(x=x_{n-1}10^{n-1}+x_{n-2}10^{n-2}+...+x_010^0.x_{-1}10^{-1}+x_{-2}10^{-2}+...\)

十进制小数转二进制小数

转换可分为整数部分和小数部分，两部分分开算，最后用一个小数点合并在一起。整数部分的转换用除2法，小数部分用乘2法。下面举例说明：

将十进制13.35转成二进制

整数部分为13，
>>> divmod(13,2)
(6, 1)
>>> divmod(6,2)
(3, 0)
>>> divmod(3,2)
(1, 1)
>>> divmod(1,2)
(0, 1)
将4个除2的余数从下到上串起来，就是 1101。

小数部分为0.35，
0.35*2 = 0.7  整数部分为0
0.7*2 = 1.4   整数部分为1，拿掉整数部分继续
0.4*2 = 0.8   整数部分为0
0.8*2 = 1.6   整数部分为1，拿掉整数部分继续
0.6*2 = 1.2   整数部分为1，拿掉整数部分继续
0.2*2 = 0.4   整数部分为0
0.4*2 = 0.8   整数部分为0
... 开始循环了，需要rounding
从上到下将整数部分串起来，就是 0.010110011001100...

最后将两部分拼在一起，得到 1101.010110011001100...

前一节介绍的二进制展开的表达方式，就是理解这个转换方法的关键。对于整数部分，每一次除2得到的余数，就是拿到\(2^0\)的系数\(x_0\)。对于小数部分，每一次乘2得到的整数部分的值，就是拿到\(2^-1\)的系数\(x_{-1}\)。

BF16（bfloat16）

深度学习促使了人们对新的浮点数格式的兴趣。通常（深度学习）算法并不需要64位，甚至32位的浮点数精度。更低的精度可以使在内存中存放更多数据成为可能，并且减少在内存中copy数据的时间。低精度浮点数的电路也会更加简单。这些好处结合在一起，带来了明显了计算速度的提升。

bfloat16，BF16格式的浮点数已经成为深度学习事实上的标准。已有一些深度学习“加速器”支持了这种格式，比如Google的TPU。Intel的处理与在未来也会支持。

BF16浮点数在格式，介于FP16和FP32之间。（FP16和FP32是 IEEE 754-2008定义的16位和32位的浮点数格式。）

|--------+------+----------+----------|
| Format | Bits | Exponent | Fraction |
|--------+------+----------+----------|
| FP32   |   32 |        8 |       23 |
| FP16   |   16 |        5 |       10 |
| BF16   |   16 |        8 |        7 |
|--------+------+----------+----------|

BF16的指数位比FP16多，跟FP32一样，不过小数位比较少。这样设计说明了设计者希望在16bits的空间中，通过降低精度（比FP16的精度还低）的方式，来获得更大的数值空间（Dynamic Range）。

浮点数的比较

由于精度的问题，浮点数的比较，我们一般认为，当两个数的差的绝对值在一个很小的范围内时就判定为相等，当不相等时，再做大小的比较。示例代码如下：

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

int main(){
    float a = 1.0;
    float b = 0.9;
    // 0.00000011920928955078
    printf("Standard Epsilon: %.20f\n", FLT_EPSILON);

    // if b+0.1 == a
    if(fabs(b+0.1-a) < FLT_EPSILON)
        printf("==\n"); // show out

    float c = 1.01;
    // not identical then ...
    if((fabs(a-c)>FLT_EPSILON) && (a<c))
        printf("a<c\n");  // show out

    return 0;
}

在代码中处理Inf和NaN

number/0.0或log(0)等计算，会产生inf。有inf的计算表达式，或者开负数根号，会产生NaN。C标准库中有几个接口，可以用来判断inf或NaN。inf可以参与比较大小，但NaN不行（返回false）。神奇的是，他俩本身都是true。

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

int main(){
    float a = 1.0/0.0;  // +inf

    if(isinf(a))
        printf("inf\n");
    if(!isfinite(a))
        printf("not finite\n");
    if(!isnormal(a))
        printf("not normal\n");
    if(!isnan(a))
        printf("not nan\n");

    if(a > 123456.0)  // inf is comparable
        printf("+inf is the biggest!\n");
    if(a && -a)
        printf("inf is true!\n");

    double b = sqrt(-1.0); // NaN
    if(!isinf(b))
        printf("not inf\n");
    if(!isfinite(b))
        printf("not finite\n");
    if(!isnormal(b))
        printf("not normal\n");
    if(isnan(b))
        printf("nan\n");

    if((b<=1.0) || (b>=1.0)) // NaN is not comparable
        printf("not shown!\n");
    if(b && -b)
        printf("But NaN itself is true!\n");

    return 0;
}

浮点数的排序

#include <stdio.h>
#include <stdlib.h>
#include <assert.h>

static int comp_float(const void *a, const void *b){
    assert(sizeof(int) == sizeof(float));
    int sign_a = (*(int*)a) >> (sizeof(int)*8-1);
    int sign_b = (*(int*)b) >> (sizeof(int)*8-1);
    unsigned int ua = *(unsigned int*)a;
    unsigned int ub = *(unsigned int*)b;
    if(sign_a && sign_b)
        return ua<ub ? 1 : ua==ub ? 0 : -1;
    if(!sign_a&&!sign_b)
        return ua<ub ? -1 : ua==ub ? 0 : 1;
    return sign_a ? -1 : 1;
}

int main(){
    float f[9] = { -1.234,
                   1.235,
                   -2.15,
                   44.99,
                   -88.99,
                   0.56,
                   -9.0,
                   1.235,
                   55.5 };

    qsort(f, sizeof(f)/sizeof(float), sizeof(float), comp_float);

    for(int i=0; i<(int)(sizeof(f)/sizeof(float)); ++i)
        printf("%f\n", f[i]);

    return 0;
}

浮点数转整数

当浮点数的整数部分在目标整数类型的范围内时，这种转换就是truncate，直接去掉小数部分。而当浮点数的整数部分的值超过了目标整数类型的范围时，这是个危险的UB。我测试的结果是，如果是目标是int，值为最小负数，如果目标是uint，值为0。

$ cat test3.c
#include <stdio.h>

int main(void) {
    double a = 1.56789;
    int b = a;     // truncated to int
    printf("%d\n", b);
    int c = 50/3;  // truncated to int
    printf("%d\n", c);
    return 0;
}

$ gcc -Wall -Wextra test3.c -o test  # no warnings
$ ./test
1
16

这是显式的类型转换，而且是narrow type conversion，但没有编译告警，这就是C语言的特性，浮点数转整数，直接truncate，去掉小数部分，只留下整数部分。

所谓的Integer Math

当整数计算中带有除法，按前述介绍，truncate后会损失精度，如果转float计算，性能堪忧。这种情况，可以尝试一些所谓的integer math的技巧，在保持整数计算的情况下，尽量保持精度。下面介绍几个简单的与浮点数计算有关的思路。

只要9*C的结果不溢出，Rearrange后计算，能够有更好的精度。只损失一份精度，而不是C份。

F = (9/5)C + 32
  = 1.8C + 32            -- but we can't have 1.8, so multiply by 10

sum = 10F = 18C + 320    -- 1.8 is now 18: now all integer operations

F = sum/10

扩大10倍后，所有计算就都是integer operation了。只需要在最后还原。如果扩大的倍数是2的幂次，最后还可以用位移的方法还原，速度更快。

1.0C + 0.5C + 0.25C + ...   OR  C + (C >> 1) + (C >> 2) + ...

《CSAPP》Data Lab，浮点数部分

这三道题不算简单，需要对浮点数的二进制表达方式有深刻理解，输入输出全都是big endian。测试通过：

$ ./btest
Score   Rating  Errors  Function
...
 4      4       0       floatScale2
 4      4       0       floatFloat2Int
 4      4       0       floatPower2

//float
/*
 * floatScale2 - Return bit-level equivalent of expression 2*f for
 *   floating point argument f.
 *   Both the argument and result are passed as unsigned int's, but
 *   they are to be interpreted as the bit-level representation of
 *   single-precision floating point values.
 *   When argument is NaN, return argument
 *   Legal ops: Any integer/unsigned operations incl. ||, &&. also if, while
 *   Max ops: 30
 *   Rating: 4
 */
unsigned floatScale2(unsigned uf) {  // assume Big Endian
    unsigned int e = (uf&0x7F800000) >> 23;
    unsigned int m = (uf&0x007FFFFF);
    unsigned int s = (uf&0x80000000);
    if(e == 0)  // denormalized
        return (uf<<1) | s;
    if((e==255) && (m!=0))  // nan
        return uf;
    e += 1;
    if(e >= 255)  // inf
        return 0x7F800000|s;
    // normalized
    return s|(e<<23)|m;
}
/*
 * floatFloat2Int - Return bit-level equivalent of expression (int) f
 *   for floating point argument f.
 *   Argument is passed as unsigned int, but
 *   it is to be interpreted as the bit-level representation of a
 *   single-precision floating point value.
 *   Anything out of range (including NaN and infinity) should return
 *   0x80000000u.
 *   Legal ops: Any integer/unsigned operations incl. ||, &&. also if, while
 *   Max ops: 30
 *   Rating: 4
 */
int floatFloat2Int(unsigned uf) {
    unsigned int e = (uf&0x7F800000) >> 23;
    int s = (uf&0x80000000);
    unsigned int fm = (uf&0x007FFFFF) | 0x00800000;

    if(e == 255)
        return 0x80000000;
    if(e <= 126)
        return 0;

    if(s)
        s = -1;
    else
        s = 1;

    if(e <= 150)  // 150 = 127+23
        return s*(fm>>(23-(e-127)));
    if(e <= 157)  // 157 = 127+23+7
        return s*(fm<<(e-150));

    return 0x80000000;
}
/*
 * floatPower2 - Return bit-level equivalent of the expression 2.0^x
 *   (2.0 raised to the power x) for any 32-bit integer x.
 *
 *   The unsigned value that is returned should have the identical bit
 *   representation as the single-precision floating-point number 2.0^x.
 *   If the result is too small to be represented as a denorm, return
 *   0. If too large, return +INF.
 *
 *   Legal ops: Any integer/unsigned operations incl. ||, &&. Also if, while
 *   Max ops: 30
 *   Rating: 4
 */
unsigned floatPower2(int x) {
    int e = 0x80;  // e of float 2

    if(x == 0)
        return 0x3F800000;  // float 1

    if(x > 0){
        e += x - 1;
        if(e >= 255)
            return 0x7F800000;  // float inf
    }

    if(x < 0){
        e += x;
        if(e <= 0){
            int rss = (~e) + 1;
            if((rss>=23) || (rss<0))  // abs(0x80000000)==0x80000000
                return 0;
            return 0x00400000 >> rss;
        }
    }

    return e << 23;
}

浮点数与x86-64平台

早期的x86芯片，通过x87协处理器支持浮点数计算，后来x87才被集成到CPU内部，感谢摩尔定律。即使是现在（2023年底），Intel的CPU依然支持x87相关的指令。x87指令支持IEEE754的单精度和双精度浮点数，同时也支持一个非标的80bits的扩展精度。因此，我们在CPU架构的资料中，总能看到一组80bits的寄存器。现在编译器基本都已经不再使用x87的指令来编译浮点数运算，因为x87的浮点数计算采用stack架构，浮点数计算过程需要大量的memory操作，速度很慢。总的来说，x87和复用寄存器的MMX（第一代SIMD技术）都已经过时了。

现代编译器基本都采用SSE指令（甚至AVX）来处理scalar浮点数的计算。值得了解的一点信息是，SSE技术有一个32bits的控制和状态寄存器，叫做MXCSR。有两个bit叫做RC，Rounding Control，它有4个值，分别对应了IEEE754标准定义的4种rounding的方式。默认的00表示银行家舍入法。