2020q3 Homework (quiz4)

contributed by < hankluo6 >

測驗 1

int hammingDistance(int x, int y) {
    return __builtin_popcount(x OP y);
}

`OP`

(c) ^

使用 XOR 會將兩個有差異的 bit 設為 1，相同的設為 0，所以回傳的值就是 hamming distance。

延伸問題

舉出不使用 gcc extension 的 C99 實作程式碼

int hammingDistance(int x, int y) {
    int v = x ^ y;
    v = v - ((v >> 1) & 0x55555555);
    v = (v & 0x33333333) + ((v >> 2) & 0x33333333);
    return ((v + (v >> 4) & 0xF0F0F0F) * 0x1010101) >> 24;
}

參考 Bit Twiddling Hacks，模擬 popcount() 的行為。

練習 LeetCode 477. Total Hamming Distance，你應當提交 (submit) 你的實作 (C 語言實作)，並確保執行結果正確且不超時

一開始使用兩個 for 迴圈，但是會有 TLE 的問題，改成對每個數字的 bit 著手，這樣可以確保經過 32 次迴圈後一定能結束。觀察在 numsSize 數字內的同一個 bit 造成的 hamming distance 為 1 的數量 * 0 的數量，這是因為每個 1 與每個 0 都可以互相配對，此種方法為將原本需要搜尋所有可能的 pair，移除兩個 bit 相同的情形，進而減少程式碼。
```
int totalHammingDistance(int* nums, int numsSize){
    int total = 0;
    int n;
    int count;
    for (;;) {
        n = 0;
        count = 0;
        for (int j = 0; j < numsSize; ++j) {
            if (nums[j] & 1)
                ++n;
            if (!nums[j])
                ++count;
            nums[j] >>= 1;
        }
        if (count == numsSize)
            break;
        total += (numsSize - n) * n;
    }
    return total;
}
```
這邊使用額外一個變數來讓迴圈提早結束，此程式可過 leetcode 測試：

Runtime: 44 ms, faster than 71.56% of C online submissions for Total Hamming Distance.
Memory Usage: 6.7 MB, less than 37.62% of C online submissions for Total Hamming Distance.

看起來還能再改進。
研讀錯誤更正碼簡介及抽象代數的實務應用，摘錄其中 Hamming Distance 和 Hamming Weight 的描述並重新闡述，舉出實際的 C 程式碼來解說;
- 搭配閱讀 Hamming codes, Part I 及 Hamming codes, Part II，你可適度摘錄影片中的素材，但應當清楚標註來源
(7,4) Hamming Code

假設我們要傳入的為 4 bit 之資料，而我們會將原始資料再加上 3 bit 的 parity bit，形成 7 bit 的資料來傳輸。

前 4 個 bit 可以對應到下圖 set 當中交集的 1, 2, 3, 4 部分，而剩餘的 3 個 bit 則是 set 外圍的部分。
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 →
parity bit 的選擇是依據每個圓圈裡面 1 的個數，只要圓圈內 1 的個數為奇數，則將該圓對應的 parity bit 設為 1，反之則設為 0，目的是讓每個圓圈裡面 1 的個數都為偶數。

接著我們就能將上述規則寫成代數形式，這邊加法為不含 carry 的二進位加法，也可以當成 XOR 運算子：

