---
tags: linux2022
---

# 2022q1 Homework2 (quiz2)
contributed by < `BensonYang1999` >

> [第 2 周測驗題](https://hackmd.io/@sysprog/linux2022-quiz2)

## 測驗 1：兩數取平均
### 解題
由於兩個無號整數相加時有可能會發生溢位，因此考慮使用以下不同方法
1. 已知兩數大小，可直接使用利用以下方式

uint32_t average(uint32_t low, uint32_t high)
{
    return low + (high - low) / 2;
}

2. 利用位移運算子可避免使用方法1中的除法，且不需預先知道兩數的大小

uint32_t average(uint32_t a, uint32_t b)
{
    return (a >> 1) + (b >> 1) + (a & b & 1);
}

其中的 `a & b & 1` 是用來補償小數進位的部份，原理為判斷兩數是否都是奇數
以下面的例子為例
`a=5, b=7 => a>>1=2, b>>1=3 => added=5`
因此需要補加上1，才會得到最終的數值6

3. 利用 bit-wise operator 實作加法器

在[你所不知道的 C 語言：數值系統篇](https://hackmd.io/@sysprog/c-numerics)有提到利用 bit-wise operator 實作加法器 `x + y = x ^ y + ( x & y ) << 1`
由加法式可以進階推導出 `(x + y) / 2 = (x & y) + ((x ^ y) >> 1`
因此可以利用上式實作兩數平均

uint32_t average(uint32_t a, uint32_t b)
{
    return (a & b) + ((a ^ b) >> 1);
}

:::info
延伸問題:

- [x] 解釋上述程式碼運作的原理
- [x] 比較上述實作在編譯器最佳化開啟的狀況，對應的組合語言輸出，並嘗試解讀 (可研讀 CS:APP 第 3 章)
- [ ] 研讀 Linux 核心原始程式碼 include/linux/average.h，探討其 Exponentially weighted moving average (EWMA) 實作，以及在 Linux 核心的應用
> 移動平均（Moving average)，又稱滾動平均值、滑動平均，在統計學中是種藉由建立整個資料集合中不同子集的一系列平均數，來分析資料點的計算方法。
移動平均通常與時間序列資料一起使用，以消除短期波動，突出長期趨勢或周期。短期和長期之間的閾值取決於應用，移動平均的參數將相應地設置。例如，它通常用於對財務數據進行技術分析，如股票價格、收益率或交易量。它也用於經濟學中研究國內生產總值、就業或其他宏觀經濟時間序列。
:::

### 編譯器最佳化後的組合語言比較

> 利用 gcc 中的 `-S` 指令可以輸出組合語言，額外加上 `-O3` 可確保編譯啟開啟最高等級的最佳化

方法1

#include <stdio.h>
#include <stdint.h>
uint32_t average(uint32_t low, uint32_t high)
{
    return low + (high - low) / 2;
}
int main()
{
    printf("%u\n", average(5, 7));
}

方法2

#include <stdio.h>
#include <stdint.h>
uint32_t average(uint32_t a, uint32_t b)
{
    return (a >> 1) + (b >> 1) + (a & b & 1);
}
int main()
{
    printf("%u\n", average(5, 7));
}

方法3

#include <stdio.h>
#include <stdint.h>
uint32_t average(uint32_t a, uint32_t b)
{
    return (a & b) + ((a ^ b) >> 1);
}
int main()
{
    printf("%u\n", average(5, 7));
}

利用diff比較方法1與方法2的差別

$ diff test1.s test2.s
10,11c10,12
< 	movl	%esi, %eax
< 	subl	%edi, %eax
---
> 	movl	%edi, %eax
> 	movl	%esi, %edx
> 	andl	%esi, %edi
12a14,16
> 	shrl	%edx
> 	andl	$1, %edi
> 	addl	%edx, %eax

* 雖然方法1的組合語言較精簡，但會用到 subl 除法運算，因此可能會花很多 clock cycles 去執行運算
* 這兩個方法的比較並不恰當，因方法1需額外加入比較兩個數字大小比較的判斷式，會需要額外的指令

利用diff比較方法2與方法3的差別

$ diff test2.s test3.s
11d10
< 	movl	%esi, %edx
12a12
> 	xorl	%esi, %eax
14,16d13
< 	shrl	%edx
< 	andl	$1, %edi
< 	addl	%edx, %eax

* 方法3使用 xor 取代移動、位移、相加、 and 指令，因此判斷可以節省執行時間
:::warning
使用 `diff -u` 命令，更好觀察。另外，方法 `1` 應該補充數值範圍的檢查程式碼。
:notes: jserv
:::

---

## 測驗 2：bitwise max
### 解題

#include <stdint.h>
uint32_t max(uint32_t a, uint32_t b)
{
    return a ^ ((a ^ b) & -(a < b));
}

考慮 `a < b` 及 `a >= b` 兩種情況
1. 當 `a < b`
& 右邊會等於 0xFFFFFFFF ，即可以省略 `& -(a < b)` 這部份
剩下的部份 `a ^ (a ^ b)` 可進一部簡化為 `b`
因此最後的算式 => `return b;`

2. 當 `a >= b`
& 右邊會等於 0 ，代表 `((a ^ b) & 0)` 會等於0
因此最後的算式 => `return a;`

:::info
延伸問題:
- [x] 解釋上述程式碼運作的原理
- [x] 針對 32 位元無號/有號整數，撰寫同樣 branchless 的實作
- [ ] Linux 核心也有若干 branchless / branch-free 的實作，例如 lib/sort.c:

/*
 * Logically, we're doing "if (i & lsbit) i -= size;", but since the
 * branch is unpredictable, it's done with a bit of clever branch-free
 * code instead.
 */
__attribute_const__ __always_inline
static size_t parent(size_t i, unsigned int lsbit, size_t size)
{
    i -= size;
    i -= size & -(i & lsbit);
    return i / 2;
}

請在 Linux 核心原始程式碼中，找出更多類似的實作手法。請善用 git log 檢索。
:::

### min/max branchless 實作
討論無號整數的 `min`，可以透過簡單修改 `max` 的程式碼達成

#include <stdint.h>
uint32_t min(uint32_t a, uint32_t b)
{
    return a ^ ((a ^ b) & -(a > b));
}

有號整數的部份，由於此運算方式可以適用於有號數，因此只需修改變數型態即可

#include <stdint.h>
int32_t max(int32_t a, int32_t b)
{
    return a ^ ((a ^ b) & -(a < b));
}
int32_t min(int32_t a, int32_t b)
{
    return a ^ ((a ^ b) & -(a > b));
}

---

## 測驗 3：GCD
### 解題
參考本題的註解

#include <stdint.h>
uint64_t gcd64(uint64_t u, uint64_t v)
{
    if (!u || !v) return u | v; // no.1
    int shift; // 紀錄共除2的幾次方
    for (shift = 0; !((u | v) & 1); shift++) {
        u /= 2, v /= 2;
    } // no.3
    while (!(u & 1))
        u /= 2; // no.4
    do {
        while (!(v & 1))
            v /= 2; // no.5
        if (u < v) {
            v -= u;
        } else {
            uint64_t t = u - v;
            u = v;
            v = t;
        }
    } while (COND);
    return RET;
}

從第14行開始沒有在註解中說明，查閱網路資料後可得知其為 [Euclid's algorithm](https://en.wikipedia.org/wiki/Greatest_common_divisor) ，其結束條件為兩數相等，因此判斷 `COND` 的空格應該填入 `v` ，也就是被減數不等於0時
而最後RET的地方需要把一開始同除的2次方給補回去，因此應該是 `u << shift`

:::info
延伸問題:

- [x] 解釋上述程式運作原理;
- [x] 在 x86_64 上透過 __builtin_ctz 改寫 GCD，分析對效能的提升;
- [ ] Linux 核心中也內建 GCD (而且還不只一種實作)，例如 lib/math/gcd.c，請解釋其實作手法和探討在核心內的應用場景。
:::

### __builtin_ctz 改寫 GCD
查閱 [GNU 說明書](https://gcc.gnu.org/onlinedocs/gcc/Other-Builtins.html)，得知本函數可以取得變數中從最低位數有幾的0，也就可以判斷其因數包含2的多少次方，因此可以用來加速原本使用 while 判斷的地方，也就是說程式可以修改為

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;
    }
    u >>= __builtin_ctz(u);
    do {
        v >>= __builtin_ctz(v);
        if (u < v) {
            v -= u;
        } else {
            uint64_t t = u - v;
            u = v;
            v = t;
        }
    } while (v);
    return u << shift;
}

仔細觀察後發現，發現第4~8行其實可以整並，改寫為

int shift = __builtin_ctz(u) < __builtin_ctz(v) ? __builtin_ctz(u) : __builtin_ctz(v);
u >>= __builtin_ctz(u);
v >>= __builtin_ctz(v);

完整程式碼

uint64_t gcd64(uint64_t u, uint64_t v)
{
    if (!u || !v) return u | v;
    int shift = __builtin_ctz(u) < __builtin_ctz(v) ? __builtin_ctz(u) : __builtin_ctz(v);
    u >>= __builtin_ctz(u);
    v >>= __builtin_ctz(v);
    do {
        v >>= __builtin_ctz(v);
        if (u < v) {
            v -= u;
        } else {
            uint64_t t = u - v;
            u = v;
            v = t;
        }
    } while (v);
    return u << shift;
}

利用傳入大數的方式，讓執行時間變長，可以明顯比較效能差異
我利用 linux 內建指令 `time` 進行計時，並帶入 uint64 的最大值 18446744073709551615 與其 -1 ，結果如下

original version:
$ time ./a.out
real	0m0.003s
user	0m0.002s
sys	0m0.000s

__builtin_ctzll version:
$ time ./a.out
real	0m0.003s
user	0m0.002s
sys	0m0.000s

兩者時間都太短了，因此改變方式，使用計算相同數字多次計時，結果如下

original version:
$ time ./a.out
real	0m1.372s
user	0m1.372s
sys	0m0.000s

__builtin_ctz version:
$ time ./a.out
real	0m1.100s
user	0m1.100s
sys	0m0.000s

大概加快 `25%`

---

## 測驗 4：bit array 使用
### 解題
依照本題的提示，考慮使用 `-bitwise` ，由於再 C 語言裡面負號的計算方式是將所有bits反轉後再加 1 ，因此可以利用此特性，與 `bitwise` 做 & 運算，及可以達成找到最低的 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 = bitset & -bitset;
            int r = __builtin_ctzll(bitset);
            out[pos++] = k * 64 + r;
            bitset ^= t;                           
        }
    }
    return pos;
}

:::info
延伸問題:

- [ ] 解釋上述程式運作原理，並舉出這樣的程式碼用在哪些真實案例中;
- [ ] 設計實驗，檢驗 ctz/clz 改寫的程式碼相較原本的實作有多少改進？應考慮到不同的 bitmap density;
- [ ] 思考進一步的改進空間;
- [ ] 閱讀 Data Structures in the Linux Kernel 並舉出 Linux 核心使用 bitmap 的案例;
:::