21、深入理解浮点运算：原理、指令与应用

深入理解浮点运算：原理、指令与应用

1. 浮点运算基础

在之前处理数值时，我们通常使用整数，但整数无法表示分数。为了表示分数，我们引入了浮点数。例如在C语言中，使用 float 和 double 数据类型来处理浮点数。

整数运算具有精确性，比如两个整数相加总能得到精确结果。而浮点数运算则存在舍入误差，结果往往是近似值。不过，浮点数也有诸多优势，它能表示极小和极大的数值，采用科学记数法，将数值分为符号、尾数和指数三部分。符号位用于标识数值正负（0为正，1为负），数值大小由 magnitude = mantissa x 2^exponent 计算得出。

计算机系统中浮点数的实现会因效率或遵循标准的不同而有所差异。像Intel 32位处理器等多数处理器遵循IEEE 754浮点标准，这种标准有助于不同计算机系统间的数据交换和高效数值软件库的编写。

浮点单元（FPU）支持三种浮点数格式，其中两种用于外部，一种用于内部。外部格式有单精度（32位）和双精度（64位），在C语言中分别对应 float 和 double ；内部格式为80位的扩展格式，FPU的所有内部寄存器都是80位，以便存储扩展格式的浮点数。

为增强处理器的数值运算能力，可使用专门硬件进行浮点运算。早期的80X87数值协处理器与80X86系列处理器配合使用，如8087协处理器为8086和8088处理器提供了高速数值处理能力，相比8086处理器上的等效软件功能，执行时间提高了约百倍。从80486处理器开始，浮点单元被集成到处理器内部，无需外部数值处理器。

2. 浮点单元组织

浮点单元提供了多个寄存器，可分为数据寄存器、控制和状态寄存器以及指针寄存器三组。由于指针寄存器涉及的编程异常处理程序超出了本文范围，这里重点介绍前两组寄存器。

2.1 数据寄存器

FPU有八个浮点寄存器用于存放浮点操作数，为浮点指令提供必要的操作数。与处理器的通用寄存器（如EAX和EBX）不同，这些寄存器被组织成寄存器栈。我们可以使用 ST0 、 ST1 等单独访问这些寄存器，但这些名称并非静态分配， ST0 指的是栈顶（TOS）寄存器，后续依次为 ST1 到 ST7 。状态寄存器中有一个3位的栈顶指针用于标识TOS寄存器。

每个数据寄存器可存储扩展精度的浮点数（80位），而非单精度（32位）或双精度（64位）格式。这是因为这些寄存器通常保存中间结果，使用扩展格式可提高最终结果的准确性。每个寄存器的状态和内容由一个2位的标记字段指示，八个寄存器共需16位标记，这些标记存储在标记寄存器中。

2.2 控制和状态寄存器

这组寄存器包括三个16位寄存器：控制寄存器、状态寄存器和标记寄存器。
- FPU控制寄存器 ：用于为程序员提供多种处理选项的控制。其最低六位包含六个浮点异常的掩码， PC 和 RC 位分别控制精度和舍入，各用两位指定四种可能的控制。舍入控制选项如下：
- 00 — 四舍五入到最接近的值
- 01 — 向下舍入
- 10 — 向上舍入
- 11 — 截断
精度控制可将内部操作精度设置为低于默认精度，以兼容早期精度较低的FPU。精度选项如下：
- 00 — 24位（单精度）
- 01 — 未使用
- 10 — 53位（双精度）
- 11 — 64位（扩展精度）
- FPU状态寄存器 ：该16位寄存器记录FPU的状态。四个条件码位（ C0 - C3 ）会根据浮点算术运算结果更新，类似于处理器的标志寄存器。其中三个位与处理器标志寄存器的对应关系如下：
| FPU标志 | CPU标志 |
| ---- | ---- |
| C0 | CF |
| C2 | PF |
| C3 | ZF |
缺失的 C1 位用于指示栈下溢/上溢情况。这些位可用于条件分支，使用时需将状态字复制到CPU标志寄存器，分两步完成：先使用 fstsw 指令将状态字存储在 AX 寄存器，再使用 sahf 指令将这些值加载到标志寄存器，之后即可使用条件跳转指令。