${\begin{matrix} x 1 + x 2 + x 3 + x 5 = 0 \\ x 1 + x 2 + x 4 + x 6 = 0 \\ x 1 + x 3 + x 4 + x 7 = 0 \end{matrix}$
將此等式寫成矩陣形式
$H$ ：

$H = (\begin{matrix} 1 & 1 & 1 & 0 & 1 & 0 & 0 \\ 1 & 1 & 0 & 1 & 0 & 1 & 0 \\ 1 & 0 & 1 & 1 & 0 & 0 & 1 \end{matrix})$
因為
$H \vec{x} = 0$ ，我們要求的矩陣
$C = (x 1, x 2, . . ., x 7)$ 就為矩陣
$H$ 的 null space。而
$r a n k (H) = 3, r a n k (C) = 7 - r a n k (H) = 4$ ，所以
$C$ 可以選擇 4 個變數，而剩餘 3 個則可以透過這些變數來決定。

我們只需透過矩陣乘法就能找出錯誤的 bit，假設傳出的資料
$\vec{y} = \vec{x} + \vec{e}$ ，其中
$\vec{x}$ 為正確的資料，
$\vec{e}$ 為錯誤向量，所以可以寫成
$H \vec{y} = H \vec{x} + H \vec{e} = H \vec{e}$ ，只要相乘的結果不為
$\vec{0}$ ，則表示有錯誤。

在習慣上，我們可以將 parity bit 放置在 2 的次方的位置上，這樣矩陣相乘的結果向量就是錯誤 bit 的位置，影片參考 Hamming codes part 2。

除了利用矩陣外，hamming code 其實可以用 XOR 來達成上面的運算，因為 parity bit 用來判斷奇偶數個 1，可以用對應的 bit 做 XOR 來取得，如果將 parity bit 放置在 2 的指數次方位置，parity bit 檢測的範圍為特定的 bit 位置，如 (7,4) hamming code 中：
- 001 位置的 parity bit 用來檢查 011, 101, 111 位置 1 的個數
- 010 位置的 parity bit 用來檢查 011, 110, 111 位置 1 的個數
- 100 位置的 parity bit 用來檢查 101, 110, 111 位置 1 的個數
其他資料的 bit 按順序放入剩餘的 bit，這樣檢查時只需對每個為 1 的 bit 的二進位位置做 XOR，結果就為錯誤 bit 的位置。這是因為錯誤的位置會剛好影響到對應的 parity bit，這樣可以有效的減少實作的程式碼。
```
void decoder(int receive)
{
    int error = 0;
    for (int i = 0; i < 32; ++i) {
        if (receive & (1 << i)) {
            error ^= i;
        }
    }
    if (!error) {
        receive ^= (1 << error);
    }
    return receive;
}
```
使用 32 bit integer 來儲存 32 bit 的資料，可以看到程式碼相當簡潔，只需將每個為 1 的 bit 的位置做 XOR，其結果就是 error 的位置。
研讀 Reed–Solomon error correction，再閱讀 Linux 核心原始程式碼的 lib/reed_solomon 目錄，抽離後者程式碼為單獨可執行的程式，作為 ECC 的示範

測驗 2

#include <stdlib.h>

typedef struct {
    AAA;
    int max_level;
} TreeAncestor;

TreeAncestor *treeAncestorCreate(int n, int *parent, int parentSize)
{
    TreeAncestor *obj = malloc(sizeof(TreeAncestor));
    obj->parent = malloc(n * sizeof(int *));

    int max_level = 32 - __builtin_clz(n) + 1;
    for (int i = 0; i < n; i++) {
        obj->parent[i] = malloc(max_level * sizeof(int));
        for (int j = 0; j < max_level; j++)
            obj->parent[i][j] = -1;
    }
    for (int i = 0; i < parentSize; i++)
        obj->parent[i][0] = parent[i];

    for (int j = 1;; j++) {
        int quit = 1;
        for (int i = 0; i < parentSize; i++) {
            obj->parent[i][j] = obj->parent[i][j + BBB] == -1
                                    ? -1
                                    : obj->parent[obj->parent[i][j - 1]][j - 1];
            if (obj->parent[i][j] != -1) quit = 0;
        }
        if (quit) break;
    }
    obj->max_level = max_level - 1;
    return obj;
}

int treeAncestorGetKthAncestor(TreeAncestor *obj, int node, int k)
{
    int max_level = obj->max_level;
    for (int i = 0; i < max_level && node != -1; ++i)
        if (k & (CCC))
            node = obj->parent[node][i];
    return node;
}

void treeAncestorFree(TreeAncestor *obj)
{
    for (int i = 0; i < obj->n; i++)
        free(obj->parent[i]);
    free(obj->parent);
    free(obj);
}

`AAA`

(b) int **parent

`BBB`

(b) (-1)

`CCC`

(f) 1 << i

AAA 是用來記憶每個 node i 的第 k 個 parent，可看到下方程式碼使用到二為陣列，所以這邊應為 **parent。

使用 dynamic programming 的技巧、Top down 的方式來解題。從樹的 root 節點開始向下，因每個 i 節點的第 j 個 parent 就是其 parent 的第 j - 1 個 parent，所以更新公式可以寫成 obj->parent[i][j] = obj->parent[obj->parent[i][j - 1]][j - 1];，但還要包括邊界條件，當第 j 個 parent 的 parent 沒有時，設定為 -1。

所以 BBB 是用來檢查是否到達邊界，故為 -1。

CCC 用來判斷 k 的 bit，運用每個 bit 來取得對應數字的 parent，再透過迴圈 loop k 的 bit，最後就能取得到第 k 個 parent node，所以 CCC 為 1 << i。

延伸問題

指出可改進的地方，並實作對應的程式碼;
- max_level 在實際上多分配一個 node 的空間，因為下方 for (int j = 1;; j++) 是用此位置來判斷是否結束，實在是多此一舉，將 max_level 改為真正的高度並設定迴圈條件 for (int j = 1; j < max_level; j++)。
- 在 treeAncestorFree() 中有用到 obj->n，但結構中沒有這項成員，將之加入結構中並新增 obj->n = n 在 treeAncestorCreate() 當中。
- Create 中第三個迴圈是以 row 為主，可能會有 locality 的問題，將陣列紀錄順序對調。
在 treeAncestorCreate 函式內部，若干記憶體被浪費，請撰寫完整測試程式碼，並透過工具分析;

在 LeetCode 1483. Kth Ancestor of a Tree Node 提交你的實作，你應該要在 C 語言項目中，執行時間領先 75% 的提交者;

將上述內容實作後並提交：

TreeAncestor *treeAncestorCreate(int n, int *parent, int parentSize)
{
    TreeAncestor *obj = malloc(sizeof(TreeAncestor));
    int max_level = 32 - __builtin_clz(n);
    obj->parent = malloc(max_level * sizeof(int *));

    for (int i = 0; i < max_level; i++) {
        obj->parent[i] = malloc(n * sizeof(int));
        for (int j = 0; j < n; j++)
            obj->parent[i][j] = -1;
    }
    for (int i = 0; i < parentSize; i++)
        obj->parent[0][i] = parent[i];

    for (int i = 1; i < max_level; i++) {
        for (int j = 0; j < parentSize; j++) {
            obj->parent[i][j] = obj->parent[i - 1][j] == -1
                                    ? -1
                                    : obj->parent[i - 1][obj->parent[i - 1][j]];
        }
    }
    obj->max_level = max_level;
    obj->n = n;
    return obj;
}

Runtime: 244 ms, faster than 98.45% of C online submissions for Kth Ancestor of a Tree Node.
Memory Usage: 62.9 MB, less than 6.50% of C online submissions for Kth Ancestor of a Tree Node.

測驗 3

#define MSG_LEN 8

static inline bool is_divisible(uint32_t n, uint64_t M) {
    return n * M <= M - 1;
}

static uint64_t M3 = UINT64_C(0xFFFFFFFFFFFFFFFF) / 3 + 1;
static uint64_t M5 = UINT64_C(0xFFFFFFFFFFFFFFFF) / 5 + 1;

int main(int argc, char *argv[]) {
    for (size_t i = 1; i <= 100; i++) {
        uint8_t div3 = is_divisible(i, M3);
        uint8_t div5 = is_divisible(i, M5);
        unsigned int length = (2 << div3) << div5;

        char fmt[MSG_LEN + 1];
        strncpy(fmt, &"FizzBuzz%u"[(MSG_LEN >> KK1) >> (KK2 << KK3)], length);
        fmt[length] = '\0';

        printf(fmt, i);
        printf("\n");
    }
    return 0;
}

is_divisible 可以參考 quiz2 的第三題與第四題。

根據題意：我們若能精準從給定輸入 i 的規律去控制 start 及 length，即可符合 FizzBuzz 題意，可以把 div3, div5 與對應的值寫下：

div3	div5	(MSG_LEN >> KK1) >> (KK2 << KK3)
0	0	8
1	0	0
0	1	4
1	1	0

分析：

可以看到對 MSG_LEN 的運算元只有 >>，所以在 div3 與 div5 皆為 0 的情形時，可得到 (KK1 == 0) && ((KK2 << KK3) == 0) 的關係。

在 div3 為 0 且 div5 為 1 時，MSG_LEN 應要右移一個 bit，此有兩種可能：

KK1 = 1, KK2 = 0, KK3 = ?
KK1 = 0, KK2 = 1, KK3 = 0

在 div3 為 1 且 div5 為 0 時，MSG_LEN 應要右移四個 bit 以上，此有七種可能：

KK1 = 3, KK2 = 0, KK3 = ?
KK1 = 2, KK2 = 1, KK3 = 1
KK1 = 2, KK2 = 2, KK3 = 0
KK1 = 1, KK2 = 3, KK3 = 0
KK1 = 0, KK2 = 1, KK3 = 2
KK1 = 0, KK2 = 2, KK3 = 1
KK1 = 0, KK2 = 4, KK3 = 0

在 div3 與 div5 皆為 1 時，與上方結果一樣。

帶數字就能得到結果。

`KK1`

(a) div5

`KK2`

(d) div3

`KK3`

(c) 2

延伸問題

評估 naive.c 和 bitwise.c 效能落差

避免 stream I/O 帶來的影響，可將 printf 更換為 sprintf

#include <stdio.h>
#include <time.h>
#include <string.h>
#include <stdint.h>
#include <stdbool.h>

#define MSG_LEN 8

static inline bool is_divisible(uint32_t n, uint64_t M) {
    return n * M <= M - 1;
}

static uint64_t M3 = UINT64_C(0xFFFFFFFFFFFFFFFF) / 3 + 1;
static uint64_t M5 = UINT64_C(0xFFFFFFFFFFFFFFFF) / 5 + 1;

void naive() {
    for (unsigned int i = 1; i < 1000000; i++) {
        char s[MSG_LEN + 1];
        if (i % 3 == 0) 1sprintf(s, "Fizz");
        if (i % 5 == 0) sprintf(s, "Buzz");
        if (i % 15 == 0) sprintf(s, "FizzBuzz");
        if ((i % 3) && (i % 5)) sprintf(s, "%u", i);
        sprintf(s, "\n");
    }
}

void bitwise() {
    for (size_t i = 1; i <= 1000000; i++) {
        uint8_t div3 = is_divisible(i, M3);
        uint8_t div5 = is_divisible(i, M5);
        unsigned int length = (2 << div3) << div5;

        char fmt[MSG_LEN + 1];
        strncpy(fmt, &"FizzBuzz%u"[(MSG_LEN >> div5) >> (div3 << 2)], length);
        fmt[length] = '\0';
    }
}

int main(int argc, char *argv[]) {
    double time_taken;
    struct timespec start, end;
    clock_gettime(CLOCK_MONOTONIC, &start);
    bitwise();
    clock_gettime(CLOCK_MONOTONIC, &end);
    time_taken = (end.tv_sec - start.tv_sec) * 1e9;
    time_taken = (time_taken + (end.tv_nsec - start.tv_nsec)) * 1e-9;
    printf("%f,", time_taken);
    clock_gettime(CLOCK_MONOTONIC, &start);
    naive();
    clock_gettime(CLOCK_MONOTONIC, &end);
    time_taken = (end.tv_sec - start.tv_sec) * 1e9;
    time_taken = (time_taken + (end.tv_nsec - start.tv_nsec)) * 1e-9;
    printf("%f\n", time_taken);

    return 0;
}

因 bitwise 少了分支與餘數運算，效能比 naive 好非常多。

分析 Faster remainders when the divisor is a constant: beating compilers and libdivide 一文的想法 (可參照同一篇網誌下方的評論)，並設計另一種 bitmask，如「可被 3 整除則末位設為 1」「可被 5 整除則倒數第二位設定為 1」，然後再改寫 bitwise.c 程式碼，試圖運用更少的指令來實作出 branchless;
- 參照 fastmod: A header file for fast 32-bit division remainders on 64-bit hardware

測驗 4

#include <stdint.h>
#define ffs32(x)  ((__builtin_ffs(x)) - 1)
size_t blockidx;
size_t divident = PAGE_QUEUES;
blockidx = offset >> ffs32(divident);
divident >>= ffs32(divident);
switch (divident) {
    CASES
}

`CASES` 至少該包含哪些數字: (複選)

(b) 3
(d) 5
(f) 7

__builtin_ffs

Returns one plus the index of the least significant 1-bit of x, or if x is zero, returns zero.

回傳從右邊數來第一個 1 的位置

#define ffs32(x) ((__builtin_ffs(x)) - 1)

回傳從右邊數來第一個 1 以前有幾個 0。
i.e. ffs32() == __builtin_ctz()

考慮到整數 PAGE_QUEUES 可能有以下數值:

(1 or 2 or 3 or 4)
(5 or 6 or 7 or 8) × (
$2^{n}$ ), n from 1 to 14

可以寫成以下：

(1 or 2 or 3 or 4)
(5 or 7) × (
$2^{n}$ ), n from 1 to 14
(3 or 4) × (
$2^{n}$ ), n from 2 to 15

可以把相同倍數的合併：

(1 or 3 or 5 or 7) × (
$2^{n}$ ), n from 0 to 15

這邊雖然某些數的 n 範圍不同，但我們可以擴充成 n 皆為相同，只要有涵蓋到原本的範圍就好。

blockidx = offset >> ffs32(divident);
divident >>= ffs32(divident);

程式碼中這兩步是把上面公式中的

2^{n}

部份移除，剩下部份就要依靠 CASE 來處理：

switch (divident) {
    case 3: blockidx /= 3; 
              break;
    case 5: blockidx /= 5; 
              break;
    case 7: blockidx /= 7; 
              break;
}

因除以 1 沒有意義，所以不必寫在 CASE 中。

延伸問題

解釋程式運算原理，可搭配 Microsoft Research 的專案 snmalloc 來解說;
- 對應的原始程式碼 src/mem/sizeclass.h
snmalloc 裡的 round_by_sizeclass() 用來計算
$⌊ \frac{x}{y} ⌋ \times y$ 的值，而其中除法的部分使用 fastdiv 來計算，並與本題一樣利用 1, 3, 5, 7 當作 special case，便能依照對應的 case 加速運算。

練習 LeetCode 51. N-Queens，應選擇 C 語言並善用上述的 __builtin_ffs，大幅降低記憶體開銷;

參考使用位元計算加速求解 N-Queen 與 RinHizakura 同學的程式碼：

int len;
char cpy[] = ".................";
int sizes[] = { 1,      1,       0,        0,        2,        10,    4,
               40,     92,      352,      724,      2680,     14200, 73712,
               365596, 2279184, 14772512, 95815104, 666090624 };

bool recursive(char*** ans, int n, int row, uint16_t l, uint16_t m, uint16_t r) {
    if (row == n) {
        ++len;
        if (len == sizes[n])
            return true;
        for (int i = 0; i < n; ++i) {
            memcpy(ans[len][i], ans[len - 1][i], n);
        }
        return false;
    }
    uint16_t full = 0xffffffff;
    uint16_t accept_cols = ~(l | m | r);
    for (int col = __builtin_ffs(accept_cols) - 1; col < n; col = __builtin_ffs(accept_cols) - 1) {
        ans[len][row][col] = 'Q';
        if (recursive(ans, n, row + 1, (l | (1 << col)) << 1, (m | (1 << col)), (r | (1 << col)) >> 1))
            return true;
        ans[len][row][col] = '.';
        accept_cols ^= (1 << col);
    }
    return false;
}

char*** solveNQueens(int n, int* returnSize, int** returnColumnSizes) {
    len = 0;
    char*** ans;
    ans = malloc(sizeof(char**) * sizes[n]);
    *returnColumnSizes = malloc(sizes[n] * sizeof(int));
    *returnSize = sizes[n];
    for (int i = 0; i < sizes[n]; ++i) {
        (*returnColumnSizes)[i] = n;
        ans[i] = malloc(sizeof(char*) * n);
        for (int j = 0; j < n; ++j) {
            ans[i][j] = malloc(sizeof(char) * (n + 1));
            memcpy(ans[i][j], cpy, n);
            ans[i][j][n] = '\0';
        }
    }
    if (!sizes[n])
        return ans;
    recursive(ans, n, 0, 0, 0, 0);
    return ans;        
}

將可放置皇后的位置 bit 設為 1，使用三個變數 l, m, r 記錄目前不能放置的位置，因為是從 row 0 迴圈到 row n，所以不用管 row 上重複的情形，而 col 則是將前一個 row 皇后的位置的左右一個 bit，加上中間的 bit 加到 bitmask 中防止選到不能選的位置。

因 bit 為 1 表示可選擇之位置，0 表示不能選擇之位置，故使用 __builtin_ffs 能快速找到可選擇的 col 位置。

Leetcode 結果：
Runtime: 4 ms, faster than 96.97% of C online submissions for N-Queens.
Memory Usage: 7.2 MB, less than 27.27% of C online submissions for N-Queens.

2020q3 Homework (quiz4)

測驗 1

OP

延伸問題

(7,4) Hamming Code

測驗 2

AAA

BBB

CCC

延伸問題

測驗 3

KK1

KK2

KK3

延伸問題

測驗 4

CASES 至少該包含哪些數字: (複選)

延伸問題

tags: linux2020

Read more

Linked List Sort

2021q1 Homework5 (sort)

2022q1 Homework5 (quiz5)

2022q1 Homework5 (quiz6)

`OP`

`AAA`

`BBB`

`CCC`

`KK1`

`KK2`

`KK3`

`CASES` 至少該包含哪些數字: (複選)

tags: `linux2020`