# 2024 Homework4 (quiz3+4) contributed by < `SimonLee0316` > ## 第 3 週測驗題 ### 測驗一 版本一使用 `log2` 來找出最高位元的1,由那個數開始往下開始計算平方根 版本二使用 `while`迴圈來找到最高位元的1,由那個數開始往下開始計算平方根 * 而這兩個版本找出最左邊的1都是個一個慢慢找,而且需要每次迴圈都做乘法運算 版本三使用 `Digit-by-digit calculation` 直式開方法來實做 透過將要求的平方根拆成 $(a_{n}+a_{n-1}+...+a_0)^2,a_m=2^m$ or $a_m=0$ 整理過後可得 \begin{align*} N^2&=a_{n}^2+[2(a_{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 \\ &=a_{n}^2+[2(P_n+a_{n-1}]a_{n-1}+[2P_{n+1})+a_{n-2}]a_{n-2}+...+[2P_{1}+a_0]a_0 \end{align*} $P_0=a_n+a_{n+1}+a_{n+2}+...++a_0 即是我們要求的平方根 N$ 觀察可以得 \begin{align*} P_m=P_{m+1}+a_m \end{align*} 從高位遞減都去計算 $N^2-P_m^2$ 即可以求得解,但是這樣的運算成本太高所以使用上一輪差值來計算 $P_m$ 上一輪差值 $X_{m+1}$ 減去 $Y_{m}$ 得到$X_m$ * $X_{m}=N^2-P_{m+1}^2=X_{m+1}-Y_m$ * $Y_m=P_m^2-p_{m+1}^2=2P_{m+1}a_m+a_m^2$ $Y_m$ 拆成 $c_m,d_m$ * $c_{m-1} = \begin{cases} c_m/2 + d_m & if \ a_m = 2^m \\ c_m/2 & if \ a_m = 0 \end{cases}$ * $d_{m-1}= d_m/4$ 演算法求 $P_0$ 從 最高位 $a_n$ 開始 * 差值一開始為 N ,$X_{n+1}=N$ * $c_n=0$ * $d_n=(2^{m})^2$ ```c int i_sqrt(int n) { if (n <= 1) /* Assume x is always positive */ return x; int c = 0; for (int d = 1UL << ((31 - __builtin_clz(x)) & ~1UL); d; d >>= 2) { int y = c + d; c >>= 1; if (n >= y) n -= y, c += d; } return c; } ``` ### 測驗二 針對正整數在相鄰敘述進行 mod 10 和 div 10 操作,如何減少運算成本? 在原始範例中 餘數和商數 都要做一次除法運算,這運算成本非常高 參考 `Hacker's Delight`,採用 bitwise 來實做除法,主要概念為把 10看成乘 $\frac{1}{10}$,而因為10無法用 $2^n$表示,所以只能找出非常接近 ${10}$的數$\frac{a}{2^n}$這裡取 $a=13,2^n=128,\frac{128}{13}\approx 9.84$, 所以 乘$\frac{1}{10}$ 現在變成 $\frac{13}{128}$ 在開始乘之前我們先將乘 13 弄成 $2^n$ 的形式才好做位移 13可以拆乘 8+4+1 可以表示成 $\frac{13}{8}=\frac{1}{8}+\frac{1}{2}+1$,因為有除與8的情況,再位移時會將數值位移掉所以要先將被位移掉的數字先存起來做完之後在加回去 * 保留後三位元 ```c d0 = q & 0b1; d1 = q & 0b11; d2 = q & 0b111; ``` * $tmp \times \frac{13}{8}\times8=tmp\times13$ ```c ((((tmp >> 3) + (tmp >> 1) + tmp) << 3) ``` * 把後面被移掉的補回去 ```c ((((tmp >> 3) + (tmp >> 1) + tmp) << 3) + d0 + d1 + d2) ``` * $tmp\div128$ 即得到商 ```c q = ((((tmp >> 3) + (tmp >> 1) + tmp) << 3) + d0 + d1 + d2) >> 7; ``` * 餘數為除數-商($q$)*除數(10) 因為除數是 10不是20的冪,所以要將$(q\times4+q)\times2=10\times{q}$ ```c r = tmp - (((q << 2) + q) << 1); ``` 完整程式碼 ```c #include <stdint.h> #include <stdio.h> int main() { for (int n = 0; n <= 19; n++) { unsigned d2, d1, d0, q, r; d0 = q & 0b1; d1 = q & 0b11; d2 = q & 0b111; q = ((((n >> 3) + (n >> 1) + n) << 3) + d0 + d1 + d2) >> 7; r = n - (((q << 2) + q) << 1); printf("q: %d r: %d\n", q, r); } return 0; } ``` 包裝成函式 ```c #include <stdint.h> void divmod_10(uint32_t in, uint32_t *div, uint32_t *mod) { uint32_t x = (in | 1) - (in >> 2); /* div = in/10 ==> div = 0.75*in/8 */ uint32_t q = (x >> 4) + x; x = q; q = (q >> 8) + x; q = (q >> 8) + x; q = (q >> 8) + x; q = (q >> 8) + x; *div = (q >> CCCC); *mod = in - ((q & ~0x7) + (*div << DDDD)); } ``` 探討程式如何包裝 首先這裡可以看成 $x-\frac{1}{4}\times{x}=\frac{3}{4}\times{x}$,(in | 1)這裡是強至讓 in 變成奇數不知道為何要這樣做 ```c uint32_t x = (in | 1) - (in >> 2); /* div = in/10 ==> div = 0.75*in/8 */ ``` 再來 q 是 $\frac{3}{64}\times{x}+\frac{3}{4}\times{x}=\frac{102}{128}\times{x}\approx0.8$ ```c uint32_t q = (x >> 4) + x; ``` 為了更接近 0.8,而做了四次逼近 ```c q = (q >> 8) + x; q = (q >> 8) + x; q = (q >> 8) + x; q = (q >> 8) + x; ``` 商數最後要除8為我們的除數是$\frac{8}{10}$,要讓它變成$\frac{1}{10}$ ```c *div = (q >> 3); ``` 餘數要保留被移掉的最後三 bit,被除數-(商數*除數)*10 ```c *mod = in - ((q & ~0x7) + (*div << 1)); ``` ### 測驗三 版本一就是一直將數字右移直到數字為0,這麼做就是為了找到最高位元是1的位元,其答案就是`log2`的答案 版本2減少一開始從最最高位開始右移,用二分搜尋法的方式找到最高位元是1的位元,其答案就是`log2`的答案 ```c static size_t ilog2(size_t i) { size_t result = 0; while (i >= 65536) { result += 16; i >>= 16; } while (i >= 256) { result += 8; i >>= 8; } while (i >= 16) { result += 4; i >>= 4; } while (i >= 2) { result += 1; i >>= 1; } return result; } ``` 使用GNU extension `__builtin_clz` ```c int ilog32(uint32_t v) { return (31 - __builtin_clz(v)); } ``` 在 [linux/log2.h](https://github.com/torvalds/linux/blob/master/include/linux/log2.h) 的 `ilog2()`為例,它可以被用來初始化全域變數 ### 測驗五 測驗五是將測驗三的運算變成 bitwise 的運算,加法可以用 or 運算代替 一開始的 `x--` 如果做最後結果會多1 , 在最後要在判斷 `x>1`防止 `x==0x3`或`x==0x2`的情況 ```c int ceil_ilog2(uint32_t x) { uint32_t r, shift; x--; r = (x > 0xFFFF) << 4; x >>= r; shift = (x > 0xFF) << 3; x >>= shift; r |= shift; shift = (x > 0xF) << 2; x >>= shift; r |= shift; shift = (x > 0x3) << 1; x >>= shift; return (r | shift | x > 1) + 1; } ``` 改進程式碼,使其得以處理 x = 0 的狀況,並仍是 branchless x = 0的情況只需要確保 `x--` 這裡不能是`0xFFFFFFFF`,所以我們要做的是確保 x=x-1是我們預期的(7-1=6,0-1=0) 所以這邊x=x-1當遇到 x 是0的情況時 * $x = \begin{cases} x + 1 & if \ x >0 \\ x & if \ x = 0 \end{cases}$ ```c x-=!!x ``` ## 第4周測驗題 ### 測驗一 population count ,世紀算術值得二進位表示中共有多少個位元是`1` ```c unsigned popcount_naive(unsigned v) { unsigned n = 0; while (v) v &= (v - 1), n = -(~n); return n; } ``` 將 `v&=(v-1)` 可以將最右邊是 1 的 bit 清成 0,藉由重複做這個步驟可以得到數字的2進制表示有多少個 1 而 `n = -(~n)` 這裡的設計更是巧妙,會等同於做了`n++`,以下為證明 $-n=\lnot{n}+1$ 模仿了2補數的步驟 $-(\lnot{n})=\lnot(\lnot{n})+1=n+1$ 上述的範例執行時間取決於數值中 1 的個數,改寫成常數時間的實做 對於 32 bit 的無號整數,popcount 可以寫成以下數學式 $popcount(x)=x- {\lfloor\frac{x}{2}\rfloor}- {\lfloor\frac{x}{4}\rfloor}-...-{\lfloor\frac{x}{2^{31}}\rfloor}$ 以下程式碼就在實做 $x- {\lfloor\frac{x}{2}\rfloor}- {\lfloor\frac{x}{4}\rfloor}-{\lfloor\frac{x}{8}\rfloor}$ ```c n = (v >> 1) & 0x77777777; v -= n; n = (n >> 1) & 0x77777777; v -= n; n = (n >> 1) & 0x77777777; v -= n; ``` 最後的 v 會是每四個 bit 的 set bit 數量 最後在將每四個 bit的 set bit 數量加起來即可 ```c v = (v + (v >> 4)) & 0x0F0F0F0F v *= 0x01010101 v=>>24 ``` [Total Hamming Distance](https://leetcode.com/problems/total-hamming-distance/) ```c int totalHammingDistance(int* nums, int numsSize) { int total = 0;; for (int i = 0;i < numsSize;i++) for (int j = 0; j < numsSize;j++) total += __builtin_popcount(nums[i] ^ nums[j]); return total >> 1; } ``` 每兩個數值算了兩次所以最後要除2 ### 測驗二 [Incomplete tic-tac-toe](https://gist.github.com/jserv/9088f6f72a35a1fcaaf1c932399a5624) 對井字遊戲做了 100萬次模擬每次選擇的路線都是隨機的,在 `play_random_game` 中,先使用 `xorshift32` 隨機選擇要走那一步,而使用 `move_masks` 來計算分數,剛開始看的時候會看不懂為何要這樣做,在判斷是否勝利的時候,勝利條件總共可以分為 8 種,分別為橫的( 3 條),直的( 3 條),斜的( 2 條),而每一條總共有三個地方可以放,總結以上所說用 32 bit 來表示的話 | bit | 31-20 | 19-8 | 7-0 | |:----:|:--------:|:--------:|:--------:| | 代表 | 三條橫的 | 三條直的 | 兩條斜的 | 範例: | 0 | 1 | 2 | 3 | |:---:|:----:|:----:|:----:| | **1** | o | | | | **2** | | | | | **3** | | | | | 0 | 1 | 2 | 3 | |:---:|:----:|:----:|:----:| | **1** | 1 | 0 | 0 | | **2** | 0 | 0 | 0 | | **3** | 0 |0 | 0 | 橫的(`0x400`):&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;直的(`0x400`): * 第一條1(1,1) 0 (1,2) 0 (1,3) &emsp;第一條1(1,1) 0 (2,1) 0 (3,1) * 第二條0(2,1) 0 (2,2) 0 (2,3) &emsp;第二條0(1,2) 0 (2,2) 0 (3,2) * 第三條0(3,1) 0 (3,2) 0 (2,3) &emsp;第三條0(1,3) 0 (2,3) 0 (3,3) 斜的(`0x40`): * 第一條1(1,1) 0 (2,2) 0 (3,3) * 第一條0(3,1) 0 (2,2) 0 (3,1) 所以我們可以得出下了(1,1)我們的值是: `0x40040040` ,每次都將值加入到自己的分數內 我們可以很直覺的看出勝利條件就是某一條線的值是7,對應到下面程式 ```c return (player_board + 0x11111111) & 0x88888888; ``` mod3 討論: 如果除數符合 $2^k\pm1$ 的形式就可以用以下提到的方法 以除數 3 為例, $1\equiv1(mod~3) 且 2\equiv-1(mod~3)$ 可以歸納成 $2^k = \begin{cases} 1(mod~3) & if \ k ~even \\ -1(mod~3) & if \ k ~odd \end{cases}$ 照以上方式來推論如果我們要知道某個數$mod~3$ 是多少的話將偶數位元是 1 的數量乘 1 + 上奇數位元數量乘 -1 我們就可以得到答案。 1. 偶數位元是 1 : `popcount(n & 0x55555555)` ,`0x5=0101`和它做 `and` 只有偶數位元是 1 才會是 1 2. 奇數位元是 1 : `popcount(n & 0xAAAAAAAA)` ,`0xA=1010`和它做 `and` 只有奇數位元是 1 才會是 1 最後將 (1) * 1 + (2) * -1 可以看成`n = popcount(n & 0x55555555) - popcount(n & 0xAAAAAAAA)` 又可以進一步化簡成為 `n = popcount(n ^ 0xAAAAAAAA) - 16` ,這裡的 n,範圍是 $-16\leq n \leq16$,而了避免負數,將數值 +39,使現在的範圍變成 $23\leq n \leq55$ 我們的目標是要作到使 $0\leq n \leq2$,所以我們這邊再做一次模同餘 3,但這次我們不用在 `popcount(n ^ 0xAAAAAAAA) - 16`,因為 n 的範圍已經備縮小到 $23\leq n \leq55$ 這意味著只需要 xor 6 個bit即可(因為 55 只需要 6 bit 就能表示),經過這次模同餘 3 之後,範圍現在縮小到 $-1\leq n \leq1$,這裡其實就已經是答案了,但是因為 $2\equiv-1(mod~3)$,我們要將 -1 變成 2 ,使用`n + ((n >> 31) & 3)`來達到,這裡如果是負數的話使用算術右移,左邊會補上 sign ,所以當你是 `-1` 做 `&3` 會得到 3,反之是 0,這就可以達到 $0\leq n \leq2$。 以下為程式碼 ```c int mod3(unsigned n) { n = popcount(n ^ 0xAAAAAAAA) + 23; n = popcount(n ^ 0x2A) - 3; return n + ((n >> 31) & 3); } ``` 再來討論另外一種變形使用 look up table 使用上面得到的結論 `popcount(n ^ 0xAAAAAAAA) - 16`,會發現它少了 -16 這個步驟,如果把 `popcount(n ^ 0xAAAAAAAA) - 16` ,他的範圍是$-16\leq n \leq16$ ,那 `popcount(n ^ 0xAAAAAAAA)`的範圍就是$0\leq n \leq32$,所以 0 實際上答案是 -16,用下面表來說明更清楚。 以下用 K 來代表 `popcount(n ^ 0xAAAAAAAA)` | $K$ | $K-16$ | $mod~3$ | | --- | ------ | ------- | | 0 | -16 | 2 | | 1 | -15 | 0 | | 2 | -14 | 1 | | 3 | -13 | 2 | | ... | ... | ... | | 32 | 16 | 1 | ### 測驗三 xtree 的主程式在 `treeint.c` 它裡面使用了 fuction hook 方式將 `treeint_xt.c` 內提供的函式定用在 `treeint.c` 內。 `treeint_xt.c` 提供了五個 xtree 的操作。 XTree 的平衡機制與 AVL tree 相似,利用 hint 來評估是否需要做 rotate,然而 xt_node 的 hint 的計算方式是,其左右節點 hint 最大值 +1 並且只有在被呼叫到 xt_update 時才會被更新。 以下討論 xtree.c 內提供了 XTree 的各種操作 * `__xt_destroy` 刪掉 XTree * `xt_first` 、 `xt_last` 找出最左最右的葉節點 * `xt_rotate_left`、`xt_rotate_right` 對節點做左旋轉和右旋轉 * `xt_balance` 計算該節點左右子點的 `hint` 差值 * `xt_update` 如果 balance 值小於 -1 ,做右旋轉 如果大於 1 ,做左旋轉 如果hint 原本是 0 ,但後來被改了,在遞迴做 `xt_update` * `__xt_insert` 插入 node 到XTree ,且同時對 XTree 做 `xt_update` * `__xt_remove` 如果要刪的節點存在右子點,則找到該節點的子樹的最左葉節點來替代要被刪的節點,之後再對該點柚子點做update 如果要刪的節點存在左子點,則找到該節點的子樹的最右葉節點來替代要被刪的節點,之後再對該點左子點做update