2015-03-09

Page Rank 迭代計算 C++

1. Problem
2. Input
3. Sample Input
4. Sample Output
5. Solution

Problem

每一個網站都有超連結 (hyperlink) 到另一個網頁，因此會貢獻一部分的 page rank 給另一個網站，假想移動是隨機的，那麼轉移機會是 1 / 出度。

把網站之間的關係轉為 Graph，問題可以轉換成馬可夫過程。光是馬可夫過程，很容易造成 page rank 失算，因為有些點並沒有出度 (連接到別的網站)，數值就會隨著迭代集中到這個孤島上 (或者說蜘蛛網)，造成其他點的 page rank 皆為 0。為了防止這一點，給予 page rank 的總和$S$，並且要求馬可夫過程，轉移只有$\beta$ 的信任度，剩餘的$1 - \beta$ 將隨機選擇全局的網站進行轉移。

Input

每組第一行輸入兩個實數$\beta$,$S$。

第二行輸入一個整數 N，表示有多少個點。接下來會有 N 行，每一行的第一個數字，表示 i-th 點會連到其他網站的個數$m_{i}$，接著同一行會有$m_{i}$ 個數字，表示連入的網站編號。

所有網站編號應為 0 … N-1。

Sample Input

Sample Output

A 0.300
B 0.405
C 2.295
   
A 1.163
B 0.644
C 1.192

Solution

以下是簡單實作的結果。迭代 100 次，就能收斂程度就相當足夠。

#include <stdio.h> 
#include <math.h>
#include <vector>
#include <string.h>
using namespace std;
// Markov process, math
const int MAXN = 64;
const int MAXIT = 100;
struct Matrix {
    double v[MAXN][MAXN];
    int row, col; // row x col
    Matrix(int n, int m, int a = 0) {
        memset(v, 0, sizeof(v));
        row = n, col = m;
        for(int i = 0; i < row && i < col; i++)
            v[i][i] = a;
    }
    Matrix operator*(const Matrix& x) const {
        Matrix ret(row, x.col);
        for(int i = 0; i < row; i++) {
            for(int k = 0; k < col; k++) {
                if (v[i][k])
                    for(int j = 0; j < x.col; j++) {
                        ret.v[i][j] += v[i][k] * x.v[k][j];
                    }
            }
        }
        return ret;
    }
    Matrix operator+(const Matrix& x) const {
        Matrix ret(row, col);
        for(int i = 0; i < row; i++) {
            for(int j = 0; j < col; j++) {
                ret.v[i][j] = v[i][j] + x.v[i][j];
            }
        }
        return ret;
    }
    Matrix operator^(const int& n) const {
        Matrix ret(row, col, 1), x = *this;
        int y = n;
        while(y) {
            if(y&1)	ret = ret * x;
            y = y>>1, x = x * x;
        }
        return ret;
    }
    Matrix powsum(int k) {
        if (k == 0) return Matrix(row, col, 0);
        Matrix vv = powsum(k/2);
        if (k&1) {
            return vv * (Matrix(row, col, 1) + vv) + vv;
        } else {
            return vv * (Matrix(row, col, 1) + vv);
        }
    }
};
#define eps 1e-6
int cmpZero(double x) {
    if (fabs(x) < eps)
        return 0;
    return x < 0 ? -1 : 1;
}
int main() {
    int N;
    double beta, S;
    while (scanf("%lf %lf", &beta, &S) == 2) {
        vector<int> g[MAXN], invg[MAXN];
        scanf("%d", &N);
        for (int i = 0; i < N; i++) {
            int M, x;
            scanf("%d", &M);
            for (int j = 0; j < M; j++) {
                scanf("%d", &x);
                g[i].push_back(x);
                invg[x].push_back(i);
            }
        }
        
        Matrix r(N, 1);
        for (int i = 0; i < N; i++)
            r.v[i][0] = 1.0 / N;
        
        for (int it = 0; it < MAXIT; it++) {
            Matrix next_r(N, 1);
            for (int i = 0; i < N; i++) {
                for (int j = 0; j < invg[i].size(); j++) {
                    int x = invg[i][j];
                    next_r.v[i][0] += beta * r.v[x][0] / g[x].size();
                }
            }
            
            double sum = 0;
            for (int i = 0; i < N; i++)
                sum += next_r.v[i][0];
            for (int i = 0; i < N; i++)
                next_r.v[i][0] += (S - sum) / N;
            r = next_r;
        }
        for (int i = 0; i < N; i++)
            printf("%c %.3lf\n", i + 'A', r.v[i][0]);
    }
    return 0;
}
/*
0.7 3
3
2 1 2
1 2
1 2
0.85 3
3
2 1 2
1 2
1 0
*/

Morris' Blog