Try   HackMD

Linux 核心專題: 重做 fibdrv

執行人: seasonwang0905
期末專題解說影片

Image Not Showing Possible Reasons
  • The image file may be corrupted
  • The server hosting the image is unavailable
  • The image path is incorrect
  • The image format is not supported
Learn More →
提問清單

  • ?

任務簡述

依據 fibdrv 作業規範,繼續投入 Linux 核心模組和相關程式的開發。

TODO: 紀錄閱讀作業說明中所有的疑惑

閱讀 fibdrv 作業規範,包含「作業說明錄影」和「Code Review 錄影」,本著「誠實面對自己」的心態,在本頁紀錄所有的疑惑,並與授課教師預約討論。

過程中,彙整 Homework3 學員的成果,挑選至少三份開發紀錄,提出值得借鏡之處,並重現相關實驗。

TODO: 回覆「自我檢查清單」

回答「自我檢查清單」的所有問題,需要附上對應的參考資料和必要的程式碼,以第一手材料 (包含自己設計的實驗) 為佳

TODO: 以 sysprog21/bignum 為範本,實作有效的大數運算

理解其中的技巧並導入到 fibdrv 中,並留意以下:

  • 在 Linux 核心模組中,可用 ktime 系列的 API;
  • 在 userspace 可用 clock_gettime 相關 API;
  • 善用統計模型,除去極端數值,過程中應詳述你的手法
  • 分別用 gnuplot 製圖,分析 Fibonacci 數列在核心計算和傳遞到 userspace 的時間開銷,單位需要用 us 或 ns (自行斟酌)

彙整及重現加速 Fibonacci 運算

參考 eecheng87 的實作,他在 fast doubling 的基礎上加入了 __builtin_clz() 來加速 Fibonacci 的運算。

Fast doubling 無 __builtin_clz() 版本:

unsigned long long double_fib1(int n)
{
    unsigned long long f1 = 0, f2 = 1, t1, t2;
    for (int i = 31; i >= 0; i--) {
        t1 = f1 * (f2 * 2 - f1);
        t2 = f2 * f2 + f1 * f1;
        f1 = t1;
        f2 = t2;
        if ((n & (1 << i))) {
            t1 = f1 + f2;
            f1 = f2;
            f2 = t1;
        }
    }
    return f1;
}

Fast doubling 有 __builtin_clz() 版本:

unsigned long long double_fib2(int n)
{
    unsigned long long f1 = 0, f2 = 1;
    int i = 31 - __builtin_clz(n);
    for (; i >= 0; i--) {
        unsigned long long t1, t2;
        t1 = f1 * (f2 * 2 - f1);
        t2 = f2 * f2 + f1 * f1;
        f1 = t1;
        f2 = t2;
        if ((n & (1 << i))) {
            t1 = f1 + f2;
            f1 = f2;
            f2 = t1;
        }
    }
    return f1;
}

在加入 __builtin_clz() 前,不論欲求的 n 值為何,迴圈必定會執行 32 ,加入後,可以減少迭代的次數,避免在 n 值很小的時候執行過多的計算

再來,參考並修改了 Hsuhaoxiang 的作法,他先計算了傳入的 n 值 (這裡為 k) 在 bitwise 中有多少個 1 並紀錄,以此切換 a 及 b 的位置

unsigned long long double_fib3(int k)
{
    int bin_k[8];
    int count = 0;
    while (k) {
        bin_k[count] = k % 2;
        k = k >> 1;
        count++;
    }
    unsigned long long int a = 0, b = 1, t1, t2, temp;
    for (int i = count - 1; i >= 0; i--) {
        t1 = a * ((2 * b) - a);
        t2 = b * b + a * a;
        a = t1;
        b = t2;
        if (bin_k[i] == 1) {
            temp = a + b;
            a = b;
            b = temp;
        }
    }
    return a;
}

這裡直接計算 k 的 1 位元個數,然後由 fast doubling 可減少一半迭代數的特性算出 count ,在效能上已和使用 __builtin_clz() 類的手法接近,改善的方式可能要從記憶體的使用或是加入大數運算

加入大數運算

因為 c 語言沒有合適的資料型態可以存儲超過

2641
264
的數值,故我們需要使用其他手段來作大數運算,參考 eecheng87 的作業中所述,他使用以下結構來存取一個動態長度的陣列

typedef struct {
    char digits[MAXDIGITS]; /* represent the number */
    int signbit;            /* 1 if positive, -1 if negative */
    int lastdigit;          /* index of high-order digit */
} bignum;

加法的實作為