状态寄存器用三位维护栈顶（TOS）信息，八个浮点寄存器被组织成循环缓冲区，TOS标识栈顶寄存器，其值会随栈的压入和弹出操作更新。最低六位给出六个异常的状态，无效操作异常可能由栈操作或算术操作引起，栈故障位可指示无效操作的原因。我们还可以使用 C1 位进一步区分栈下溢（ C1 = 0 ）和上溢（ C1 = 1 ）。溢出和下溢异常在数值过大或过小时发生，通常在执行浮点算术指令时出现。精度异常表示操作结果无法精确表示，如表示1/3这样的分数时。除零异常类似于处理器产生的除零错误异常，非正规异常在算术指令对非正规操作数进行操作时产生。
- 标记寄存器 ：该寄存器存储数据寄存器的状态和内容信息。每个寄存器用两位表示以下信息：
- 00 — 有效
- 01 — 零
- 10 — 特殊（无效、无穷大或非正规）
- 11 — 空
最低两位用于 ST0 寄存器，后续依次类推。该标记字段可识别相关寄存器是否为空，若不为空，则标识其内容为有效数字、零或特殊值（如无穷大）。

3. 浮点指令

FPU提供了多种浮点指令，用于数据移动、算术运算、比较和超越运算等，还有加载常用常量（如π）和处理器控制字的指令。下面介绍一些常见指令。

3.1 数据移动指令

数据移动由加载和存储两种类型的指令支持。
- 加载指令 ：通用加载指令格式为 fld src ，将 src 压入FPU栈，即递减TOS指针并将 src 存储在 ST0 。 src 操作数可以在寄存器或内存中，若在内存中，可以是单精度（32位）、双精度（64位）或扩展（80位）浮点数。由于寄存器以扩展格式存储数字，单精度或双精度数在存储到 ST0 之前会转换为扩展格式。此外，还有一些将常量压入栈的指令，如下表所示：
| 指令 | 描述 |
| ---- | ---- |
| fldz | 将 +0.0 压入栈 |
| fld1 | 将 +1.0 压入栈 |
| fldpi | 将 π 压入栈 |
| fldl21 | 将 log₂10 压入栈 |
| fldl2e | 将 log₂e 压入栈 |
| fldlg2 | 将 log₁₀2 压入栈 |
| fldln2 | 将 logₑ2 压入栈 |
加载整数可使用 fild 指令， src 操作数必须是位于内存中的16位或32位整数，该指令将整数转换为扩展格式并压入栈。
- 存储指令 ：格式为 fst dest ，将栈顶值存储在 dest 。 dest 可以是FPU寄存器或内存，内存操作数可以是单精度、双精度或扩展浮点数。若 dest 是单精度或双精度操作数，寄存器值会转换为目标格式。需要注意的是，该指令不会从栈中移除值，仅复制其值。若要复制并弹出栈顶值，可使用 fstp dest 指令。还有整数版本的存储指令 fist dest ，将 ST0 中的值转换为有符号整数并存储在内存中的 dest 位置，转换时使用 RC （舍入控制）字段。其弹出版本 fistp dest 执行类似转换并弹出栈顶值。

3.2 算术指令

