# F05: review
contributed by < `afcidk` >
* [F05: review](https://hackmd.io/s/S1-wcjavN)
## W2 測驗 2
population count 簡稱 popcount 或叫 sideways sum，是計算數值的二進位表示中，有多少位元是 1，在一些場合下很有用，例如計算 0-1 稀疏矩陣 (sparse matrix)或 bit array 中非 0 元素個數、計算兩個字串的 Hamming distance。Intel 在 2008 年 11 月 Nehalem 架構的處理器 Core i7 引入 SSE4.2 指令集，其中就有 CRC32 和 POPCNT 指令，POPCNT 可處理 16-bit, 32-bit, 64-bit 整數。
對應到 C 程式的實作:
```clike
unsigned popcnt_naive(unsigned v) {
unsigned n = 0;
while (v)
v &= (v - 1), n = -(~n);
return n;
}
```
可改寫為以下常數時間的實作:
```clike
unsigned popcnt(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;
}
```
### 原理
我們先看到
```clike
n = (v >> 1) & 0x77777777;
v -= n;
n = (n >> 1) & 0x77777777;
v -= n;
n = (n >> 1) & 0x77777777;
v -= n;
```
由 0x7 的二進位表示 0b0111 可以猜測這邊以 4 bit 為一組取出各組的 lsb 3 bit。
展開這3個部份會得到的 v 如下表示
$\begin{align}
v&=(2^3a_3+2^2a_2+2^1a_1+2^0a_0)-(2^2a_3+2^1a_2+2^0a_1)-(2^1a_3+2^0a_2)-(2^0a_3) \\
&= a_3+a_2+a_1+a_0
\end{align}$
這裡的 $a_n$ 是以 4 bit 為一組的，到了這個步驟我們已經把每 4 bit 的 population count 整理好了。v 現在表示為 $A_0A_1A_2A_3A_4A_5A_6A_7A_8$。
> $A_n$=$a_0a_1a_2a_3$
```clike
v = (v + (v >> 4)) & 0x0F0F0F0F;
```
類似上面的步驟，只是這次我們以 **4** 個 4 bit 為單位，這裡我們將偶數單位的 population count 加到奇數單位上 ($A_1=A_1+ A_2$, $A_3=A_2+ A_3$...) 再把偶數單位清空，這時候我們得到以 8 bit 為單位的 $v=B_1B_2B_3B_4$
```
A0 A1 A2 A3 A4 A5 A6 A7 <-- v
+) A0 A1 A2 A3 A4 A5 A6| A7 <-- v>>4
&) 0 F 0 F 0 F 0 F <-- mask
-------------------------------
B1 B2 B3 B4
```
```clike
v *= 0x01010101;
```
對 v 乘上 0x01010101，我們用直式來看
```
B1 B2 B3 B4
x 1 1 1 1
--------------------
| B1 B2 B3 B4
B1 | B2 B3 B4
B1 B2 | B3 B4
B1 B2 B3 | B4
----------------------------
will be cut|32 bit unsigned integer
```
發現在 第 8\*3=24 bit 的地方會得到 $B_1+B_2+B_3+B_4$，因為我們使用的是 32 bit 無號整數，而 $B$ 以 8 bit 為一單位，後面的數字會被截斷因此不需要做特別處理。
取得 population count 為 `v>>24`
### 應用場景
## W2 測驗 3
考慮以下程式碼:
```clike=
#include <stdio.h>
#define cons(x, y) (struct llist[]){{x,y}}
struct llist { int val; struct llist *next; };
int main() {
struct llist *list = cons(9, cons(5, cons(4, cons(7, NULL))));
struct llist *p = list;
for (; p; p = p->next)
printf("%d", p->val);
printf("\n");
return 0;
}
```
我們需要補完的是 cons 的巨集，我們看到在第 7 行的地方被包了 4 層 `cons`。因為我們直接把展開後的 `cons(.., cons(.., ....))` 指派給 list，所以可以知道展開後的巨集應該是 compound literal 或是有需要的話再加上 designated initializer。
> As a GNU extension, GCC allows **initialization of objects with static storage** duration by compound literals (which is not possible in ISO C99, because the initializer is not a constant).
如果我們想要在 compound literal 中配置空間給 structure 中的 `*next` 的話，我們必須要在前面加上 `(struct llist[])` 告訴 compiler 這個地方需要配置空間。
:::info
找不太到相關資料為什麼在巨集的地方需要兩個大括號，只用一個括號會出現 -Wmissing-braces 警告
:::
## W3(上) 測驗 2
考慮以下計算平方根的程式:
```clike
#include <stdint.h>
/* The algorithm uses the property that computing x² is trivial compared to
* sqrt. It will search the biggest x so that x² <= v, assuming we compute
* sqrt(v).
*/
int32_t sqrt_I2I(int32_t v) {
uint32_t r = v; // r = v - x²
uint32_t b = 0x40000000; // a²
uint32_t q = 0; // 2ax
while (b > 0) {
uint32_t t = q + b; // t = 2ax + a²
q >>= 1; // if a' = a/2, then q' = q/2
if (r >= t) { // if (v - x²) >= 2ax + a²
r -= t; // r' = (v - x²) - (2ax + a²)
q += b; // if x' = (x + a) then ax' = ax + a²,
// thus q' = q' + b
}
b >>= 2; // if a' = a/2, then b' = b / 4
}
return q;
}
```
該函式利用 $(x + a)^2 = x^2 + 2ax + a^2$ 漸進找出平方根。不過上面的程式僅能處理整數，我們引入 fixed point 來表達浮點數，改寫程式碼如下:
```clike=
/* A fixed point is a 32 bit value with the comma between
* the bits 15 and 16, where bit 0 is the less significant bit of the value.
* FIXME: Larger numbers might cause overflow and rounding error.
*/
typedef int32_t fixed;
fixed sqrt_I2F(int32_t v) {
if (!v) return 0;
uint32_t r = v;
uint32_t b = 0x40000000;
uint32_t q = 0;
while (b > 0) {
uint32_t t = q + b;
if (r >= t) {
r -= t, q = t + b;
}
r<<=1, b>>=1;
}
if (r >= q) ++q;
return q;
}
```
### 原理
我們先理解 `sqrt_I2I` 的運作方式。
給定一個 32 bit 無號數 v，回傳一個無條件設去小數點的整數。
我們設定兩個變數 a, x。a 的初始值是 0x40000000，往下收斂到 0。x 的初始值是 0，會慢慢變大，直到 $(a+x)^2 \gt v$，這時候我們就可以得到 $x \approx sqrt(v)$
在 while 迴圈中，我們每個 iteration 都會將範圍除以二減少搜尋範圍，當發現 $v-x^2 \ge a^2+2ax$ 的時候，代表 $(a+x)^2 \le v$ 有成立，這時候我們更新 a 以及 x 的值，但是還不能停下來，因為我們希望 $x^2$ 最逼近 $v$。
當 $a$ 變成 0 的時候，代表我們沒辦法繼續逼近了，於是回傳最後的 x 值。
再來以 `sqrt_I2I` 的概念擴展成 fixed point 版本，從程式碼上面看可以發現其實我們不需要做太多更改，唯一的差異就是在每次該如何收斂 a 以及 x。
看到關鍵的收斂方式
```clike
r <<= 1;
b >>= 1;
```
r 每次會向左位移，相對來說我們可以視為小數點往右位移了。這樣有點像是二進位版本的[十分逼近法](https://medium.com/@pakls/%E5%8D%81%E5%88%86%E9%80%BC%E8%BF%91%E6%B3%95%E9%80%9F%E7%AE%97%E6%B3%95-46c7a29f4ff6)，我們每次都會往下一個小數點逼近。
### 修正 overflow
可能發生 overflow 的地方在第 17 行的 `r <<= 1` 中，我們每次將 r 向右位移為了得到更精準的小數。因為整數的表示位數有限，向又位移可能就會導致 overflow 的發生。
上面的程式碼沒有修正這個問題, 舉例來說，當我們輸入 2000000000 時，會得到 -23923.929688 (轉換回浮點數後)。
第1個解法是我們可以先將 `fixed` type 改為 `uint32_t`，因為我們不預期應該得到負數，如果得到負數就是從 overflow 得到的，又因為我們把他視為有號整數，所以才會出錯。
第2個方式是我們可以將 r 的型態從 `uint32_t` 改為 `uint64_t`，輸入的型態依然使用 `int32_t`，透過這樣的方式，我們可以確保 r 向右位移不會產生 overflow，因為演算法中，我們最多位移 30 次，加上 32 bit 的整數並不會產生 overflow。
> 30 = log2(0x40000000)
> 30 + 32 = 62 bits < 64 bits
實際在程式碼驗證時，我們可以在 `sqrt_I2F` 中加上 `assert(r < r<<1);` 並讓函式的參數跑過整個 32 bit 大小。
Sign in with Wallet
Connect another wallet