2016-01-22

[讀檔操作] 有限記憶體找中位數

1. 背景
2. 題目描述
3. 輸入格式
4. 輸出格式
5. 範例輸入
6. 範例輸出
7. 提示
1. 7.1. p10067-in.dat
8. Solution
1. 8.1. 二分搜尋
2. 8.2. 分桶算法

背景

現在的系統架構中，內存 (Memory，台灣翻譯：記憶體，大陸翻譯：內存) 的大小仍無法完全像硬碟機一樣。如檔案大小 64GB，內存只有 4GB，處理檔案無法全部 In Memory，當然最近的硬體技術也逐漸朝著 All In Memory 的方式進行加速。回過頭來，在內存不足的情況下，來練習如何處理檔案大小遠大於內存的情況吧。

題目描述

給予一個 binary file 的檔案名稱，裡面用 4 個位元組為一個 signed 32-bits 整數，有數個直到 EOF，保證恰好有奇數個，請輸出中位數為何？

輸入格式

標準輸入串流只有一行一個字串，表示檔案名稱 $F$。檔案大小最多 16 MB，而你的程序被限制最多使用 4 MB 的內存。

每個整數 $x$，限制條件 $0 \le x \le 10^9$

輸出格式

輸出一行整數 $Median$ 為 binary file 的中位數。

範例輸入

1	p10067-in.dat

範例輸出

97537111

提示

二分搜尋，分桶

download p10067-in.dat

p10067-in.dat

1
2
3

10825098
97537111
208681850

Solution

這一個問題很久以前見過，對岸通常叫做「海量資料找中位數」約束在記憶體不夠的情況下，如何高效率找到中位數，操作時允許重複讀檔。

遲遲沒有 Online Judge 進行記憶體用量檢測以及開放讀寫檔案，如今有機會擔任台灣大學計算機程式設計的課程助教，架了一個允許亂來的 Online Judge - Judge Girl (中文翻譯：批改娘系統)，讀寫檔案也不會被封鎖！洽詢網站，目前只針對課堂學生開放。

回到問題，題目給訂 16MB 的整數序列，最多 $n = 4 \times 10^6$，由於記憶體不足沒辦法直接宣告陣列大小為 $n$ 的整數陣列，又由於數字可能會重複，就沒辦法利用 bitmask 進行壓縮，只能倚靠不斷地讀檔案進行計算。若以數值範圍 $[0, V]$，大致放分成兩種策略

二分搜尋 (AC, 250ms)，效率 $\mathcal{O}(V \log V)$，開檔次數 $\mathcal{O}(\log V)$。
分桶算法 (AC, 70ms)，效率 $\mathcal{O}(N)$，開檔次數 $\mathcal{O}(1)$。

讀取檔案時，在有限記憶體空間，配置一塊作為緩衝區，一次讀入一坨子序列，一個一個讀取效率非常差，別輕忽從磁碟讀取資料的速度。剩餘空間再進行紀錄用途。從效能比較 二分搜尋 慢於 分桶算法。

二分搜尋

二分答案 $x$，接著線性掃描序列，計算有多少個整數小於等於 $x$，藉此找到中位數。

#include <stdio.h>
int fsize(FILE *fp) {
    int prev = ftell(fp);
    fseek(fp, 0L, SEEK_END);
    int size = ftell(fp);
    fseek(fp, prev, SEEK_SET);
    return size;
}
#define MAXELE (3<<20)/4
int countLess(FILE *fp, int ans) {
    static int A[MAXELE];
    int ret = 0, x, n = 0;
    fseek(fp, 0, SEEK_SET);
    while ((n = fread(A, sizeof(int), MAXELE, fp))) {
        for (int i = 0; i < n; i++)
            ret += A[i] <= ans;
    }
    return ret;
}
int main() {
    char fName[128];
    scanf("%s", fName);
    FILE *fin = fopen(fName, "rb");
    
    int n = fsize(fin) / sizeof(int);
    
    int m = n/2 + 1;
    int l = 0, r = 1000000000, mid, ret = 0;
    
    while (l <= r) {
        mid = (l + r)/2;
        int tn = countLess(fin, mid);
        if (tn < m)
            l = mid + 1, ret = mid+1;
        else
            r = mid - 1;
    }
    printf("%d\n", ret);
    return 0;
}

分桶算法

分桶算法分成兩次掃瞄，假設內存分配最多 $I$ 個 32bits 整數 ($I \ll N$)。

由於記憶體限制大小關係，將數值範圍分成 $I$ 堆，如 [0, V/I-1]、[V/I, 2V/I-1]、… 等，計算區間範圍內的個數，讀取一次檔案，可以得到中位數介於 [iV/I, (i+1)V/I-1]。
接著再劃分一次，直到區間長度為 1，中位數便可得到。

#include <stdio.h>
#include <assert.h>
#include <string.h>
int fsize(FILE *fp) {
    int prev = ftell(fp);
    fseek(fp, 0L, SEEK_END);
    int size = ftell(fp);
    fseek(fp, prev, SEEK_SET);
    return size;
}
#define MAXELE ((1<<20)/4)
static int A[MAXELE];
static int B[MAXELE];
int countLess(FILE *fp, int ans) {
    int ret = 0, x, n = 0;
    fseek(fp, 0, SEEK_SET);
    while ((n = fread(A, sizeof(int), MAXELE, fp))) {
        for (int i = 0; i < n; i++)
            ret += A[i] <= ans;
    }
    return ret;
}
int main() {
    char fName[128];
    scanf("%s", fName);
    FILE *fin = fopen(fName, "rb");
    
    int n = fsize(fin) / sizeof(int);
    
    int m = n/2 + 1;
    int l = 0, r = 1000000000, mid, ret = 0;
    
    int WIDTH = r, SHIFT;
    
    for (SHIFT = 1; (1<<SHIFT) < MAXELE; SHIFT++);
    WIDTH = r / (1<<SHIFT);
    
    memset(B, 0, sizeof(B));
    fseek(fin, 0, SEEK_SET);
    while ((n = fread(A, sizeof(int), MAXELE, fin))) {
        for (int i = 0; i < n; i++) {
            assert((A[i]>>SHIFT) < MAXELE);
            B[A[i]>>SHIFT]++;
        }
    }
    
    int BASE = 0, SELECT;
    for (int i = 0, sum = m; i < MAXELE; i++) {
        sum -= B[i];
        if (sum <= 0) {
            BASE = i * (1<<SHIFT);
            SELECT = sum + B[i];
            break;
        }
    }
    
    memset(B, 0, sizeof(B));
    fseek(fin, 0, SEEK_SET);
    while ((n = fread(A, sizeof(int), MAXELE, fin))) {
        for (int i = 0; i < n; i++) {
            if (A[i] - BASE < MAXELE && A[i] >= BASE)
                B[A[i] - BASE]++;
        }
    }
    for (int i = 0, sum = SELECT; i < MAXELE; i++) {
        sum -= B[i];
        if (sum <= 0) {
            printf("%d\n", i + BASE);
            break;
        }
    }
    return 0;
}

Morris' Blog