Try   HackMD

2023q1 Homework3 (fibdrv)

contributed by < ericlai1021 >

開發環境

$ 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):                  8
  On-line CPU(s) list:   0-7
Vendor ID:               GenuineIntel
  Model name:            Intel(R) Core(TM) i5-8250U CPU @ 1.60GHz
    CPU family:          6
    Model:               142
    Thread(s) per core:  2
    Core(s) per socket:  4
    Socket(s):           1
    Stepping:            10
    CPU max MHz:         3400.0000
    CPU min MHz:         400.0000
    BogoMIPS:            3600.00
Virtualization features: 
  Virtualization:        VT-x
Caches (sum of all):     
  L1d:                   128 KiB (4 instances)
  L1i:                   128 KiB (4 instances)
  L2:                    1 MiB (4 instances)
  L3:                    6 MiB (1 instance)
NUMA:                    
  NUMA node(s):          1
  NUMA node0 CPU(s):     0-7

$ uname -r
5.19.0-35-generic

前期準備

  • 取得原始程式碼後,進行編譯及測試
  • 觀察產生的 fibdrv.ko 核心模組
$ modinfo fibdrv.ko

得到以下輸出

filename:       /home/eric/linux2023/fibdrv/fibdrv.ko
version:        0.1
description:    Fibonacci engine driver
author:         National Cheng Kung University, Taiwan
license:        Dual MIT/GPL
srcversion:     9A01E3671A116ADA9F2BB0A
depends:        
retpoline:      Y
name:           fibdrv
vermagic:       5.19.0-35-generic SMP preempt mod_unload modversions
  • 觀察 fibdrv.ko 核心模組在 Linux 核心掛載後的行為
    模組載入核心:
$ sudo insmod fibdrv.ko
  • 查看 /dev 目錄下可以發現的確有一個名為 fibonacci 的 character device
$ ls -l /dev/fibonacci
crw------- 1 root root 509, 0  三  16 13:44 /dev/fibonacci

繼續觀察

$ cat /sys/class/fibonacci/fibonacci/dev
509:0

可以注意到 509 這個數字,不同電腦設備會產生不同數字,試著對照 fibdrv.c ,發現關鍵在 alloc_chrdev_region() , 此函式會向 linux kernel 申請設備號。

$ cat /sys/module/fibdrv/version
0.1

這和 fibdrv.c 透過 MODULE_VERSION 所指定的版本號碼相同。

$ lsmod | grep fibdrv
fibdrv                 16384  0

$ cat /sys/module/fibdrv/refcnt
0

這兩道命令的輸出都是 0,意味著目前的 reference counting

  • 掛載 fibdrv.ko 核心模組後,執行 client 並觀察輸出
$ sudo ./cilent
Reading from /dev/fibonacci at offset 0, returned the sequence 0.
Reading from /dev/fibonacci at offset 1, returned the sequence 1.
Reading from /dev/fibonacci at offset 2, returned the sequence 1.
...
Reading from /dev/fibonacci at offset 92, returned the sequence 7540113804746346429.
Reading from /dev/fibonacci at offset 93, returned the sequence 7540113804746346429.
Reading from /dev/fibonacci at offset 94, returned the sequence 7540113804746346429.

可以發現在 offset 大於 92 之後都會得到同樣的值,猜測是因為 overflow ,於是先去查看 client.c ,基本上該檔案只是對 /dev/fibonacci 此裝置做開檔、讀檔、寫檔的動作,所以 fibonacci sequence 的計算都在 fibdrv.c 檔案內,馬上去查看,發現以下程式碼

/* MAX_LENGTH is set to 92 because
 * ssize_t can't fit the number > 92
 */
#define MAX_LENGTH 92

再往下看到 fib_device_lseek() ,此函式中會判斷欲找尋的位置是否大於 MAX_LENGTH ,若大於則回傳 MAX_LENGTH ,因此要讀取大於 F(92) 的數值都會讀到 F(92)。

Linux 效能分析

參考 作業說明 (三) 當中的設定

限定 CPU 給特定的程式使用

使用 isolcpus 這個 Linux 核心起始參數,修改 /etc/default/grub 設定檔,在 GRUB_CMDLINE_LINUX_DEFAULT 欄位加入 isolcpus=0 使開機時保留第 0 個 CPU 核給那些被 taskset 特別指定的行程使用。

