# 2024q1 Homework4 (quiz3+4)
contributed by < [ChenFuhuangKye](https://github.com/kyehuang) >
## 2024q1 第 3 週測驗題
### 測驗 `1`
#### 根據 [Digit-by-digit calculation](https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Digit-by-digit_calculation) 提出直式開方法的變形
假設欲求的平方根的 $N^2$ 拆成
\begin{aligned}
N^2 = (a_n + a_{n-1} + a_{n-2} + ... + a_0)^2 \\ a_m = 2^m\quad or \quad a_m = 0
\end{aligned}
將其乘開
\begin{aligned}
N^2 = a_n^2 + [2P_n+a_{n-1}]a_{n-1} + [2P_{n-1}+a_{n-2}]a_{n-2} + ... + [2P_1+a_0]a_0 \\
P_m = a_n+a_{n-1}+a_{n-2}+...+a_m
\end{aligned}
而 $P_0 = a_n+a_{n-1}+a_{n-2}+...+a_0$ 即為所求的平方根 $N^2$
可以知道
\begin{aligned}
P_m = P_{m+1}+a_m
\end{aligned}
我們需要判斷 $a_m$ 是 $2^m$ 或是 $0$。
$a_m$ 可以透過試驗 $P_m^2 <= N^2$ 是否成立,然而這個做法成本太高。
另一種想法就是從上一輪的差值 $X_{m+1}$ 減去 $Y_m$ 得到 $X_m$
* $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$
接著將 $Y_m$ 拆分成 $c_m$ 及 $d_m$
* $c_m=2P_{m+1}a_m = P_{m+1}2^{m+1}$
* $d_m=a_m^2=(2^m)^2$
* $Y_m = \begin{cases} c_m+d_m & if \ a_m = 2^m \\ 0 & if \ a_m = 0 \end{cases}$
可以發現 $c_0=P_02^1 = 2P_0$ , 因此求出 $c_0/2$ 即為所求的平方根 $N^2$。
而
* $c_{m-1}=P_m2^m=(P_{m+1}+a_m)2^m=P_{m+1}2^m+a_m2^m = \begin{cases} c_m/2+dm & if \ a_m =2^m\\ c_m/2 & if \ a_m=0 \end{cases}$
* $d_{m-1} = d_m/4$
可以知道每一輪 $c$ 需要除以2 且 $d$ 需要除以4。
因此程式碼可以改寫成
```c
int i_sqrt(int x)
{
if (x <= 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 (x >= y)
x -= y, c += d;
}
return c;
}
```
__builtin_clz(x) 是求出 x 的最高位之前的 0 的個數,跟 ffs (find first bit in word) 功能一樣,因此可以將程式碼改成。
```diff
- for (int d = 1UL << ((31 - __builtin_clz(x)) & ~1UL); d; d >>= 2) {
+ for (int d = 1UL << ((31 - __ffs(x)) & ~1UL); d; d >>= 2) {
```
### 測驗 `2`
參照 [Hacker's Delight](https://web.archive.org/web/20180517023231/http://www.hackersdelight.org/divcMore.pdf) 以及 [公式推導](https://hackmd.io/@sysprog/linux2024-quiz3#%E6%B8%AC%E9%A9%97-2) 可以知道只要找到一個除數 $x$ 且 $9.55<=x<=10$ 。假設 $n=ab$ ( a 是十位數、 b 是個位數), $n$ 除以 $x$ 的商即為 $n$ 除以 10 的商,如此一來可以透過 bitwise operation 取得商跟餘數。
此外也可以透過 $0.8n$ 除以 $8$ 得到 $0.1n$ , 而 $0.1n$ 既為 div 10 的結果,因此找到一個 $x$ , 使得 $x*n\cong0.8n$。
```c
uint32_t q = (x >> 4) + x;
```
$q=\frac{3}{4*16}n + \frac{3}{4}n = \frac{51}{64}n =0.796875n$ , $x=0.796875$
透過 4 輪的逼近。
```c
q = (q >> 8) + x;
```
| 第幾輪 | x |
|:------:|:------------------ |
| 1 | 0.79998779296875 |
| 2 | 0.7999999523162842 |
| 3 | 0.7999999998137355 |
| 4 | 0.7999999999992724 |
我們可以得到 $q\cong0.8n$ ,將 $q$ 除以 $8$ 可以得到 $div$。
$div = \frac{q}{8} \cong \frac{0.8n}{8} \cong 0.1n$
```c
(q & ~0x7)
```
清除 q 的後三位元,相當於 div 乘 8,再加上 div 乘 2 , 可以算出 $mod=n-(8*div+2*div)=n-10*div$。
因此程式可以改寫成
```diff
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));
+ *div = (q >> 3);
+ *mod = in - ((q & ~0x7) + (*div << 1));
}
```
```c
uint32_t x = (in | 1) - (in >> 2);
```
要避免出現 20 的倍數,因為 20 乘上 x 會出現接近 16 的值 15。 當 $q=15$ , 20 的 div 為 1 。 所以將每一個數字轉換成奇數,避免出現 20 的倍數。
### 測驗 `3`
版本一的作法透過不斷把目標整數 i 除以 2 直到其變為 0 , 同時使用一個計數器 log 來記錄右移的次數, 這個過程是計算 $log_2(i)$ 的整數部分。因此目標整數出現第一個 1 的位數,對應 $log_2(i)$ 的整數部分。
版本二的做法是確認目標整數 i 第一個 1 的位數。 result 是一個計數器計算目前右移的次數,假設有一個整數 x 且 x 是 2 的冪,如果 i 大於 x 時,將 result 加上 $log2(x)$ 並把 i 右移 $log2(x)$ 次。
```diff
static size_t ilog2(size_t i)
{
size_t result = 0;
- while (i >= AAAA) {
+ while (i >= 65536) {
result += 16;
i >>= 16;
}
- while (i >= BBBB) {
+ while (i >= 256) {
result += 8;
i >>= 8;
}
- while (i >= CCCC) {
+ while (i >= 16) {
result += 4;
i >>= 4;
}
while (i >= 2) {
result += 1;
i >>= 1;
}
return result;
}
```
版本三直接利用 `__builtin_clz` 直接取得目標整數 i 第一個 1 的位數,為了處理 i 為 0 的情況, 將 i 與 1 進行 `i | 1` ,確保 i 不為 0 。
### 測驗 `4`
程式碼的資料型態是 `unsigned long` ,然而平均值會是小數,這邊是利用定點數,自行定義小數點的位置,如此一來就可以利用 `bitwise` 來計算小數。
程式碼中透過 `ewma_read` 函數讀取 `EWMA` , 發現讀取時會先把 S 右移 `factor` 次 , 因此在 `ewma_add` 所有新加入的值都會先左移 `factor` 次 。
因此在 `ewma_add` 中 `EWMA` 的公式可以看成
$f*S_t =\begin{cases} f*Y_0 & t=0\\ f*\alpha Y_t+(1-\alpha)*S_{t-1} & t>0\end{cases}$
>$f$ 為 $2^{factor}$ 、 $\alpha$ 為 $2^{-weight}$。
當 $t=0$ 時 , 將 `val` 左移 `factor` 次即可。
當 $t>0$ 時 , 這邊先將 $\alpha$ 提出去 $S_t=f*\alpha Y_t+(1-\alpha)*S_{t-1} = \alpha [f*Y_t + (\frac{1}{\alpha} -1)S_{t-1}] = \alpha [f*Y_t + (\frac{1}{\alpha}*S_{t-1} -S_{t-1})]$
* $\alpha=2^{-weight}$ , 代表右移`weight` 次
* $\frac1\alpha={2^{weight}}$ , 代表左移 `weight` 次
```diff
struct ewma *ewma_add(struct ewma *avg, unsigned long val)
{
avg->internal = avg->internal
- ? (((avg->internal << EEEE) - avg->internal) +
- (val << FFFF)) >> avg->weight
+ ? (((avg->internal << avg->weight) - avg->internal) +
+ (val << avg->factor)) >> avg->weight
: (val << avg->factor);
return avg;
}
```
### 測驗 `5`
程式碼中 `r` 、 `shift` 和 `x>1` 可以視為 `ilog2` 計算 $log_2(i)$ 的整數部分的結果 , [ceil](https://man7.org/linux/man-pages/man3/ceil.3.html) 為找到不小於 x 的整數值,如何利用 `ilog2` 的結果,計算出 $[log_2(i)]$ 需要考慮以下兩點。
第一點,如果 $log_2(i)$ 有小數 , 將 `ilog2(i)` 的結果加 1 即為答案。
第二點,i 為 2 的冪也就是 $log_2(i)$ 只有整數,不能直接將 `ilog2` 的結果加 1 ,因此透過改變 i 的值,先將 i 減 1 ,再將`ilog2(i)` 的結果加 1 即為答案。
## 2024q1 第 4 週測驗題
### 測驗 `1`
假設有三個數字 $A$ $B$ $C$ 。
$total = (A\bigoplus A+A\bigoplus B+A\bigoplus C)+(B\bigoplus A+B\bigoplus B+B\bigoplus C)+(C\bigoplus A+C\bigoplus B+C\bigoplus C)$
因為 $A\bigoplus A =0$
$total = (A\bigoplus B+A\bigoplus C)+(B\bigoplus A+B\bigoplus C)+(C\bigoplus A+C\bigoplus B)$
$total = 2*(A\bigoplus B+A\bigoplus C+B\bigoplus C)$
所以須將 total 除以2。
```diff
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 >> AAAA;
+ return total >> 1;
}
```
### 測驗 `2`
`popcount(n ^ 0xAAAAAAAA) - 16` 的結果介於 -16~16 ,為了避免負數,選擇加上足夠大的 3 的倍數 39 , 改成 `n = popcount(n ^ 0xAAAAAAAA) + 23` ,此時的 n 介於 23~55 ,需要再進行一次化簡。 55 出現第一個 1 的位數為 5 。 因此本次的 m 要改成 `0x2AA` , `n = popcount(n ^ 0x2A) - 3` 的結果介於 -2~2 , 當結果是負數需要加 3 。因此 `return n + ((n >> 31) & 3)` 有誤,應改成 `return n + (((int)n >> 31) & 3)`。
#### 棋盤解析
在棋盤中,有 8 種不同的連線方式可以構成勝利條件,這包括 3 行、 3 列以及 2 個對角線的連線,因此每 4 個位元 (nibble) 可以作為一的單位來判斷棋盤上是否成功形成連線。
```c
static const uint32_t move_masks[9] = {
0x40040040, 0x20004000, 0x10000404,
0x04020000, 0x02002022, 0x01000200,
0x00410001, 0x00201000, 0x00100110,
};
```
假設第一行成功形成連線,也就是第 0 個、第 1 個以及第 2 個元素做 or 運算,可以得到 `0x70044444` ,假設第三行成功形成連線,也就是第 6 個、第 7 個以及第 8 個元素做 or 運算,可以得到 `0x00711111`。
因此可以看出如果有一組 nibble 的和為 7 代表連線成功,所以將每一組 nibble 加 1 ,跟 8 作 and 運算時就不會是 0。
```diff
- return (player_board + BBBB) & 0x88888888;
+ return (player_board + 0x11111111) & 0x88888888;
```
[Hacker's Delight](https://web.archive.org/web/20180517023231/http://www.hackersdelight.org/divcMore.pdf) 中有提及 mod 7 的作法。
```c
int remu7(unsigned n) {
static char table[75] = {0,1,2,3,4,5,6, 0,1,2,3,4,5,6,
0,1,2,3,4,5,6, 0,1,2,3,4,5,6, 0,1,2,3,4,5,6,
0,1,2,3,4,5,6, 0,1,2,3,4,5,6, 0,1,2,3,4,5,6,
0,1,2,3,4,5,6, 0,1,2,3,4,5,6, 0,1,2,3,4};
n = (n >> 15) + (n & 0x7FFF); // Max 0x27FFE.
n = (n >> 9) + (n & 0x001FF); // Max 0x33D.
n = (n >> 6) + (n & 0x0003F); // Max 0x4A.
return table[n];
}
```
可以發現 `n = (n >> 15) + (n & 0x7FFF)` 為作一次化簡的步驟,但詳細的推導過程還在釐清。
```diff
- x = (x >> CCCC) + (x & UINT32_C(0x7FFF));
+ x = (x >> 15) + (x & UINT32_C(0x7FFF));
```
### 測驗 `3`