2024q1 Homework4 (quiz3+4)

contributed by < youjiaw >

第 3 週測驗題

測驗 1

Digit-by-digit calculation 提及了一種計算平方根的方式，其原理是將欲求平方根的

$N^2$ 拆成：

$$N^2=(a_n+a_{n-1}+a_{n-2}+...+a_0{)}^2,a_m=2^{m}or{a}_m=0$$
並根據

$(x+y{）}^2=x^2+2xy+y^2$ 將其擴展為：

$$N^2=a_n^2+[2a_n+a_{n-1}]a_{n-1}+[2(a_n+a_{n-1})+a_{n-2}]a_{n-2}+...+[2(\sum _{i=1}^n{a}_i)+a_0]a_0$$

以 n=2 舉例，

$N^2=(a_2+a_1+a_0{)}^2$ 可以用下方的矩陣表示

$$\left[\begin{array}{ccc}{a}_0{a}_0& a_1{a}_0& a_2{a}_0\\ a_0{a}_1& a_1{a}_1& a_2{a}_1\\ a_0{a}_2& a_1{a}_2& a_2{a}_2\end{array}\right]$$
將算式調整後可得

$$N^2=a_2^2+[2a_2+a_1]a_1+[2(a_2+a_1)+a_0]a_0$$
分別代表以下三個部分

$$\left[\begin{array}{ccc}& & \\ & & ⋮\\ & ...& a_2{a}_2\end{array}\right]$$

$$\left[\begin{array}{ccc}\ddots & ⋮& ⋮\\ ...& a_1{a}_1& a_2{a}_1\\ ...& a_1{a}_2& \ddots \end{array}\right]$$

$$\left[\begin{array}{ccc}{a}_0{a}_0& a_1{a}_0& a_2{a}_0\\ a_0{a}_1& \ddots & ⋮\\ a_0{a}_2& ⋮& \end{array}\right]$$

接著令

$P_m=\sum _{i=m}^n{a}_i$，則

$P_0=a_n+a_{n-1}+...+a_0$ 就會是我們要求的平方根 N，因此我們需要找出每個

$a_m$ 的值。

想檢測

$a_m$ 的值是

$2^m$ 還是 0，可以計算

$P_m^2$ 與

$N^2$ 的差值，若

$P_m^2\le N^2$ 就代表

$a_m$ 為

$2^m$，否則為0，但是此方法每輪的運算成本太高了。

若是令

• $X_m=N^2-P_m^2=X_{m+1}-Y_m$
• $Y_m=P_m^2-P_{m+1}^2=(P_{m+1}+a_m{)}^2-P_{m+1}^2=2P_{m+1}{a}_m+a_m^2$

就可以透過紀錄

$P_{m+1}$ 來取代

$P_m^2$ 的計算。

最後將

$Y_m$ 拆為

$c_m$、

$d_m$，可得

• $c_m=2P_{m+1}{a}_m$
• $d_m=a_m^2$
• $Y_m=\{\begin{array}{ll}{c}_m+d_m& \text{if }{a}_m=2^m\\ 0& \text{if }{a}_m=0\end{array}$

並藉由位元運算，由

$c_m$、

$d_m$ 推出下一輪的

$c_{m-1}$、

$d_{m-1}$。

• $c_{m-1}=P_m{2}^m=(P_{m+1}+a_m)2^m=P_{m+1}{2}^m+a_m{2}^m=\{\begin{array}{ll}{c}_m/2+d_m& \text{if }{a}_m=2^m\\ c_m/2& \text{if }{a}_m=0\end{array}$
• $d_{m-1}={\displaystyle \frac{d_m}{4}}$

最後求得的

$c_{-1}$ 即為所求的

$P_0=N$。

程式碼如下：

int i_sqrt(int x)
{
    if (x <= 1) /* Assume x is always positive */
        return x;

    int z = 0;
    for (int m = 1UL << ((31 - __builtin_clz(x)) & ~1UL); m; m >>= 2) {
        int b = z + m;
        z >>= 1;
        if (x >= b)
            x -= b, z += m;               
    }
    return z;
}

在確認 x 的範圍後，會使用三個變數，z、m 與 b，分別對應到先前的

$c_m$、

$d_m$ 與

$Y_m$。
在每輪迴圈中，無論

$a_m$ 值為何，

$c_m$ 都會除以 2，

$d_m$ 則是除以 4，因此 AAAA 為 2，BBBB 為 1。

但是我還沒理解為什麼要做 & ~1UL，相當於限制左移位元數為偶數，目前猜測是因為

$d_m=a_m^2=2^{2m}$，所以 m 僅可能為 2 的偶數次方。

嘗試用 ffs / fls 取代 __builtin_clz

ffs 會回傳最低有效位元的索引，從 1 開始計算，若 x 為 0 則會回傳 0。因此可以不斷將 x 右移 ffs(x) 個位元，直到 x 為 0 時就停止，最後累計的位移次數減一就會是最高有效位元的索引。

因此加入以下程式碼

int tmp = x, msb = 0;
while (tmp) {
    int index = ffs(tmp);
    tmp >>= index;
    msb += index;
}

並更改迴圈條件為

-    for (int m = 1UL << ((31 - __builtin_clz(x)) & ~1UL); m; m >>= 2) {
+    for (int m = 1UL << ((msb - 1) & ~1UL); m; m >>= 2) {

而 fls 是直接回傳最高有效位元的索引

，因此更改如下

-    for (int m = 1UL << ((31 - __builtin_clz(x)) & ~1UL); m; m >>= 2) {
+    for (int m = 1UL << ((fls(x) - 1) & ~1UL); m; m >>= 2) {

測驗 3

此程式碼目的是取得輸入值 i 在二進位表示中最高有效位元的索引。

版本一是將 i （

$2^n$）不斷右移一個位元直到為 0，此時右移次數即為 n，迴圈執行次數也為 n。
版本二則是依照 i 的大小，提供四種右移的方式，只要 i 大於

$2^k$，就一次右移 k 個位元，讓迴圈的執行次數可以小於 n。

static size_t ilog2(size_t i)
{
    size_t result = 0;
    while (i >= 2^16) {
        result += 16;
        i >>= 16;
    }
    while (i >= 2^8) {
        result += 8;
        i >>= 8;
    }
    while (i >= 2^4) {
        result += 4;
        i >>= 4;
    }
    while (i >= 2) {
        result += 1;
        i >>= 1;
    }
    return result;
}

版本三使用了 __builtin_clz 改寫，但因為輸入值不能為 0，所以先與 1 做 | 運算。

return (31 - __builtin_clz(v | 1));

測驗 5

此程式碼以測驗三做延伸，計算

$\lceil lo{g}_2(x)\rceil$。

若 (x > 0xFFFF) 回傳 1，則 r 的值為 10000，因此會將 x 右移 16 個位元。

r = (x > 0xFFFF) << 4;                                   
x >>= r;

對應到測驗三的

while (i >= 2^16) {
    result += 16;
    i >>= 16;
}

而因為每次左移的位元數都不同，代表 r 和 shift 不會有相同位元都為 1 的情況，因此計算 result 的方式可以用 r |= shift 代替。

shift = (x > 0xFF) << 3;
x >>= shift;
r |= shift;
shift = (x > 0xF) << 2;
x >>= shift;
r |= shift;

最後再判斷 x 是否大於 1。

shift = (x > 0x3) << 1;
x >>= shift;
return (r | shift | x > 1) + 1;

改進程式碼

最直覺的想法就是當 x 為 0 時減 0，其他情況減 1。目前採用下方的方法。

-    x--;
+    x -= (x > 0)

第 4 週測驗題

測驗 1