int add_bignum(bignum *a, bignum *b, bignum *c)
{
    ..
    if (a->lastdigit < b->lastdigit)
        return add_bignum(b, a, c);

    int k = c->lastdigit = a->lastdigit + 1;

    c->digits[k--] = '\0';
    carry = 0;
    n_carry = 0;
    for (i = b->lastdigit - 1, j = a->lastdigit - 1; i >= 0; i--, j--) {
        carry = b->digits[i] - '0' + a->digits[j] - '0' + carry;
        c->digits[k--] = (carry % 10) + '0';
        carry = carry / 10;
        if (carry)
            n_carry++;
    }
    for (; j >= 0; j--) {
        carry = a->digits[j] - '0' + carry;
        c->digits[k--] = (carry % 10) + '0';
        carry = carry / 10;
        if (carry)
            n_carry++;
    }
    if (carry)
        c->digits[k] = carry + '0';
    else {
        char string[MAXDIGITS];
        strlcpy(string, &c->digits[1], MAXDIGITS);
        strlcpy(c->digits, string, MAXDIGITS);
        c->lastdigit = c->lastdigit - k - 1;
    }
    return n_carry;
}

上面省去了判斷加負數的部份,加負數可以替代成減法的形式。我們需要以字串來紀錄大數,使用時再將其一個一個輸出成大數的樣子

再來是減法的部份

int subtract_bignum(bignum *a, bignum *b, bignum *c)
{
    // int borrow; /* has anything been borrowed? */
    // int v; /* placeholder digit */
    register int i, j, op = 0; /* counter */
    int n_borrow;
    int temp;
    c->signbit = PLUS;

    if ((a->signbit == MINUS) || (b->signbit == MINUS))
    {
        b->signbit = -1 * b->signbit;
        n_borrow = add_bignum(a, b, c);
        b->signbit = -1 * b->signbit;
        return n_borrow;
    }
    ..
        
    ..
    int k = c->lastdigit = MAX(a->lastdigit, b->lastdigit);
    n_borrow = 0;
    c->digits[k--] = '\0';
    for (i = a->lastdigit - 1, j = b->lastdigit - 1; j >= 0; i--, j--) {
        temp = a->digits[i] - '0' - (b->digits[j] - '0' + op);

        if (temp < 0) {
            temp += 10;
            op = 1;
            n_borrow++;
        } else
            op = 0;
        c->digits[k--] = temp + '0';
    }
    while (op) {
        temp = a->digits[i--] - op - '0';
        if (temp < 0) {
            temp += 10;
            op = 1;
            n_borrow++;
        } else
            op = 0;
        c->digits[k--] = temp + '0';
    }
    for (; i >= 0; i--)
        c->digits[k--] = a->digits[i];
    for (i = 0; !(c->digits[i] - '0'); i++)
        ;
    c->lastdigit = c->lastdigit - i;
    if (i == a->lastdigit)
        strlcpy(c->digits, "0", MAXDIGITS);
    else {
        char string[MAXDIGITS];
        strlcpy(string, &c->digits[i], MAXDIGITS);
        strlcpy(c->digits, string, MAXDIGITS);
    }
    return n_borrow;
}

減法中有很多借位的算法,目前還在嘗試理解

這邊有個疑問:加法、減法輸出 n_carryn_borrow 進位或是借位目的是什麼,用在哪裡?

最後是乘法

void multiply_bignum(bignum *a, bignum *b, bignum *c)
{
    // long int n_d;
    register long int i, j, k = 0;
    short int num1[MAXDIGITS], num2[MAXDIGITS], of = 0, res[MAXDIGITS] = {0};
    // n_d = (a->lastdigit < b->lastdigit) ? b->lastdigit : a->lastdigit;
    // n_d++;
    for (i = 0, j = a->lastdigit - 1; i < a->lastdigit; i++, j--)
        num1[i] = a->digits[j] - 48;
    for (i = 0, j = b->lastdigit - 1; i < b->lastdigit; j--, i++)
        num2[i] = b->digits[j] - 48;
    res[0] = 0;
    for (j = 0; j < b->lastdigit; j++) {
        for (i = 0, k = j; i < a->lastdigit || of; k++, i++) {
            if (i < a->lastdigit)
                res[k] += num1[i] * num2[j] + of;
            else
                res[k] += of;
            of = res[k] / 10;
            res[k] = res[k] % 10;
        }
    }
    for (i = k - 1, j = 0; i >= 0; i--, j++)
        c->digits[j] = res[i] + 48;
    c->digits[j] = '\0';
    c->lastdigit = k;
    c->signbit = a->signbit * b->signbit;
}

因為字串沒有所謂的相乘,所以只能把數字擷取出來再做乘法,最後把結果存到 c.->digits