• 加法指令 ：基本加法指令格式为 fadd src ，将内存中（ src 处）的浮点数与 ST0 中的数相加，并将结果存储回 ST0 。 src 处的值可以是单精度或双精度数，该指令不会弹出栈。双操作数版本 fadd dest, src 允许指定目标寄存器， src 和 dest 必须是FPU寄存器，同样不弹出栈，若要弹出可使用 faddp dest, src 指令。还可以使用 fiadd src 指令添加整数， src 是位于内存中的16位或32位整数。
• 减法指令 ：减法指令格式与加法类似。 fsub src 执行 ST0 = ST0 - src 操作。双操作数版本 fsub dest, src 执行 dest = dest - src ，也有弹出版本 fsubp dest, src 。由于减法不满足交换律，还有反向减法操作 fsubr src ，执行 ST0 = src - ST0 ，同样有双操作数和弹出版本。若要减去整数，可使用 fisub 进行标准减法，或 fisubr 进行反向减法，整数必须位于内存中。
• 乘法指令 ：乘法指令有多种版本，内存操作数版本 fmul src 将内存中（ src 处）的32位或64位浮点数与 ST0 中的值相乘，并将结果存储在 ST0 。双操作数版本 fmul dest, src 执行 dest = dest * src ，有弹出版本 fmulp dest, src ，还有一种特殊的无操作数弹出版本 fmulp ，将 ST0 和 ST1 相乘。使用 fimul src 指令可将 ST0 的内容与内存中存储的整数相乘， src 处的值可以是32位或64位整数。
• 除法指令 ：内存版本的除法指令 fdiv src 将 ST0 的内容除以 src 并将结果存储在 ST0 ， src 操作数可以是内存中的单精度或双精度浮点值。双操作数版本 fdiv dest, src 执行 dest = dest / src ，弹出版本使用 fdivp 。若要将 ST0 除以整数，可使用 fidiv 指令。与减法指令类似，除法指令也有反向版本，如 fdivr src 执行 ST0 = src / ST0 。

3.3 比较指令

比较指令 fcom src 用于比较 ST0 中的值与 src ，并设置FPU标志。 src 操作数可以在内存或寄存器中， C1 位用于指示栈溢出/下溢情况，其他三位（ C0 、 C2 和 C3 ）用于指示关系如下：
| 关系 | C3 C2 C0 |
| ---- | ---- |
| ST0 > src | 0 0 0 |
| ST0 = src | 1 0 0 |
| ST0 < src | 0 0 1 |
| 不可比较 | 1 1 1 |
若指令中未给出操作数，则比较栈顶两个值（即 ST0 与 ST1 ），有弹出版本 fcomp 。还有双弹出版本 fcompp ，比较 ST0 与 ST1 ，更新FPU标志并弹出两个值。使用 ficom src 可将栈顶与内存中的整数值进行比较， src 可以是16位或32位整数，也有弹出版本 ficomp 。比较栈顶与零可使用 ftst 指令，该指令无操作数，比较栈顶值与0.0并更新FPU标志。 fxam 指令用于检查 ST0 中的数值类型，返回其符号到 C1 标志位（0为正，1为负），并在其余三位（ C0 、 C2 和 C3 ）返回以下信息：
| 类型 | C3 | C2 | C0 |
| ---- | ---- | ---- | ---- |
| 不支持 | 0 | 0 | 0 |
| NaN | 0 | 0 | 1 |
| 正常 | 0 | 1 | 0 |
| 无穷大 | 0 | 1 | 1 |
| 零 | 1 | 0 | 0 |
| 空 | 1 | 0 | 1 |
| 非正规 | 1 | 1 | 0 |

3.4 其他指令

• fchs ：改变 ST0 中数值的符号。
• fabs ：将 ST0 中的值替换为其绝对值。
• fldcw src ：将内存中 src 处的16位值加载到FPU控制字寄存器。
• fstcw dest ：存储控制字，执行后四个标志位（ C0 - C3 ）未定义。
• fstsw dest ：存储状态字， dest 可以是16位内存位置或 AX 寄存器，结合 sahf 指令可使用条件跳转指令，执行后四个标志位（ C0 - C3 ）未定义。

4. 第一个浮点程序：数组求和

下面通过一个示例展示如何使用浮点指令计算双精度数组的和。

C程序（Program 22.1）

#include <stdio.h>
#define SIZE 10

int main(void)
{
    double value[SIZE];
    int i;
    extern double array_fsum(double*, int);
    printf("Input %d array values:\n", SIZE);
    for (i = 0; i < SIZE; i++)
        scanf("%lf", &value[i]);
    printf("Array sum = %lf\n", array_fsum(value, SIZE));
    return 0;
}

