---
tags: 進階電腦系統理論與實作
---
# 2020q3 Homework2 (quiz2)
contributed by < `RinHizakura` >
> [2020q3 第 2 週測驗題](https://hackmd.io/@sysprog/2020-quiz2)
> [GitHub](https://github.com/RinHizakura/NCKU-System-Software-Quiz/tree/master/quiz2)
## 測驗 `1`
### 程式原理
對於單個字元 1 byte = 8 bits 的 mask 是 0x80，所以推廣到 8 個字元很好理解，答案 `MMM` = `(d) 0x8080808080808080`。
下面也舉例會使錯誤的選項出現判斷錯誤的測資(可以用 8 位元 `is_ascii` 和 64 位元 `is_ascii` 的輸出異同來比較):
* `( a ) 0x0080008000800080`
我們知道，對於 1 byte 的 mask 是 0x80，因此對於此選項，相當於是字串中可以先以每 8 個一數形成的 `payload` 中 index 1、3、5、7 的字元會被錯誤的檢查，舉例來說:
```cpp
char a[] = "12345678";
a[7] = 200;
```
相當於是
| String | | '7' | '6' | '5' | '4' | '3' | '2' | '1' |
|:---------- |:--- | --- | --- | --- | --- | --- |:--- |:--- |
| ASCII(hex) | C8 | 37 | 36 | 35 | 34 | 33 | 32 | 31 |
| Mask(hex) | 00 | 80 | 00 | 80 | 00 | 80 | 00 | 80 |
| & | 00 | 00 | 00 | 00 | 00 | 00 | 00 | 00 |
可以看到 `C8` 沒有對應到 0x80 的 mask，而 00 這個 mask 等於是不做任何檢查，於是此選項的 mask 會出現判斷錯誤。
所以會導致這個選項出現錯誤的測資為:
```cpp
char a[] = "12345678"; /* 或任意合法的 ASCII 形成的字串 */
a[7] = 200;/* 或者調整 a[i] = {128 至 255 的任意整數}
* i 為 1、3、5、7 等奇數 index
* 且非 8 個一數後剩餘的字元 */
```
* `( b ) 0x80 00 80 00 80 00 80 00`
情況類似 (a)，這次是 `payload` 中 index 0、2、4、6 的字元會被錯誤的檢查。
所以會導致這個選項出現錯誤的測資為:
```cpp
char a[] = "12345678"; /* 或任意合法的 ASCII 形成的字串 */
a[6] = 200;/* 或者調整 a[i] = {128 至 255 的任意整數}
* i 為 0、2、4、6 等偶數 index
* 且非 8 個一數後剩餘的字元 */
```
* `( c ) 0x0808080808080808`
相當於是對每個字元都用 08 這個 mask 做檢查，也就是說，只有字元對應的二進位數字中，第 4 個 bit 不為 1 時才被當成 ASCII。
所以會導致這個選項把 ASCII 判斷成非 ASCII 錯誤的測資例如下面的 `a`、`b` 字串，`c` 字串則會正確判斷:
```cpp
char a[] = "HHHHHHHH"; /* 因為 H = 0100 1000 */
char b[] = "OOOOOOOO"; /* 因為 O = 0100 1111 */
^
這裡是 1
char c[] = "PPPPPPPP"; /* 因為 P = 0101 0000 */
^
這裡是 0
```
另一方面，也會把非 ASCII 判斷成 ASCII ，例如以下測資:
```cpp
char a[] = "PPPPPPPP"; /* 或任意第四個 bit 不為 1 的 ASCII 形成的字串 */
a[6] = 200;
```
* `( e ) 0x8888888888888888`
相當於是對每個字元都用 88 這個 mask 做檢查，也就是說，只有字元對應的二進位數字中，第 4、第 8 個 bit 都不為 1 時才被當成 ASCII。
所以 ( c ) 選項中那些會把 ASCII 判斷成非 ASCII 錯誤的測資在這裡同樣會有錯
```cpp
char a[] = "HHHHHHHH"; /* 因為 H = 0100 1000 */
char b[] = "OOOOOOOO"; /* 因為 O = 0100 1111 */
^ ^
這裡是 1
```
不過這裡不會把非 ASCII 判斷成 ASCII ，例如以下測資可以得到預期結果:
```cpp
char a[] = "PPPPPPPP"; /* 或任意第四個 bit 不為 1 的 ASCII 形成的字串 */
a[6] = 200; /* 因為 200 = 1100 1000 */
^ ^
這裡是 1
```
* `( f ) 0xFFFFFFFFFFFFFFFF`
這個 mask 會把所有對應二進位非 0(`'\0'`) 的字元都判斷成非 ASCII...
### Why memcpy？
> 答案是: 因為一個字元是 8 bits，而我們想要一次對 64 個 bits 做處理，所以我們需要 `memcpy` 來將 8 bits * 8 個字元轉為一個 64 bits 的 `payload` 儲存，並用其來與 `0x8080808080808080` 做 and 運算。
雖然上面的敘述沒有錯，但是這裡我更想關注的細節是：「可不可以用 casting 的方法來取得 payload 呢？」，也就是把
```c=
uint64_t payload;
memcpy(&payload, str + i, 8);
```
改成
```c=
uint64_t *payload = (uint64_t *) str;
```
乍看之下，似乎兩者沒有甚麼分別，而且真的去執行的話也有機會得到同樣的結果。事實上，如果放在 intel x86 的系統下，因為硬體可以支援 unalign 的記憶體操作，我們可能就會誤認程式的正確性。但是在 ARM 等不支援 [unaligned memory access](https://heyrick.eu/armwiki/Unaligned_data_access) 的架構上，這樣不合法的操作可能直接導致程式無法執行或者產生非預期的結果。
我們可以打開 [AddressSanitizer](https://en.wikipedia.org/wiki/AddressSanitizer) 來檢查是否有 unaligned memory access:
```
gcc -fsanitize=alignment *.c
```
然後透過類似以下的呼叫來測試使用 `memcpy` 的版本，和使用 casting 的版本。
```c=
char a[] = "ABCDEFGHIJKL";
is_ascii(a + 1, strlen(a) - 1 ));
```
在 casting 的版本中，應該可以得到類似以下的訊息。
```
runtime error: load of misaligned address 0x7ffc3959cad1 for type 'uint64_t', which requires 8 byte alignment
```
那麼，為甚麼這裡 `memcpy` 不會有問題呢？讓我們來瞧瞧 glibc 中的 [memcpy.c](https://github.com/lattera/glibc/blob/master/string/memcpy.c) 是如何實作的。
`memcpy` 會先檢查需要 copy 的 bytes 數量，如果數量很少(OP_T_THRES 除 i386 架構是 8 以外都設定為 16)，就直接逐 bytes 複製，而我們每次 `memcpy` 只要求複製 8 bytes，在這種情況下，alignment 自然就不會有影響。
```c=
if (len >= OP_T_THRES)
{
...
}
/* There are just a few bytes to copy. Use byte memory operations. */
BYTE_COPY_FWD (dstp, srcp, len);
```
對於更長的複製，實際上也無需擔心 alignment 的問題，在以下的程式中可以看到，`memcpy` 會先調整 len 讓 `dstp` aligned，而對於 `srcp` 的 alignment 的處理則在 `WORD_COPY_FWD` 的實作內部，這裡就暫時不細究。
```c=
/* Copy just a few bytes to make DSTP aligned. */
len -= (-dstp) % OPSIZ;
BYTE_COPY_FWD (dstp, srcp, (-dstp) % OPSIZ);
/* Copy whole pages from SRCP to DSTP by virtual address manipulation,
as much as possible. */
PAGE_COPY_FWD_MAYBE (dstp, srcp, len, len);
/* Copy from SRCP to DSTP taking advantage of the known alignment of
DSTP. Number of bytes remaining is put in the third argument,
i.e. in LEN. This number may vary from machine to machine. */
WORD_COPY_FWD (dstp, srcp, len, len);
/* Fall out and copy the tail. */
```
因此，`memcpy` 的使用不但是用來得到 64 bits 的 payload，更可以避免因為對記憶體的 unaligned memory access 產生錯誤。
### 進階應用
#### level 1
延伸上面的技巧，實現：「給予一個已知長度的字串，檢測裡頭是否包含有效的英文大小寫字母」。
事實上，此實作的答案就在測驗 `5` 中，為求程式的可讀性，這裡也順便把 `PACKED_BYTE` 搬過來使用。
關於 `A`、`Z` 等和 `lower_mask` 、 `upper mask` 的操作原理請直接參考測驗 `5`。總而言之，如果字元是小寫的英文子母，則 `lower_mask` 對應的 8 bits 就會是 `0x80`，大寫的字母亦然，`upper_mask` 對應的 8 bits 會是 `0x80`，因此，如果 `str` 的內容僅由大小寫字母組成，則 `lower_mask | upper_mask` 會得到 `0x8080808080808080`，我們可以藉此來判斷一個已知長度的字串，是否包含有效的英文大小寫字母。
```c=
bool is_alpha(const char str[], size_t size)
{
if (size == 0)
return false;
int i = 0;
while ((i + 8) <= size) {
uint64_t payload;
memcpy(&payload, str + i, 8);
uint64_t A = payload + PACKED_BYTE(128 - 'A' );
uint64_t Z = payload + PACKED_BYTE(128 - 'Z' - 1);
uint64_t a = payload + PACKED_BYTE(128 - 'a' );
uint64_t z = payload + PACKED_BYTE(128 - 'z' - 1);
uint64_t lower_mask = (a ^ z) & PACKED_BYTE(0x80);
uint64_t upper_mask = (A ^ Z) & PACKED_BYTE(0x80);
if ((lower_mask | upper_mask ) ^ PACKED_BYTE(0x80)){
return false;
}
i += 8;
}
while (i < size) {
if ( (str[i] < 65 || (str[i] > 90 && str[i] < 97) || str[i] > 122))
return false;
i++;
}
return true;
}
```
#### level 2
Intel® Streaming SIMD Extensions 4 (Intel® SSE4) 是一系列的 SIMD 指令集，其中在 [SSE4.2](https://en.wikipedia.org/wiki/SSE4#SSE4.2) 添加的 STTNI（String and Text New Instructions)，可以在單指令中對 16 bytes 的字節進行字元比對或搜尋等，加速文檔的解析。
參照 [Intel Intrinsics Guide 4.2](https://software.intel.com/sites/landingpage/IntrinsicsGuide/#techs=SSE4_2)，首先嘗試通過 STTNI 指令改寫與上面的 level 1 同等的功能。
```c=
/* 編譯時要加上 -msse4，也要 include 下面兩個 library */
#include <emmintrin.h>
#include <nmmintrin.h>
uint64_t get_m128i_mask(__m128i var)
{
// Here we assume var's last half bits are not needed
uint64_t val[2];
memcpy(val, &var, sizeof(val));
return val[0];
}
bool is_alpha_sse4(const char str[], size_t size)
{
if (size == 0)
return false;
const __m128i ranges = _mm_setr_epi8('a', 'z', 'A', 'Z', 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
const __m128i diff = _mm_set1_epi8(0xFF);
const uint8_t mode =
_SIDD_UBYTE_OPS |
_SIDD_CMP_RANGES |
_SIDD_UNIT_MASK;
unsigned int i = 0;
while ((i + 8) <= size) {
uint64_t payload;
memcpy(&payload, str + i, 8);
uint64_t *payload_ptr = &payload;
__m128i *mem = (__m128i *)(payload_ptr);
const __m128i chunk = _mm_loadu_si128(mem);
const __m128i mask = _mm_cmpestrm(ranges, 4, chunk, 8, mode);
if(get_m128i_mask(_mm_xor_si128(mask, diff)))
return false;
i += 8;
}
while (i < size) {
if ( (str[i] < 65 || (str[i] > 90 && str[i] < 97) || str[i] > 122))
return false;
i++;
}
return true;
}
```
* [`_mm_setr_epi8`](https://software.intel.com/sites/landingpage/IntrinsicsGuide/#text=_mm_setr_epi8&expand=4995) 可以設定一個 128 bits integer 上的值
* [`_mm_set1_epi8`](https://software.intel.com/sites/landingpage/IntrinsicsGuide/#text=_mm_set1_epi8&expand=4966) 可以想成是 128 bits 的 `PACKED_BYTE`，所以這裡是把 diff 設成一個 0xFF....FF 的 __m128i
* 迴圈的概念與前面的概念相同，先提取 8 bytes 到 payload，再把一個指向 payload 的 pointer 轉型成 __m128i 的 pointer
:::warning
1. 事實上，因為 __m128i 可以容納 128 bits，其實我們可以一次提取 16 bytes，在 `str` 較長時會有更佳的效能
2. 這裡先用另一個空間儲存再將指標轉型是為了符合 alignment 和指令 API 的 input 型別，因為還不熟悉 API 的使用所以搞得有點複雜，研究之後有更佳的寫法會再調整
:::
* `_mm_loadu_si128` 取出 chunk 後，透過 `_mm_cmpestrm` 去做比對
* `_mm_cmpestrm (__m128i a, int la, __m128i b, int lb, const int imm8)` 會比較長度 `la` 的 `a` 和長度 `lb` 的 `b`，而比較的模式會根據 imm8 的實際設定去調整，可以看到這裡我們使用的3個參數是：
* _SIDD_UBYTE_OPS: 表示 `__m128i` 中的字元是以 8 bits 為一單位的
* _SIDD_CMP_RANGES: 用範圍的模式來比較(*For each character in a, determine if b[0] <= c <= b[1] or b[1] <= c <= b[2]...*)
* _SIDD_UNIT_MASK: 使回傳的 mask 是 byte/word mask
* 如果符合，得到的 mask 就會是 0xFF...FF，所以可以透過 XOR `diff` 去檢查，使用 `get_m128i_mask` 去得到 XOR 後的結果
延伸到標點符號的比對的話，例如我們假設可以接受的符號是：`{space}`、 `！`、 `,` 、`-`、 `.` 、`？` ，我們只需要對上面的程式稍做修改：
* ranges
```c=
const __m128i ranges = _mm_setr_epi8('a', 'z', 'A', 'Z', 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
```
改成
```c=
#define SPACE 32
const __m128i ranges = _mm_setr_epi8('a', 'z', 'A', 'Z', '?', '?', ',', '.', \
SPACE, '!', 0, 0, 0, 0, 0, 0);
```
:::warning
:warning:注意這裡的字元順序是有意義的！要使判斷符合前面提到 `_SIDD_CMP_RANGES` 的條件。
:::
* `_mm_cmpestrm` 要比較的 element 數量增加
```c=
_mm_cmpestrm(ranges, 10, chunk, 8, mode);
```
* 最後逐字處理的部份也做對應的調整(有點醜XD 其實最好是寫成 macro 或者 enum 可讀性會好一點)
```c=
while (i < size) {
/* 1. lower
* 2. upper
* 3. ?
* 4. , - .
* 5. {space} ! */
if ( !((str[i] >= 65 && str[i] <= 90)\
|| (str[i] >= 97 && str[i] <= 122)\
|| str[i] == 63 \
|| (str[i] >= 44 && str[i] <= 46) \
|| (str[i] >= 32 && str[i] <=33)))
return false;
i++;
}
```
詳細請參照 Github 連結。
## 測驗 `2`
### 程式原理
我們可以通過觀察每個字元的編碼來推理。
| 字元 | b7 | b6 | b5 | b4 | b3 | b2 | b1 | b0 |
| ----- | --- | --- | --- | --- | --- | --- |:--- |:--- |
| `'0'` | 0 | 0 | 1 | 1 | 0 | 0 | 0 | 0 |
| `'9'` | 0 | 0 | 1 | 1 | 1 | 0 | 0 | 1 |
| `'A'` | 0 | 1 | 0 | 0 | 0 | 0 | 0 | 1 |
| `'F'` | 0 | 1 | 0 | 0 | 0 | 1 | 1 | 0 |
| `'a'` | 0 | 1 | 1 | 0 | 0 | 0 | 0 | 1 |
| `'f'` | 0 | 1 | 1 | 0 | 0 | 1 | 1 | 0 |
可以首先觀察出幾個特性：
* `'0'` - `'9'` 的後 4 位(b3 - b0)事實上就是要轉換成的數值
* b6 是判別字元是 `'0'` - `'9'` 系列，還是其他的準則
再從 return 的結果 `(in + shift) & 0xf` 去推論的話：
* return 的結果就是 `(in + shift)` 的後四位元
* MASK 是 `0b01000000`，也就是 `(c) 0x40`
* `'0'` - `'9'`： `letter = 0b00000000`
* 其他： `letter = 0b01000000`
* `'a'` / `'A'` 的後 4 bit 是 0b0001，所以要 shift 到目標的 0b1010 需要 0b1001。
* `AAA = (a) 3`
* `BBB = (a) 6`
順便補上測試正確的程式碼：
```c=
int main(){
char a[] = "123456789AaBbCcDdEeFf";
for(int i = 0; i < strlen(a);i++)
printf("%c %d\n", a[i], hexchar2val(a[i]));
return 0;
}
```
```
1 1
2 2
3 3
4 4
5 5
6 6
7 7
8 8
9 9
A 10
a 10
B 11
b 11
C 12
c 12
D 13
d 13
E 14
e 14
F 15
f 15
```
### Hexspeak
擴充前面的功能，使可以接受如 `0xDEADBEAF` 等字串，並且輸出對應的十進位數值。
```c=
uint64_t hexspeak(char str[])
{
int len = strlen(str);
assert(len <= 18);
assert(str[0] == '0' && (str[1] == 'x' || str[1] == 'X'));
uint8_t tmp[8] = {
0,
};
if (!is_bigendian()) {
unsigned int count = 0;
for (int i = len - 1; i >= 2; i--) {
if ((len - i) & 1)
tmp[count] = hexchar2val(str[i]);
else
tmp[count++] += hexchar2val(str[i]) << 4;
}
} else {
unsigned int count = 7;
for (int i = len - 1; i >= 2; i--) {
if ((len - i) & 1)
tmp[count] = hexchar2val(str[i]);
else
tmp[count--] += hexchar2val(str[i]) << 4;
}
}
uint64_t *ans = (uint64_t*)tmp;
return *ans;
}
```
* 首先檢查輸入字串的長度，這裡我限定得是 `0x` 或者 `0X` 開頭，並且字串長度不能大於 18(也就是只能轉成 (18 - 2)* 4 = 64 個位元)
* 先用一個 tmp 陣列，把轉換後的結果存在裡面，這裡需要注意： 因為 tmp 每個單元是 8 bits，而 str 中每個字元對應的數字是 4 bits 的表示，因此例如連續的兩個字元 `0xAB` 在 tmp 中要存成 `16*0xA + 0xB`
* 為減少計算量，我透過 casting 的技巧把 `tmp[8]` 轉成一個 64 位元的整數，但是因為 endian 的關係，需要對應硬體的 endian 正確的擺放 `tmp` 中的位置
> 可將上方程式碼的第 12 到第 28 行改為透過巨集來區隔實作:
> ```cpp
> #if __BYTE_ORDER == __LITTLE_ENDIAN
> ...
> #else /* Big endian */
> ...
> #endif
> ```
> 記得要 `#include <endian.h>`, 後者在 POSIX 相容環境通常存在，判斷 endianness 的程式碼應該儘量在編譯時期確認。
> :notes: jserv
> :ok_hand: 已經修改至 GitHub 上！
:::warning
延伸問題： 不知道實務上會不會儘量避免寫出根據 endian 而有差異的程式呢？或者根據硬體的架構這樣的檢查和對應的實作可能是有瑕疵的？
:::
## 測驗 `3`
### 程式原理
從函式的結果是輸出餘數來看，也就是求 $n\mod d = n - (n/d) \times d$。不難想像 `n - quotient * D` 的 `quotient` 就是 `n / D` 的整數部份。而從題目裡的提示
$$
\dfrac{n}{d} = n \times (\dfrac{2^N / d}{2^N})
$$
不難看出 `quotient = ((__uint128_t) M * n) >> 64;` 暗示著右移 64 bits 就是除以 $2^N$ ，所以 N = 64。因此 `quotient = (M * n) >> 64` 就是 $n / d$ 的整數部份，所以 `M` 是 ${(2^{64} - 1) / d} + 1$，`XXX` 的答案是 `(f) 0xFFFFFFFFFFFFFFFF`
### jemalloc 的快速除法
在 [include/jemalloc/internal/div.h](https://github.com/jemalloc/jemalloc/blob/dev/include/jemalloc/internal/div.h) 和 [src/div.c](https://github.com/jemalloc/jemalloc/blob/dev/src/div.c) 也有相關的技巧運用。在 div.c 的註解中，有比較詳細數學推導：
* 假設 n = q * d 皆為整數，而我們的目標是求出 q = n / d
* 在上述前提下，對於式子(注意這裡所有都是的除法除至小數點以下，而非取整)
$$
\left\lfloor{\left\lceil\frac{2^k}{d}\right\rceil \times \frac{n}{2^k}} \right\rfloor
$$
可以寫成
$$
= \left\lfloor{\frac{2^k+r}{d} \times \frac{n}{2^k}} \right\rfloor
$$
其中 k 為任意整數，而 r = $-(2^k)\mod d$
* 展開上述的式子，可以得到:
$$
= \left\lfloor{\frac{2^k}{d} \times \frac{n}{2^k} + \frac{r}{d} \times \frac{n}{2^k}} \right\rfloor \\
= \left\lfloor{\frac{n}{d} + \frac{r}{d} \times \frac{n}{2^k}} \right\rfloor
$$
* 如前所述，由於前提為 $\frac{n}{d}$ 為整數，因此上面的式子可以寫成：
$$
\frac{n}{d} + \left\lfloor{\frac{r}{d} \times \frac{n}{2^k}} \right\rfloor
$$
也就是說，如果可以滿足
$$
\frac{r}{d} \times \frac{n}{2^k} < 1
$$
則 n / d 可以被改寫成 $\left\lfloor{\left\lceil\frac{2^k}{d}\right\rceil \times \frac{n}{2^k}} \right\rfloor$
* 因為 r < d 且 r / d < 1 總是滿足，為使 $n / 2^k$ < 1，對於不會超出 $2^{32}$ 的 n，可以令 k = 32 來滿足此條件
#### `div_init`
```cpp
void
div_init(div_info_t *div_info, size_t d) {
/* Nonsensical. */
assert(d != 0);
/*
* This would make the value of magic too high to fit into a uint32_t
* (we would want magic = 2^32 exactly). This would mess with code gen
* on 32-bit machines.
*/
assert(d != 1);
uint64_t two_to_k = ((uint64_t)1 << 32);
uint32_t magic = (uint32_t)(two_to_k / d);
/*
* We want magic = ceil(2^k / d), but C gives us floor. We have to
* increment it unless the result was exact (i.e. unless d is a power of
* two).
*/
if (two_to_k % d != 0) {
magic++;
}
div_info->magic = magic;
#ifdef JEMALLOC_DEBUG
div_info->d = d;
#endif
}
```
對數學的推導有一定的了解後，整個快速除法的部份並不難理解，`div_init` 首先初始化除數相關的資料：
* 除數不允許是 `0` (不能除 0) 或 `1` (因為 r / d < 1 會不成立)
* 計算出 magic number = $\lceil \frac{2^{32}}{d} \rceil$。可以看到計算的方式是先求出 $\lfloor 2^{32} / d \rfloor$，再判斷 $2^{32}$ 是不是可以被 d 整除，如果不行，對 magic 加一(原因請讀註解)。這其實也呼應原題中對 `M` 的計算，雖然計算的思路不同，但同樣是在取 $\lceil \frac{2^{32}}{d} \rceil$
* magic number 的運用會取代除 d 的運算，在 DEBUG 下則仍會在 struct 中保留 d 的資訊
#### `div_compute`
```c=
static inline size_t
div_compute(div_info_t *div_info, size_t n) {
assert(n <= (uint32_t)-1);
/*
* This generates, e.g. mov; imul; shr on x86-64. On a 32-bit machine,
* the compilers I tried were all smart enough to turn this into the
* appropriate "get the high 32 bits of the result of a multiply" (e.g.
* mul; mov edx eax; on x86, umull on arm, etc.).
*/
size_t i = ((uint64_t)n * (uint64_t)div_info->magic) >> 32;
#ifdef JEMALLOC_DEBUG
assert(i * div_info->d == n);
#endif
return i;
}
```
n / d 的計算轉為呼叫 `div_compute` 來實現，計算被替換成:
$$
\frac{n \times \lceil \frac{2^{32}}{d} \rceil}{2^{32}}
$$
### 快速除法
要抽離 jemalloc 的快速除法，直接將程式中的的資料結構拉出來使用即可：
```cpp
typedef struct div_info_s div_info_t;
struct div_info_s {
uint32_t magic;
};
void div_init(div_info_t *div_info, size_t d) {
assert(d != 0);
assert(d != 1);
uint64_t two_to_k = ((uint64_t)1 << 32);
uint32_t magic = (uint32_t)(two_to_k / d);
if (two_to_k % d != 0) {
magic++;
}
div_info->magic = magic;
}
static inline size_t div_compute(div_info_t *div_info, size_t n) {
assert(n <= (uint32_t)-1);
size_t i = ((uint64_t)n * (uint64_t)div_info->magic) >> 32;
return i;
}
```
為了實驗是否實作是否有顯著的效益，結合統計方法(詳細請參考 [perf.py](https://github.com/RinHizakura/NCKU-System-Software-Quiz2/blob/master/scripts/perf.py))，透過以下的實驗來驗證，並視覺化結果。
```cpp
static uint64_t MAX = ((uint64_t)(1) << 32) - 1;
int main()
{
struct timespec tt1, tt2;
clock_gettime(CLOCK_MONOTONIC, &tt1);
for (uint64_t i = 2; i < 5000; i++) {
div_info_t DIV;
div_init(&DIV, i);
uint64_t num = rand() % MAX;
clock_gettime(CLOCK_MONOTONIC, &tt1);
size_t ans1 = div_compute(&DIV, num);
clock_gettime(CLOCK_MONOTONIC, &tt2);
long long time1 = (long long) (tt2.tv_sec * 1e9 + tt2.tv_nsec) -
(long long) (tt1.tv_sec * 1e9 + tt1.tv_nsec);
clock_gettime(CLOCK_MONOTONIC, &tt1);
size_t ans2 = num / i;
clock_gettime(CLOCK_MONOTONIC, &tt2);
long long time2 = (long long) (tt2.tv_sec * 1e9 + tt2.tv_nsec) -
(long long) (tt1.tv_sec * 1e9 + tt1.tv_nsec);
printf("%ld, %lld, %lld\n", i, time1, time2);
}
}
```
嘗試了各種最佳化參數，結果各自如下，其中左圖為所有除數的實驗結果，右圖則排除最左側的資料(divisor = 2、3)，放大整個圖的後面部份以利更詳細的觀察：
* `-O0`
* `-O1`
* `-O2`
觀察這個結果，有許多問題可以思考，後續要補充相關的實驗來驗證：
* 在沒有最佳化的情況下，fastdiv 是有明顯的優於一般的除法的，但當最佳化的程度越高，可以看到兩者的表現越來越接近，是否編譯器背後也透過同樣技巧加速？或者有其他可能可以最佳化之處被編譯器處理？例如因為被除數是隨機產生的，如果被除數小於除數，那麼直接最佳化成回傳 0 即可，不需要真正去運算。
* 為了避免固定的被除數可能會被編譯器最佳化掉，我刻意讓其為隨機更動，讓同一個除數下使用相同的被除數以減少誤差，但是不固定的被除數是否有其他的問題呢？有沒有更多實驗的細節需要被考量？
### 最佳化的觀察與實驗改進
因為我只想比對關鍵的除法部份，在儘量不影響產生出的 assembly code 有差別的情況下，透過下面的化簡程式來產生 assembly code:
```c=
#include <stdlib.h>
#include <stdint.h>
static uint64_t MAX = ((uint64_t)(1) << 32) - 1;
int main(){
for (uint64_t i = 2; i < 5000; i++) {
uint64_t num = rand() % MAX;
size_t ans = num / i;
}
return 0;
}
```
來比較一下 `-O0` 和 `-O2` 的 assembly code 的關鍵迴圈部份:
* `-O0` :
```c=
.L3:
call rand@PLT
cltq
movq MAX(%rip), %rcx
movl $0, %edx
divq %rcx
movq %rdx, -16(%rbp)
movq -16(%rbp), %rax
movl $0, %edx
divq -24(%rbp)
movq %rax, -8(%rbp)
addq $1, -24(%rbp)
.L2:
cmpq $4999, -24(%rbp)
jbe .L3
```
> divq S 是把 %rdx:%rax 除以 S，然後商會存在 stored in %rax，餘數會存在 stored in %rdx
第一個 divq 是把 `rand() % MAX` 的結果算出放在 `-16(%rbp)`，第二個 divq 則是求出 ans 的值。總而言之，可以看到，在 `-O0` 的狀況下，除法的指令被使用。
* `-O2` :
```c=
main:
.LFB16:
.cfi_startproc
endbr64
pushq %rbx
.cfi_def_cfa_offset 16
.cfi_offset 3, -16
movl $4998, %ebx
.p2align 4,,10
.p2align 3
.L2:
call rand@PLT
subq $1, %rbx
jne .L2
```
在 `-O2` 下觀察 assembly 後不禁恍然大悟，因為我計算出的結果 `ans` 根本沒有要使用，所以就根本沒必要真的去計算囉，所以程式就只是乖乖的去數完 for 迴圈而已XD。如果實際去使用該 `ans` (例如 printf)，事實上就會使用到除法指令了。
所以當我更動原始程式，額外去印出計算的結果，再去觀察 -O2 下時間的差異時，就可以看出其中的分別，算是給自己的愚蠢上了一課了:D
### 修正錯誤的程式理解
根據方鈺學同學的指教(相當感謝 :pray:)，考慮到正確性來說，這裡的測資其實是有問題的。雖然我認為實驗上的時間效率仍然是足以採信的，不過老師也說了：
> 既然都可能不正確了，還談什麼效率呢？
> :notes: jserv
因此，這裡我們調整可以執行出正確結果的測資，再次進行實驗。
```cpp
int main()
{
struct timespec tt1, tt2;
clock_gettime(CLOCK_MONOTONIC, &tt1);
for (uint64_t i = 2; i < 5000; i++) {
div_info_t DIV;
div_init(&DIV, i);
uint64_t num = rand() % 10000;
num *= i;
clock_gettime(CLOCK_MONOTONIC, &tt1);
size_t ans1 = div_compute(&DIV, num);
clock_gettime(CLOCK_MONOTONIC, &tt2);
long long time1 = (long long) (tt2.tv_sec * 1e9 + tt2.tv_nsec) -
(long long) (tt1.tv_sec * 1e9 + tt1.tv_nsec);
clock_gettime(CLOCK_MONOTONIC, &tt1);
size_t ans2 = num / i;
clock_gettime(CLOCK_MONOTONIC, &tt2);
long long time2 = (long long) (tt2.tv_sec * 1e9 + tt2.tv_nsec) -
(long long) (tt1.tv_sec * 1e9 + tt1.tv_nsec);
printf("%ld, %lld, %lld, %ld, %ld\n", i, time1, time2, ans1, ans2);
}
}
```
為了符合 div 的假設，因此需讓被除數 n 可以整除 d，才可以得到正確的計算結果，而兩者的時間效率差異則如下圖所示:
## 測驗 `4`
### 程式原理
前面我們也提到 M 對應到 $\lceil \frac{2^{32}}{d} \rceil$，因此 n * M 可以想成是把數字位移到剩下 n / D 的小數部份，而 M - 1 則是 1 / D 的小數部份。要判斷 n 是否可以整除 D，思路是判斷 n/D 的小數部份是否小於 1/D (因為小數部份只能是 0、1/D、2/D...D-1/D)，因此答案是 `YYY` = `(c) M - 1`
> 延伸閱讀： [Faster remainders when the divisor is a constant: beating compilers and libdivide](https://lemire.me/blog/2019/02/08/faster-remainders-when-the-divisor-is-a-constant-beating-compilers-and-libdivide/)
> 等你研讀 jemalloc 原始程式碼後，應該會有新的啟發
> :notes: jserv
## 測驗 `5`
### 程式原理
* `PACKED_BYTE(b)` 的作用是取出 `b` 的後 8 個位元(假設是 0xAB)，並且變成 0xABABABABABABABAB
* 因此 `VV1` 就是我們在測驗 `1` 中有用到的技巧 = `(e) 0x80`。
* 不難猜測到，變數 `A` 和 `Z` 的用途就是判斷 chunk 中各字元的實際內容是否是 `A` 到 `Z` 之間，先思考一下如果是對單個字元 c，判斷式就是：
```c=
if(c - 'A' >= 0 && c - 'Z' <= 0)
```
* 再思考一下，`128 - 'A'` 是從 `'A'` 到 128 的距離，也就是說，如果 `c + 128 - 'A'` 的第一個(最左邊)的 bit 為 1，表示 c >= 'A'。
* 同理 `128 - 'Z' - 1` 是從 `'Z'` 到 127 的距離，如果 `c + 128 - 'Z' - 1` 的第一個(最左邊)的 bit 是 1，表示 c > `'Z'` ，如果是 0，表示 c <= `'Z'`。
* 因此如果 c 同時滿足 >= `'A'` 且 <= `'Z'`，`c + 128 - 'A'` XOR `c + 128 - 'Z' + 1` 的結果會是 0b10000000，否則就是 0b00000000。
* 又如果 s 是 `'A'` ~ `'Z'` 之間的字元，最終對 chunk 處理的 mask 是 0b00100000(對應 strtolower 的 `1 << 5`)，也就是前面的結果 `0b10000000 >> 2` 所得到的 mask，不然就是 0b00000000。
因此:
* `VV2` = `(b) 0`
* `VV3` = `(c) -1`
* `VV4` = `(e) 0x80`
順便附上一個可以很好檢驗正確與否的測資：
```c=
int main()
{
char str[] =
"?@ABYZ[\\^";
int n = strlen(str);
strlower_vector(str, n);
if(strcmp(str,"?@abyz[\\^")==0)
printf("yes!!\n");
puts(str);
return 0;
}
```
### `char str[]` or `char *str`？
如果我們把前述程式的 `char str[]` 至換成 `char *str` 會看到 `Segmentation fault (core dumped)` 的錯誤訊息，為甚麼呢？
在規格書 **J.2 Undefined behavior** 有一項是：
> — The program attempts to modify a string literal (6.4.5).
**6.4.5 String literals** 中的第 7、8 點中則提到：
> In translation phase 7, a byte or code of value zero is appended to each multibyte character sequence that results from a string literal or literals. The multibyte character sequence is then used to initialize an array of static storage duration and length just sufficient to contain the sequence.
> It is unspecified whether these arrays are distinct provided their elements have the appropriate values. If the program attempts to modify such an array, the behavior is undefined.
string literal 是靜態配置的，嘗試去更動其內容是 undefined behavior。
而 `char str[]` 和 `char *str` 在語言中的行為，根據 C 語言規格書 **6.7.9 Initialization** 的第 32 點：
> EXAMPLE 8 The declaration
>
> **char s[] = "abc", t[3] = "abc";**
> defines ‘‘plain’’ char array objects s and t whose elements are initialized with character string literals. This declaration is identical to
>
> **char s[] = { 'a', 'b', 'c', '\0' },**
> **t[] = { 'a', 'b', 'c' };**
>
> The contents of the arrays are modifiable. On the other hand, the declaration
>
> **char \*p = "abc";**
>
> defines p with type ‘‘pointer to char’’ and initializes it to point to an object with type ‘‘array of char’’ with length 4 whose elements are initialized with a character string literal. If an attempt is made to use p to modify the contents of the array, the behavior is undefined
`char str[]` 的行為是用 string literal 建立一個 char array object 複本，所以物件會被配置在 stack 中，其內容是可以被改變的;而 `char *str` 則是建立一個直接指向 string literal 的指標，如前所述，嘗試去改寫 string literal 的內容是 undefined behavior。
## 測驗 `6`
### 程式原理
#### 由淺入深
讓我們把題目從 「有一個數字只會出現一次，其他數字都會出現 2 次」的解題思考起。
所以我們就可以對應每個 bit 的出現次數用一個 bit 去紀錄：首先初始為 0，如果第 n 個 bit 出現第一次，設為 1，再出現一次，設為 0。如此一來，紀錄的 bit 為 1 就表示對應的 bit 只出現過一次，那些 bit 編碼起來的數字就是我們想要的答案。
然而該如何實現呢？我們可以把 XOR 當成「兩個 bit 相加並且不進位」的運算，則在上面的題意下，我們可以透過 XOR 的特性，用一個整數的 n 個 bit 對應第 n 個 bit 出現的次數，每當該 bit 出現，就把那個 bit XOR 1，則最後該 interger 上的數字就會是那一個只出現一次的數字。
#### 進階思考
推廣到「出現三次」的話該怎麼做呢？我們可以延續前面的想法，但現在要紀錄 bit 出現的狀態就需要 3 種。要對此紀錄，我們可以用兩個整數的 n 個 bit 對應第 n 個 bit 出現的次數。首先初始為 00，如果第 n 個 bit 出現第一次，設為 01，第二次，設為 10，第三次設回 00。如此一來，狀態為 01 的 bit 編碼起來的數字就是我們想要的答案。
那麼這該如何實作呢？我們可以觀察下面的狀態變化表來思考。
| input | higher lower | higher' lower' |
| ----- | ------------ |:-------------- |
| 0 | 00 | 00 |
| 0 | 01 | 01 |
| 0 | 10 | 10 |
| 1 | 00 | 01 |
| 1 | 01 | 10 |
| 1 | 10 | 00 |
我們可以把 bit 的變化彙整到程式邏輯：
* 如果第 n 個 bit 為 1，則先對 lower 做 XOR 1
* 觀察上面的表格，如果當前狀態是 10 的話，lower bit XOR 1 不是預想的結果(狀態是 10 的時候，lower bit 是 0 變成 0)。要調整成 0 的話，可以判斷現在的狀態是不是 10，也就是 higher bit，把 higher 取成補數 0 做 mask
* 如果第 n 個 bit 為 1，則先對 higher 做 XOR 1
* 觀察上面的表格，如果當前狀態是 00 的話，higher bit XOR 1 不是預想的結果(狀態是 00 的時候，higher bit 是 0 變成 0)。要調整成 0 的話，可以判斷現在的狀態是不是 00，也就是 lower bit，把 lower 取成補數 0 做 mask(注意原本狀態是 00 的話，到 higher 的處理時 lower 已經變成 1)
因此
`KKK` = `(c) ~higher`
`JJJ` = `(d) ~lower`
### 不同的解題思路
我們可以換個角度思考此題的解法。由於數字只會有出現三次和一次的兩種，如果我們可以設計一種神奇(?)的加法運算: 對於每個 bit ，如果出現 3 次 1 時就變成 0，則經過這個神奇的加法後，得到的結果就是我們的答案。舉例來說 [6, 6, 1, 6, 7, 7, 7]:
| | bit2 | bit0 | bit1 |
| ------------ | ---------------------------- | --------------------------- | ---------------------------- |
| 6 | <font color="#f00">1 </font> | <font color="#f00">1</font> | 0 |
| 6 | <font color="#f00">1 </font> | <font color="#f00">1</font> | 0 |
| 1 | 0 | 0 | <font color="#f00">1 </font> |
| 6 | <font color="#f00">1 </font> | <font color="#f00">1</font> | 0 |
| 7 | <font color="#ff0">1 </font> | <font color="#ff0">1</font> | <font color="#f00">1 </font> |
| 7 | <font color="#ff0">1 </font> | <font color="#ff0">1</font> | <font color="#f00">1 </font> |
| 7 | <font color="#ff0">1 </font> | <font color="#ff0">1</font> | 1 |
| magic sum(?) | 0 | 0 | 1 |
因此我們可以透過每個 bit 分開計算有幾個 1，得到正確答案，有幾點需要注意的：
* 因為數字只會出現 3 次或 1 次，所以 bit % 3 理論上只該等於 0 或 1
* 要寫成 `(unsigned int)(bit % 3) << i` 是因為 int << 31 是 undefined behavior
```c=
int singleNumber(int* nums, int numsSize){
int bit_count = 0;
for (int i = 0; i < 32 ; i++){
int bit = 0;
for(int idx = 0; idx < numsSize; idx++){
bit += (nums[idx] >> i) & 1;
}
bit_count |= (unsigned int)(bit % 3) << i;
}
return bit_count;
}
```
附上在 leetcode 中執行正確且通過的證明：
Sign in with Wallet
Connect another wallet