GRUB_CMDLINE_LINUX_DEFAULT="quiet splash isolcpus=0"

修改完輸入 sudo update-grub 以更新設定,重新開機後確認設定是否生效。

$ cat /sys/devices/system/cpu/isolated
0

$ taskset -cp 1
pid 1's current affinity list: 1-7

可以發現第 0 個 CPU 確實被保留下來,再使用 taskset 查看 PID=1 的行程的處理器親和性,此時已經不包含第 0 個 CPU。

時間量測

使用 ktimer 相關的 API 來量測時間。 改寫 fibdrv.c 當中的 fib_read()fib_write() 兩個函式。
首先宣告一個 ktime_t 的變數:

static ktime_t kt;

針對 fib_read() 的修改如下:

static long long fib_time_proxy(long long k)
{
    kt = ktime_get();
    long long result = fib_sequence(k);
    kt = ktime_sub(ktime_get(), kt);

    return result;
}

/* calculate the fibonacci number at given offset */
static ssize_t fib_read(struct file *file,
                        char *buf,
                        size_t size,
                        loff_t *offset)
{
    return (ssize_t) fib_time_proxy(*offset);
}

針對 fib_write() 的修改如下:

/* write operation is skipped */
static ssize_t fib_write(struct file *file,
                         const char *buf,
                         size_t size,
                         loff_t *offset)
{
    return ktime_to_ns(kt);
}

我們可發現到 write 在此核心模組暫時沒作用,所以可以用來輸出上一次執行 fib_sequence() 的時間。
client.c 中先用 lseek(fd, i, SEEK_SET) 將檔案的讀寫位置指向特定位置 (即第 i 個位置),接著再用 read() 取得 Fib(i) 的數值,最後再使用 write() 取得計算 Fib(i) 的時間,相對應的程式碼如下:

for (int i = 0; i <= offset; i++) {
    lseek(fd, i, SEEK_SET);
    sz = read(fd, buf, 1);
    sz = write(fd, write_buf, strlen(write_buf));
}

改進費氏數列

Fast Doubling

參考作業說明 (一)(四) 以迭代方式實作 bottom-up Fast Doubling。

使用 Fast Doubling 可以得出以下公式

F(2k)=F(k)[2F(k+1)F(k)]F(2k+1)=F(k+1)2+F(k)2

對應的虛擬碼如下:

Fast_Fib(n)
    a = 0; b = 1;       // m = 0
    for i = (number of binary digit in n) to 1
        t1 = a*(2*b - a);
        t2 = b^2 + a^2;
        a = t1; b = t2; // m *= 2
        if (current binary digit == 1)
            t1 = a + b; // m++
            a = b; b = t1;
    return a;

搭配 bitwise 操作求得 number of binary digit in n 以及判斷 current binary digit 是否為 1。

  • 首先討論如何判斷當前位元是否為 1,想法很簡單,就是將 1 往左位移 i 個位元再跟 n 做 AND 運算即可求得 n 的第 i 位元是否為 1。
    對應的程式碼如下:
...
if (n & (1 << i)) {
    ...
}
...
  • 接著來討論如何求得 n 的二進制位數,這邊可以選擇是否搭配 clz / ctz 一類的指令來實作。
    • 使用 clz 指令取得 n 的 leading 0-bits 數量
      因為 n 的型態為 long long ,所以首先將 sizeof(long long) 向左位移 3 個位元即可求得 long long 所需的位元數,接著再減去 n 的 leading 0-bits 數量就可以求得 n 的二進制位數。
      對應的程式碼如下:
      ​​​​​​​​/* calculate the number of binary digits in k using clz */
      ​​​​​​​​int nbits = (sizeof(n) << 3) - __builtin_clzll(n);
      
    • 相反的,如果不使用 clz 指令的話,我的作法是簡單的透過一個迴圈產生跟 clz 指令同樣的結果
      想法是從 n 的 MSB 往下逐一判斷並搭配變數 count 表示 n 的 leading 0-bits 數,若當前位元為 1 則跳出迴圈,否則讓 count++ , 最後 (sizeof(long long) << 3) - count 即為 n 的 leading 0-bits 數量。
      對應的程式碼如下:
      ​​​​​​​​/* calculate the number of binary digits in k without using clz */
      ​​​​​​​​int nbits = 0;
      ​​​​​​​​int llsize = sizeof(n) << 3;
      ​​​​​​​​
      ​​​​​​​​for (int i = llsize - 1; i >= 0; i--) {
      ​​​​​​​​    if (n >> i)
      ​​​​​​​​        break;
      
      ​​​​​​​​    nbits++;
      ​​​​​​​​}
      ​​​​​​​​nbits = llsize - nbits;
      

