---
tags: linux2022
---
# 2022q1 Homework3 (fibdrv)
contributed by < `ibat10clw` >
> [fibdrv](https://github.com/ibat10clw/fibdrv)
## 費氏數列
### 版本一 (DP)
為理解 client 和 fibdrv 之間的互動，我的第一個版本只是將 `f[i] = f[i-1] + f[i-2]` 以字元陣列儲存，並使用字元加法做改寫，同時根據作業說明中 clinet.c 的提示使用 `copy_to_user` 將 fibdrv 計算結果的字元陣列取出。
```c
static void add(char *a, char *b, char *out)
{
int len_a = strlen(a);
int len_b = strlen(b);
int i = 0;
int carry = 0;
int diff = len_a - len_b;
for (i = len_b - 1; i >= 0; i--) {
int sum = (a[i + diff] - '0') + (b[i] - '0') + carry;
out[i + diff] = '0' + sum % 10;
carry = sum / 10;
}
for (i = diff - 1; i >= 0; i--) {
int sum = (a[i] - '0') + carry;
out[i] = '0' + sum % 10;
carry = sum / 10;
}
if (carry) {
if (diff == 0)
len_a++;
for (i = len_a; i >= 0; i--) {
out[i + 1] = out[i];
}
out[len_a + 1] = 0;
out[0] = '0' + carry;
}
out[len_a] = 0;
}
```
加法部份是將傳入的字元陣列以較短的為準對齊，依序從最右邊相加到比較短陣列結束後，同時需要紀錄當前的相加是否需要進位，若有進位需求的話，則需要把結果累計到下一次的計算，再把剩下的數字填入 buf。
```c
static long my_fibo(int k, char *buf)
{
char f3[SIZE];
char f2[SIZE] = "1";
char f1[SIZE] = "0";
for (int i = 2; i <= k; ++i) {
add(f2, f1, f3);
strncpy(f1, f2, strlen(f2) + 1);
strncpy(f2, f3, strlen(f3) + 1);
}
if (k == 0) {
copy_to_user(buf, f1, strlen(f1) + 1);
return 0;
}
if (k == 1) {
copy_to_user(buf, f2, strlen(f2) + 1);
return 1;
}
copy_to_user(buf, f3, strlen(f3) + 1);
return k;
}
```
關於費氏數列的計算是以 DP 的方式依序將 `f(n-1)` 與 `f(n-2)` 相加後把結果存到 `f3` ，之後分別把 `f2` 和 `f3` 複製給 `f1` 和 `f2`以便進行下一輪的計算，而 buf 是由 `bif_read` 中的參數傳入的，因此最後只要將結果放入 `copy_to_user` 即可在 client 中將其印出。
### 版本二 (Fast Doubling)
若要使用 fast doubling 的方法計算費氏數列還需要引入乘法與減法的操作
```c
void sub(char *a, char *b, char *out)
{
int len_a = strlen(a);
int len_b = strlen(b);
int i = 0;
int borrow = 0;
int idx = len_a - len_b;
for (i = len_b - 1; i >= 0; i--) {
int diff = (a[i + idx] - '0') - (b[i] - '0') - borrow;
if (diff < 0) {
diff += 10;
borrow = 1;
} else
borrow = 0;
out[i + idx] = '0' + diff;
}
for (i = idx - 1; i >= 0; i--) {
int diff = (a[i] - '0') - borrow;
if (diff < 0) {
diff += 10;
borrow = 1;
} else
borrow = 0;
out[i] = '0' + diff;
}
out[len_a] = 0;
i = 0;
while (out[i] == '0')
i++;
if (i != 0) {
for (int j = 0; j < len_a; ++j)
out[j] = out[j + i];
}
}
```
減法與加法的想法是類似的，不過從需要檢查進位變成需要檢查借位，同樣的以較短的輸入為準對齊後從右邊開始計算回來，每個回合將對應的位數相減後判斷是否需要借位，借助著直式減法的思路，當需要借位時可以將當前的差 +10 後把前一個位數 -1 ，透過這種方法便可以達成以字元鎮列為儲存單位的減法
最後的 while 迴圈是判斷是否有 leading zero ，若有的話就將有效位數全部向前挪動
```c
void mul(char *a, char *b, char *out)
{
int len_a = strlen(a);
int len_b = strlen(b);
int N = len_a + len_b;
int i;
memset(out, '0', N);
out[N] = 0;
for (i = len_b - 1; i >= 0; i--) {
for (int j = len_a - 1; j >= 0; j--) {
int digit1 = a[j] - '0';
int digit2 = b[i] - '0';
int pos = 1 + i + j;
int carry = out[pos] - '0';
int multiplication = digit1 * digit2 + carry;
out[pos] = (multiplication % 10) + '0';
out[pos - 1] += (multiplication / 10);
}
}
i = 0;
while (out[i] == '0')
i++;
if (i != 0) {
for (int j = 0; j < N; ++j)
out[j] = out[j + i];
}
}
```
字元的乘法與前面的流程稍有不同，但一樣是模擬直式的計算流程取得結果
首先把兩字元陣列的長度相加，並初始化為 0 以儲存計算結果
接著兩層的 for 迴圈分別看過所有存在陣列內的數字，兩兩計算相乘的結果，同時每次計算時還得先找出要將結果填在 out 的哪個位置，若有進位的話則在將其填在 pos - 1 的地方，當所有的數字都看過後即可取得結果
最底下的 while 迴圈同樣是去除存在陣列中的 leading zero
```c=
static long fast_fib(int k, char *buf)
{
if (k == 0) {
copy_to_user(buf, "0", 2);
return k;
}
if (k == 1) {
copy_to_user(buf, "1", 2);
return k;
}
int stack[100];
int top = 0;
while (k > 1) {
stack[top++] = k;
k /= 2;
}
char f1[SIZE] = "1";
char f2[SIZE] = "1";
char tmp1[SIZE];
char tmp2[SIZE];
char res1[SIZE];
while (top--) {
int n = stack[top];
// f1 = f1 * (2 * f2 - f1)
mul("2", f2, tmp1);
sub(tmp1, f1, tmp2);
mul(tmp2, f1, res1);
// f2 = f1 ** 2 + f2 ** 2
mul(f1, f1, tmp1);
mul(f2, f2, tmp2);
add(tmp2, tmp1, f2);
if (n % 2 == 1) {
add(f2, res1, tmp1);
strncpy(f1, f2, strlen(f2) + 1);
strncpy(f2, tmp1, strlen(tmp1) + 1);
} else {
strncpy(f1, res1, strlen(res1) + 1);
}
}
copy_to_user(buf, f1, strlen(f1) + 1);
return k;
}
```
在使用 fast doubling 的方法時需要根據當前的 fibonacci number 是偶數項或奇數項來決定下次的計算方法，因此這邊使用了一個 stack 來紀錄項數的奇偶，以便決定後續的計算方式。
當項數為奇數的時候，就用原始定義的 `fib(n) = fib(n-1) + fib(n-2)` 把數列往前推進一項，等待全部計算完畢後同樣使用 copy_to_user 將結果傳出。
### 效能評測
根據作業說明的提示，可以使用 `ktime` 相關的 API 在 kernel 中紀錄計算 fibonacci 的時間，再搭配還沒有被使用上的 `fib_write` 從 client 中獲得紀錄的時間
```c
static long long fib_time_proxy(long long k, char *buf)
{
kt = ktime_get();
long long result = fast_fib(k, buf);
kt = ktime_sub(ktime_get(), kt);
return result;
}
```
在呼叫 fast_fib 的前後先計算時間差，之後在從 client 呼叫 write 得到結果
```c
for (int i = 0; i <= offset; i++) {
lseek(fd, i, SEEK_SET);
clock_gettime(CLOCK_REALTIME, &t1);
sz = read(fd, buf, sizeof(buf));
clock_gettime(CLOCK_REALTIME, &t2);
sz = write(fd, write_buf, strlen(write_buf));
printf("%d %lld %ld %lld\n", i, sz, t2.tv_nsec - t1.tv_nsec, (t2.tv_nsec - t1.tv_nsec) - sz);//k, u, k2u
}
```
此外在 client 中同樣在 read 前後各紀錄一次時間，就可以計算執行時間
而 `copy_to_user` 的開銷則可以透過總執行時間扣掉在 kernel 的時間，空缺的那一塊就是 `copy_to_user` 的開銷
## 大數運算
上面的費氏數列實作雖然能計算出正確的答案，但由於使用字串的方式模擬 10 進位的運算，每次運算相當於只能處理 1 個位數，以現在 CPU 的規格來說並沒有完全發揮他的效能，因此先全面改寫運算的規則後再重新實作費氏數列。
### 用於儲存大數的資料結構
```c
typedef struct _bn {
unsigned int number[BNSIZE];
unsigned int size;
} bn;
```
`unsigned int` 共可以儲存 32 bit 的資料，使用一個陣列儲存運算的結果，並且根據 size 來決定目前有效資料的長度，LSB 存於 index 0 的位置
示意圖如下:
若要解釋以上資料則從 index 0 開使將數字以 2 進位的方式連接，可得以下結果:
`00000000 00000000 000000001 11000011`
如此便可以將 2 進位的格式在輸出時轉換為 10 進位
### 將大數以十進位的方式輸出
這個部份參考了 [KYG-yaya573142](https://hackmd.io/@KYWeng/rkGdultSU#%E5%A4%A7%E6%95%B8%E9%81%8B%E7%AE%97) 的實作
```c=
static char *bn_to_string(const bn *src)
{
// log10(x) = log2(x) / log2(10) ~= log2(x) / 3.322
size_t len = (8 * sizeof(int) * src->size) / 3 + 2;
char *s = kmalloc(len, GFP_KERNEL);
char *p = s;
memset(s, '0', len - 1);
s[len - 1] = '\0';
for (int i = src->size - 1; i >= 0; i--) {
for (unsigned long d = 1U << 31; d; d >>= 1) {
// double logical not can make the result to 1
int carry = !!(d & src->number[i]);
for (int j = len - 2; j >= 0; j--) {
s[j] += s[j] - '0' + carry;
carry = (s[j] > '9');
if (carry)
s[j] -= 10;
}
}
}
// skip leading zero
while (p[0] == '0' && p[1] != '\0') {
p++;
}
memmove(s, p, strlen(p) + 1);
return s;
}
```
一般將 2 進位的數字轉換為 10 進位時，可以借助下列公式
$\displaystyle\sum_{b=1}^{N} {2^N}\times b$ `for b = lsb to msb`
範例如下: `00000000 00000000 000000001 11000011`
根據 bit=1 所在的位置可以得出該 2 進位數的 10 進位表示為
$2^8+2^7+2^6+2^1+2^0$
由於輸入的 bn 可能非常的大，因此需要直接將其轉為字元陣列的形式保存，演算法描述如下
1. 從 bn 的 MSB 開始讀過所有的有效位元(根據size取捨)
* 對應程式碼第 10 行的 for 迴圈
3. 由於 bn 的每個欄位是 32 bit，因此使用 32 bit 的無號數對每個 bit 做 bitwise and 判斷是否為1，若出現的是 1 就將字元陣列 +1
* 對應程式碼 11 行的 for 迴圈
5. 從初始化字元陣列的尾巴開始把每個位數扣掉 ASCII 的偏移後與自己相加，相當於乘以 2 的效果，由於是從 MSB 開始做轉換，當第一個 bit 1 出現到轉換結束執行的次數恰好能將其恢復成 power of 2，同時也需要判斷字元陣列表示的十進位數是否需要進位
* 對應程式碼 14 行的 for 迴圈
如此便能夠在輸出時將 bn 轉化為十進位表示的字串，並且透過 `copy_to_user` 將結果傳到 userspace
### 大數加法
bn 的儲存格式在選擇時使用了 unsigned integer，unsigned integer 在 overflow 發生的時候會將結果 `mod` 該型態的最大值後再做保存，根據這個特性只需要判斷是否發生 overflow，若發生的話就將 bn 的下一個欄位 +1 來做進位，如此便能完成大數加法
```c=
#define MAXSLOT 0xffffffffUL+1UL
static void add_bn(bn *a, bn *b, bn *c)
{
unsigned long carry = 0;
int i;
c->size = a->size;
for (i = 0; i < b->size; ++i) {
unsigned int sum = a->number[i] + b->number[i] + carry;
carry =
((unsigned long) a->number[i] + b->number[i] + carry) >= MAXSLOT;
c->number[i] = sum;
}
for (; i < a->size; ++i) {
unsigned int sum = a->number[i] + carry;
carry = ((unsigned long) a->number[i] + carry) >= MAXSLOT;
c->number[i] = sum;
}
if (carry) {
c->size++;
c->number[i]++;
}
}
```
其中 MAXSLOT 定義為 `UINT_MAX+1` 是為了 `branchless` 的實作，只要當前要相加的 bn 和 carry 加在一起後判斷是否會超過 `MAXSLOT` 若會的話則代表 overflow 發生了，這時就要將 carry 設為 1，這部份的操作則是透過 logical operator 來達成
Sign in with Wallet
Connect another wallet