# 2020q3 Homework3 (quiz3) contributed by < `Veternal1226` > ###### tags: `sysprog2020` --- ### 測驗 `1` 依據 ISO/IEC 9899:TC2 (即 C99) 標準的 6.5.7 章節 (第 84 到第 85 頁): >“The result of E1 >> E2 is E1 right-shifted E2 bit positions. … If E1 has a signed type and a negative value, the resulting value is implementation-defined.” 在 C 語言中,對有號型態的物件進行算術右移 (arithmetic right shift,以下簡稱 ASR 或 asr) 歸屬於 unspecified behavior,包含 Cppcheck 在內的靜態分析器會嘗試找出 C 程式中這項 unspecified behavior。考慮到程式碼的相容性,我們想實作一套可跨越平台的 ASR,假設在二補數系統,無論標的環境的 shift 究竟輸出為何,我們的 ASR 總可做到: $−5≫^{arith}1=−3$ 上述 $−3$ 正是 $\tfrac{−5}{2}$ 的輸出,並向 $−∞$ 進位 (rounding)。 針對有號整數的跨平台 ASR 實作如下: ```cpp=1 int asr_i(signed int m, unsigned int n) { const int logical = (((int) -1) OP1) > 0; unsigned int fixu = -(logical & (OP2)); int fix = *(int *) &fixu; return (m >> n) | (fix ^ (fix >> n)); } ``` 請補完程式碼。 作答區 OP1 = (b) >> 1 OP2 = (c\) m < 0 提示: 透過 gcc/clang 編譯程式碼時,可加上 -fsanitize=undefined 編譯參數,在執行時期若遇到未定義行為,會得到以下錯誤訊息: >runtime error: load of misaligned address 0x7ffd8a89be8f for type ‘int’, which requires 4 byte alignment 補充資訊: Arm Assembly Language Programming: Arithmetic Shift Operations #### 解題想法 //TODO :::success 延伸問題: 1. 解釋上述程式運作原理,應提供數學證明和 Graphviz 製作的圖解; 2. 練習實作其他資料寬度的 ASR,可參照 C 語言:前置處理器應用篇 撰寫通用的巨集以強化程式碼的共用程度; ::: //TODO --- ### 測驗 `2` 在 C 語言:數值系統篇 中,我們介紹過 power of 2 (漢字可寫作「二的冪」)。 >女星楊冪的名字裡,也有一個「冪」字。且,她取這個名字,就有「次方」的意義,因為她一家三口都姓楊,所以是「楊的三次方」。 GCC extension 其中兩個是 ctz 和 clz: >Built-in Function: `int __builtin_ctz (unsigned int x)` > - Returns the number of trailing 0-bits in x, starting at the least significant bit position. If x is 0, the result is undefined. >Built-in Function: `int __builtin_clz (unsigned int x)` > - Returns the number of leading 0-bits in x, starting at the most significant bit position. If x is 0, the result is undefined. 我們可用 `__builtin_ctz` 來實作 LeetCode 342. Power of Four,假設 int 為 32 位元寬度。實作程式碼如下: ```cpp=1 bool isPowerOfFour(int num) { return num > 0 && (num & (num - 1))==0 && !(__builtin_ctz(num) OPQ); } ``` 請補完程式碼。 作答區 OPQ = (f) & 0x1 #### 解題想法 先觀察4的冪次 |Dec|Bin|`__builtin_ctz(num)`| |-|-|-| |4|0000 0100|3| |16|0001 0000|5| |64|0100 0000|7| 推測要找到 3, 5, 7 的共通點,也就是第 0 bit 為 1 ,但又得排除`__builtin_ctz(num)` 為 1 的情況,此條件正好為第二項判斷式 `(num & (num - 1))==0`,原判斷式應為`__builtin_ctz(num) & 0x1 == 1`,可改寫成`!(__builtin_ctz(num) & 0x1)`,因此這題選擇 ==(f) & 0x1==。 :::success 延伸問題: 1. 解釋上述程式運作原理; ::: `num > 0` :判斷輸入的數是否為正數,若為負數或 0 則 false。 `(num & (num - 1))==0` :判斷輸入的數是否為偶數,若為奇數則 false。 `!(__builtin_ctz(num) & 0x1)` :依照 `__builtin_ctz(num)` 回傳結果判斷是否最右 bit 為 1 ,若回傳結果為0,代表 num 為 4 的冪次方。 :::success 延伸問題: 3. 練習 LeetCode 1009. Complement of Base 10 Integer 和 41. First Missing Positive,應善用 clz; ::: #### LeetCode 1009. Complement of Base 10 Integer 想法: 先做輸入為 0 的例外處理:直接回傳 1。 設計一個 mask,此 mask 用來選擇有效的位元數。 >例如:$(5)_{10}=(101)_{2}$ , mask 應為 $111$。 >例如:$(10)_{10}=(1010)_{2}$ , mask 應為 $1111$。 題目要求將輸入值 N 在有效的位元數內做反向,因此將輸入值 N 與此 mask 做 XOR 運算,得到的值即為答案(回傳值)。 程式碼: ```cpp=1 int bitwiseComplement(int N){ if (!N) return 1; int mask = 0xffffffff >> __builtin_clz(N); return N ^ mask; } ``` 結果: ![](https://i.imgur.com/NGgikr2.png) --- #### LeetCode 41. First Missing Positive WA紀錄: 第一次:忘記處理空陣列,result未初始化。 第二次:忘記處理完全排序好的陣列(如:[1][1,2,3]等等),結果的產生方式要與因為缺失而跳出的處理方式不同,因此不使用result存結果再回傳,改為發現缺失時直接回傳。 原本的錯誤寫法: ```cpp=17 for(i=0;i<numsSize;i++) { if(nums[i] != i+1) { result=i+1; break; } } return result; ``` 想法: 遍歷整個陣列,將陣列中的元素交換到其值減一的位置 (例如 2 應該交換到 nums[1] 中) ,最後再檢查陣列中的index是否與其值相同 (程式第 18 行) 。 交換前,需先檢查值是否大於 0 且小於 numsSize , ```cpp=9 if(nums[i]<numsSize && nums[i]>0) ``` 若要交換的值已經在正確位置或要交換的目標位置已經有正確值,則不進行交換, ```cpp=11 if(pos != i && nums[i] != nums[pos]) ``` 因此能避開因相同數值造成的無窮迴圈,若成功交換後要記得重新檢查換到該位置的值是否能再次換到正確位置。 ```cpp=13 i-=1 ``` 交換部分寫成函式形式以增加可讀性。 ```cpp=1 void swap(int *a, int *b) { int t=*a; *a=*b; *b=t; return; } ``` 完整程式碼: ```cpp=1 void swap(int *a, int *b) { int t=*a; *a=*b; *b=t; return; } int firstMissingPositive(int* nums, int numsSize) { int i,pos; for(i=0;i<numsSize;i++) { if(nums[i]<numsSize && nums[i]>0) { //only swap vaild value int pos = nums[i] - 1; if(pos != i && nums[i] != nums[pos]) { //avoid same value swap(&nums[i],&nums[pos]); i-=1; //recheck swaped value } } } for(i=0;i<numsSize;i++) { if(nums[i] != i+1) { return i+1; } } return numsSize+1; } ``` 結果: ![](https://i.imgur.com/gk6OyFH.png) :::success 延伸問題: 2. 改寫上述程式碼,提供等價功能的不同實作,儘量降低 branch 的數量; 4. 研讀 2017 年修課學生報告,理解 clz 的實作方式,並舉出 Exponential Golomb coding 的案例說明,需要有對應的 C 原始程式碼和測試報告; >x-compressor 是個以 Golomb-Rice coding 為基礎的資料壓縮器,實作中也用到 clz/ctz 指令,可參見 Selecting the Golomb Parameter in Rice Coding。 ::: // TODO --- ### 測驗 `3` 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 程式的實作: ``` unsigned popcnt_naive(unsigned v) { unsigned n = 0; while (v) v &= (v - 1), n = -(~n); return n; } ``` 函式 `popcnt_naive()` 利用不斷清除 LSB 直到輸入的數值 `v` 為 0。 v &= (v - 1) 的作用是將 LSB 設為 0 (**reset LSB**) 舉例來說,假設輸入數值為 20: ```cpp 0001 0100 # 20 ; LSB in bit position 2 0001 0011 # 20 - 1 0001 0000 # 20 & (20 - 1) ``` >類似的操作還有 x & -x,將 x 的 LSB 取出 (**isolate LSB**) n = -(~n) 等同 n++,因為在二補數系統中, $−n= ∼n+1$ $−(∼n)=n+1$ 因此 `popcnt_naive()` 的執行時間取決於輸入數值中 1 (set bit) 的個數。可改寫為以下常數時間複雜度的實作: ```cpp=1 unsigned popcnt_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; } ``` 對於一個 32 bit 的無號整數,popcount 可以寫成以下數學式: $popcnt(x)=x−⌊\frac{x}{2}⌋−⌊\tfrac{x}{4}⌋−...−⌊\tfrac{x}{2^{31}}⌋$ 假設 $x=b_{31}...b_3b_2b_1b_0$,先看看 $x[3:0]$ 4 個位元,用以上公式可以計算得: $(2^3b_3 + 2^2b_2+2^1b_1+2^0b_0)−(2^2b_3+2^1b_2+2^0b_1)−(2^1b_3+2^0b_2)−2^0b_3$ >$⌊\tfrac{x}{2}⌋$ 相當於 C 表達式中 `x >> 1` 稍微改寫可得到: $(2^3−2^2−2^1−2^0)b_3+(2^2−2^1−2^0)b_2+(2^1−2^0)b_1+2^0b_0$ 因此 popcnt 的一般式可改寫為: $popcnt(x)=\displaystyle\sum_{n=0}^{31}(2^n−\sum_{i=0}^{n-1}2^i)b^n=\sum_{n=0}^{31}b^n$ 因為 $2^n−\displaystyle\sum_{i=0}^{n-1}2^i=1$,只要對應的 $b^n$ 為 1,這個 bit 就會在 popcnt 的總和中加一,剛好對應 `popcnt_naive()`,因此映證上述數學式確實可計算出 population count。 且一個 32 位元無號整數最多有 32 個 1 (set bit),剛好可用一個 byte 表示,所以可分成幾個區塊平行計算,最後再全部加總到一個 byte 中,進行避免檢查 32 次。 `popcnt_branchless` 實作一開始以**每 4 個位元 (nibble) 為一個單位**計算 1 的個數,利用最初的公式計算 $x−⌊\dfrac{x}{2}⌋−⌊\dfrac{x}{4}⌋−⌊\dfrac{x}{8}⌋$ 關鍵的程式碼,逐行解釋: ``` n = (v >> 1) & 0x77777777; v -= n; n = (n >> 1) & 0x77777777; v -= n; n = (n >> 1) & 0x77777777; v -= n; ``` 1. `n = (v >> 1) & 0x77777777` : 將輸入數值 v 除以 2,得到 $⌊\dfrac{v}{2}⌋$ ``` b_31 b_30 b_29 b_28 ... b7 b6 b5 b4 b3 b2 b1 b0 // v 0 b_31 b_30 b_29 ... 0 b7 b6 b4 0 b3 b2 b1 // (v >> 1) & 0x77777777 ``` 2. `v -= n` : 計算結果相當於 $v−⌊\dfrac{v}{2}⌋$ 3. `n = (n >> 1) & 0x77777777` : 再對 n 除以 2,得到 $⌊\dfrac{v}{4}⌋$ 4. `v -= n` : 計算出 $v−⌊\dfrac{v}{2}⌋−⌊\dfrac{v}{4}⌋$ 5. 和 6. 重複同樣動作 6. 最後這段結束後計算出 $v−⌊\dfrac{v}{2}⌋−⌊\dfrac{v}{4}⌋−⌊\dfrac{v}{8}⌋$ ,得到每 4 個位元為一個單位中 set bit 的個數 **`v = (v + (v >> 4)) & 0x0F0F0F0F`** : 將每 4 個位元中 set bit 的個數加到 byte 中: 1. 假設 $B_n$ 代表第 n 個 nibble (4 位元) 中的數值 ``` B7 B6 B5 B4 B3 B2 B1 B0 // v 0 B7 B6 B5 B4 B3 B2 B1 // (v >> 4) ``` 2. 加總可得到: ``` // (v + (v >> 4)) B7 (B7+B6) (B6+B5) (B5+B4) (B4+B3) (B3+B2) (B2+B1) (B1+B0) ``` 3. 最後使用 0x0F0F0F0F 做 mask 可得 ``` // (v + (v >> 4)) & 0x0F0F0F0F 0 (B7+B6) 0 (B5+B4) 0 (B3+B2) 0 (B1+B0) ``` **`v *= 0x01010101`** : 在最後一道敘述中,將 `v` 乘上 `0x01010101`。寫成直式如下 ``` 0 A6 0 A4 0 A2 0 A0 x 0 1 0 1 0 1 0 1 --------------------------------------------------- 0 A6 0 A4 0 A2 0 A0 0 A6 0 A4 0 A2 0 A0 0 0 A6 0 A4 0 A2 0 A0 0 0 A6 0 A4 0 A2 0 A0 0 --------------------------------------------------- ↑_______________________A6+A4+A2+A0 ``` 我們可發現期望輸出就在原本 $A_6$ 的位置 ($2^7$),因此將 `v` 右移 24 bits 即為所求,剩下的位數會 overflow ,右移後不影響結果。 - 假設 $A=B_7+B_6$, $B=B_5+B_4$, $C=B_3+B_2$, $D=B_1+B_0$, 根據分配律可得: ``` v * 0x01010101 = (A + B + C + D) (B + C + D) (C + D) (D) |<-- 1 byte -->|<-- 1 byte -->|<-- 1 byte -->|<-- 1 byte -->| ``` **`return v >> 24`** : 最後得到的結果會放在 Most significant byte 中,因此向右位移 24 位元,即為所求的 popcount 值。 GCC 提供對應的內建函式: >`__builtin_popcount(x)`: 計算 x 的二進位表示中,總共有幾個 1 使用示範: ```cpp=1 int x = 5328; // 00000000000000000001010011010000 printf("%d\n", __builtin_popcount(x)); // 5 ``` 針對 LeetCode 1342. Number of Steps to Reduce a Number to Zero,我們可運用 gcc 內建函式實作,假設 `int` 寬度為 32 位元,程式碼列表如下: ```cpp=1 int numberOfSteps (int num) { return num ? AAA + 31 - BBB : 0; } ``` 請補完程式碼。 作答區 AAA = (a) __builtin_popcount(num) BBB = (b) __builtin_clz(num) #### 解題想法 numberOfSteps的結果可分為: 需右移 n 次與 反向位元 (1變0) m 次,也就是共 n + m 次,因為選項的 `__builtin_popcount(num)` 與 `__builtin_clz(num)` 的回傳值不為負數。 因為推測AAA為增加次數,BBB為減少次數,31表示把最左側的 1 右移是最右側需要的次數。 若輸入值的二進制最右側有 x 個 0 ,表示可以少 x 次右移,所以 ==BBB = (b) `__builtin_clz(num)`==。 若輸入值的二進制共有 y 個 1 ,表示需反向 y 次,所以 ==AAA = (a) `__builtin_popcount(num)`==。 :::success 延伸問題: 1. 解釋上述程式運作原理; ::: ```cpp=3 num ? __builtin_popcount(num) + 31 - __builtin_clz(num) : 0; ``` 若輸入的值為 0 則直接回傳 0 不做處理;若不為 0 則計算 numberOfSteps 結果,算法為 `__builtin_popcount(num) + 31 - __builtin_clz(num)` ,計算輸入值的二進制有幾個 bit 為 1 與最多須向右位移幾次,兩者相加即為符合要求的回傳值。 :::success 2. 改寫上述程式碼,提供等價功能的不同實作並解說; >提示: Bit Twiddling Hacks 中提及許多 bitwise operation 的使用,如 bit inverse、abs 等等,是極好的參考材料。 3. 避免用到編譯器的擴充功能,只用 C99 特徵及 bitwise 運算改寫 LeetCode 1342. Number of Steps to Reduce a Number to Zero,實作出 branchless 解法,儘量縮減程式碼行數和分支數量; ::: //TODO --- ### 測驗 `4` 考慮以下 64-bit GCD (greatest common divisor, 最大公因數) 求值函式: ```cpp=1 #include <stdint.h> uint64_t gcd64(uint64_t u, uint64_t v) { if (!u || !v) return u | v; while (v) { uint64_t t = v; v = u % v; u = t; } return u; } ``` 改寫為以下等價實作: ```cpp=1 #include <stdint.h> uint64_t gcd64(uint64_t u, uint64_t v) { if (!u || !v) return u | v; int shift; for (shift = 0; !((u | v) & 1); shift++) { u /= 2, v /= 2; } while (!(u & 1)) u /= 2; do { while (!(v & 1)) v /= 2; if (u < v) { v -= u; } else { uint64_t t = u - v; u = v; v = t; } } while (XXX); return YYY; } ``` 請補完程式碼。 作答區 XXX = (b) v YYY = (e) u << shift #### 解題想法 :::success 延伸問題: 1.解釋上述程式運作原理; ::: ```cpp=3 if (!u || !v) return u | v; ``` 例外處理,若 u 或 v 有一個為 0 ,則直接回傳另一個數。 ```cpp=5 for (shift = 0; !((u | v) & 1); shift++) { u /= 2, v /= 2; } ``` 若 u 跟 v 均為偶數 (`!((u | v) & 1)`),則將兩者同時除以 2 直到其中一個變成奇數,並記錄除的次數於 shift。 接下來進行輾轉相除法: ```cpp=8 while (!(u & 1)) u /= 2; ``` 若 u 為偶數,則把 u 除 2 直到不能再除。 ```cpp=10 do { while (!(v & 1)) v /= 2; if (u < v) { v -= u; } else { uint64_t t = u - v; u = v; v = t; } } while (XXX); ``` ```cpp=11 while (!(v & 1)) v /= 2; ``` 若 v 為偶數,則把 v 除 2 直到不能再除。 ```cpp=13 if (u < v) { v -= u; } else { uint64_t t = u - v; u = v; v = t; } ``` 先確認 u 與 v 誰是大數誰是小數(比大小),一次迴圈後的結果變成`v = 大數 - 小數`,`u = 小數` ```cpp=20 } while (XXX); ``` 依照輾轉相除法定義,這邊應該要判斷 `大數 - 小數` 是否為零,因此 ==XXX = (b) v==。 ```cpp=21 return YYY; ``` 最後回傳的結果為兩數的最大公因數,也就是結束迴圈後的 `小數`乘以 $2^{shift}$ 等價於選項中的 ==YYY = (e) u << shift== :::success 延伸問題: 2.在 x86_64 上透過 `__builtin_ctz` 改寫 GCD,分析對效能的提升; ::: --- ### 測驗 `5` 在影像處理中,bit array (也稱 `bitset`) 廣泛使用,考慮以下程式碼: ```cpp=1 #include <stddef.h> size_t naive(uint64_t *bitmap, size_t bitmapsize, uint32_t *out) { size_t pos = 0; for (size_t k = 0; k < bitmapsize; ++k) { uint64_t bitset = bitmap[k]; size_t p = k * 64; for (int i = 0; i < 64; i++) { if ((bitset >> i) & 0x1) out[pos++] = p + i; } } return pos; } ``` 可用 ctz 改寫為效率更高且等價的程式碼: ```cpp=1 #include <stddef.h> size_t improved(uint64_t *bitmap, size_t bitmapsize, uint32_t *out) { size_t pos = 0; uint64_t bitset; for (size_t k = 0; k < bitmapsize; ++k) { bitset = bitmap[k]; while (bitset != 0) { uint64_t t = KKK; int r = __builtin_ctzll(bitset); out[pos++] = k * 64 + r; bitset ^= t; } } return pos; } ``` 請補完程式碼。 作答區 KKK = (e) bitset & -bitset #### 解題想法 :::success 1. 解釋上述程式運作原理 ::: 從題目要求與ctz版程式碼可得知ctz版的8~13行迴圈的功能如下: `uint64_t t = KKK;` t 為題目要求的mask。 `int r = __builtin_ctzll(bitset);` 找出該二進制數字中,最右邊的 1 的位置。 `out[pos++] = k * 64 + r;` 將該位置紀錄於out陣列中。 `bitset ^= t; ` 對 bitset 做 XOR mask。 由上述行為推測,此mask作用應該是要對某個位置的 bit 做反向,因為第 10 行的 `__builtin_ctzll` 是找出該二進制中最右邊(LSB)的 1,因此此mask應設計為找出最右邊(LSB)的 1,因此只有 ==(e) bitset & -bitset== 符合要求。 :::success 延伸問題: 1. 舉出這樣的程式碼用在哪些真實案例中; ::: //TODO :::success 2. 設計實驗,檢驗 clz 改寫的程式碼相較原本的實作有多少改進?應考慮到不同的 bitmap density; 3. 思考進一步的改進空間; ::: //TODO