實驗結果 (尚未排除干擾效能分析的因素)

  • 比較使用遞迴方法與 Fast Doubling (搭配 clz 指令)計算費氏數列的時間

可以發現 Fast Doubling 在計算 Fib(0) 時會花費大量時間,猜測是跟使用 clz 指令有關,目前正在做實驗研究,若我們將 y 軸範圍縮小到 2500 ns 會得到以下結果

即便尚未排除干擾效能的因素,仍可以看到 Fast Doubling 明顯優於遞迴方式

  • 目前對於 Fast Doubling 計算 Fib(0) 花費大量時間的解法為在 fib_sequence() 的一開始加上提前中止的條件處理
    ​​​​if (k <= 0)
    ​​​​    return k;
    

排除干擾效能分析的因素

  • 抑制 address space layout randomization (ASLR)

    ​​​​$ sudo sh -c "echo 0 > /proc/sys/kernel/randomize_va_space"
    
  • 設定 scaling_governor 為 performance。準備以下 shell script,檔名為 performance.sh:

    ​​​​for i in /sys/devices/system/cpu/cpu*/cpufreq/scaling_governor
    ​​​​do
    ​​​​    echo performance > ${i}
    ​​​​done
    

    之後再用 $ sudo sh performance.sh 執行

  • 針對 Intel 處理器,關閉 turbo mode

    ​​​​$ sudo sh -c "echo 1 > /sys/devices/system/cpu/intel_pstate/no_turbo"
    
  • 關閉 SMT (Intel 稱它為 Hyper-Threading)

    ​​​​$ sudo sh -c "echo off > /sys/devices/system/cpu/smt/control"
    

實驗結果 (排除干擾效能分析的因素)

  • 比較使用遞迴方法與 Fast Doubling (搭配 clz 指令)計算費氏數列的時間

可以看出來兩種方式的曲線都明顯平滑許多

  • 比較在 Fast Doubling 上有無使用 clz 指令的差別

可以看到使用 clz 指令對計算費氏數列有顯著的改善。

計算 Fib(92) 之後的 Fibonacci 數

參考作業說明 (四) 初步支援大數運算

支援大數運算結構體

為了能正確計算 92 項以後的費氏數列,引入長度可變動的數值表示法,藉由動態配置不同大小記憶體來呈現更大範圍的整數,結構體的定義如下
值得注意的地方是,我透過觀察發現即便是使用 fast doubling 在計算費氏數列的過程中,所有會用到自定義結構體 bn 宣告的變數 (e.g., a, b, t1, t2 ) 皆為非負整數,所以我就沒有定義 sign 來表示數值的正負。

/* 
 * number[size - 1] = msb, number[0] = lsb
 * sign = 1 for negative number
 */
typedef struct _bn {
    unsigned int *number;
    unsigned int size;
} bn;
  • number : 指向無號整數陣列的開頭
  • size : 配置給 number array 的大小,單位為 sizeof(unsigned int)

iterative 實作方法

初步嘗試使用大數運算來實作迭代作法的費氏數列計算,首先會需要加入以下功能

  1. bn_to_string ,將 bn 型態變數轉換成字串
  2. bn_add ,將兩個 bn 型態變數相加
  3. bn_cpy ,將一個 bn 型態變數複製到另一個 bn 型態變數(即做 assign 的動作)

在講解以上功能之前要先介紹支援 bn 型態的記憶體配置方法以及釋放方法,這邊為了能更加彈性的使用 bn 結構體,因此皆採用動態配置記憶體的方式來使用 bn 結構體

bn_alloc

