# 2017q3 Homework1 (clz) ###### tags: `sysprogram` contributed by <`st9007a`> ## [Github](https://github.com/st9007a/clz-tests) ### Reviewed by `amikai` - [726ac7e](https://github.com/st9007a/clz-tests/commit/726ac7e50c35684aa63144192b5ebff22b775de8) 此 commit 做的是對整個架構的大修改, 建議加上修改的動機, 讓其他協作者 (如果很多人開發) 了解為什麼要有這些修改 - 8 bit, 16bit, 32bit 三張圖, 沒有對三張圖的比較做詳細的解釋, 講述了數據不穩定, 但沒有說緣由 以上為小弟簡短的建議, 如有冒犯請見諒 ## 開發環境 ``` Linux 4.4.0-92-generic Ubuntu 16.04.1 LTS L1d cache: 32K L1i cache: 32K L2 cache: 256K L3 cache: 6144K ``` ## 針對 32-bits 數值對各版本之 clz 比較效能 利用 fork 過來的程式，可以快速的實測各版本的 clz 對 32-bits 整數的效能比較 ![](https://i.imgur.com/gGTUDwm.png) 從上圖可以知道 byte-shift 跟 binary search 所使用的 cycle 數平均是比較少的，然後 cycle 數 400 到 500 之間也有點分布，試著隨便測其他區間： ![](https://i.imgur.com/isOLevl.png) ![](https://i.imgur.com/Zlpgy4m.png) 分布狀況看起來都差不多，不過一部分的數據會跳到 400 至 500 cycles 之間，這是由於作業系統multi task，如果有其他 task 要優先值行的話，作業系統會先中斷現在執行的 task，等完成其他 task 在跳回來繼續執行，所以花費的 cycle 數會比較高 ## 解析原始碼 執行完效能比較後，來看看 `clz.h` 裡面是怎麼實作個版本的 clz 演算法 ### recursive ```c static const int mask[] = {0, 8, 12, 14}; static const int magic[] = {2, 1, 0, 0}; //static inline __attribute((always_inline)) unsigned clz2(uint32_t x,int c) { // *INDENT-OFF* if (!x && !c) return 32; uint32_t upper = (x >> (16 >> c)); uint32_t lower = (x & (0xFFFF>>mask[c])); if (c == 3) return upper ? magic[upper] : 2 + magic[lower]; return upper ? clz2(upper, c + 1) : (16 >> (c)) + clz2(lower, c + 1); // *INDENT-ON* } ``` `c` 這個變數紀錄每次遞迴的層數，每次遞迴會將 `x` 分成前半部( upper )跟後半部( lower )，如果前半部不為零，則回傳前半部的 clz，若前半部為零，則回傳前半部的 bit 數量加上後半部的 clz，而倒數第二行的停止條件則是第四次遞迴時，前半部跟後半部只剩 2 個 bits，所以直接回傳查表對應的 clz ### iteration ```c static inline __attribute((always_inline)) unsigned clz(uint32_t x) { // *INDENT-OFF* int n = 32, c = 16; do { uint32_t y = x >> c; if (y) { n -= c; x = y; } c >>= 1; } while (c); return (n - x); // *INDENT-ON* } ``` 迴圈裡面的 `y` 用來記錄 `x` 的前半部，`n` 用來記錄 clz，如果前半部不為零就用前半部去做下一個迴圈，同時 `n` 減去後半部的 bit 的數量( 因為至少後半部的 bit 不會算進 clz )，迴圈做到最後當 `c` 變成 0 跳出迴圈時 `x` 還剩下一個 bit， `n` 直接減去最後一個 bit 的值就能得到 clz　 ### byte shift ```c static inline __attribute((always_inline)) unsigned clz(uint32_t x) { // *INDENT-OFF* if (x == 0) return 32; int n = 1; if ((x >> 16) == 0) { n += 16; x <<= 16; } if ((x >> 24) == 0) { n += 8; x <<= 8; } if ((x >> 28) == 0) { n += 4; x <<= 4; } if ((x >> 30) == 0) { n += 2; x <<= 2; } n = n - (x >> 31); return n; // *INDENT-ON* } ``` `n` 用來記錄 clz，利用 if 來分別確認 `x` 的前 16 , 8 , 4 , 2 個 bits 是否為零，如果為零則將他們加入 `n`，並且左移 `x` 確保下個 if 能檢查到正確的位置，如果不為零則直接做下一個 if 檢查，最後 `x` 剩下 最前面的 1 bit 沒有被檢查到，做法跟 iteration 一樣，直接減去 bit 的值就是 clz ### binary search ```c static inline __attribute((always_inline)) unsigned clz(uint32_t x) { // *INDENT-OFF* if (x == 0) return 32; int n = 0; if (x <= 0x0000FFFF) { n += 16; x <<= 16; } if (x <= 0x00FFFFFF) { n += 8; x <<= 8; } if (x <= 0x0FFFFFFF) { n += 4; x <<= 4; } if (x <= 0x3FFFFFFF) { n += 2; x <<= 2; } if (x <= 0x7FFFFFFF) { n += 1; x <<= 1; } return n; // *INDENT-ON* } ``` 原理跟 byte shift 一樣，只是改成用比大小的方式來確認前 16 , 8 , 4 , 2 個 bits 是否為零，最後一個 if 取代 byte shift 版本的 `n = n - (x >> 31);` 檢查最前面的 bit ### harely ```c static inline __attribute((always_inline)) unsigned clz(uint32_t x) { // *INDENT-OFF* static uint8_t const Table[] = { 32, 31, 0, 16, 0, 30, 3, 0, 15, 0, 0, 0, 29, 10, 2, 0, 0, 0, 12, 14, 21, 0, 19, 0, 0, 28, 0, 25, 0, 9, 1, 0, 17, 0, 4, 0, 0, 0, 11, 0, 13, 22, 20, 0, 26, 0, 0, 18, 5, 0, 0, 23, 0, 27, 0, 6, 0, 24, 7, 0, 8, 0, 0, 0 }; // *INDENT-ON* /* Propagate leftmost 1-bit to the right */ x = x | (x >> 1); x = x | (x >> 2); x = x | (x >> 4); x = x | (x >> 8); x = x | (x >> 16); /* x = x * 0x6EB14F9 */ x = (x << 3) - x; /* Multiply by 7. */ x = (x << 8) - x; /* Multiply by 255. */ x = (x << 8) - x; /* Again. */ x = (x << 8) - x; /* Again. */ return Table[(x >> 26)]; } ``` 從演算法名稱可以預期到` 0x6EB14F9` 是 De Bruijn 數列，一開始前五個 or 運算將 `x` 的最高有效位往右的 bit 都變成 1，接著可以從作業共筆知道存在著 De Bruijn 數列 M 使得 $(a \times M) >> n$ 出來的結果彼此不衝突，其中 $a = 0...000111...1_2$，在這裡 $M = 6EB14F9_{16} , n = 26$，接著 hash 出來的結果直接去查表就能得到 clz ## 修改程式碼 由於作業要求有提到 C11 的 `_Generic`，研究了一下把程式改成支援 8-bit, 16-bit, 32-bit, 64-bit unsigned 整數，首先將標頭檔套用 `_Generic` ```c #ifndef CLZ_H #define CLZ_H #include <stdint.h> #if defined(recursive) #define clz(x) \ _Generic((x), \ uint8_t: clz_u8, \ uint16_t: clz_u16, \ uint32_t: clz_u32, \ uint64_t: clz_u64 \ )(x, 0) unsigned clz_u8(uint8_t x, int c); unsigned clz_u16(uint16_t x, int c); unsigned clz_u32(uint32_t x, int c); unsigned clz_u64(uint64_t x, int c); #else #define clz(x) \ _Generic((x), \ uint8_t: clz_u8, \ uint16_t: clz_u16, \ uint32_t: clz_u32, \ uint64_t: clz_u64 \ )(x) unsigned clz_u8(uint8_t x); unsigned clz_u16(uint16_t x); unsigned clz_u32(uint32_t x); unsigned clz_u64(uint64_t x); #endif #endif ``` 接著分別在 `iteration.c`, `recursive.c`, `byte.c` , `binary.c` 根據原本的函式修改出 8-bit, 16-bit, 64-bit的版本，由於全部程式碼太多行，且只是按 32-bit 改出其它版本，所以就不附上程式碼，harely 還沒實作，需先找到對應的 De Bruijn 數列 改好程式碼後，把 Makefile 中的 `-std=gnu99` 改成 `-std=c11`，然後來編譯看看 ``` cc -Wall -std=c11 -g -DDEBUG -O0 -o iteration \ -MMD -MF .iteration.d \ -Diteration main.c main.c: In function ‘get_cycles’: main.c:14:5: error: ‘asm’ undeclared (first use in this function) asm volatile ("CPUID\n\t" ^ main.c:14:5: note: each undeclared identifier is reported only once for each function it appears in main.c:14:9: error: expected ‘;’ before ‘volatile’ asm volatile ("CPUID\n\t" ^ main.c: In function ‘get_cycles_end’: main.c:24:5: error: ‘asm’ undeclared (first use in this function) asm volatile("RDTSCP\n\t" ^ main.c:24:9: error: expected ‘;’ before ‘volatile’ asm volatile("RDTSCP\n\t" ^ Makefile:30: recipe for target 'iteration' failed make: *** [iteration] Error 1 ``` 會發現他不認得 `asm` 這個保留字，查了一下原因，知道了 `asm` 是 C 語言的 extension，要使用 extension 必須把 gcc 的 `-std=c11` 改成 `-std=gnu11`，接著再試著編譯看看 ``` make: Warning: File 'Makefile' has modification time 0.76 s in the future cc -Wall -std=gnu11 -g -DDEBUG -O0 -o iteration \ -MMD -MF .iteration.d \ -Diteration main.c /tmp/ccU6IkJY.o: In function `main': /home/st9007a/project/sysprog/clz-tests/main.c:107: undefined reference to `clz_u32' collect2: error: ld returned 1 exit status Makefile:30: recipe for target 'iteration' failed make: *** [iteration] Error 1 ``` 發現它竟然認不得 `clz_u32`，所以打開 `Makefile` 看看是什麼狀況 ``` %: $(SRCS_common) %.c clz.h $(CC) $(CFLAGS_common) -o $@ \ -MMD -MF .$@.d \ -D$(shell echo $(subst .o,,$@)) $(SRCS_common) ``` 原來實作 clz function 的檔案沒有被一起編譯，所以改成如下 ``` %: $(SRCS_common) %.c clz.h $(CC) $(CFLAGS_common) -o $@ \ -MMD -MF .$@.d \ -D$(shell echo $(subst .o,,$@)) $(SRCS_common) $@.c ``` 在最後面加入 `$@.c`，把實作 clz function 的檔案一起編譯 ``` $ make cc -Wall -std=gnu11 -g -DDEBUG -O0 -o iteration \ -MMD -MF .iteration.d \ -Diteration main.c iteration.c cc -Wall -std=gnu11 -g -DDEBUG -O0 -o binary \ -MMD -MF .binary.d \ -Dbinary main.c binary.c cc -Wall -std=gnu11 -g -DDEBUG -O0 -o byte \ -MMD -MF .byte.d \ -Dbyte main.c byte.c cc -Wall -std=gnu11 -g -DDEBUG -O0 -o recursive \ -MMD -MF .recursive.d \ -Drecursive main.c recursive.c ``` 大功告成！來用原本的測試測試看看如何 ![](https://i.imgur.com/f9uEymk.png) 跟先前的版本差不多，看來效能沒有因為我的修改被改爛，接著來測看看不同的bit width，首先把 `main.c` 做一點小修改 ```c #if defined(u8) uint8_t i; #elif defined(u16) uint16_t i; #elif defined(u32) uint32_t i; #else uint64_t i; #endif ``` 把要測試的數字根據不同 bit width 宣告成不同型態，`Makefile` 也加入對應 flag 來自動化 ``` bit_width ?= 32 from ?= 67100000 to ?= 67110000 CFLAGS_common += -Du$(bit_width) run: $(EXEC) for method in $(EXEC); do\ taskset -c 1 ./$$method $(from) $(to); \ done ``` ### 8-bit 結果 ![](https://i.imgur.com/OBftzix.png) ### 16-bit 結果 ![](https://i.imgur.com/xxWFqIN.png) ### 64-bit 結果 ![](https://i.imgur.com/1PosoUA.png) 可以看到效能是差不多的，8-bit 由於數字只有 256 個，所以只看 100 到 200 的範圍，其數據比較不穩定 ## 應用場合 ### Division Algorithm 參考 [A Fast Hi Precision Fixed Point Divide](http://me.henri.net/fp-div.html) 中，對於硬體實作的除法演算法要做 32 個 iterations 花費 97 個 cycles，程式如下： ```asm @ Entry r0: numerator (lo) must be signed positive @ r2: deniminator (den) must be non-zero and signed negative idiv: lo .req r0; hi .req r1; den .req r2 mov hi, #0 @ hi = 0 adds lo, lo, lo .rept 32 @ repeat 32 times adcs hi, den, hi, lsl #1 subcc hi, hi, den adcs lo, lo, lo .endr mov pc, lr @ return @ Exit r0: quotient (lo) @ r1: remainder (hi) ``` 在被除數與除數有 leading zero 的狀況下，這個演算法前幾個( 被除數與除數共有的leading zero 數量 ) iterations 就等同於只單純的左移 `hi` 跟 `lo` 這兩個暫存器，因此可以預先計算好 clz 然後先移動 `hi` , `lo` 來減少 iterations 的數量，參考資料中提到，在 ARM v5 架構下使用 clz 指令可以減少 9 到 18 個 cycles ### Integer Square Root 整數平方根的算法有很多種，這裡指的是牛頓法( Newton's Method )，牛頓法的概念是數字 $a$ 的平方根可以透過給定一個初始值 $g_0$ 並反覆計算以下公式來找到： $g_{n+1} = (g_n + a / g_n)\times 0.5$ 在計算過程中 $g$ 的值會趨近於 $a$ 的平方根，因此可以知道 $g_0$ 的值會影響運算時間，在 [Hacker's Delight](https://books.google.com.tw/books?id=VicPJYM0I5QC&printsec=frontcover&hl=zh-TW#v=onepage&q&f=false) 281頁中提到，應找大於等於 $a$ 之平方根的 $2 ^ n$ 整數作為 $g_0$，尋找的方法就可以透過 clz 去算出來 附上程式碼： ```clike int isqrt(unsigned x) { int s, g0, g1; if (x <= 1) return x; s = 16 - clz(x - 1) / 2; g0 = 1 << s; g1 = (g0 + x >> s) >> 1; while (g1 < g0) { g0 = g1; g1 = (g0 + x / g0) >> 1; } return g0; } ``` 關鍵在於 `s = 16 - clz(x - 1) / 2` 這個部分， `s` 的目的是算出最接近的數字符合 $2^s >= \sqrt{x}$，可以將程式碼換成數學式如下： $s = 16 - \lfloor\ \frac{1}{2}(31 - \lfloor log_2 (x - 1) \rfloor)\ \rfloor >= 16 - \frac{1}{2}(31 - log_2 (x - 1))$ $= \frac{1}{2} + log_2 \sqrt{x - 1}$ 所以 $2 ^ s >= \sqrt{2} + \sqrt{x - 1} >= \sqrt{2} + \sqrt{x} - 1 > \sqrt{x}$ ## 參考資料 [ktvexe 同學共筆](https://hackmd.io/s/BJBZn6Q6) [baimao8437 同學共筆](https://hackmd.io/s/BkZOVSIcx) [wrong clock cycle measurements with rdtsc](https://stackoverflow.com/questions/19941588/wrong-clock-cycle-measurements-with-rdtsc) [C11 - Generic Selections](http://www.robertgamble.net/2012/01/c11-generic-selections.html) [stackoverflow -- Writing function definition in header files in C++](https://stackoverflow.com/questions/453372/writing-function-definition-in-header-files-in-c) [what “inline \_\_attribute\_\_((always_inline))” means in the function?](https://stackoverflow.com/questions/22767523/what-inline-attribute-always-inline-means-in-the-function) [C 語言的各整版本](https://www.crifan.com/summary_c_language_version_c89_amd1_c99_c11/) [A Fast Hi Precision Fixed Point Divide](http://me.henri.net/fp-div.html) [ARM Keil](http://www.keil.com/support/man/docs/armclang_asm/armclang_asm_pge1425890373464.htm) [carolc0708 同學共筆](https://hackmd.io/s/H1EkbFDT)
Last changed by  
 
Uno Wu
760
Published on HackMD