有了以上大數運算,就可以運用到 Fibonacci 的數列上了

Fast doubling 演算法支援大數運算

bignum a, b;
bignum big_two;
int_to_bignum(0, &a);
int_to_bignum(1, &b);
int_to_bignum(2, &big_two);
for (int i = 31 - __builtin_clz(k); i >= 0; i--) {
    bignum t1, t2;
    bignum tmp1, tmp2;
    multiply_bignum(&b, &big_two, &tmp1);
    (void) subtract_bignum(&tmp1, &a, &tmp2);
    multiply_bignum(&a, &tmp2, &t1);

    multiply_bignum(&a, &a, &tmp1);
    multiply_bignum(&b, &b, &tmp2);
    (void) add_bignum(&tmp1, &tmp2, &t2);
    copy(&a, &t1);
    copy(&b, &t2);
    if ((k & (1 << i)) > 0) {
        (void) add_bignum(&a, &b, &t1);
        copy(&a, &b);
        copy(&b, &t1);
    }
}
return a;

輸入第 100 、200 項時,看看輸出結果為何

354224848179261915075                         \\ n = 100
280571172992510140037611932413038677189525    \\ n = 200

比對 The Fibonacci numbers 發現相同,大數運算成功執行。

研讀上述費氏數列相關材料(包含論文),摘錄關鍵手法,並思考 clz / ctz 一類的指令對 Fibonacci 數運算的幫助。請列出關鍵程式碼並解說

James L. Holloway 發表於 1988 年的碩士論文 〈Algorithms for Computing Fibonacci Numbers Quickly〉針對多種 Fibonacci 數列之演算法進行模擬及分析,包括遞迴、重複相加、Binet's Formula (公式解)、Vorobev等式、二項式係數法等等逐項討論,以實際成果來展示各演算法在 Fibonacci 數列之項數 n 改變的情況下,其執行時間及效率究竟如何,最後驗證其論文結果。

  • 遞迴 (natural recursive)

fib(n)={0ifn=01ifn=1fib(n1)+fib(n2)ifn2

在介紹函式遞迴的時,常見的例子即計算 Fibonacci 數,因重複呼叫函式的關係,演算法時間複雜度來到

O(λ1n),其中
λ1>1.5
,隨著 n 值增加,執行時間會以指數倍成長,故在 n 值較大時,幾乎不會選用此演算法。

  • 重複相加 (repeated addition)

原理與前述相同,只不過使用了 for 迴圈來取代函式遞迴的動作,效率卻有顯著的提昇,隨著 n 值增加,執行時間則呈線性增長。另外,論文中提到的 K-Fibonacci 數列也有相關的演算法,不過這裡暫不討論此部份。

透過 gnuplot 製圖:

橫軸代表 Fibonacci 數列第 n 項,即

fib(n),縱軸則是其執行時間,不難看出僅僅是 n = 10 ,遞迴演算法的執行時間已來到最初的 6 倍了,其效率可說是驚人的低落。重複相加的演算法可以作為 gnuplot 的練習及效能的參考。

接著,列出改善計算 Fibonacci 數列的演算法

F(2n)=2F(n)×F(n+1)F(n)2F(2n+1)=F(n)2+F(n+1)2

根據上述關係,著手撰寫程式碼

unsigned long long fast_doubling(int n)
{
    if (n <= 2)
        return !!n;
    
    unsigned long long f1 = 1,f2 = 1, t1, t2, temp;
    int count = n / 2;
    for (int i = 0; i < count; i++)
    {
        t1 = 2 * f1 * f2 - f1 * f1;
        t2 = f1 * f1 + f2 * f2;
        temp = f2;
        f2 = f1 + temp;
        f1 = temp;
    }

    if (n & 1)
        return t2;
    return t1;
}

寫完後,馬上來測試效能如何,結果如下圖所示

藍色的線是參考上述提及的 eecheng87 的 fast doubling 無硬體加速版本的演算法,綠色則是我實作的部份,可以發現綠色的執行時間與 n 呈線性成長,效能只有比重複相加的演算法好一些

提問:在 fast doubling 演算法中,最多只需要做

n/2 次的迭代即可算完每項
F(n)
,即使我設定 n 不會超過 32 次,輸出結果是正確的,但執行時間仍就無法降低呢?是否是因為多了 f2 = f1 + temp 這一項呢?又或者應該要支援大數運算後再重新比較效能看看?

你用 perf 一類的工具去觀察執行的表現,會發現裡頭存在耗時的指令

Image Not Showing Possible Reasons
  • The image file may be corrupted
  • The server hosting the image is unavailable
  • The image path is incorrect
  • The image format is not supported
Learn More →
jserv