# 2024q1 Homework4 (quiz3+4)
contributed by < [`Wufangni`]("https://github.com/Wufangni") >
## 第三週測驗
### 測驗 1
`i_sqrt` 為找出平方根的函式,下面為版本一的作法:
``` c
#include <math.h>
int i_sqrt(int N)
{
int msb = (int) log2(N);
int a = 1 << msb;
int result = 0;
while (a != 0) {
if ((result + a) * (result + a) <= N)
result += a;
a >>= 1;
}
return result;
}
```
使用 `log2` 找出最高位元數 `msb` ,接著利用 `if ((result + a) * (result + a) <= N)` 判別本次迴圈的結果取平方是否小於原數,若不是則持續往下找,直到找出平方可小於原數的結果即為他的平方根。
版本二使用迴圈方式替代 `log2` 函式找出最高位元 `msb`:
``` c
int i_sqrt(int N)
{
int msb = 0;
int n = N;
while (n > 1) {
n >>= 1;
msb++;
}
int a = 1 << msb;
int result = 0;
while (a != 0) {
if ((result + a) * (result + a) <= N)
result += a;
a >>= 1;
}
return result;
}
```
``` c
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;
}
```
### 測驗 2
``` c
carry = tmp / 10;
tmp = tmp % 10;
```
因 `mod` 的計算量較大,根據除法原理利用已計算出來的 `tmp / 10` 可將 `%` 改寫成下面形式:
``` c
carry = tmp / 10;
tmp = tmp - carry * 10;
```
接著著手改寫 `divmod_10` 方法:
``` 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 >> 3);
*mod = in - ((q & ~0x7) + (*div << 1));
}
```
`(in | 1) - (in >> 2)` 可看成 $in-({in \over 4})$ ,化簡開可得 $0.75in$,因此 `q` 的計算結果為 ${3 \over 4}in \cdot {1 \over 16} + {3 \over 4}in$, 再整理後可得 ${102 \over 128}in$,做了四次 `q = (q >> 8) + x;` 讓 `q` 的結果更趨近於 $0.8in$ , 最後可計算 `div` 結果,${8 \over 10}in \cdot {1 \over 8}={1 \over 10}in$。
`mod` 計算則利用 `(q & ~0x7)` 方式,將 `~0x7` 看待成 `0x11...1000` 與 `q` 做 `&` 計算,相等於 `div` 乘上 $8$ ,再與 `(*div << 1)` 相加可得結果為 $div \cdot 8 + div \cdot 2 = div \cdot 10$,根據上述提到的 `mod` 簡化算法 `tmp = tmp - carry * 10` 代入數值,`mod` 即為 $in - div \cdot 10$ 結果。
### 測驗 3
定義一個 `ilog2` 用來找出以 2 為底的對數,每次迴圈中 i 的二進位值都會右移一格,用 `log` 作為計數器計算有幾個 $2^n, n\ge0$ 。
``` c
int ilog2(int i)
{
int log = -1;
while (i) {
i >>= 1;
log++;
}
return log;
}
```
上述方法因每個 `bit` 都需掃過一次做計算,因此改良成下面版本,在迴圈中新增判別式判斷 i 是否大於某個 $2^n, n\ge0$ ,若達成條件則一次移動 n 格 `bit` ,以減少總迴圈次數。
``` 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 中提供的 `__builtin_clz` ,因輸入不能為零所以設置 `v | 1` 作為輸入。
``` c
int ilog32(uint32_t v)
{
return (31 - __builtin_clz(v | 1));
}
```
### 測驗 4
預先定義 `ewma` 結構, `internal` 為當前 `ewma` 的值, `factor` 為縮放因子,`weight` 用來決定衰減的 $\alpha$ 大小。
``` c
struct ewma {
unsigned long internal;
unsigned long factor;
unsigned long weight;
};
```
初始化 `ewma` 結構體,因 `factor` 及 `weight` 必須為 2 的冪,使用 `is_power_of_2` 排除可能報錯的值,接著使用測驗 3 提及的 `ilog2` 取出 2 的次方。
``` c
void ewma_init(struct ewma *avg, unsigned long factor, unsigned long weight)
{
if (!is_power_of_2(weight) || !is_power_of_2(factor))
assert(0 && "weight and factor have to be a power of two!");
avg->weight = ilog2(weight);
avg->factor = ilog2(factor);
avg->internal = 0;
}
```
創建 `*ewma_add` 更新 `ewma` ,添加新資料 `val` ,若 `avg->internal` 為空值代表 `val` 為 $Y_0$,直接填入 `val << avg->factor` ,其餘情況則按照 $\alpha Y_t \ + \ (1-\alpha)\cdot S_{t-1}$ 更新 $S_{t}$。
$\alpha$ 以$\cfrac{1}{2^{avg->weight}}$ 表示,$\alpha Y_t$ 則可以以 `(val << avg->factor)) >> avg->weight` 實現。
$(1 -\alpha)$ 可寫成 $\cfrac{2^{avg->weight}-1} {2^{avg->weight}}$,乘回 `avg->internal` 展開可得 $[(2^{avg->weight})(avg\to internal)-(avg\to internal)]\ /\ 2^{avg->weight}$,前半部分可以 `((avg->internal << avg->weight) - avg->internal)` 實現,最後因 `(val << avg->factor)) >> avg->weight` 也需 `>> avg->weight` 因此將此移至算式尾端當分母。
``` c
struct ewma *ewma_add(struct ewma *avg, unsigned long val)
{
avg->internal = avg->internal
? (((avg->internal << avg->weight) - avg->internal) +
(val << avg->factor)) >> avg->weight
: (val << avg->factor);
return avg;
}
```
## 第四周測驗題
### 測驗 1
計算數值的二進位表示法中有被設置的 `bit` 數量, `popcount_naive` 為其中一種方式:
``` c
unsigned popcount_naive(unsigned v)
{
unsigned n = 0;
while (v)
v &= (v - 1), n = -(~n);
return n;
}
```
利用 `v &= (v - 1)` 重置數值中 `LSB` 位元,並使用 `n = -(~n)` 方式對 `n` 加一,直到所有 `LSB` 位元重置為止。
以另一種方式做改寫:
``` c
unsigned popcount_branchless(unsigned v)
{
unsigned n;
n = (v >> 1) & 0x77777777;
v -= n;
n = (n >> 1) & 0x77777777;
v -= n;
n = (n >> 1) & 0x77777777;
v -= n;
v = (v + (v >> 4)) & 0x0F0F0F0F;
v *= 0x01010101;
return v >> 24;
}
```
以 32bit 舉例,population count 可以下列數學式表示:
$popcount(x)=x-\lfloor \cfrac{x}{2} \rfloor - \lfloor \cfrac{x}{4}\rfloor-...-\lfloor\cfrac{x}{32} \rfloor$
為了避免每個 bit 都做一次計算,可以以每 4 個 bit 同時做計算,即:$popcount(x)=x-\lfloor \cfrac{x}{2} \rfloor - \lfloor \cfrac{x}{4}\rfloor-\lfloor\cfrac{x}{8} \rfloor$
`n = (v >> 1) & 0x77777777` 得出$\lfloor \cfrac{x}{2} \rfloor$,持續累加三次得到以 4 bits 為單位的 `set bit` 數量。
``` c
v = (v + (v >> 4)) & 0x0F0F0F0F;
v *= 0x01010101;
return v >> 24;
```
利用一次位移一個 4 bits 單位做相加乘上 `0x0F0F0F0F` 過濾多餘的內容,可得以下數列, $Bn$ 為第 n 個單位 `set bit` 數量。
```c
0 (B7+B6) 0 (B5+B4) 0 (B3+B2) 0 (B1+B0)
```
最後乘上 `0x01010101` 並取出在原 A6 位置的結果 (B7+B6+B5+B4+B3+B2+B1+B0),即為總 `set bit` 數目。
``` c
int hammingDistance(int x, int y)
{
return __builtin_popcount(x ^ y);
}
```
Hamming distance 為兩數列不相同位元之數目,透過 GCC 內建函式庫中提供 `__builtin_popcount(x)` 可計算 `set bit` 數量,可利用 $XOR$ 特性(相同為0,不同為1)計算 $A\oplus B$ ,得出數列 A 與數列 B 相異位元的數列,即為 `A ^ B`。
``` 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;
}
```
### 測驗 3
從 `xtree.h` 中先觀察 `xt_node` 結構:
``` c
struct xt_node {
short hint;
struct xt_node *parent;
struct xt_node *left, *right;
};
```
`*parent` 作為指向父節點的指標,`*left` 及 `*right` 分別指向左右子樹的根節點,`hint` 則作為判斷是否處於平衡的變數。
利用 `treeint_xt.c` 中提供對 `xt_tree` 的基本操作,完成 `__xt_remove` 及 `xt_insert` 等節點操作。
``` c
if (xt_right(del)) {
struct xt_node *least = xt_first(xt_right(del));
if (del == *root)
*root = least;
xt_replace_right(del, least);
xt_update(root, xt_right(least));
return;
}
```