平方根函數演進史：從二分法到Quake III神奇常數0x5f3759df

引言：一個函數引發的性能探索

sqrt() 是科學計算、圖形學、物理引擎中調用頻率最高的函數之一。本文追溯其算法演進：二分法 → 牛頓迭代法 → 快速逆平方根，揭示那個傳奇常數0x5f3759df背後的數學美學與工程智慧。

第一階段：二分法 — 直觀但緩慢

float SqrtByBisection(float n) {
    if(n < 0) return n;
    float low = 0, up = n, mid, last;
    do {
        mid = (low + up) / 2;
        if(mid * mid > n) up = mid;
        else low = mid;
    } while(fabs(mid - last) > EPS);
    return mid;
}

性能：收斂線性，需~30次迭代達到float精度，比庫函數慢數百倍。

第二階段：牛頓迭代法 — 二次收斂

求解 f(x) = x² - a = 0，切線迭代公式：

x_{n+1} = (x_n + a/x_n) / 2

float SqrtByNewton(float a) {
    float x = a;
    do {
        float last = x;
        x = (x + a/x) / 2;
    } while(fabs(x - last) > EPS);
    return x;
}

性能：誤差指數級減少，~5次迭代達到精度，但仍比庫函數慢數十倍。

第三階段：快速逆平方根 — 位級黑客藝術

核心代碼 (Quake III Arena)

float Q_rsqrt(float number) {
    long i;
    float x2 = number * 0.5F, y = number;
    const float threehalfs = 1.5F;
    
    i = *(long*)&y;                       // 1. 浮點位模式 → 整數
    i = 0x5f3759df - (i >> 1);            // 2. 神奇常數 + 移位 = 初始猜測
    y = *(float*)&i;                      // 3. 整數 → 浮點位模式
    y = y * (threehalfs - x2 * y * y);    // 4. 一次牛頓迭代修正
    return y;                             // 返回 1/√number
}

IEEE 754位模式洞察

Float32：值 = (-1)^符號 × 2^(指數-127) × (1.尾數)

對數域近似：log₂(1/√x) ≈ -½ log₂(x)

整數解讀：I_x ≈ 2²³(log₂x + 127 - 0.5)

目標：I_y ≈ 2²³(127 - ½(log₂x + 127 - 0.5)) = 0x5f3759df - I_x/2

常數演進史

性能對比 (現代x86-64, -O3)

std::sqrt      : 1.00x (基準, 硬件指令SQRTSD)
Newton (3 iter): 2.3x slower
Q_rsqrt (1 iter): 1.8x slower (但輸出1/√x)
rsqrtss指令    : 0.7x faster (近似, 需修正)

關鍵結論：現代CPU硬件指令SQRTSD/RSQRTSS已全面超越軟件算法。快速逆平方根的歷史意義在於無硬件平方根指令的早期3D時代（RISC、嵌入式、GPU早期），以極低成本換取可用精度。

現代最佳實踐

• 通用高精度：直接用std::sqrt / sqrtf，編譯器會發射硬件指令
• 圖形學/物理引擎：常需1/√x（歸一化），用rsqrtss + 1次牛頓修正
• 嵌入式無FPU：定點牛頓法或查表+插值
• 向量化：AVX2/AVX-512 _mm256_rsqrt_ps 批量處理

總結

從二分法到0x5f3759df，見證了「數學洞察 + 硬件特性利用 + 工程權衡」的經典範例。雖然硬件進步使其退出性能舞台，但位級操作思維、IEEE 754深度理解、近似計算藝術，永遠是系統程序員的核心內功。

1. 上一篇文章： Bash實現C代碼行數統計：從粗糙腳本到工業級工具的演進之路
2. 下一篇文章： C語言通用HashMap實現：從設計原理到生產級代碼

關於C語言

C# 3.0 匿名方法實現與應用 (1)：從 C# 2.0 匿名方法到隱式型別

C# 2.0 匿名方法回顧 回呼方法名直接傳遞，C

位元組序大小端完整解析：原理、檢測程式、網路傳輸與跨平台相容性

位元組序大小端完整解析：原理、檢測程式、網路傳輸與跨平台相容

C 語言內存分配全景圖：從虛擬地址空間到 malloc 實現原理

C 語言內存分配全景圖：從虛擬地址空間到 malloc 實現

C# 設計模式無廢話系列十九：觀察者模式（三）—— 代碼執行解析與擴充性展示

代碼執行結果 代碼深度解析 本例採用 C# 事件機