先配置空間給 bn 結構體,再配置 sizeof(unsigned int) * size 的大小給 number 指標變數,最後將所有數值初始化為 0
參考教材 kmalloc vs. vmalloc 及網路材料 difference between kmalloc and vmalloc 後選擇使用 kmalloc 來配置記憶體是因為 kmalloc 會配置實體連續的記憶體空間,接著就可以簡單的使用 memset 來初始化此空間,但當沒有足夠的實體記憶體空間時 kmalloc 就會回傳錯誤,通常是回傳 NULL

/*
 * allocate a bn structure with the given size
 * the value is initialized to 0
 */
bn *bn_alloc(size_t size)
{
    bn *new = (bn *) kmalloc(sizeof(bn), GFP_KERNEL);
    new->number =
        (unsigned int *) kmalloc(sizeof(unsigned int) * size, GFP_KERNEL);
    // QUESTION? physical or virtual contiguous memory allocation
    memset(new->number, 0, sizeof(unsigned int) * size);
    new->size = size;

    return new;
}

bn_free

依序將 number 以及 bn 結構體釋放

/*
 * free entire bn data structure
 * return 0 on success, -1 on error
 */
int bn_free(bn *src)
{
    if (src == NULL)
        return -1;
    kfree(src->number);
    kfree(src);
    return 0;
}

bn_to_string

由於大數沒辦法直接以數值的形式列出,這裡改用字串來呈現,轉換的概念很簡單,就是將大數除以 10 取餘數取得個位數,接著再將大數更新為商數,不斷重複此步驟直到大數為 0

memset(s, '0', len - 1);
s[len - 1] = '\0';

int cur = len - 2;
while (!bn_is_zero(src)) {
    unsigned long long remainder = 0U;
    // src->number[size - 1] is msb
    for (int i = src->size - 1; i >= 0; i--) {
        unsigned int quotient = 0U;
        // walk through every bit of bn
        for (unsigned int d = 1U << 31; d > 0; d >>= 1) {
            if (src->number[i] & d) {
                remainder = (remainder << 1) | 1U;
            } else {
                remainder = (remainder << 1) | 0U;
            }
            quotient = (quotient << 1) | (remainder / 10U);
            remainder %= 10;
        }
        src->number[i] = quotient;
    }
    s[cur] += remainder;
    cur--;
}

bn_add

實作加法的概念很直觀,就是從 lsb 開始相加,將進位的部份加到下一輪的相加結果

/* c = a + b */
void bn_add(const bn *a, const bn *b, bn *c)
{
    bn_resize(c, MAX(a->size, b->size) + 1);

    unsigned long long carry = 0;
    for (int i = 0; i < c->size; i++) {
        unsigned int tmp1 = (i < a->size) ? a->number[i] : 0;
        unsigned int tmp2 = (i < b->size) ? b->number[i] : 0;
        carry += (unsigned long long) tmp1 + tmp2;

        c->number[i] = carry & 0x00000000FFFFFFFF;
        carry >>= 32;
    }

    if (!c->number[c->size - 1] && c->size > 1)
        bn_resize(c, c->size - 1);
}

bn_cpy

使用 memcpy 將一個 bn 結構體的數值部份複製到另一個 bn 結構體的數值部份

/*
 * copy the value from src to dest
 * return 0 on success, -1 on error
 */
int bn_cpy(bn *dest, const bn *src)
{
    if (bn_resize(dest, src->size) < 0)
        return -1;
    memcpy(dest->number, src->number, src->size * sizeof(unsigned int));
    return 0;
}

fast doubling 實作方法

使用大數運算實作 fast doubling 會需要額外實作以下功能

  1. bn_sub ,將兩個 bn 型態變數相減
  2. bn_lshift ,將 bn 型態變數向左位移
  3. bn_mult ,將兩個 bn 型態變數相乘

bn_sub

實作減法的概念類似加法,一樣從 lsb 開始相減,若相減結果為負數,則向下一位借位

