# 2023q1 Homework3 (fibdrv)
contributed by <`gaberplaysgame`>
## 開發環境

$ gcc --version
gcc (Ubuntu 11.3.0-1ubuntu1~22.04) 11.3.0
$ lscpu
Architecture:            x86_64
  CPU op-mode(s):        32-bit, 64-bit
  Address sizes:         39 bits physical, 48 bits virtual
  Byte Order:            Little Endian
CPU(s):                  12
  On-line CPU(s) list:   0-11
Vendor ID:               GenuineIntel
  Model name:            11th Gen Intel(R) Core(TM) i5-11400 @ 2.60GHz
    CPU family:          6
    Model:               167
    Thread(s) per core:  2
    Core(s) per socket:  6
    Socket(s):           1
    Stepping:            1
    CPU max MHz:         4400.0000
    CPU min MHz:         800.0000
    BogoMIPS:            5184.00

## 研讀費氏數列相關材料並摘錄關鍵手法
### 手法
#### 遞迴法
利用 $fib(n+2)=fib(n)+fib(n+1), fib(0)=0, fib(1)=1$ 的定義最容易得到的方法，雖然只要進行 $n$ 次遞迴，但實際時間複雜度並不等於 $O(n)$，因而費時。

「什麼是費氏數列」系列影片有用到一些方法來改進 python 內的費氏數列運算。
- [PART 2](https://www.youtube.com/watch?v=TA0Dpx0LOeY) 有提到 python 內使用 lru_cache 來存已經計算過的數值的方法。
- [PART 3](https://www.youtube.com/watch?v=WyDn6wiwW74) 使用 tail recursive 的方法，將尾數兩個的數值傳給下一層的函數來避免重複計算。
- [PART 4](https://www.youtube.com/watch?v=iCSNHH45EeI) 使用迭代的方式來計算，程式碼與[費氏數列分析](https://hackmd.io/@sysprog/fibonacci-number)的基本程式碼一致，雖然能夠處理大數運算，但根據分析所言利用 timeit 去測時可以發現執行時間隨數量級而變長。

#### Binet’s Formula （公式解）
費氏數列（以下用 $fib(n)$ 代替）有以下公式解：$fib(n) = \frac{1}{\sqrt{5}}(\frac{1+\sqrt{5}}{2})^n-\frac{1}{\sqrt{5}}(\frac{1-\sqrt{5}}{2})^n$

雖然可以將遞迴呼叫次數降到 $O(1)$，但是由於浮點數精度的問題，在進行計算 $\sqrt5$ 時會導致誤差放大。引用費氏數列分析的文獻，可得知公式解其實沒有比遞迴法來的快：
>$√5$ 可以用牛頓法來算，假設每個 $n$ bit 除法需要 $n^2$ 個 bit 運算，但每一輪只需要留下前面 $Nn+logn$ 個 bit. 整體上說來需要的 bit 運算數中最大宗的大約是 $logn×(Nn+logn)2$, 再加上一些較不顯著的運算 —— **其實沒有明顯地比反覆加法快**。

#### Vorobev 等式
當 $n>=1$ 時，存在特性 $f(n+k)=f(n−1)×f(k)+f(n)×f(k+1)$，若將 $n=k$ 代入並經過一些換算可得 $f(2n)=2×f(n+1)×f(n)−(f(n))^2$，遞迴次數為 $O(logn)$，但執行時間則並非上值。

[Fast Doubling](https://chunminchang.github.io/blog/post/calculating-fibonacci-numbers-by-fast-doubling) 原則上利用這個等式求值。
- $f(2n)=2×f(n+1)×f(n)−(f(n))^2＝f(n)[2f(n+1)-f(n)]$
- $f(2n+1)=f(n+1)^2+f(n)^2$

#### 配合 GNU 函式庫
GNU 的 [GMP](https://gmplib.org/) 與 [MPFR](https://www.mpfr.org/) 函式庫可以將記憶體拿來當作快取，以進行接近於無限位的小數運算。經過運算對比可以發現再大數上公式解不會比 Vorobev 等式還來的快

### [clz / ctz](https://en.wikipedia.org/wiki/Find_first_set) 一類的指令對 Fibonacci 數運算的幫助

參考作業描述中關於費波那契數列以 Fast Doubling 實作的虛擬碼，新建測試檔案 fibonacci.c，內包含 `fib()` 與 `fast_fib()`：

uint64_t fib(uint64_t n)
{
    if (!n)
        return 0;
    else if (n <= 2)
        return 1;
    return fib(n-1) + fib(n-2);
}

uint64_t fast_fib(uint64_t n)
{
    if (!n)
        return 0;
    else if (n <= 2)
        return 1;
    uint64_t a = 0, b = 1;
    
    for (int i = 63 - __builtin_clzll(n), j = 1UL << i; 
         i >= 0; 
         --i, j >>= 1
    ) {
        uint64_t t1 = a * (2 * b - a);  //fib(2a)
        uint64_t t2 = b * b + a * a;    //fib(2a+1)
        
        if (n & j) {
            a = t2;
            b = t1 + t2;
        } else {
            a = t1;
            b = t2;
        }
    }

    return a;
}

`fast_fib()` 支援無號 64 位元的費波那契數運算，比使用遞迴法的 `fib()` 速度更快且可以在短時間內達成 $n > 50$ 的運算。（尚未支持超過 $2^{64}$ 的大數運算）

這邊利用 quiz2 提到的 `__builtin_clz` 變體函式來去快速取得有效位元個數，並利用變數 `j` 來對 `n` 做 bitmask 來判斷當前的位元為 1 或 0。

if-else 區塊透過進一步的 bitwise 可以改寫為：

int64_t mask = !(n & j) - 1; // 結果會是 0 (0x00000000) 或 -1 (0xFFFFFFFF)
a = (t2 & mask) + (t1 & ~mask);
b = t2 + (t1 & mask);

## 改進大數運算
參考 [Small/Short String Optimization](https://github.com/AdrianHuang/fibdrv/tree/big_number) 進行改寫，利用字串的方式來儲存數值。若使用長度為 128 的字串可使能正確計算的數值來到 $fib(500)$ 。

void bn_add (char *n1, char *n2, char *n) {
    int it, digit = 0, carry = 0;
    int size_n1 = strlen(n1), size_n2 = strlen(n2);
    
    // 這邊預設 n-1 的長度 >= n-2，若沒有預設需要做額外處理
    for (it = 0; it < size_n2; it++) { 
        digit = (n1[it] - '0') + (n2[it] - '0') + carry;
        n[it] = digit % 10 + '0';
        carry = digit / 10;
    }
    for (it = size_n2; it < size_n1; it++) {
        digit = (n1[it] - '0') + carry;
        n[it] = digit % 10 + '0';
        carry = digit / 10;
    }

    if (carry)
        n[it] = carry + '0';
    n[it+1] = '\0';
}

## 效能測試
參考 [yanjiew 的筆記](https://hackmd.io/@yanjiew/linux2023q1-fibdrv#%E6%95%88%E8%83%BD%E5%88%86%E6%9E%90%E5%B7%A5%E5%85%B7%E8%A8%AD%E5%AE%9A%E8%88%87%E6%B8%AC%E8%A9%A6)，利用 static 變數儲存 ktime 數值並利用 write 回傳，並改寫 client.c 使其輸出結果至 out.txt。

static inline long long clock_gettime_ns()
{
    struct timespec t;
    clock_gettime(CLOCK_REALTIME, &t);
    return t.tv_sec * 1e9 + t.tv_nsec;
}

for (int i = 0; i <= offset; i++) {
    lseek(fd, i, SEEK_SET);
    long long temp = clock_gettime_ns();
    sz = read(fd, buf, sizeof(buf));
    long long utime = clock_gettime_ns() - temp;
    buf[sz] = 0;
    printf("Reading from " FIB_DEV
           " at offset %d, returned the sequence "
           "%s.\n",
           i, buf);
    long long ktime = write(fd, write_buf, sizeof(write_buf));
    printf("The kernel time from fib function is %lld.\n", kt);
}

> make
> sudo insmod ./fibdrv.ko
> sudo ./client > out.txt
> gnuplot scripts/gnuplot.gp
> eog plot.png

利用 GNUplot 繪製折線圖可以得到下圖，`ktime` 的數值大致隨著 `offset` 變多而往上升，且有一些極端案例存在。

將範圍放大到 500 可以得到下圖，但重複測試可以發現每次繪製的圖表落差巨大，需要排除干擾效能分析的因素：

### 將行程固定在特定的 CPU 中執行
修改 /etc/default/grub 檔案，修改為 `GRUB_CMDLINE_LINUX_DEFAULT="quiet splash isolcpus=11"`，指定 CPU 11 不會被排程器排程。並利用 `sudo update-grub` 儲存設定。

> taskset -cp 1
pid 1 目前的近似者清單：0-10

執行上述指令可發現修改成功。程式改以 `sudo taskset -c 11 ./client > out.txt` 執行（有多次執行），可以看到折線變得稍微平緩，但在 offset = 420 左右仍有極端峰值。

### 排除干擾效能分析的因素
參考[說明文件](https://hackmd.io/@sysprog/linux2023-fibdrv/%252F%2540sysprog%252Flinux2023-fibdrv-c)的提示，將所有指令皆執行過後再次測量時間，可以發現原有的極端峰值被移除，但關閉 SMT 後導致原有設定的 CPU 11 不能使用，再次檢查 taskset -cp 1 得知近似者僅有 0~5，故使用 CPU 5 執行。

> sudo sh -c "echo 0 > /proc/sys/kernel/randomize_va_space"
> sudo sh scripts/performance.sh
> sudo sh -c "echo 1 > /sys/devices/system/cpu/intel_pstate/no_turbo"
> sudo sh -c "echo off > /sys/devices/system/cpu/smt/control"
> taskset -cp 1
pid 1 目前的近似者清單：0-5


## 加速大數運算
### 紀錄認知與疑惑
遞迴法的 Fast Doubling 為使用 Top-down 實作，對遞迴的步驟進行剖析可發現有重複計算，為改進效能而改使用 Bottom-up Fast Doubling，即上方（`fast_fib()`）所設計的版本。

範例所使用的 bn 結構體以矩陣的方式實作，配合 size 與 sign 兩數值來紀錄矩陣的長度與數的符號。
- 加減法實作較為簡單，利用矩陣的相同欄的值相加再計算進位。
- 乘法實作中使用長乘法，以並改以 unsigned long long 儲存計算結果避免乘法溢出。
- bn 結構體可以被看成 $a=\sum\limits_{x = 0}^{size-1}{2^{32x}}number_{i}$，故乘法的最大 $size = a_{size}+b_{size} + 1$
- 利用雙層 for 迴圈來計算乘法數值，並用 `bn_mult_add` 來對 $number_{i+j}$ 與 $number_{i+j+1}$ 進行結果相加。
- **疑惑**：`bn_mult_add` 的計算值到 $i = c_{size}-1$ 才結束運算，但根據下方 `x >>= 32` 與 unsigned long long 為 64 位元可知執行兩次就會結束運算，若加上 carry 最多三次，那中止條件為何不改為 `i < offset+3`，是為了避免執行 `c->number[i]` 時出現 segfault 嗎？
- bn 結構體轉換成字串由於會需要經過二進位至十進位的轉換，反而比較複雜，這邊列舉一些疑惑的點：
- 長度 len 以 `8*sizeof(int)*size`取得，為何不直接利用 32*size 而要經過 sizeof(int) 運算？
- 原本的 len 運算是位元數 / 3.322，這邊除以 3 可以保證 len 一定比實際長度大，且後面會運用到 `memmove` 來去處理最後產出長度。
- carry 是得出的特定位元，只會是 1 或 0。由於十進位需要 $carry\times 2^{nth}$ 轉換才能得出正確數值，但程式碼中沒有看見，只有看見對於溢位的操作運算，所以看不懂是如何轉成十進位的。

## 自我檢查清單
- [x] 研讀上述 ==Linux 效能分析的提示== 描述，在自己的實體電腦運作 GNU/Linux，做好必要的設定和準備工作 $\to$ 從中也該理解為何不希望在虛擬機器中進行實驗;
- [x] 研讀上述費氏數列相關材料 (包含論文)，摘錄關鍵手法，並思考 [clz / ctz](https://en.wikipedia.org/wiki/Find_first_set) 一類的指令對 Fibonacci 數運算的幫助。請列出關鍵程式碼並解說
- [x] 複習 C 語言 [數值系統](https://hackmd.io/@sysprog/c-numerics) 和 [bitwise operation](https://hackmd.io/@sysprog/c-bitwise)，思考 Fibonacci 數快速計算演算法的實作中如何減少乘法運算的成本;
- [x] 學習指出針對大數運算的加速運算和縮減記憶體操作成本的舉措，紀錄你的認知和疑惑
- [ ] 注意到 `fibdrv.c` 存在著 `DEFINE_MUTEX`, `mutex_trylock`, `mutex_init`, `mutex_unlock`, `mutex_destroy` 等字樣，什麼場景中會需要呢？撰寫多執行緒的 userspace 程式來測試，觀察 Linux 核心模組若沒用到 mutex，到底會發生什麼問題。嘗試撰寫使用 [POSIX Thread](https://en.wikipedia.org/wiki/POSIX_Threads) 的程式碼來確認。 $\to$ 搭配閱讀〈[並行和多執行緒程式設計](https://hackmd.io/@sysprog/concurrency)〉
- [x] 詳閱第 2 和第 3 週所有教材及作業描述 [(一)](/@sysprog/linux2023-fibdrv-a), [(二)](/@sysprog/linux2023-fibdrv-b), [(三)](/@sysprog/linux2023-fibdrv-c), [(四)](/@sysprog/linux2023-fibdrv-d), [(五)](/@sysprog/linux2023-fibdrv-e)，記錄你的啟發和疑惑
- [x] 在 GitHub 上 fork [fibdrv](https://github.com/sysprog21/fibdrv)，目標是修正 Fibonacci 數計算並改善 `fibdrv` 核心模組的計算效率，過程中需要量化執行時間，可在 Linux 核心模組和使用者空間 (userspace) 中測量
- [x] 原本的程式碼只能列出到 $Fibonacci(100)$ 而且==部分輸出是錯誤的數值==，請修改程式碼，列出後方更多數值 (注意: 檢查正確性和數值系統的使用)
- [x] 務必研讀上述 ==Linux 效能分析的提示== 的描述，降低效能分析過程中的干擾;
- [x] 在 Linux 核心模組中，可用 ktime 系列的 API;
- [x] 在 userspace 可用 [clock_gettime](https://linux.die.net/man/2/clock_gettime) 相關 API;
- [ ] 善用統計模型，除去極端數值，過程中應詳述你的手法
- [x] 分別[用 gnuplot 製圖](https://hackmd.io/@sysprog/Skwp-alOg)，分析 Fibonacci 數列在核心計算和傳遞到 userspace 的時間開銷，單位需要用 us 或 ns (自行斟酌)
- [ ] 嘗試解讀上述實驗的時間分佈，特別是隨著 Fibonacci 數列增長後，對於 Linux 核心的影響。可修改 VFS 函式，允許透過寫入 `/dev/fibonacci` 裝置檔案來變更不同的 Fibonacci 求值的實作程式碼。
- [ ] 用 printk 固然方便，但其執行時間易受各種因素而不穩定，除了透過統計來觀察，也可改用核心的 sysfs 介面，後者可讓我們自行宣告 show 和 store 介面來提供讀寫操作，可參見 [Sample kobject implementation](https://elixir.bootlin.com/linux/latest/source/samples/kobject/kobject-example.c) (注意: 切換到對應的 Linux 核心版本)。
- [ ] 逐步縮減 Fibonacci 的執行時間，過程中要充分量化
- [ ] 嘗試研讀 [sysprog21/bignum](https://github.com/sysprog21/bignum) 的實作，理解其中的技巧並與你的實作進行效能比較，探討個中值得取鏡之處。
- [ ] 當 $Fib(n)$ 隨著 $n$ 越來越大時，由 Linux 核心傳遞字串形式的十進位就成為效能瓶頸
* 儘量降低由核心傳遞到應用程式的資料量
* 確保應用程式不需要過多的計算量才可重現 Fibonacci 數列
* BONUS: 研讀〈[Integer Encoding Algorithm 筆記](https://kkc.github.io/2021/01/30/integer-compression-note/)〉所提及的策略，在 Linux 核心模組先對數值進行編碼和壓縮，再由應用程式予以顯示 $Fib(n)$ 數值