汇编语言程序（Program 22.2）

; This procedure receives an array pointer and its size
; via the stack. It computes the array sum and returns
; it via ST0.
segment .text
global array_fsum
array_fsum:
    enter 0,0
    mov EDX, [EBP+8]    ; copy array pointer
    mov ECX, [EBP+12]   ; copy array size
    fldz                ; ST0 = 0 (sum is in ST0)
add_loop:
    jecxz done          ; exit loop if index is zero
    dec ECX             ; update the array index
    fadd qword[EDX+ECX*8] ; ST0 = ST0 + array_element
    jmp add_loop
done:
    leave
    ret

这个程序中，C程序负责用户界面，请求用户输入数组值，然后调用汇编语言过程 array_fsum 计算数组的和。汇编语言过程将数组指针复制到 EDX ，数组大小复制到 ECX ，使用 fldz 指令将 ST0 初始化为零，通过循环使用 fadd 指令累加数组元素，最终结果存储在 ST0 中并返回。

通过以上内容，我们对浮点运算有了更深入的了解，包括浮点单元的组织、各种浮点指令的使用以及如何编写简单的浮点程序。后续我们将通过更多示例进一步展示浮点指令的应用。

5. 浮点指令应用示例

5.1 示例一：求解二次方程的根

在这个示例中，我们要解决二次方程 $ax^2 + bx + c = 0$ 的根的问题。二次方程的两个根定义如下：
$root1 = \frac{-b + \sqrt{b^2 - 4ac}}{2a}$
$root2 = \frac{-b - \sqrt{b^2 - 4ac}}{2a}$

当 $b^2 > 4ac$ 时，根为实数；否则为虚数。

C程序（Program 22.3）

#include <stdio.h>

int main(void)
{
    double a, b, c, root1, root2;
    extern int quad_roots(double, double, double, double*, double*);
    printf("Enter quad constants a, b, c: ");
    scanf("%lf %lf %lf", &a, &b, &c);
    if (quad_roots(a, b, c, &root1, &root2))
        printf("Root1 = %lf and root2 = %lf\n", root1, root2);
    else
        printf("There are no real roots.\n");
    return 0;
}

汇编语言程序（Program 22.4）

; It receives three constants a, b, c and pointers to two
; roots via the stack. It computes the two real roots if
; they exist and returns them in root1 & root2. In this
; case, EAX = 1. If no real roots exist, EAX = 0.
%def a qword[EBP+8]
%def b qword[EBP+16]
%def c qword[EBP+24]
%def root1 dword[EBP+32]
%def root2 dword[EBP+36]

segment .text
global quad_roots
quad_roots:
    enter 0, 0
    fid a           ; a
    fadd ,ST0       ; a,a
    fid a           ; a,a,a
    fid c           ; c,a,a,a
    fmulp           ; a*c,a,a
    fadd            ; 2a*c,a
    fadd            ; 4a*c
    fchs            ; -4a*c
    fid b           ; b,-4a*c
    fid b           ; b,b,-4a*c
    fmulp           ; b*b,-4a*c
    faddp           ; b*b - 4a*c
    ftst            ; compare (b*b - 4a*c) with 0
    fstsw AX        ; store status word in AX
    sahf            ; move AX to flags register
    jb no_real_roots ; jump if (b*b - 4a*c) < 0

    fsqrt           ; sqrt(b*b - 4a*c)
    fid b           ; b,sqrt(b*b - 4a*c)
    fchs            ; -b,sqrt(b*b - 4a*c)
    fadd            ; -b + sqrt(b*b - 4a*c),sqrt(b*b - 4a*c)
    fid a           ; a,-b + sqrt(b*b - 4a*c),sqrt(b*b - 4a*c)
    fadd            ; 2a,-b + sqrt(b*b - 4a*c),sqrt(b*b - 4a*c)
    fdiv            ; (-b + sqrt(b*b - 4a*c))/2a,sqrt(b*b - 4a*c)
    mov EAX, root1  ;
    fstp qword[EAX] ; store root1

    fchs            ; -sqrt(b*b - 4a*c)
    fid b           ; b,-sqrt(b*b - 4a*c)
    fsubp           ; -b - sqrt(b*b - 4a*c)
    fid a           ; a,-b - sqrt(b*b - 4a*c)
    fadd            ; 2a,-b - sqrt(b*b - 4a*c)
    fdivrp          ; (-b - sqrt(b*b - 4a*c))/2a
    mov EAX, root2  ;
    fstp qword[EAX] ; store root2
    mov EAX, 1      ; real roots exist
    jmp done