/* c = a - b */
void bn_sub(const bn *a, const bn *b, bn *c)
{
    bn_resize(c, MAX(a->size, b->size));

    long long borrow = 0;
    for (int i = 0; i < c->size; i++) {
        unsigned int tmp1 = (i < a->size) ? a->number[i] : 0;
        unsigned int tmp2 = (i < b->size) ? b->number[i] : 0;

        borrow = (long long) tmp1 - tmp2 - borrow;
        if (borrow < 0) {
            c->number[i] = borrow + (1LL << 32);
            borrow = 1;
        } else {
            c->number[i] = borrow;
            borrow = 0;
        }
    }

    int clz_ints = bn_clz(c) / 32;
    if (clz_ints == c->size)
        --clz_ints;

    bn_resize(c, c->size - clz_ints);
}

bn_shift

bn 結構體左移 shift 個 bit ,因為左移超過 32 bits 的處理相對複雜,因此就先做出比較簡單的版本

/* dest = src << shift */
void bn_lshift(bn *dest, bn *src, size_t shift)
{
    size_t z = bn_clz(src);
    shift %= 32;  // only handle shift within 32 bits
    if (!shift)
        return;

    if (shift > z)
        bn_resize(src, src->size + 1);

    bn_cpy(dest, src);
    /* bit shift */
    for (int i = src->size - 1; i > 0; i--)
        dest->number[i] =
            src->number[i] << shift | src->number[i - 1] >> (32 - shift);
    dest->number[0] <<= shift;
}

bn_mult

實作大數乘法的概念類似硬體乘法器 (如下圖),從 multiplier 的 lsb 開始判斷,若當前 bit 為 1 ,則將 multiplicand 加至 product ,再將 product 往左位移 1 bit ,重複上述直到 multiplier 的判斷到達 msb 為止

/* c = a * b */

bn *tmp1 = bn_alloc(1);
bn_cpy(tmp1, a);
if (b->number[0] & 1U)
    bn_add(c, tmp1, c);

for (int i = 1; i < bn_digit(b); i++) {
    bn_lshift(tmp1, tmp1, 1);
    if (b->number[i / 32] & (1U << (i % 32)))
        bn_add(c, tmp1, c);
}
bn_free(tmp1);

int clz_ints = bn_clz(c) / 32;
if (clz_ints == c->size)
    --clz_ints;

bn_resize(c, c->size - clz_ints);

正確計算
F(92)
以後的數值

使用以上實作出來的大數 API 來正確計算 92 項以後的費氏數列,程式邏輯皆參照原始的實作,首先展示迭代演算法

void bn_fib(bn *dest, unsigned int n)
{
    bn_resize(dest, 1);
    if (n < 2) {  // Fib(0) = 0, Fib(1) = 1
        dest->number[0] = n;
        return;
    }
    bn *a = bn_alloc(1);
    bn *b = bn_alloc(1);
    dest->number[0] = 1;
    for (unsigned int i = 1; i < n; i++) {
        bn_cpy(b, dest);
        bn_add(dest, a, dest);
        bn_cpy(a, b);
    }
    bn_free(a);
    bn_free(b);
}

接著是 fast doubling 演算法的實作

void bn_fib_fast(bn *dest, unsigned int n)
{
    bn_resize(dest, 1);
    if (n < 2) {  // Fib(0) = 0, Fib(1) = 1
        dest->number[0] = n;
        return;
    }

    int nbits = (sizeof(n) << 3) - __builtin_clz(n);
    bn *a = bn_alloc(1);
    a->number[0] = 0U;
    bn *b = bn_alloc(1);
    b->number[0] = 1U;

    for (int i = nbits - 1; i >= 0; i--) {
        bn *t1 = bn_alloc(1);
        bn_lshift(t1, b, 1);
        bn_sub(t1, a, t1);
        bn_mult(a, t1, t1);

        bn *t2 = bn_alloc(1);
        bn_mult(b, b, b);
        bn_mult(a, a, a);
        bn_add(a, b, t2);

        bn_cpy(a, t1);
        bn_cpy(b, t2);

        if (n & (1U << i)) {
            bn_add(a, b, t1);
            bn_cpy(a, b);
            bn_cpy(b, t1);
        }

        bn_free(t1);
        bn_free(t2);
    }
    bn_cpy(dest, a);

    bn_free(a);
    bn_free(b);
}

使用 /scripts/verify.py 來驗證,至少能夠正確計算至第 3000 項

改善大數運算

後續改進皆更新在 2023期末專題