no_real_roots:
    sub EAX, EAX    ; EAX = 0 (no real roots)

done:
    leave
    ret

在这个程序中，C程序负责用户界面，请求用户输入常数 $a$、$b$ 和 $c$，然后将这三个值以及两个指向根的指针传递给汇编语言过程 quad_roots 。汇编语言过程接收五个参数，通过一系列浮点指令计算 $b^2 - 4ac$ 的值，使用 ftst 指令检查其是否为负数，若为负数则跳转到 no_real_roots 标记处，返回 0 表示没有实数根；否则计算两个根并存储在相应的指针位置，返回 1 表示有实数根。

5.2 示例二：数组求和的内联汇编版本

在这个示例中，我们使用内联汇编方法重写 array_fsum 过程。使用内联汇编时，需要使用 AT&T 语法，对于浮点指令，使用以下后缀明确操作数大小：
- s ：单精度
- l ：双精度
- t ：扩展精度

内联汇编代码（Program 22.5）

double array_fsum(double* value, int size)
{
    double sum;
    asm("fldz; "
        "add_loop: jecxz done; "
        "decl %%ecx; "
        "faddl (%%ebx,%%ecx,8); "
        "jmp add_loop; "
        "done: "
        : "=t"(sum)
        : "b"(value), "c"(size)
        : "cc");
    return(sum);
}

在这个代码中，我们使用 asm 关键字嵌入汇编代码。通过 =t 输出说明符将变量 sum 映射到浮点寄存器，将 value 映射到 EBX 寄存器， size 映射到 ECX 寄存器。代码逻辑与之前的汇编语言过程类似，先将 ST0 初始化为零，然后通过循环使用 faddl 指令累加数组元素，最终结果存储在 sum 变量中并返回。

6. 总结

我们对浮点单元的组织进行了简要介绍，重点关注了 FPU 提供的寄存器。FPU 有八个浮点数据寄存器，它们被组织成栈结构。浮点指令包括多种算术和非算术指令，我们详细讨论了其中一些指令，如数据移动、算术运算、比较等指令。最后，通过几个示例展示了这些浮点指令的实际应用，包括计算数组的和、求解二次方程的根等。

通过深入学习浮点运算的原理、指令和应用，我们可以更好地处理涉及浮点数的计算任务，提高程序的数值计算能力和准确性。在实际编程中，根据具体需求合理选择和使用浮点指令，能够优化程序的性能和精度。

下面用 mermaid 流程图展示求解二次方程根的流程：

graph TD;
    A[开始] --> B[输入 a, b, c];
    B --> C[计算 b^2 - 4ac];
    C --> D{b^2 - 4ac < 0?};
    D -- 是 --> E[无实数根，返回 0];
    D -- 否 --> F[计算 sqrt(b^2 - 4ac)];
    F --> G[计算 (-b + sqrt(b^2 - 4ac))/2a];
    G --> H[存储 root1];
    F --> I[计算 (-b - sqrt(b^2 - 4ac))/2a];
    I --> J[存储 root2];
    H --> K[返回 1];
    J --> K;
    E --> L[结束];
    K --> L;

通过这个流程图，我们可以更清晰地看到求解二次方程根的程序逻辑，根据判别式的值决定是否有实数根，并进行相应的计算和结果返回。